#include #include "RooRealVar.h" #include "RooDataSet.h" #include #include "RooFormulaVar.h" #include "RooAddPdf.h" #include "RooExponential.h" #include "TPad.h" #include "TCanvas.h" #include "TAxis.h" #include "RooPlot.h" #include "TTree.h" #include "TCut.h" #include "TMath.h" #include "TFile.h" #include "TCanvas.h" #include "RooHist.h" #include "RooGenericPdf.h" #include "RooTruthModel.h" #include "RooGaussModel.h" #include "RooDecay.h" #include "RooProdPdf.h" #include "RooEffProd.h" using namespace std; using namespace RooFit; void massfit(){ ofstream gttt0,gttt1,gttt,gttt3,gttpure,gttt4,gttpure1,gttst0,gttst1; gttst0.open("subsample_data_result_ch0_v04_2500_st0-try1.txt", ios::out | ios::app); gttt0.open("subsample_data_result_ch0_v04_2500_try1.txt", ios::out | ios::app); gttt.open("subsample_data_result_v04_2500_try1.txt", ios::out | ios::app); gttt3.open("subsample_data_result_failed_ch0_v04_2500_try1.txt", ios::out | ios::app); gttpure.open("subsample_data_result_ch0_v04_2500_pure_try1.txt", ios::out | ios::app); int size,unsucess0=0,unsucess1=0; cout<<"size for fitting "<>size; for(int fc=0;fcGet("events"); int nentries_ = tree->GetEntries(); cout<<"\nSubsample Number "< total entries in data tree = " << nentries_ << endl; double init=5.17; double fin=5.5; double bin=2*(fin-init)*100; RooRealVar* m=new RooRealVar("m","M(B^{+})[GeV/c^{2}]",init,fin); RooRealVar* treco=new RooRealVar("treco","t_{reco}[ps]",1.4,12.5) ; RooRealVar* trecoe=new RooRealVar("trecoe","terr[ps]",0.009,0.30); RooRealVar* t1=new RooRealVar("t1","t",1.63, 0.0, 10.0) ; RooRealVar* t=new RooRealVar("t","t",1.637, 0.0, 10.0) ; RooRealVar* chan=new RooRealVar("chan","",-3,5); RooDataSet* data1=new RooDataSet("data1","ntple candidates", RooArgSet(*m,*treco,*trecoe,*chan) , Import(*tree)); TCut c1="chan==0"; TCut c2="chan==1"; TCut c3="m>5.4&&m<5.5"; RooDataSet* datach0=(RooDataSet*)data1->reduce(c1); RooDataSet* datach1=(RooDataSet*)data1->reduce(c2); RooDataSet* sbdata=(RooDataSet*)data1->reduce(c3); ////-------------------------------------------- ///..........Data_Bmass..... ... ////------------------------------------------- RooRealVar* lambda=new RooRealVar("lambda","lambda",4.441,-34,34); RooExponential* Expo=new RooExponential("Expo","Exponential background",*m,*lambda); RooRealVar* a0=new RooRealVar("a0","",1); RooRealVar* a1=new RooRealVar("a1","",-10,10); RooRealVar* a2=new RooRealVar("a2","a2",0.4,0,20); RooPolynomial* poly=new RooPolynomial("poly","",*m,RooArgList(*a0,*a2)); RooRealVar* B1=new RooRealVar("B1", "B_1comb", 0.3758, 0.0 ,1); RooFormulaVar* B2=new RooFormulaVar("B2", "B_2comb", "1.-@0", RooArgList(*B1)); RooBernstein* comb=new RooBernstein("comb", "mass_comb", *m, RooArgSet(*B1,*B2)); RooRealVar* mean=new RooRealVar("mean", "mean",5.27,5.24, 5.30); RooRealVar* sigma1=new RooRealVar("sigma1","sigma1 of Gaussian",0.0259,0., 1.0); RooRealVar* sigma2=new RooRealVar("sigma2","sigma2 of Gaussian",0.0986, 0., 1.0); RooGaussian* gausian1=new RooGaussian("gausian1","gaus1", *m, *mean, *sigma1); RooGaussian* gausian2=new RooGaussian("gausian2","gaus2", *m, *mean, *sigma2); RooRealVar* fg1=new RooRealVar("fg1","fg1", 0.893, 0., 1.); RooAddPdf* mass1=new RooAddPdf("mass1","mass1", RooArgList(*gausian1,*gausian2), *fg1); RooRealVar* nbkg=new RooRealVar("nbkg","Background fraction",13,0.0,100); RooRealVar* Nsig=new RooRealVar("Nsig","signal fraction",59,0,100); RooAddPdf* model=new RooAddPdf("model","",RooArgSet(*mass1,*Expo),RooArgSet(*Nsig,*nbkg)); mean->setVal(5.279); // 1D fit to Background Decay time distribution RooRealVar* meant=new RooRealVar("meant", "mean", 0); RooRealVar* sigmat=new RooRealVar("sigmat","sigma",1.0); RooGaussModel* gm=new RooGaussModel("gm","gm", *treco,*meant,*sigmat,*trecoe); ////...background... RooRealVar*fracbkg1=new RooRealVar("fracbkg1","",0.25,0,1); RooRealVar* TauBkg3=new RooRealVar("TauBkg3","",1.50,0,3); RooTruthModel* CTM =new RooTruthModel("CTM","",*treco); RooDecay* Bck_Ctau3 =new RooDecay("Bck_Ctau3","",*treco,*TauBkg3,*CTM,RooDecay::SingleSided); // Decay time and mass fit 2D with the Decay time uncertainty TFile* EffFile = TFile::Open("BsFitEffWS_ancut_totaleff_difint_v04_test1.root"); RooWorkspace* wEff = (RooWorkspace*) EffFile->Get("EffWorkSpace"); RooFormulaVar* effForch0 =(RooFormulaVar*) wEff->function("effForch1"); RooFormulaVar* effForch1 =(RooFormulaVar*) wEff->function("effFor"); RooDecay* decay_sig=new RooDecay("decay_sig","decay",*treco,*t,*gm,RooDecay::SingleSided) ; RooEffProd* CtEffSig =new RooEffProd("CtEffSig","",*decay_sig,*effForch0); // treco->setRange("partreco",1.4,12.5); RooProdPdf* SigTot=new RooProdPdf ("SigTot","",RooArgList(*mass1,*CtEffSig) ); RooProdPdf* BkgTot=new RooProdPdf ("BkgTot","",RooArgList(*Expo,*Bck_Ctau3) ); RooAddPdf* TotPdf=new RooAddPdf("TotPdf"," Signal + Bkg Pdf",RooArgList(*SigTot,*BkgTot),RooArgList(*Nsig,*nbkg)); t->setVal(1.5823); RooMsgService::instance().setGlobalKillBelow(RooFit::FATAL); RooMsgService::instance().setSilentMode(kTRUE); RooMsgService::instance().setStreamStatus(1,false); RooMsgService::instance().getStream(1).removeTopic(Integration) ; RooMsgService::instance().getStream(1).removeTopic(Minimization) ; RooMsgService::instance().getStream(1).removeTopic(Fitting) ; RooMsgService::instance().getStream(1).removeTopic(NumIntegration) ; RooMsgService::instance().getStream(1).removeTopic(Optimization) ; RooMsgService::instance().getStream(1).removeTopic(ObjectHandling) ; RooMsgService::instance().getStream(1).removeTopic(Eval) ; RooMsgService::instance().Print() ; RooFitResult* fitres=TotPdf->fitTo(*datach0,Extended(kTRUE),ConditionalObservables(*trecoe),Save(kTRUE),NumCPU(7)); fitres->Print("v"); double nsigch0=Nsig->getVal(); double nbkgch0=nbkg->getVal(); double nbkgch0err=nbkg->getError(); double men0=mean->getVal(); double Tau0=t->getVal(); double tauer=t->getError(); int st0=fitres->status(); int qual0=fitres->covQual(); cout<<"tau "<nbkgch0err){ gttpure<Delete(); f->Close(); TotPdf->Delete(); datach0->Delete(); data1->Delete(); } gttt0.close(); gttt3.close(); gttpure.close(); gttst0.close(); }