RooMCStudy

Roofit and RooStats
1

I am using RooMCStudy to run over 5000 toy sets

mcs->generateAndFit(5000,nFullAllEvents,kTRUE,"data/toymc_%04d.dat");

I’m checking how often the null hypothesis fluctuations will fake signal. Some rare fluctuations push fake signal up to 40 events based on the plot nX v.s. numGenerated. I’d like to see this distribution. I understand you can save the generated toy sets to data files, which I do above. I then grab the file which has more than 40 events with:

    for(Int_t i=1;i<5000;i++){
      if((mcs->fitParDataSet().get(i)->getRealValue("numGamma")) > 40){                                  
        cout << "i: " << i << " " << mcs->fitParDataSet().get(i)->getRealValue("numGamma") << endl;      
      }                                                                                                  
    }

This gives the iteration of the file. But when I try and fit the data in this dat file I don’t get the +40 events from the fit that I expect.

Is there some way to save the results of these 5000 fits as images?
If not, is there a better way for me to pinpoint these extreme distributions and see what the resulting fit looks like?

Thanks,
james

2

Hi James,

What you do is in principle the correct way, and you should be able to reproduce the fit outside of RooMCStudy if you take the corresponding input dataset. Do you make sure that the initial conditions for your
’manual’ fit are the same as used in RooMCStudy and that all fit options are indeed the same? (RooMCStudy uses the parameters of the original input p.d.f. as starting point for each fit)

On the other hand, if you just want to see what the p.d.f looks like for a given toymc fit result, you can do it as follows

RooAbsPdf* myPdf ; // your model (also used as input to RooMCStudy)

RooArgSet* params = myPdf->getParameters(obs) ; // where obs is your observable

// Set parameter values to that from a given row in the fitpar dataset
*params = *mcs->fitParDataSet().get(i)

// Now plot the p.d.f.
RooPlot* frame = obs.frame() ;
myPdf->plotOn(frame) ;
frame->Draw()

Wouter