Hi all,
since I have multiple channels in my analysis, with different binnings etc, I’d like to put them all together in a for loop. I’m encountering a segmentation violation which brought me to just test a much simpler case, i.e. I’ve edited $ROOTSYS/tutorials/roofit/rf501_simultaneouspdf.C to be like what you see below (different RooRealVar for each channel, and use a map<string,RooDataSet*> for the data), and in this way I find different fit results from the nominal case.
This tells either I’m doing something wrong or there is a bug somewhere. Any help ?
With this implementation, I also wonder about how ROOT knows which RooRealVar in the RooArgSet refers to which pdf/dataset.
void rf501_simultaneouspdf()
{
// C r e a t e m o d e l f o r p h y s i c s s a m p l e
// -------------------------------------------------------------
// Create observables
RooRealVar x("x","x",-8,8) ;
RooRealVar x1("x1","x1",-8,8) ; // NEW
// Construct signal pdf
RooRealVar mean("mean","mean",0,-8,8) ;
RooRealVar sigma("sigma","sigma",0.3,0.1,10) ;
RooGaussian gx("gx","gx",x,mean,sigma) ;
// Construct background pdf
RooRealVar a0("a0","a0",-0.1,-1,1) ;
RooRealVar a1("a1","a1",0.004,-1,1) ;
RooChebychev px("px","px",x,RooArgSet(a0,a1)) ;
// Construct composite pdf
RooRealVar f("f","f",0.2,0.,1.) ;
RooAddPdf model("model","model",RooArgList(gx,px),f) ;
// C r e a t e m o d e l f o r c o n t r o l s a m p l e
// --------------------------------------------------------------
RooRealVar sigma1("sigma1","sigma1",0.3,0.1,10) ; // NEW
// Construct signal pdf.
// NOTE that sigma is shared with the signal sample model
RooRealVar mean_ctl("mean_ctl","mean_ctl",-3,-8,8) ;
RooGaussian gx_ctl("gx_ctl","gx_ctl",x1,mean_ctl,sigma1) ; // NEW
// Construct the background pdf
RooRealVar a0_ctl("a0_ctl","a0_ctl",-0.1,-1,1) ;
RooRealVar a1_ctl("a1_ctl","a1_ctl",0.5,-0.1,1) ;
RooChebychev px_ctl("px_ctl","px_ctl",x1,RooArgSet(a0_ctl,a1_ctl)) ; // NEW
// Construct the composite model
RooRealVar f_ctl("f_ctl","f_ctl",0.5,0.,1.) ;
RooAddPdf model_ctl("model_ctl","model_ctl",RooArgList(gx_ctl,px_ctl),f_ctl) ;
// G e n e r a t e e v e n t s f o r b o t h s a m p l e s
// ---------------------------------------------------------------
// Generate 1000 events in x and y from model
RooDataSet *data = model.generate(RooArgSet(x),100) ;
RooDataSet *data_ctl = model_ctl.generate(RooArgSet(x1),2000) ; // NEW
// C r e a t e i n d e x c a t e g o r y a n d j o i n s a m p l e s
// ---------------------------------------------------------------------------
// Define category to distinguish physics and control samples events
RooCategory sample("sample","sample") ;
sample.defineType("physics") ;
sample.defineType("control") ;
// NEW starts
std::map<std::string,RooDataSet*> dataset;
dataset["physics"] = data;
dataset["control"] = data_ctl;
RooArgSet Observables;
Observables.add(x);
Observables.add(x1);
// NEW ends
// Construct combined dataset in (x,sample)
RooDataSet combData("combData","combined data",Observables,Index(sample),Import(dataset)) ; // NEW
// RooDataSet combData("combData","combined data",x,Index(sample),Import("physics",*data),Import("control",*data_ctl)) ;
// C o n s t r u c t a s i m u l t a n e o u s p d f i n ( x , s a m p l e )
// -----------------------------------------------------------------------------------
// Construct a simultaneous pdf using category sample as index
RooSimultaneous simPdf("simPdf","simultaneous pdf",sample) ;
// Associate model with the physics state and model_ctl with the control state
simPdf.addPdf(model,"physics") ;
simPdf.addPdf(model_ctl,"control") ;
// P e r f o r m a s i m u l t a n e o u s f i t
// ---------------------------------------------------
// Perform simultaneous fit of model to data and model_ctl to data_ctl
simPdf.fitTo(combData) ;
}