29 using namespace RooFit ;
30 using namespace RooStats ;
37 TString card_name(
"rs603_card_WsMaker.rs");
38 ofstream ofile(card_name);
39 ofile <<
"// The simplest card for combination\n\n";
40 ofile <<
"gauss1 = Gaussian(x[0,100],mean1[50,0,100],4);\n";
41 ofile <<
"flat1 = Polynomial(x,0);\n";
42 ofile <<
"sb_model1 = SUM(nsig1[120,0,300]*gauss1 , nbkg1[100,0,1000]*flat1);\n\n";
43 ofile <<
"echo In the middle!;\n\n";
44 ofile <<
"gauss2 = Gaussian(x,mean2[80,0,100],5);\n";
45 ofile <<
"flat2 = Polynomial(x,0);\n";
46 ofile <<
"sb_model2 = SUM(nsig2[90,0,400]*gauss2 , nbkg2[80,0,1000]*flat2);\n\n";
47 ofile <<
"echo At the end!;\n";
50 TString card_name2(
"rs603_card.rs");
51 ofstream ofile2(card_name2);
52 ofile2 <<
"// The simplest card for combination\n\n";
53 ofile2 <<
"gauss1 = Gaussian(x[0,100],mean1[50,0,100],4);\n";
54 ofile2 <<
"flat1 = Polynomial(x,0);\n";
55 ofile2 <<
"sb_model1 = SUM(nsig1[120,0,300]*gauss1 , nbkg1[100,0,1000]*flat1);\n\n";
56 ofile2 <<
"echo In the middle!;\n\n";
57 ofile2 <<
"gauss2 = Gaussian(x,mean2[80,0,100],5);\n";
58 ofile2 <<
"flat2 = Polynomial(x,0);\n";
59 ofile2 <<
"sb_model2 = SUM(nsig2[90,0,400]*gauss2 , nbkg2[80,0,1000]*flat2);\n\n";
60 ofile2 <<
"#include rs603_included_card.rs;\n\n";
61 ofile2 <<
"echo At the end!;\n";
64 TString card_name3(
"rs603_included_card.rs");
65 ofstream ofile3(card_name3);
66 ofile3 <<
"echo Now reading the included file!;\n\n";
67 ofile3 <<
"echo Including datasets in a Workspace in a Root file...;\n";
68 ofile3 <<
"data1 = import(rs603_infile.root,\n";
69 ofile3 <<
" rs603_ws,\n";
70 ofile3 <<
" data1);\n\n";
71 ofile3 <<
"data2 = import(rs603_infile.root,\n";
72 ofile3 <<
" rs603_ws,\n";
73 ofile3 <<
" data2);\n";
79 "rs603_card_WsMaker.rs",
98 TFile
outfile(
"rs603_infile.root",
"RECREATE");
102 cout <<
"-------------------------------------------------------------------\n"
103 <<
" Rootfile and Workspace prepared \n"
104 <<
"-------------------------------------------------------------------\n";
107 HLFactory hlf_2(
"HLFactoryElaborateExample",
112 pdf1 = hlf_2.
GetWs()->
pdf(
"sb_model1");
113 pdf2 = hlf_2.
GetWs()->
pdf(
"sb_model2");
115 hlf_2.
AddChannel(
"model1",
"sb_model1",
"flat1",
"data1");
116 hlf_2.
AddChannel(
"model2",
"sb_model2",
"flat2",
"data2");
135 gROOT->SetStyle(
"Plain");
virtual Int_t Write(const char *name=0, Int_t option=0, Int_t bufsize=0)
Write this object to the current directory.
void SetName(const char *name)
Change the name of this dataset into the given name.
int AddChannel(const char *label, const char *SigBkgPdfName, const char *BkgPdfName=0, const char *datasetName=0)
Add channel for the combination.
virtual Bool_t setIndex(Int_t index, Bool_t printError=kTRUE)
Set value by specifying the index code of the desired state.
RooDataSet * GetTotDataSet()
Get the combined dataset.
RooCmdArg ProjWData(const RooAbsData &projData, Bool_t binData=kFALSE)
RooCmdArg Extended(Bool_t flag=kTRUE)
HLFactory is an High Level model Factory allows you to describe your models in a configuration file (...
RooAbsPdf * pdf(const char *name) const
Retrieve p.d.f (RooAbsPdf) with given name. A null pointer is returned if not found.
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
Plot dataset on specified frame.
RooRealVar represents a fundamental (non-derived) real valued object.
void rs603_HLFactoryElaborateExample()
virtual RooPlot * plotOn(RooPlot *frame, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none(), const RooCmdArg &arg9=RooCmdArg::none(), const RooCmdArg &arg10=RooCmdArg::none()) const
Plot (project) PDF on specified frame.
RooPlot * frame(const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none()) const
RooDataSet is a container class to hold unbinned data.
RooCategory represents a fundamental (non-derived) discrete value object.
RooRealVar * var(const char *name) const
Retrieve real-valued variable (RooRealVar) with given name. A null pointer is returned if not found...
A RooPlot is a plot frame and a container for graphics objects within that frame. ...
RooCategory * GetTotCategory()
Get the combined dataset.
RooWorkspace * GetWs()
Get the RooWorkspace containing the models and variables.
RooAbsPdf is the abstract interface for all probability density functions The class provides hybrid a...
Bool_t import(const RooAbsArg &arg, const RooCmdArg &arg1=RooCmdArg(), const RooCmdArg &arg2=RooCmdArg(), const RooCmdArg &arg3=RooCmdArg(), const RooCmdArg &arg4=RooCmdArg(), const RooCmdArg &arg5=RooCmdArg(), const RooCmdArg &arg6=RooCmdArg(), const RooCmdArg &arg7=RooCmdArg(), const RooCmdArg &arg8=RooCmdArg(), const RooCmdArg &arg9=RooCmdArg())
Import a RooAbsArg object, e.g.
RooAbsArg * arg(const char *name) const
Return RooAbsArg with given name. A null pointer is returned if none is found.
RooDataSet * generate(const RooArgSet &whatVars, Int_t nEvents, const RooCmdArg &arg1, const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none())
Generate a new dataset containing the specified variables with events sampled from our distribution...
virtual RooFitResult * fitTo(RooAbsData &data, const RooCmdArg &arg1=RooCmdArg::none(), const RooCmdArg &arg2=RooCmdArg::none(), const RooCmdArg &arg3=RooCmdArg::none(), const RooCmdArg &arg4=RooCmdArg::none(), const RooCmdArg &arg5=RooCmdArg::none(), const RooCmdArg &arg6=RooCmdArg::none(), const RooCmdArg &arg7=RooCmdArg::none(), const RooCmdArg &arg8=RooCmdArg::none())
Fit PDF to given dataset.
RooAbsPdf * GetTotSigBkgPdf()
Get the combined signal plus background pdf.
RooCmdArg Slice(const RooArgSet &sliceSet)
The RooWorkspace is a persistable container for RooFit projects.
virtual void Draw(Option_t *options=0)
Draw this plot and all of the elements it contains.