Pull histogram with multiple ranges

Roofit and RooStats
1

Hi,

I’m having a problem when drawing a pull histogram after a one dimensional fit with two ranges.
I define the “left” and “right” ranges as follows:

mass.setRange("left",377,540); mass.setRange("right",551,642);

And then I fit my RooAddPdf in the sum of the two ranges, as follows:

The fit works perfectly fine. However, when I plot the pull histogram:

RooHist* hpull = massplot->pullHist("data_hist","tot_pdf"); RooPlot* pull = mass.frame(); pull->addPlotable(hpull,"P");

only pulls in the “right” range are plotted. If I swap “left” and “right” in the fitting line, i.e.:

only pulls in the “left” range are plotted.
Is there a way to solve this problem?

Cheers,
Lorenzo

2

Hi,

Can you please post your full running script, so we can reproduce the problem ?
Thank you

Lorenzo

3

Hi Lorenzo,

here is the code:

[code]#include
#include
#include
#include
#include <TFile.h>
#include <TROOT.h>
#include <TH1.h>
#include <TTree.h>
#include <TPaveText.h>
#include <TAxis.h>
#include <TCanvas.h>
#include <TSystem.h>
#include <TStyle.h>
#include <TString.h>
#include <RooRealVar.h>
#include <RooDataSet.h>
#include <RooDataHist.h>
#include <RooGaussian.h>
#include <RooChebychev.h>
#include <RooPolynomial.h>
#include <RooCBShape.h>
#include <RooExponential.h>
#include <RooAddPdf.h>
#include <RooGenericPdf.h>
#include <RooProdPdf.h>
#include <RooEffProd.h>
#include <RooFitResult.h>
#include <RooPlot.h>
#include <RooHist.h>
#include <RooMCStudy.h>
#include <RooRealIntegral.h>
#include <TIterator.h>

