| 1 |
/****
|
| 2 |
|
| 3 |
Off peak status (RestrictToMassPeak) :
|
| 4 |
|
| 5 |
x Necessary adaptations identified
|
| 6 |
x Started working on necessary adaptations
|
| 7 |
x Necessary adaptations implemented
|
| 8 |
x Necessary adaptations tested
|
| 9 |
|
| 10 |
DONE!
|
| 11 |
|
| 12 |
****/
|
| 13 |
|
| 14 |
#include <iostream>
|
| 15 |
#include <vector>
|
| 16 |
#include <sys/stat.h>
|
| 17 |
|
| 18 |
#include <TCut.h>
|
| 19 |
#include <TROOT.h>
|
| 20 |
#include <TCanvas.h>
|
| 21 |
#include <TMath.h>
|
| 22 |
#include <TColor.h>
|
| 23 |
#include <TPaveText.h>
|
| 24 |
#include <TRandom.h>
|
| 25 |
#include <TH1.h>
|
| 26 |
#include <TH2.h>
|
| 27 |
#include <TF1.h>
|
| 28 |
#include <TSQLResult.h>
|
| 29 |
#ifndef GeneralToolBoxLoaded
|
| 30 |
#include "GeneralToolBox.C"
|
| 31 |
#endif
|
| 32 |
|
| 33 |
using namespace std;
|
| 34 |
|
| 35 |
using namespace PlottingSetup;
|
| 36 |
|
| 37 |
//--------------------------------------------- below: ttbar
|
| 38 |
|
| 39 |
void ttbar_limit_shapes_for_systematic_effect(TFile *limfile, string identifier, string mcjzb, string datajzb, int JES,vector<float> binning, TCanvas *limcan) {
|
| 40 |
dout << "Creatig shape templates ... ";
|
| 41 |
if(identifier!="") dout << "for systematic called "<<identifier;
|
| 42 |
dout << endl;
|
| 43 |
int dataormc=mc;//this is only for tests - for real life you want dataormc=data !!!
|
| 44 |
if(dataormc!=data) write_warning("limit_shapes_for_systematic_effect","WATCH OUT! Not using data for limits!!!! this is ok for tests, but not ok for anything official!");
|
| 45 |
|
| 46 |
TCut limitnJetcut;
|
| 47 |
if(JES==noJES) limitnJetcut=cutnJets;
|
| 48 |
else {
|
| 49 |
if(JES==JESdown) limitnJetcut=cutnJetsJESdown;
|
| 50 |
if(JES==JESup) limitnJetcut=cutnJetsJESup;
|
| 51 |
}
|
| 52 |
TH1F *ZOSSFP = allsamples.Draw("ZOSSFP",datajzb,binning, "JZB4limits", "events",cutmass&&cutOSSF&&limitnJetcut&&basiccut,dataormc,luminosity);
|
| 53 |
TH1F *ZOSOFP = allsamples.Draw("ZOSOFP",datajzb,binning, "JZB4limits", "events",cutmass&&cutOSOF&&limitnJetcut&&basiccut,dataormc,luminosity);
|
| 54 |
TH1F *ZOSSFN = allsamples.Draw("ZOSSFN","-"+datajzb,binning, "JZB4limits", "events",cutmass&&cutOSSF&&limitnJetcut&&basiccut,dataormc,luminosity);
|
| 55 |
TH1F *ZOSOFN = allsamples.Draw("ZOSOFN","-"+datajzb,binning, "JZB4limits", "events",cutmass&&cutOSOF&&limitnJetcut&&basiccut,dataormc,luminosity);
|
| 56 |
|
| 57 |
TH1F *LZOSSFP = allsamples.Draw("LZOSSFP",mcjzb,binning, "JZB4limits", "events",cutmass&&cutOSSF&&limitnJetcut&&basiccut,mc,luminosity,allsamples.FindSample("TTJets"));
|
| 58 |
TH1F *LZOSOFP = allsamples.Draw("LZOSOFP",mcjzb,binning, "JZB4limits", "events",cutmass&&cutOSOF&&limitnJetcut&&basiccut,mc,luminosity,allsamples.FindSample("TTJets"));
|
| 59 |
TH1F *LZOSSFN = allsamples.Draw("LZOSSFN","-"+mcjzb,binning, "JZB4limits", "events",cutmass&&cutOSSF&&limitnJetcut&&basiccut,mc,luminosity,allsamples.FindSample("TTJets"));
|
| 60 |
TH1F *LZOSOFN = allsamples.Draw("LZOSOFN","-"+mcjzb,binning, "JZB4limits", "events",cutmass&&cutOSOF&&limitnJetcut&&basiccut,mc,luminosity,allsamples.FindSample("TTJets"));
|
| 61 |
|
| 62 |
string obsname="data_obs";
|
| 63 |
string predname="background";
|
| 64 |
string signalname="signal";
|
| 65 |
if(identifier!="") {
|
| 66 |
obsname=("data_"+identifier);
|
| 67 |
predname=("background_"+identifier);
|
| 68 |
signalname="signal_"+identifier;
|
| 69 |
}
|
| 70 |
|
| 71 |
TH1F *obs = (TH1F*)ZOSSFP->Clone();
|
| 72 |
obs->SetName(obsname.c_str());
|
| 73 |
obs->Write();
|
| 74 |
TH1F *pred = (TH1F*)ZOSSFN->Clone();
|
| 75 |
pred->Add(ZOSOFN,-1.0);
|
| 76 |
pred->SetName(predname.c_str());
|
| 77 |
pred->Write();
|
| 78 |
|
| 79 |
TH1F *Lobs = new TH1F("Lobs","Lobs",binning.size()-1,&binning[0]);
|
| 80 |
TH1F *Lpred = new TH1F("Lpred","Lpred",binning.size()-1,&binning[0]);
|
| 81 |
Lobs->Add(LZOSSFP);
|
| 82 |
Lpred->Add(LZOSSFN);
|
| 83 |
TH1F *signal = (TH1F*)Lobs->Clone();
|
| 84 |
signal->Add(Lpred,-1);
|
| 85 |
signal->SetName(signalname.c_str());
|
| 86 |
signal->Write();
|
| 87 |
|
| 88 |
delete Lobs;
|
| 89 |
delete Lpred;
|
| 90 |
|
| 91 |
delete ZOSSFP;
|
| 92 |
delete ZOSOFP;
|
| 93 |
delete ZOSSFN;
|
| 94 |
delete ZOSOFN;
|
| 95 |
|
| 96 |
delete LZOSSFP;
|
| 97 |
delete LZOSOFP;
|
| 98 |
delete LZOSSFN;
|
| 99 |
delete LZOSOFN;
|
| 100 |
|
| 101 |
}
|
| 102 |
|
| 103 |
void prepare_ttbar_datacard(TFile *f) {
|
| 104 |
TH1F *dataob = (TH1F*)f->Get("data_obs");
|
| 105 |
TH1F *signal = (TH1F*)f->Get("signal");
|
| 106 |
TH1F *background = (TH1F*)f->Get("background");
|
| 107 |
|
| 108 |
ofstream datacard;
|
| 109 |
ensure_directory_exists(get_directory()+"/limits");
|
| 110 |
datacard.open ((get_directory()+"/limits/ttbardatacard.txt").c_str());
|
| 111 |
datacard << "Writing this to a file.\n";
|
| 112 |
datacard << "imax 1\n";
|
| 113 |
datacard << "jmax 1\n";
|
| 114 |
datacard << "kmax *\n";
|
| 115 |
datacard << "---------------\n";
|
| 116 |
datacard << "shapes * * ttbarlimitfile.root $PROCESS $PROCESS_$SYSTEMATIC\n";
|
| 117 |
datacard << "---------------\n";
|
| 118 |
datacard << "bin 1\n";
|
| 119 |
datacard << "observation "<<dataob->Integral()<<"\n";
|
| 120 |
datacard << "------------------------------\n";
|
| 121 |
datacard << "bin 1 1\n";
|
| 122 |
datacard << "process signal background\n";
|
| 123 |
datacard << "process 0 1\n";
|
| 124 |
datacard << "rate "<<signal->Integral()<<" "<<background->Integral()<<"\n";
|
| 125 |
datacard << "--------------------------------\n";
|
| 126 |
datacard << "lumi lnN 1.10 1.0\n";
|
| 127 |
datacard << "bgnorm lnN 1.00 1.4 uncertainty on our prediction (40%)\n";
|
| 128 |
datacard << "JES shape 1 1 uncertainty on background shape and normalization\n";
|
| 129 |
datacard << "peak shape 1 1 uncertainty on signal resolution. Assume the histogram is a 2 sigma shift, \n";
|
| 130 |
datacard << "# so divide the unit gaussian by 2 before doing the interpolation\n";
|
| 131 |
datacard.close();
|
| 132 |
|
| 133 |
|
| 134 |
|
| 135 |
}
|
| 136 |
|
| 137 |
void prepare_ttbar_limits(string mcjzb, string datajzb, float jzbpeakerrordata, float jzbpeakerrormc, vector<float> jzbbins) {
|
| 138 |
ensure_directory_exists(get_directory()+"/limits");
|
| 139 |
TFile *limfile = new TFile((get_directory()+"/limits/ttbarlimitfile.root").c_str(),"RECREATE");
|
| 140 |
TCanvas *limcan = new TCanvas("limcan","Canvas for calculating limits");
|
| 141 |
ttbar_limit_shapes_for_systematic_effect(limfile,"",mcjzb,datajzb,noJES,jzbbins,limcan);
|
| 142 |
ttbar_limit_shapes_for_systematic_effect(limfile,"peakUp",newjzbexpression(mcjzb,jzbpeakerrormc),newjzbexpression(datajzb,jzbpeakerrordata),noJES,jzbbins,limcan);
|
| 143 |
ttbar_limit_shapes_for_systematic_effect(limfile,"peakDown",newjzbexpression(mcjzb,-jzbpeakerrormc),newjzbexpression(datajzb,-jzbpeakerrordata),noJES,jzbbins,limcan);
|
| 144 |
ttbar_limit_shapes_for_systematic_effect(limfile,"JESUp",mcjzb,datajzb,JESup,jzbbins,limcan);
|
| 145 |
ttbar_limit_shapes_for_systematic_effect(limfile,"JESDown",mcjzb,datajzb,JESdown,jzbbins,limcan);
|
| 146 |
|
| 147 |
prepare_ttbar_datacard(limfile);
|
| 148 |
|
| 149 |
write_info(__FUNCTION__,"limitfile.root and datacard.txt have been generated. You can now use them to calculate limits!");
|
| 150 |
limfile->Close();
|
| 151 |
|
| 152 |
}
|
| 153 |
|
| 154 |
void Poisson_ratio_plot(int is_data,vector<float> binning, string jzb, TCanvas *can, float high=-9999, float jzbpearerrorMC=-999, float jzbpeakerrorData=-999) {
|
| 155 |
can->SetLogy(0);
|
| 156 |
cout << "About to prepare ratio plots with Poisson errors and obs/pred plots, in the new style!" << endl;
|
| 157 |
TH1F *pred = new TH1F("pred","",binning.size()-1,&binning[0]);
|
| 158 |
TH1F *obs = new TH1F("obs","",binning.size()-1,&binning[0]);
|
| 159 |
for(int ibin=0;ibin<binning.size()-1;ibin++) {
|
| 160 |
// cout << "________________________________________________________________________" << endl;
|
| 161 |
cout << "Ratio plot: Working on bin " << ibin+1 << " of " << binning.size()-1 << " for step " << is_data+1 << " of 3" << endl;
|
| 162 |
float binerr=0;
|
| 163 |
vector<float> contentanderror=get_result_between_two_fixed_jzb_values(true,binning[ibin],binning[ibin+1], jzb,jzb, is_data,jzbpearerrorMC, jzbpeakerrorData, can,false, false,false);
|
| 164 |
//[0] Bpred [1] Bpred uncert [2] Observed [3] Observed error
|
| 165 |
// cout << "Bpred : " << contentanderror[0] << " +/- " << contentanderror[1] << endl;
|
| 166 |
// cout << "Obs : " << contentanderror[2] << " +/- " << contentanderror[3] << endl;
|
| 167 |
pred->SetBinContent(ibin+1,contentanderror[0]);
|
| 168 |
pred->SetBinError(ibin+1,contentanderror[1]);
|
| 169 |
obs->SetBinContent(ibin+1,contentanderror[2]);
|
| 170 |
obs->SetBinError(ibin+1,contentanderror[3]);
|
| 171 |
}
|
| 172 |
TH1F *JZBratioforfitting=(TH1F*)obs->Clone("JZBratioforfitting");
|
| 173 |
JZBratioforfitting->Divide(pred);
|
| 174 |
TGraphAsymmErrors *JZBratio = histRatio(obs,pred,is_data,binning,true);//the last argument basically forces the function to give us the real error, instead of the stat one
|
| 175 |
|
| 176 |
JZBratio->SetTitle("");
|
| 177 |
|
| 178 |
TF1 *pol0 = new TF1("pol0","[0]",0,1000);
|
| 179 |
TF1 *pol0d = new TF1("pol0","[0]",0,1000);
|
| 180 |
JZBratioforfitting->Fit(pol0,"Q0R","",0,30);
|
| 181 |
pol0d->SetParameter(0,pol0->GetParameter(0));
|
| 182 |
|
| 183 |
can->cd();
|
| 184 |
JZBratio->GetYaxis()->SetTitle("Observed/Predicted");
|
| 185 |
JZBratio->GetXaxis()->SetTitle("JZB [GeV]");
|
| 186 |
if ( high>0 ) JZBratio->GetXaxis()->SetRangeUser(0.0,high);
|
| 187 |
JZBratio->GetYaxis()->SetNdivisions(519);
|
| 188 |
JZBratio->GetYaxis()->SetRangeUser(0.0,4.0);
|
| 189 |
JZBratio->GetYaxis()->CenterTitle();
|
| 190 |
JZBratio->GetXaxis()->CenterTitle();
|
| 191 |
JZBratio->SetMarkerSize(DataMarkerSize);
|
| 192 |
JZBratio->Draw("AP");
|
| 193 |
/////----------------------------
|
| 194 |
TPaveText *writeresult = new TPaveText(0.15,0.78,0.49,0.91,"blNDC");
|
| 195 |
writeresult->SetFillStyle(4000);
|
| 196 |
writeresult->SetFillColor(kWhite);
|
| 197 |
writeresult->SetTextFont(42);
|
| 198 |
writeresult->SetTextSize(0.03);
|
| 199 |
writeresult->SetTextAlign(12);
|
| 200 |
ostringstream jzb_agreement_data_text;
|
| 201 |
jzb_agreement_data_text<< setprecision(2) << "mean =" << pol0->GetParameter(0) << " #pm " << setprecision(1) << pol0->GetParError(0);
|
| 202 |
if(is_data==1) fitresultconstdata=pol0->GetParameter(0);// data
|
| 203 |
if(is_data==0) fitresultconstmc=pol0->GetParameter(0); // monte carlo, no signal
|
| 204 |
|
| 205 |
writeresult->AddText(jzb_agreement_data_text.str().c_str());
|
| 206 |
// writeresult->Draw("same");
|
| 207 |
// pol0d->Draw("same");
|
| 208 |
TF1 *topline = new TF1("","1.5",0,1000);
|
| 209 |
TF1 *bottomline = new TF1("","0.5",0,1000);
|
| 210 |
topline->SetLineColor(kBlue);
|
| 211 |
topline->SetLineStyle(2);
|
| 212 |
bottomline->SetLineColor(kBlue);
|
| 213 |
bottomline->SetLineStyle(2);
|
| 214 |
// topline->Draw("same");
|
| 215 |
// bottomline->Draw("same");
|
| 216 |
TF1 *oneline = new TF1("","1.0",0,1000);
|
| 217 |
oneline->SetLineColor(kBlue);
|
| 218 |
oneline->SetLineStyle(1);
|
| 219 |
oneline->Draw("same");
|
| 220 |
if(is_data==1) DrawPrelim();
|
| 221 |
else DrawMCPrelim();
|
| 222 |
TLegend *leg = new TLegend(0.5,0.55,0.94,0.89);
|
| 223 |
leg->SetTextFont(42);
|
| 224 |
leg->SetTextSize(0.04);
|
| 225 |
|
| 226 |
TString MCtitle("MC ");
|
| 227 |
if (is_data==1) MCtitle = "";
|
| 228 |
|
| 229 |
leg->SetFillStyle(4000);
|
| 230 |
leg->SetFillColor(kWhite);
|
| 231 |
leg->SetTextFont(42);
|
| 232 |
leg->AddEntry(JZBratio,MCtitle+"obs / "+MCtitle+"pred","p");
|
| 233 |
leg->AddEntry(oneline,"ratio = 1","l");
|
| 234 |
leg->AddEntry(pol0d,"fit in [0,30] GeV","l");
|
| 235 |
leg->AddEntry(bottomline,"#pm50% envelope","l");
|
| 236 |
// leg->Draw("same");
|
| 237 |
if(is_data==1) CompleteSave(can, "precise_jzb_ratio_data");
|
| 238 |
if(is_data==0) CompleteSave(can, "precise_jzb_ratio_mc");
|
| 239 |
if(is_data==2) CompleteSave(can, "precise_jzb_ratio_mc_BandS");//special case, MC with signal!
|
| 240 |
|
| 241 |
can->SetLogy(1);
|
| 242 |
if(is_data==1) DrawPrelim();
|
| 243 |
else DrawMCPrelim();
|
| 244 |
TLegend *secondleg = make_legend("Observed vs Predicted");
|
| 245 |
secondleg->AddEntry(pred,"Predicted (stat&syst band)","f");
|
| 246 |
secondleg->AddEntry(obs,"Observed","p");
|
| 247 |
secondleg->SetY1(0.7);
|
| 248 |
secondleg->SetX1(0.4);
|
| 249 |
|
| 250 |
pred->GetXaxis()->SetTitle("JZB [GeV]");
|
| 251 |
pred->GetXaxis()->CenterTitle();
|
| 252 |
pred->GetYaxis()->SetTitle("events");
|
| 253 |
pred->GetYaxis()->CenterTitle();
|
| 254 |
pred->SetFillColor(kGreen);
|
| 255 |
pred->SetMarkerColor(kWhite);
|
| 256 |
pred->SetMarkerSize(0.000001);
|
| 257 |
pred->SetLineColor(pred->GetFillColor());
|
| 258 |
pred->Draw("e3");
|
| 259 |
DrawPrelim();
|
| 260 |
obs->Draw("e1x0,same");
|
| 261 |
secondleg->Draw("same");
|
| 262 |
|
| 263 |
if(is_data==1) CompleteSave(can, "Exp/precise_obs_pred_data");
|
| 264 |
if(is_data==0) CompleteSave(can, "Exp/precise_obs_pred_mc");
|
| 265 |
if(is_data==2) CompleteSave(can, "Exp/precise_obs_pred_mc_BandS");//special case, MC with signal!
|
| 266 |
|
| 267 |
delete obs;
|
| 268 |
delete pred;
|
| 269 |
}
|
| 270 |
|
| 271 |
void Poisson_ratio_plots(string mcjzb,string datajzb,vector<float> ratio_binning,float jzbpearerrorMC,float jzbpeakerrorData) {
|
| 272 |
TCanvas *globalc = new TCanvas("globalc","Ratio Plot Canvas");
|
| 273 |
globalc->SetLogy(0);
|
| 274 |
|
| 275 |
Poisson_ratio_plot(data,ratio_binning,datajzb,globalc, jzbHigh,jzbpearerrorMC,jzbpeakerrorData);
|
| 276 |
Poisson_ratio_plot(mc,ratio_binning,mcjzb,globalc, jzbHigh,jzbpearerrorMC,jzbpeakerrorData);
|
| 277 |
Poisson_ratio_plot(mcwithsignal,ratio_binning,mcjzb,globalc, jzbHigh,jzbpearerrorMC,jzbpeakerrorData);
|
| 278 |
}
|
| 279 |
|
| 280 |
|
| 281 |
|
| 282 |
|
| 283 |
|
| 284 |
TF1* do_new_logpar_fit_to_plot(TH1F *osof) {
|
| 285 |
TCanvas *logpar_fit_can = new TCanvas("logpar_fit_can","Fit canvas for LogPar");
|
| 286 |
TF1 *logparfunc = new TF1("logparfunc",LogParabola,60,200,3);
|
| 287 |
TF1 *logparfunc2 = new TF1("logparfunc2",LogParabola,0,(osof->GetXaxis()->GetBinLowEdge(osof->GetNbinsX())+osof->GetXaxis()->GetBinWidth(osof->GetNbinsX())),3);
|
| 288 |
TF1 *logparfuncN = new TF1("logparfuncN",LogParabolaN,0,(osof->GetXaxis()->GetBinLowEdge(osof->GetNbinsX())+osof->GetXaxis()->GetBinWidth(osof->GetNbinsX())),3);
|
| 289 |
TF1 *logparfuncP = new TF1("logparfuncP",LogParabolaP,0,(osof->GetXaxis()->GetBinLowEdge(osof->GetNbinsX())+osof->GetXaxis()->GetBinWidth(osof->GetNbinsX())),3);
|
| 290 |
osof->SetMinimum(0);
|
| 291 |
osof->Fit(logparfunc,"QR");
|
| 292 |
osof->Draw();
|
| 293 |
logparfunc->SetLineWidth(2);
|
| 294 |
logparfunc2->SetParameters(logparfunc->GetParameters());
|
| 295 |
logparfuncN->SetParameters(logparfunc->GetParameters());
|
| 296 |
logparfuncP->SetParameters(logparfunc->GetParameters());
|
| 297 |
stringstream fitinfo;
|
| 298 |
fitinfo << "#Chi^{2} / NDF : " << logparfunc->GetChisquare() << " / " << logparfunc->GetNDF();
|
| 299 |
TText *writefitinfo = write_text(0.8,0.8,fitinfo.str());
|
| 300 |
writefitinfo->SetTextSize(0.03);
|
| 301 |
DrawPrelim();
|
| 302 |
writefitinfo->Draw();
|
| 303 |
logparfunc->Draw("same");
|
| 304 |
logparfunc2->Draw("same");
|
| 305 |
logparfuncN->SetLineStyle(2);
|
| 306 |
logparfuncP->SetLineStyle(2);
|
| 307 |
logparfuncN->Draw("same");
|
| 308 |
logparfuncP->Draw("same");
|
| 309 |
CompleteSave(logpar_fit_can,"Bpred_ttbar_LogPar_Fit_To_OSOF");
|
| 310 |
delete logpar_fit_can;
|
| 311 |
return logparfunc2;
|
| 312 |
}
|
| 313 |
|
| 314 |
vector<TF1*> do_new_extended_fit_to_plot(TH1F *prediction, TH1F *ossf, TH1F *osof) {
|
| 315 |
/* there are mainly two background contributions: Z+Jets (a) and ttbar (b). So:
|
| 316 |
a) The case is clear - we take the OSSF prediction - OSOF prediction, and fit a crystal ball function to it. We then extract the CB parameters.
|
| 317 |
b) For ttbar, we use the OSOF distribution and look at the [10,100] GeV JZB range, and fit our log parabola. We then extract the LP parameters.
|
| 318 |
Once we have these two components, we use the combined parameters to get the final function and we're done.
|
| 319 |
*/
|
| 320 |
//Step 1: take the OSSF prediction - OSOF prediction, and fit a crystal ball function to it
|
| 321 |
TH1F *step1cb = (TH1F*)ossf->Clone("step1cb");
|
| 322 |
step1cb->Add(osof,-1);
|
| 323 |
vector<TF1*> functions = do_cb_fit_to_plot(step1cb,PlottingSetup::JZBPeakWidthData);
|
| 324 |
TF1 *zjetscrystalball = functions[1];
|
| 325 |
|
| 326 |
//Step 2: use the OSOF distribution and look at the [10,100] GeV JZB range, and fit our log parabola
|
| 327 |
TH1F *ttbarprediction=(TH1F*)prediction->Clone("ttbarprediction");
|
| 328 |
ttbarprediction->Add(ossf,-1);//without the Z+Jets estimate, this is really just the ttbar estimate!
|
| 329 |
TF1 *ttbarlogpar = do_new_logpar_fit_to_plot(ttbarprediction);
|
| 330 |
TF1 *ttbarlogparP = new TF1("ttbarlogparP",LogParabolaP,0,(prediction->GetXaxis()->GetBinLowEdge(prediction->GetNbinsX())+prediction->GetXaxis()->GetBinWidth(prediction->GetNbinsX())),8);
|
| 331 |
TF1 *ttbarlogparN = new TF1("ttbarlogparN",LogParabolaN,0,(prediction->GetXaxis()->GetBinLowEdge(prediction->GetNbinsX())+prediction->GetXaxis()->GetBinWidth(prediction->GetNbinsX())),8);
|
| 332 |
|
| 333 |
//and now fuse the two!
|
| 334 |
TF1 *kmlp = new TF1("kmlp", CrystalBallPlusLogParabola, 0,(prediction->GetXaxis()->GetBinLowEdge(prediction->GetNbinsX())+prediction->GetXaxis()->GetBinWidth(prediction->GetNbinsX())),8);
|
| 335 |
TF1 *kmlpP= new TF1("kmlpP",CrystalBallPlusLogParabolaP,0,(prediction->GetXaxis()->GetBinLowEdge(prediction->GetNbinsX())+prediction->GetXaxis()->GetBinWidth(prediction->GetNbinsX())),8);
|
| 336 |
TF1 *kmlpN= new TF1("kmlpN",CrystalBallPlusLogParabolaN,0,(prediction->GetXaxis()->GetBinLowEdge(prediction->GetNbinsX())+prediction->GetXaxis()->GetBinWidth(prediction->GetNbinsX())),8);
|
| 337 |
double kmlp_pars[10];
|
| 338 |
for(int i=0;i<5;i++) kmlp_pars[i]=zjetscrystalball->GetParameter(i);
|
| 339 |
for(int i=0;i<3;i++) kmlp_pars[5+i]=ttbarlogpar->GetParameter(i);
|
| 340 |
ttbarlogparP->SetParameters(ttbarlogpar->GetParameters());
|
| 341 |
ttbarlogparN->SetParameters(ttbarlogpar->GetParameters());
|
| 342 |
kmlp->SetParameters(kmlp_pars);
|
| 343 |
prediction->Fit(kmlp);
|
| 344 |
|
| 345 |
kmlpP->SetParameters(kmlp->GetParameters());
|
| 346 |
kmlpN->SetParameters(kmlp->GetParameters());
|
| 347 |
|
| 348 |
// now that we're done, let's save all of this so we can have a look at it afterwards.
|
| 349 |
TCanvas *can = new TCanvas("can","Prediction Fit Canvas");
|
| 350 |
can->SetLogy(1);
|
| 351 |
prediction->SetMarkerColor(kRed);
|
| 352 |
prediction->Draw();
|
| 353 |
|
| 354 |
kmlp->SetLineColor(TColor::GetColor("#04B404"));
|
| 355 |
kmlpP->SetLineColor(TColor::GetColor("#04B404"));
|
| 356 |
kmlpN->SetLineColor(TColor::GetColor("#04B404"));
|
| 357 |
kmlp->Draw("same");
|
| 358 |
kmlpN->SetLineStyle(2);
|
| 359 |
kmlpP->SetLineStyle(2);
|
| 360 |
kmlpN->Draw("same");
|
| 361 |
kmlpP->Draw("same");
|
| 362 |
|
| 363 |
ttbarlogpar->SetLineColor(TColor::GetColor("#CC2EFA"));//purple
|
| 364 |
ttbarlogpar->Draw("same");
|
| 365 |
ttbarlogparP->SetLineColor(TColor::GetColor("#CC2EFA"));//purple
|
| 366 |
ttbarlogparN->SetLineColor(TColor::GetColor("#CC2EFA"));//purple
|
| 367 |
ttbarlogparP->SetLineStyle(2);
|
| 368 |
ttbarlogparN->SetLineStyle(2);
|
| 369 |
ttbarlogparP->Draw("same");
|
| 370 |
ttbarlogparN->Draw("same");
|
| 371 |
|
| 372 |
functions[0]->Draw("same");functions[1]->Draw("same");functions[2]->Draw("same");
|
| 373 |
|
| 374 |
TLegend *analyticalBpredLEG = make_legend("",0.5,0.55);
|
| 375 |
analyticalBpredLEG->AddEntry(prediction,"predicted","p");
|
| 376 |
analyticalBpredLEG->AddEntry(functions[1],"Crystal Ball fit","l");
|
| 377 |
analyticalBpredLEG->AddEntry(functions[0],"1#sigma Crystal Ball fit","l");
|
| 378 |
analyticalBpredLEG->AddEntry(ttbarlogparN,"TTbar fit","l");
|
| 379 |
analyticalBpredLEG->AddEntry(ttbarlogpar,"1#sigma TTbar fit","l");
|
| 380 |
analyticalBpredLEG->AddEntry(kmlp,"Combined function","l");
|
| 381 |
analyticalBpredLEG->AddEntry(kmlpN,"1#sigma combined function","l");
|
| 382 |
analyticalBpredLEG->Draw("same");
|
| 383 |
|
| 384 |
CompleteSave(can,"Bpred_MC_Analytical_Function_Composition");
|
| 385 |
delete can;
|
| 386 |
|
| 387 |
//and finally: prep return functions
|
| 388 |
vector<TF1*> return_functions;
|
| 389 |
return_functions.push_back(kmlpN);
|
| 390 |
return_functions.push_back(kmlp);
|
| 391 |
return_functions.push_back(kmlpP);
|
| 392 |
|
| 393 |
return_functions.push_back(ttbarlogparN);
|
| 394 |
return_functions.push_back(ttbarlogpar);
|
| 395 |
return_functions.push_back(ttbarlogparP);
|
| 396 |
|
| 397 |
return_functions.push_back(functions[0]);
|
| 398 |
return_functions.push_back(functions[1]);
|
| 399 |
return_functions.push_back(functions[2]);
|
| 400 |
|
| 401 |
return return_functions;
|
| 402 |
}
|
| 403 |
void do_new_prediction_plot(string jzb, TCanvas *globalcanvas, float sigma, float high, bool is_data, bool overlay_signal = false )
|
| 404 |
{
|
| 405 |
int nbins=100;
|
| 406 |
if(is_data) nbins=50;
|
| 407 |
float low=0;
|
| 408 |
float hi=500;
|
| 409 |
|
| 410 |
|
| 411 |
TH1F *RcorrJZBeemm = allsamples.Draw("RcorrJZBeemm",jzb.c_str(),nbins,low,hi, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,is_data, luminosity,allsamples.FindSample("TTJet"));
|
| 412 |
TH1F *LcorrJZBeemm = allsamples.Draw("LcorrJZBeemm",("-"+jzb).c_str(),nbins,low,hi, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,is_data, luminosity,allsamples.FindSample("TTJet"));
|
| 413 |
TH1F *RcorrJZBem = allsamples.Draw("RcorrJZBem",jzb.c_str(),nbins,low,hi, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,is_data, luminosity,allsamples.FindSample("TTJet"));
|
| 414 |
TH1F *LcorrJZBem = allsamples.Draw("LcorrJZBem",("-"+jzb).c_str(),nbins,low,hi, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,is_data, luminosity,allsamples.FindSample("TTJet"));
|
| 415 |
|
| 416 |
TH1F *RcorrJZBSBem;
|
| 417 |
TH1F *LcorrJZBSBem;
|
| 418 |
TH1F *RcorrJZBSBeemm;
|
| 419 |
TH1F *LcorrJZBSBeemm;
|
| 420 |
|
| 421 |
TH1F *lm4RcorrJZBeemm = allsamples.Draw("lm4RcorrJZBeemm",jzb.c_str(),nbins,low,hi, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,is_data, luminosity,allsamples.FindSample("LM4"));
|
| 422 |
|
| 423 |
if(RestrictToMassPeak||!PlottingSetup::UseSidebandsForcJZB) {
|
| 424 |
RcorrJZBSBem = allsamples.Draw("RcorrJZBSBem",jzb.c_str(),nbins,low,hi, "JZB [GeV]", "events", sidebandcut&&cutOSOF&&cutnJets,is_data, luminosity,allsamples.FindSample("TTJet"));
|
| 425 |
LcorrJZBSBem = allsamples.Draw("LcorrJZBSBem",("-"+jzb).c_str(),nbins,low,hi, "JZB [GeV]", "events", sidebandcut&&cutOSOF&&cutnJets,is_data, luminosity,allsamples.FindSample("TTJet"));
|
| 426 |
RcorrJZBSBeemm = allsamples.Draw("RcorrJZBSBeemm",jzb.c_str(),nbins,low,hi, "JZB [GeV]", "events", sidebandcut&&cutOSSF&&cutnJets,is_data, luminosity,allsamples.FindSample("TTJet"));
|
| 427 |
LcorrJZBSBeemm = allsamples.Draw("LcorrJZBSBeemm",("-"+jzb).c_str(),nbins,low,hi, "JZB [GeV]", "events", sidebandcut&&cutOSSF&&cutnJets,is_data, luminosity,allsamples.FindSample("TTJet"));
|
| 428 |
}
|
| 429 |
|
| 430 |
TH1F *conta = allsamples.Draw("conta",jzb.c_str(),2*nbins,-hi,hi, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,is_data, luminosity,allsamples.FindSample("TTJet"));
|
| 431 |
|
| 432 |
TFile *ttbarfile=new TFile("ttbarfile.root","RECREATE");
|
| 433 |
conta->Write();
|
| 434 |
ttbarfile->Close();
|
| 435 |
|
| 436 |
TH1F *Bpred = (TH1F*)LcorrJZBeemm->Clone("Bpred");
|
| 437 |
if(RestrictToMassPeak) {
|
| 438 |
Bpred->Add(RcorrJZBem,1.0/3.);
|
| 439 |
Bpred->Add(LcorrJZBem,-1.0/3.);
|
| 440 |
Bpred->Add(RcorrJZBSBem,1.0/3.);
|
| 441 |
Bpred->Add(LcorrJZBSBem,-1.0/3.);
|
| 442 |
Bpred->Add(RcorrJZBSBeemm,1.0/3.);
|
| 443 |
Bpred->Add(LcorrJZBSBeemm,-1.0/3.);
|
| 444 |
} else {
|
| 445 |
Bpred->Add(RcorrJZBem,1.0);
|
| 446 |
Bpred->Add(LcorrJZBem,-1.0);
|
| 447 |
}
|
| 448 |
|
| 449 |
globalcanvas->cd();
|
| 450 |
globalcanvas->SetLogy(1);
|
| 451 |
|
| 452 |
RcorrJZBeemm->SetMarkerStyle(20);
|
| 453 |
RcorrJZBeemm->GetXaxis()->SetRangeUser(0,high);
|
| 454 |
RcorrJZBeemm->Draw("e1x0");
|
| 455 |
|
| 456 |
Bpred->SetLineColor(kRed);
|
| 457 |
Bpred->SetStats(0);
|
| 458 |
Bpred->Draw("hist,same");
|
| 459 |
|
| 460 |
if ( overlay_signal ) {
|
| 461 |
lm4RcorrJZBeemm->SetLineColor(TColor::GetColor("#088A08"));
|
| 462 |
lm4RcorrJZBeemm->Draw("histo,same");
|
| 463 |
}
|
| 464 |
RcorrJZBeemm->SetMarkerSize(DataMarkerSize);
|
| 465 |
|
| 466 |
TLegend *legBpred = make_legend("",0.6,0.55);
|
| 467 |
if(is_data)
|
| 468 |
{
|
| 469 |
vector<TF1*> analytical_function = do_new_extended_fit_to_plot(Bpred,LcorrJZBeemm,LcorrJZBem);
|
| 470 |
globalcanvas->cd();
|
| 471 |
analytical_function[0]->Draw("same"); analytical_function[1]->Draw("same");analytical_function[2]->Draw("same");
|
| 472 |
legBpred->AddEntry(RcorrJZBeemm,"observed","p");
|
| 473 |
legBpred->AddEntry(Bpred,"predicted","l");
|
| 474 |
legBpred->AddEntry(analytical_function[1],"predicted fit","l");
|
| 475 |
legBpred->AddEntry(analytical_function[2],"stat. uncert.","l");
|
| 476 |
if ( overlay_signal ) legBpred->AddEntry(lm4RcorrJZBeemm,"LM4","l");
|
| 477 |
legBpred->Draw();
|
| 478 |
CompleteSave(globalcanvas,"Bpred_Data");
|
| 479 |
}
|
| 480 |
else {
|
| 481 |
vector<TF1*> analytical_function = do_new_extended_fit_to_plot(Bpred,LcorrJZBeemm,LcorrJZBem);
|
| 482 |
globalcanvas->cd();
|
| 483 |
analytical_function[0]->Draw("same"); analytical_function[1]->Draw("same");analytical_function[2]->Draw("same");
|
| 484 |
legBpred->AddEntry(RcorrJZBeemm,"MC observed","p");
|
| 485 |
legBpred->AddEntry(Bpred,"MC predicted","l");
|
| 486 |
legBpred->AddEntry(analytical_function[1],"predicted fit","l");
|
| 487 |
legBpred->AddEntry(analytical_function[2],"stat. uncert.","l");
|
| 488 |
if ( overlay_signal ) legBpred->AddEntry(lm4RcorrJZBeemm,"LM4","l");
|
| 489 |
legBpred->Draw();
|
| 490 |
CompleteSave(globalcanvas,"Bpred_MC");
|
| 491 |
}
|
| 492 |
|
| 493 |
TH1F *Bpredem = (TH1F*)LcorrJZBeemm->Clone("Bpredem");
|
| 494 |
Bpredem->Add(RcorrJZBem);
|
| 495 |
Bpredem->Add(LcorrJZBem,-1);
|
| 496 |
TH1F *BpredSBem = (TH1F*)LcorrJZBeemm->Clone("BpredSBem");
|
| 497 |
BpredSBem->Add(RcorrJZBSBem);
|
| 498 |
Bpred->Add(LcorrJZBSBem,-1);
|
| 499 |
TH1F *BpredSBeemm = (TH1F*)LcorrJZBeemm->Clone("BpredSBeemm");
|
| 500 |
BpredSBeemm->Add(RcorrJZBSBeemm);
|
| 501 |
BpredSBeemm->Add(LcorrJZBSBeemm,-1.0);
|
| 502 |
globalcanvas->cd();
|
| 503 |
globalcanvas->SetLogy(1);
|
| 504 |
|
| 505 |
RcorrJZBeemm->SetMarkerStyle(20);
|
| 506 |
RcorrJZBeemm->GetXaxis()->SetRangeUser(0,high);
|
| 507 |
RcorrJZBeemm->Draw("e1x0");
|
| 508 |
RcorrJZBeemm->SetMarkerSize(DataMarkerSize);
|
| 509 |
|
| 510 |
Bpredem->SetLineColor(kRed+1);
|
| 511 |
Bpredem->SetStats(0);
|
| 512 |
Bpredem->Draw("hist,same");
|
| 513 |
|
| 514 |
BpredSBem->SetLineColor(kGreen+2);//TColor::GetColor("#0B6138"));
|
| 515 |
BpredSBem->SetLineStyle(2);
|
| 516 |
BpredSBem->Draw("hist,same");
|
| 517 |
|
| 518 |
BpredSBeemm->SetLineColor(kBlue+1);
|
| 519 |
BpredSBeemm->SetLineStyle(3);
|
| 520 |
BpredSBeemm->Draw("hist,same");
|
| 521 |
|
| 522 |
TLegend *legBpredc = make_legend("",0.6,0.55);
|
| 523 |
if(is_data)
|
| 524 |
{
|
| 525 |
legBpredc->AddEntry(RcorrJZBeemm,"observed","p");
|
| 526 |
legBpredc->AddEntry(Bpredem,"OFZP","l");
|
| 527 |
legBpredc->AddEntry(BpredSBem,"OFSB","l");
|
| 528 |
legBpredc->AddEntry(BpredSBeemm,"SFSB","l");
|
| 529 |
legBpredc->Draw();
|
| 530 |
CompleteSave(globalcanvas,"Bpred_Data_comparison");
|
| 531 |
}
|
| 532 |
else {
|
| 533 |
legBpredc->AddEntry(RcorrJZBeemm,"MC observed","p");
|
| 534 |
legBpredc->AddEntry(Bpredem,"MC OFZP","l");
|
| 535 |
legBpredc->AddEntry(BpredSBem,"MC OFSB","l");
|
| 536 |
legBpredc->AddEntry(BpredSBeemm,"MC SFSB","l");
|
| 537 |
legBpredc->Draw();
|
| 538 |
legBpredc->Draw();
|
| 539 |
CompleteSave(globalcanvas,"Bpred_MC_comparison_ttbar");
|
| 540 |
}
|
| 541 |
delete RcorrJZBeemm;
|
| 542 |
delete LcorrJZBeemm;
|
| 543 |
delete RcorrJZBem;
|
| 544 |
delete LcorrJZBem;
|
| 545 |
delete RcorrJZBSBem;
|
| 546 |
delete LcorrJZBSBem;
|
| 547 |
delete RcorrJZBSBeemm;
|
| 548 |
delete LcorrJZBSBeemm;
|
| 549 |
delete lm4RcorrJZBeemm;
|
| 550 |
}
|
| 551 |
|
| 552 |
void do_new_prediction_plots(string mcjzb, string datajzb, float DataSigma, float MCSigma, bool overlay_signal ) {
|
| 553 |
TCanvas *globalcanvas = new TCanvas("globalcanvas","Prediction Canvas");
|
| 554 |
// do_new_prediction_plot(datajzb,globalcanvas,DataSigma,jzbHigh ,data,overlay_signal);
|
| 555 |
do_new_prediction_plot(mcjzb,globalcanvas,MCSigma,jzbHigh ,mc,overlay_signal);
|
| 556 |
}
|
| 557 |
|
