ROOT  6.07/01
Reference Guide
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
StandardHypoTestDemo.C
Go to the documentation of this file.
1 // Standard tutorial macro for hypothesis test (for computing the discovery significance) using all
2 // RooStats hypotheiss tests calculators and test statistics
3 //
4 //
5 //Author: L. Moneta
6 //
7 // Usage:
8 //
9 // root>.L StandardHypoTestDemo.C
10 // root> StandardHypoTestDemo("fileName","workspace name","S+B modelconfig name","B model name","data set name",calculator type, test statistic type, // number of toys)
11 //
12 //
13 // type = 0 Freq calculator
14 // type = 1 Hybrid calculator
15 // type = 2 Asymptotic calculator
16 // type = 3 Asymptotic calculator using nominal Asimov data sets (not using fitted parameter values but nominal ones)
17 //
18 // testStatType = 0 LEP
19 // = 1 Tevatron
20 // = 2 Profile Likelihood
21 // = 3 Profile Likelihood one sided (i.e. = 0 if mu_hat < 0)
22 //
23 
24 
25 #include "TFile.h"
26 #include "RooWorkspace.h"
27 #include "RooAbsPdf.h"
28 #include "RooRealVar.h"
29 #include "RooDataSet.h"
30 #include "RooStats/ModelConfig.h"
31 #include "RooRandom.h"
32 #include "TGraphErrors.h"
33 #include "TGraphAsymmErrors.h"
34 #include "TCanvas.h"
35 #include "TLine.h"
36 #include "TSystem.h"
37 #include "TROOT.h"
38 
42 #include "RooStats/ToyMCSampler.h"
43 #include "RooStats/HypoTestPlot.h"
44 
50 
54 
55 using namespace RooFit;
56 using namespace RooStats;
57 
58 
59 bool noSystematics = false; // force all systematics to be off (i.e. set all nuisance parameters as constat
60 double nToysRatio = 4; // ratio Ntoys Null/ntoys ALT
61 double poiValue = -1; // change poi snapshot value for S+B model (needed for expected p0 values)
62 int printLevel=0;
63 bool generateBinned = false; // for binned generation
64 
65 void StandardHypoTestDemo(const char* infile = "",
66  const char* workspaceName = "combined",
67  const char* modelSBName = "ModelConfig",
68  const char* modelBName = "",
69  const char* dataName = "obsData",
70  int calcType = 0, // 0 freq 1 hybrid, 2 asymptotic
71  int testStatType = 3, // 0 LEP, 1 TeV, 2 LHC, 3 LHC - one sided
72  int ntoys = 5000,
73  bool useNC = false,
74  const char * nuisPriorName = 0)
75 {
76 
77 /*
78 
79  Other Parameter to pass in tutorial
80  apart from standard for filename, ws, modelconfig and data
81 
82  type = 0 Freq calculator
83  type = 1 Hybrid calculator
84  type = 2 Asymptotic calculator
85 
86  testStatType = 0 LEP
87  = 1 Tevatron
88  = 2 Profile Likelihood
89  = 3 Profile Likelihood one sided (i.e. = 0 if mu < mu_hat)
90 
91  ntoys: number of toys to use
92 
93  useNumberCounting: set to true when using number counting events
94 
95  nuisPriorName: name of prior for the nnuisance. This is often expressed as constraint term in the global model
96  It is needed only when using the HybridCalculator (type=1)
97  If not given by default the prior pdf from ModelConfig is used.
98 
99  extra options are available as global paramwters of the macro. They major ones are:
100 
101  generateBinned generate binned data sets for toys (default is false) - be careful not to activate with
102  a too large (>=3) number of observables
103  nToyRatio ratio of S+B/B toys (default is 2)
104  printLevel
105 
106 */
107 
108  // disable - can cause some problems
109  //ToyMCSampler::SetAlwaysUseMultiGen(true);
110 
111  SimpleLikelihoodRatioTestStat::SetAlwaysReuseNLL(true);
112  ProfileLikelihoodTestStat::SetAlwaysReuseNLL(true);
113  RatioOfProfiledLikelihoodsTestStat::SetAlwaysReuseNLL(true);
114 
115  //RooRandom::randomGenerator()->SetSeed(0);
116 
117  // to change minimizers
118  // ROOT::Math::MinimizerOptions::SetDefaultStrategy(0);
119  // ROOT::Math::MinimizerOptions::SetDefaultMinimizer("Minuit2");
120  // ROOT::Math::MinimizerOptions::SetDefaultTolerance(1);
121 
122  /////////////////////////////////////////////////////////////
123  // First part is just to access a user-defined file
124  // or create the standard example file if it doesn't exist
125  ////////////////////////////////////////////////////////////
126  const char* filename = "";
127  if (!strcmp(infile,"")) {
128  filename = "results/example_combined_GaussExample_model.root";
129  bool fileExist = !gSystem->AccessPathName(filename); // note opposite return code
130  // if file does not exists generate with histfactory
131  if (!fileExist) {
132 #ifdef _WIN32
133  cout << "HistFactory file cannot be generated on Windows - exit" << endl;
134  return;
135 #endif
136  // Normally this would be run on the command line
137  cout <<"will run standard hist2workspace example"<<endl;
138  gROOT->ProcessLine(".! prepareHistFactory .");
139  gROOT->ProcessLine(".! hist2workspace config/example.xml");
140  cout <<"\n\n---------------------"<<endl;
141  cout <<"Done creating example input"<<endl;
142  cout <<"---------------------\n\n"<<endl;
143  }
144 
145  }
146  else
147  filename = infile;
148 
149  // Try to open the file
150  TFile *file = TFile::Open(filename);
151 
152  // if input file was specified byt not found, quit
153  if(!file ){
154  cout <<"StandardRooStatsDemoMacro: Input file " << filename << " is not found" << endl;
155  return;
156  }
157 
158 
159  /////////////////////////////////////////////////////////////
160  // Tutorial starts here
161  ////////////////////////////////////////////////////////////
162 
163  // get the workspace out of the file
164  RooWorkspace* w = (RooWorkspace*) file->Get(workspaceName);
165  if(!w){
166  cout <<"workspace not found" << endl;
167  return;
168  }
169  w->Print();
170 
171  // get the modelConfig out of the file
172  ModelConfig* sbModel = (ModelConfig*) w->obj(modelSBName);
173 
174 
175  // get the modelConfig out of the file
176  RooAbsData* data = w->data(dataName);
177 
178  // make sure ingredients are found
179  if(!data || !sbModel){
180  w->Print();
181  cout << "data or ModelConfig was not found" <<endl;
182  return;
183  }
184  // make b model
185  ModelConfig* bModel = (ModelConfig*) w->obj(modelBName);
186 
187 
188  // case of no systematics
189  // remove nuisance parameters from model
190  if (noSystematics) {
191  const RooArgSet * nuisPar = sbModel->GetNuisanceParameters();
192  if (nuisPar && nuisPar->getSize() > 0) {
193  std::cout << "StandardHypoTestInvDemo" << " - Switch off all systematics by setting them constant to their initial values" << std::endl;
194  RooStats::SetAllConstant(*nuisPar);
195  }
196  if (bModel) {
197  const RooArgSet * bnuisPar = bModel->GetNuisanceParameters();
198  if (bnuisPar)
199  RooStats::SetAllConstant(*bnuisPar);
200  }
201  }
202 
203 
204  if (!bModel ) {
205  Info("StandardHypoTestInvDemo","The background model %s does not exist",modelBName);
206  Info("StandardHypoTestInvDemo","Copy it from ModelConfig %s and set POI to zero",modelSBName);
207  bModel = (ModelConfig*) sbModel->Clone();
208  bModel->SetName(TString(modelSBName)+TString("B_only"));
209  RooRealVar * var = dynamic_cast<RooRealVar*>(bModel->GetParametersOfInterest()->first());
210  if (!var) return;
211  double oldval = var->getVal();
212  var->setVal(0);
213  //bModel->SetSnapshot( RooArgSet(*var, *w->var("lumi")) );
214  bModel->SetSnapshot( RooArgSet(*var) );
215  var->setVal(oldval);
216  }
217 
218  if (!sbModel->GetSnapshot() || poiValue > 0) {
219  Info("StandardHypoTestDemo","Model %s has no snapshot - make one using model poi",modelSBName);
220  RooRealVar * var = dynamic_cast<RooRealVar*>(sbModel->GetParametersOfInterest()->first());
221  if (!var) return;
222  double oldval = var->getVal();
223  if (poiValue > 0) var->setVal(poiValue);
224  //sbModel->SetSnapshot( RooArgSet(*var, *w->var("lumi") ) );
225  sbModel->SetSnapshot( RooArgSet(*var) );
226  if (poiValue > 0) var->setVal(oldval);
227  //sbModel->SetSnapshot( *sbModel->GetParametersOfInterest() );
228  }
229 
230 
231 
232 
233 
234  // part 1, hypothesis testing
235  SimpleLikelihoodRatioTestStat * slrts = new SimpleLikelihoodRatioTestStat(*bModel->GetPdf(), *sbModel->GetPdf());
236  // null parameters must includes snapshot of poi plus the nuisance values
237  RooArgSet nullParams(*bModel->GetSnapshot());
238  if (bModel->GetNuisanceParameters()) nullParams.add(*bModel->GetNuisanceParameters());
239 
240  slrts->SetNullParameters(nullParams);
241  RooArgSet altParams(*sbModel->GetSnapshot());
242  if (sbModel->GetNuisanceParameters()) altParams.add(*sbModel->GetNuisanceParameters());
243  slrts->SetAltParameters(altParams);
244 
245 
247 
248 
250  ropl = new RatioOfProfiledLikelihoodsTestStat(*bModel->GetPdf(), *sbModel->GetPdf(), sbModel->GetSnapshot());
251  ropl->SetSubtractMLE(false);
252 
253  if (testStatType == 3) profll->SetOneSidedDiscovery(1);
254  profll->SetPrintLevel(printLevel);
255 
256  // profll.SetReuseNLL(mOptimize);
257  // slrts.SetReuseNLL(mOptimize);
258  // ropl.SetReuseNLL(mOptimize);
259 
260  AsymptoticCalculator::SetPrintLevel(printLevel);
261 
262  HypoTestCalculatorGeneric * hypoCalc = 0;
263  // note here Null is B and Alt is S+B
264  if (calcType == 0) hypoCalc = new FrequentistCalculator(*data, *sbModel, *bModel);
265  else if (calcType == 1) hypoCalc= new HybridCalculator(*data, *sbModel, *bModel);
266  else if (calcType == 2) hypoCalc= new AsymptoticCalculator(*data, *sbModel, *bModel);
267 
268  if (calcType == 0)
269  ((FrequentistCalculator*)hypoCalc)->SetToys(ntoys, ntoys/nToysRatio);
270  if (calcType == 1)
271  ((HybridCalculator*)hypoCalc)->SetToys(ntoys, ntoys/nToysRatio);
272  if (calcType == 2 ) {
273  if (testStatType == 3) ((AsymptoticCalculator*) hypoCalc)->SetOneSidedDiscovery(true);
274  if (testStatType != 2 && testStatType != 3)
275  Warning("StandardHypoTestDemo","Only the PL test statistic can be used with AsymptoticCalculator - use by default a two-sided PL");
276 
277 
278  }
279 
280 
281  // check for nuisance prior pdf in case of nuisance parameters
282  if (calcType == 1 && (bModel->GetNuisanceParameters() || sbModel->GetNuisanceParameters() )) {
283  RooAbsPdf * nuisPdf = 0;
284  if (nuisPriorName) nuisPdf = w->pdf(nuisPriorName);
285  // use prior defined first in bModel (then in SbModel)
286  if (!nuisPdf) {
287  Info("StandardHypoTestDemo","No nuisance pdf given for the HybridCalculator - try to deduce pdf from the model");
288  if (bModel->GetPdf() && bModel->GetObservables() )
289  nuisPdf = RooStats::MakeNuisancePdf(*bModel,"nuisancePdf_bmodel");
290  else
291  nuisPdf = RooStats::MakeNuisancePdf(*sbModel,"nuisancePdf_sbmodel");
292  }
293  if (!nuisPdf ) {
294  if (bModel->GetPriorPdf()) {
295  nuisPdf = bModel->GetPriorPdf();
296  Info("StandardHypoTestDemo","No nuisance pdf given - try to use %s that is defined as a prior pdf in the B model",nuisPdf->GetName());
297  }
298  else {
299  Error("StandardHypoTestDemo","Cannnot run Hybrid calculator because no prior on the nuisance parameter is specified or can be derived");
300  return;
301  }
302  }
303  assert(nuisPdf);
304  Info("StandardHypoTestDemo","Using as nuisance Pdf ... " );
305  nuisPdf->Print();
306 
307  const RooArgSet * nuisParams = (bModel->GetNuisanceParameters() ) ? bModel->GetNuisanceParameters() : sbModel->GetNuisanceParameters();
308  RooArgSet * np = nuisPdf->getObservables(*nuisParams);
309  if (np->getSize() == 0) {
310  Warning("StandardHypoTestDemo","Prior nuisance does not depend on nuisance parameters. They will be smeared in their full range");
311  }
312  delete np;
313 
314  ((HybridCalculator*)hypoCalc)->ForcePriorNuisanceAlt(*nuisPdf);
315  ((HybridCalculator*)hypoCalc)->ForcePriorNuisanceNull(*nuisPdf);
316  }
317 
318  // hypoCalc->ForcePriorNuisanceAlt(*sbModel->GetPriorPdf());
319  // hypoCalc->ForcePriorNuisanceNull(*bModel->GetPriorPdf());
320 
321  ToyMCSampler * sampler = (ToyMCSampler *)hypoCalc->GetTestStatSampler();
322 
323  if (sampler && (calcType == 0 || calcType == 1) ) {
324 
325  // look if pdf is number counting or extended
326  if (sbModel->GetPdf()->canBeExtended() ) {
327  if (useNC) Warning("StandardHypoTestDemo","Pdf is extended: but number counting flag is set: ignore it ");
328  }
329  else {
330  // for not extended pdf
331  if (!useNC) {
332  int nEvents = data->numEntries();
333  Info("StandardHypoTestDemo","Pdf is not extended: number of events to generate taken from observed data set is %d",nEvents);
334  sampler->SetNEventsPerToy(nEvents);
335  }
336  else {
337  Info("StandardHypoTestDemo","using a number counting pdf");
338  sampler->SetNEventsPerToy(1);
339  }
340  }
341 
342  if (data->isWeighted() && !generateBinned) {
343  Info("StandardHypoTestDemo","Data set is weighted, nentries = %d and sum of weights = %8.1f but toy generation is unbinned - it would be faster to set generateBinned to true\n",data->numEntries(), data->sumEntries());
344  }
346 
347 
348  // set the test statistic
349  if (testStatType == 0) sampler->SetTestStatistic(slrts);
350  if (testStatType == 1) sampler->SetTestStatistic(ropl);
351  if (testStatType == 2 || testStatType == 3) sampler->SetTestStatistic(profll);
352 
353  }
354 
355  HypoTestResult * htr = hypoCalc->GetHypoTest();
356  htr->SetPValueIsRightTail(true);
357  htr->SetBackgroundAsAlt(false);
358  htr->Print(); // how to get meaningfull CLs at this point?
359 
360  delete sampler;
361  delete slrts;
362  delete ropl;
363  delete profll;
364 
365  if (calcType != 2) {
366  HypoTestPlot * plot = new HypoTestPlot(*htr,100);
367  plot->SetLogYaxis(true);
368  plot->Draw();
369  }
370  else {
371  std::cout << "Asymptotic results " << std::endl;
372 
373  }
374 
375  // look at expected significances
376  // found median of S+B distribution
377  if (calcType != 2) {
378 
379  SamplingDistribution * altDist = htr->GetAltDistribution();
380  HypoTestResult htExp("Expected Result");
381  htExp.Append(htr);
382  // find quantiles in alt (S+B) distribution
383  double p[5];
384  double q[5];
385  for (int i = 0; i < 5; ++i) {
386  double sig = -2 + i;
387  p[i] = ROOT::Math::normal_cdf(sig,1);
388  }
389  std::vector<double> values = altDist->GetSamplingDistribution();
390  TMath::Quantiles( values.size(), 5, &values[0], q, p, false);
391 
392  for (int i = 0; i < 5; ++i) {
393  htExp.SetTestStatisticData( q[i] );
394  double sig = -2 + i;
395  std::cout << " Expected p -value and significance at " << sig << " sigma = "
396  << htExp.NullPValue() << " significance " << htExp.Significance() << " sigma " << std::endl;
397 
398  }
399  }
400  else {
401  // case of asymptotic calculator
402  for (int i = 0; i < 5; ++i) {
403  double sig = -2 + i;
404  // sigma is inverted here
405  double pval = AsymptoticCalculator::GetExpectedPValues( htr->NullPValue(), htr->AlternatePValue(), -sig, false);
406  std::cout << " Expected p -value and significance at " << sig << " sigma = "
407  << pval << " significance " << ROOT::Math::normal_quantile_c(pval,1) << " sigma " << std::endl;
408 
409  }
410  }
411 
412 }
413 
virtual Double_t sumEntries() const =0
virtual Bool_t AccessPathName(const char *path, EAccessMode mode=kFileExists)
Returns FALSE if one can access a file using the specified access mode.
Definition: TSystem.cxx:1213
void SetBackgroundAsAlt(Bool_t l=kTRUE)
ModelConfig is a simple class that holds configuration information specifying how a model should be u...
Definition: ModelConfig.h:52
RooAbsPdf * GetPdf() const
get model PDF (return NULL if pdf has not been specified or does not exist)
Definition: ModelConfig.h:244
RooAbsPdf * MakeNuisancePdf(RooAbsPdf &pdf, const RooArgSet &observables, const char *name)
RooAbsData * data(const char *name) const
Retrieve dataset (binned or unbinned) with given name. A null pointer is returned if not found...
double normal_quantile_c(double z, double sigma)
Inverse ( ) of the cumulative distribution function of the upper tail of the normal (Gaussian) distri...
#define assert(cond)
Definition: unittest.h:542
virtual void SetName(const char *name)
Change (i.e.
Definition: TNamed.cxx:128
RooArgSet * getObservables(const RooArgSet &set, Bool_t valueOnly=kTRUE) const
Definition: RooAbsArg.h:194
std::vector< double > values
Definition: TwoHistoFit2D.C:32
HypoTestResult is a base class for results from hypothesis tests.
static const char * filename()
virtual void SetNEventsPerToy(const Int_t nevents)
Definition: ToyMCSampler.h:177
#define gROOT
Definition: TROOT.h:344
Basic string class.
Definition: TString.h:137
This class provides the plots for the result of a study performed with any of the HypoTestCalculatorG...
Definition: HypoTestPlot.h:37
double poiValue
SamplingDistribution * GetAltDistribution(void) const
void SetLogYaxis(Bool_t ly)
changes plot to log scale on y axis
void SetAltParameters(const RooArgSet &altParameters)
const std::vector< Double_t > & GetSamplingDistribution() const
Get test statistics values.
const RooArgSet * GetNuisanceParameters() const
get RooArgSet containing the nuisance parameters (return NULL if not existing)
Definition: ModelConfig.h:250
static TFile * Open(const char *name, Option_t *option="", const char *ftitle="", Int_t compress=1, Int_t netopt=0)
Create / open a file.
Definition: TFile.cxx:3851
bool generateBinned
virtual void SetTestStatistic(TestStatistic *testStatistic, unsigned int i)
Definition: ToyMCSampler.h:214
RooAbsArg * first() const
void plot(TString fname="data.root", TString var0="var0", TString var1="var1")
Definition: createData.C:17
virtual ModelConfig * Clone(const char *name="") const
clone
Definition: ModelConfig.h:76
Common base class for the Hypothesis Test Calculators.
virtual void Print(Option_t *options=0) const
Print TNamed name and title.
Definition: RooAbsArg.h:227
virtual HypoTestResult * GetHypoTest() const
inherited methods from HypoTestCalculator interface
TestStatSampler * GetTestStatSampler(void) const
void SetNullParameters(const RooArgSet &nullParameters)
void Quantiles(Int_t n, Int_t nprob, Double_t *x, Double_t *quantiles, Double_t *prob, Bool_t isSorted=kTRUE, Int_t *index=0, Int_t type=7)
Computes sample quantiles, corresponding to the given probabilities Parameters: x -the data sample n ...
Definition: TMath.cxx:1173
double normal_cdf(double x, double sigma=1, double x0=0)
Cumulative distribution function of the normal (Gaussian) distribution (lower tail).
void Print(const Option_t *="") const
RooAbsPdf * pdf(const char *name) const
Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found.
void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.
void Info(const char *location, const char *msgfmt,...)
bool noSystematics
tuple np
Definition: multifit.py:30
virtual Bool_t isWeighted() const
Definition: RooAbsData.h:92
Double_t getVal(const RooArgSet *set=0) const
Definition: RooAbsReal.h:64
virtual Double_t NullPValue() const
Return p-value for null hypothesis.
void Error(const char *location, const char *msgfmt,...)
RooRealVar represents a fundamental (non-derived) real valued object.
Definition: RooRealVar.h:37
tuple infile
Definition: mrt.py:15
virtual void setVal(Double_t value)
Set value of variable to 'value'.
Definition: RooRealVar.cxx:203
const int nEvents
Definition: testRooFit.cxx:42
R__EXTERN TSystem * gSystem
Definition: TSystem.h:545
virtual Int_t numEntries() const
Definition: RooAbsData.cxx:293
tuple w
Definition: qtexample.py:51
virtual const char * GetName() const
Returns name of object.
Definition: TNamed.h:51
RooAbsData is the common abstract base class for binned and unbinned datasets.
Definition: RooAbsData.h:37
ProfileLikelihoodTestStat is an implementation of the TestStatistic interface that calculates the pro...
ToyMCSampler is an implementation of the TestStatSampler interface.
Definition: ToyMCSampler.h:99
virtual Double_t Significance() const
familiar name for the Null p-value in terms of 1-sided Gaussian significance
void Warning(const char *location, const char *msgfmt,...)
void SetTestStatisticData(const Double_t tsd)
This class simply holds a sampling distribution of some test statistic.
TestStatistic that returns the ratio of profiled likelihoods.
virtual void Append(const HypoTestResult *other)
add values from another HypoTestResult
This class implements the Hypothesis test calculation using an hybrid (frequentist/bayesian) procedur...
const RooArgSet * GetObservables() const
get RooArgSet for observables (return NULL if not existing)
Definition: ModelConfig.h:259
void SetPValueIsRightTail(Bool_t pr)
void StandardHypoTestDemo(const char *infile="", const char *workspaceName="combined", const char *modelSBName="ModelConfig", const char *modelBName="", const char *dataName="obsData", int calcType=0, int testStatType=3, int ntoys=5000, bool useNC=false, const char *nuisPriorName=0)
tuple file
Definition: fildir.py:20
Bool_t canBeExtended() const
Definition: RooAbsPdf.h:216
bool SetAllConstant(const RooAbsCollection &coll, bool constant=true)
Definition: RooStatsUtils.h:93
const RooArgSet * GetSnapshot() const
get RooArgSet for parameters for a particular hypothesis (return NULL if not existing) ...
TestStatistic class that returns -log(L[null] / L[alt]) where L is the likelihood.
TObject * obj(const char *name) const
Return any type of object (RooAbsArg, RooAbsData or generic object) with given name) ...
virtual Double_t AlternatePValue() const
Return p-value for alternate hypothesis.
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Definition: RooAbsPdf.h:41
void Print(Option_t *opts=0) const
Print contents of the workspace.
Hypothesis Test Calculator using a full frequentist procedure for sampling the test statistic distrib...
double nToysRatio
const RooArgSet * GetParametersOfInterest() const
get RooArgSet containing the parameter of interest (return NULL if not existing)
Definition: ModelConfig.h:247
void SetGenerateBinned(bool binned=true)
Definition: ToyMCSampler.h:233
int printLevel
Int_t getSize() const
RooAbsPdf * GetPriorPdf() const
get parameters prior pdf (return NULL if not existing)
Definition: ModelConfig.h:256
float * q
Definition: THbookFile.cxx:87
virtual void SetSnapshot(const RooArgSet &set)
set parameter values for a particular hypothesis if using a common PDF by saving a snapshot in the wo...
Hypothesis Test Calculator based on the asymptotic formulae for the profile likelihood ratio...
virtual Bool_t add(const RooAbsArg &var, Bool_t silent=kFALSE)
Add element to non-owning set.
Definition: RooArgSet.cxx:448
The RooWorkspace is a persistable container for RooFit projects.
Definition: RooWorkspace.h:42