void fit(string fileName){

using namespace RooFit;

//Open data file
TFile* inputFile = TFile::Open(fileName.c_str());
TTree* inputTree;
RooRealVar M_OS("M_OS","M_OS",500,630,"MeV/c^{2}");
M_OS.setBins(130);
TList* keyList = inputFile->GetListOfKeys();
char treeName[50];
sprintf (treeName, "%s", keyList->At(0)->GetName()); 
inputTree = (TTree*)inputFile->Get(treeName);

//Open Ds2eta shape file
TFile* histfile = TFile::Open("Ds2eta.root");
TH1D *hist = (TH1D*) histfile->Get("hetaDs");

RooDataHist *pdfh = new RooDataHist("pdfh", "pdfh", RooArgList(M_OS), hist);

//Set ranges
M_OS.setRange("sx",500,544);
M_OS.setRange("dx",552,630);

//Import data from TTree and create a RooDataHist
RooDataSet *data = new RooDataSet("data","data",RooArgList(M_OS),Import(*inputTree),CutRange("sx,dx"));
RooDataHist *datah = new RooDataHist("datah","datah",RooArgList(M_OS),*data);

Int_t nEntries = data->numEntries();
cout << "Number of entries = " << nEntries << endl;

//The mass plot
RooPlot *massplot = M_OS.frame();
datah->plotOn(massplot,Name("datahistogram"));


// **********************THE MASS FIT********************************
// **********************THE PARAMETERS***************************

RooRealVar numBkg("numBkg","numBkg", 2.0e+05, 1.4321e+05, 2.8321e+05); 
RooRealVar bkg_p0("bkg_p0","bkg_p0",-0.1,-1,1);
RooRealVar bkg_p1("bkg_p1","bkg_p1",5e-03,-1,1);
RooRealVar bkg_p2("bkg_p2","bkg_p2", 0 , -1, 1 );
RooRealVar bkg_p3("bkg_p3","bkg_p3",0,-1,1);
RooRealVar bkg_p4("bkg_p4","bkg_p4",0,-1,1);
RooRealVar bkg_p5("bkg_p5","bkg_p5",-0.0006,-1,1);

RooRealVar numHist("numHist","numHist", 1700 );


// **********************THE PDFS*********************************

RooHistPdf* histpdf = new RooHistPdf("histpdf","histpdf",M_OS,*pdfh,5);

RooChebychev* bkg = new RooChebychev("bkg","bkg", M_OS, RooArgList(bkg_p0, bkg_p1, bkg_p2, bkg_p3));

RooFormulaVar* step = new RooFormulaVar("step", "step", "M_OS < 544 || M_OS > 552", RooArgList(M_OS));
RooFormulaVar* step2 = new RooFormulaVar("step2", "step2", "M_OS < 544 || M_OS > 552", RooArgList(M_OS));

RooEffProd* blind_bkg = new RooEffProd("blind_bkg", "blind_bkg", *bkg, *step);
RooEffProd* blind_eta = new RooEffProd("blind_eta", "blind_eta", *histpdf, *step2);

RooAddPdf* blind = new RooAddPdf("blind", "blind", RooArgList(*bkg,*histpdf), RooArgList(numBkg,numHist));


// **********************THE FIT*********************************

RooFitResult *mass_result;

mass_result = blind->fitTo(*datah,Extended(),Save(),Range("sx,dx"));

mass_result->Print();

//Check fit results
RooArgList list = mass_result->floatParsFinal();
for (Int_t i = 0; i < list.getSize(); i++)
{
	RooRealVar* temp = (RooRealVar*)list.at(i);
if ((temp->getVal() + 3.*temp->getError()) > temp->getMax())
{
	cout << "Parameter " << temp->GetName() << " exceeds range at 3 sigma!" << endl;
}
else if ((temp->getVal() - 3.*temp->getError()) < temp->getMin())
{
	cout << "Parameter " << temp->GetName() << " exceeds range at 3 sigma!" << endl;
}
}

//We update the mass plot from earlier to show the fit to the data. 
blind->plotOn(massplot, Components(RooArgSet(*bkg)), LineStyle(5), LineColor(kBlack), Name("fit"));
blind->plotOn(massplot, Components(RooArgSet(*histpdf)), LineStyle(5), LineColor(kBlack), Name("fit2"));
blind->plotOn(massplot, LineColor(kBlue), Name("blindtot"));

massplot->GetXaxis()->SetTitle("#it{m}(#it{#pi}^{ #plus}#it{#pi}^{ #minus}) [MeV/#it{c}^{2}]");
massplot->GetYaxis()->SetTitle("Candidates / (1 MeV/#it{c}^{2})");
massplot->SetTitle("");

RooHist* hresid = massplot->pullHist("datahistogram","blindtot");
RooPlot* resid = M_OS.frame();
resid->addPlotable(hresid,"P"); 
resid->SetTitle("");
resid->GetXaxis()->SetTitle("#it{m}(#it{#pi}^{ #plus}#it{#pi}^{ #minus}) [MeV/#it{c}^{2}]");
gStyle->SetPadLeftMargin(0.1);

TCanvas *c = new TCanvas("datac","datac",1440,900);
TPad *pad1 = new TPad("pad1", "The pad 80 of the height",0.0,0.4,1.0,1.0);
TPad *pad2 = new TPad("pad2", "The pad 20 of the height",0.0,0.0,1.0,0.4);
pad1->Draw();
pad2->Draw();
pad1->cd();
massplot->Draw();
pad2->cd();
resid->Draw();

}

void fitSpectrum(){

vector<string>  files;

//Set files
files.push_back("Spectrum_fit.root");

for (int i = 0; i < files.size(); i++)
{
	cout << "Input file = " << files.at(i) << endl;
fit(files.at(i));
}

}[/code]

Find attached also the necessary root files to run the script.
I have tested it with ROOT 6.04/02 and RooFit v3.60 and here is the result:

test.pdf (36.6 KB)

Thank you for your time.

Cheers,
Lorenzo
Spectrum_fit.root (622 KB)
Ds2eta.root (8.02 KB)

4

Hi,

Thank you for attaching the files. I can reproduce your problem now. I’ll investigate it

Best

Lorenzo

5

Hi,

You need to explicitly retrieve the two different curve which are drawn when fitting the separate ranges.
Attached you find your modified script which uses RooHist::MakePullHist instead of RooPlot::pullHist

Best Regards

LorenzorootForum_21406.C (5.7 KB)

6

Hi,

Thank you, now everything is working properly.

Cheers,
Lorenzo