30 |
|
dout << "Info : The lepton requirement is currently set to " << (const char*) leptoncut << endl; |
31 |
|
} |
32 |
|
|
33 |
< |
void find_peaks(float &MCPeak,float &MCPeakError, float &DataPeak, float &DataPeakError, float &MCSigma, float &DataSigma, stringstream &result, bool doPUreweighting = true ) |
33 |
> |
|
34 |
> |
|
35 |
> |
void find_one_peak_combination(TCut specialcut, float &MCPeak,float &MCPeakError, float &DataPeak, float &DataPeakError, float &MCSigma, float &DataSigma, stringstream &result, bool doPUreweighting = true, string saveas="") |
36 |
|
{ |
37 |
|
// Temporarily switch off PU reweighting, if asked |
38 |
|
TCut weightbackup=cutWeight; |
39 |
|
if ( !doPUreweighting ) cutWeight="1.0"; |
40 |
< |
|
40 |
> |
|
41 |
> |
int nbins=100; |
42 |
> |
if(PlottingSetup::DoBTag) nbins=25; |
43 |
> |
|
44 |
|
TCanvas *tempcan = new TCanvas("tempcan","Temporary canvas for peak finding preparations"); |
45 |
< |
TH1F *rawJZBeemmMC = allsamples.Draw("rawJZBeemmMC",jzbvariablemc,100,-50,50, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,mc, luminosity); |
46 |
< |
TH1F *rawJZBeemmData = allsamples.Draw("rawJZBeemmData",jzbvariabledata,100, -50,50, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,data, luminosity); |
47 |
< |
TH1F *rawJZBemMC = allsamples.Draw("rawJZBemMC",jzbvariablemc,100,-50,50, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,mc, luminosity); |
48 |
< |
TH1F *rawJZBemData = allsamples.Draw("rawJZBemData",jzbvariabledata,100, -50,50, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets,data, luminosity); |
45 |
> |
TH1F *rawJZBeemmMC = allsamples.Draw("rawJZBeemmMC",jzbvariablemc,nbins,-50,50, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets&&specialcut,mc, luminosity); |
46 |
> |
TH1F *rawJZBeemmData = allsamples.Draw("rawJZBeemmData",jzbvariabledata,nbins, -50,50, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets&&specialcut,data, luminosity); |
47 |
> |
TH1F *rawJZBemMC = allsamples.Draw("rawJZBemMC",jzbvariablemc,nbins,-50,50, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets&&specialcut,mc, luminosity); |
48 |
> |
TH1F *rawJZBemData = allsamples.Draw("rawJZBemData",jzbvariabledata,nbins, -50,50, "JZB [GeV]", "events", cutmass&&cutOSOF&&cutnJets&&specialcut,data, luminosity); |
49 |
|
TH1F *rawttbarjzbeemmMC; |
50 |
|
|
51 |
|
if(method==doKM) { |
52 |
|
//we only need this histo for the KM fitting... |
53 |
< |
rawttbarjzbeemmMC = allsamples.Draw("rawttbarjzbeemmMC",jzbvariablemc,100, -50,50, "JZB [GeV]", "events",cutmass&&cutOSSF&&cutnJets,mc,luminosity,allsamples.FindSample("TTJet")); |
54 |
< |
MCPeak=find_peak(rawJZBeemmMC, rawttbarjzbeemmMC, -40, 40, mc, MCPeakError,MCSigma,method); |
55 |
< |
DataPeak=find_peak(rawJZBeemmData, rawJZBeemmData, -40, 40, data, DataPeakError,DataSigma,method); |
53 |
> |
rawttbarjzbeemmMC = allsamples.Draw("rawttbarjzbeemmMC",jzbvariablemc,nbins, -50,50, "JZB [GeV]", "events",cutmass&&cutOSSF&&cutnJets&&specialcut,mc,luminosity,allsamples.FindSample("TTJet")); |
54 |
> |
MCPeak=find_peak(rawJZBeemmMC, rawttbarjzbeemmMC, -40, 40, mc, MCPeakError,MCSigma,method,saveas); |
55 |
> |
DataPeak=find_peak(rawJZBeemmData, rawJZBeemmData, -40, 40, data, DataPeakError,DataSigma,method,saveas); |
56 |
> |
delete rawttbarjzbeemmMC; |
57 |
|
} |
58 |
|
else { |
59 |
|
TH1F *reducedMC = (TH1F*)rawJZBeemmMC->Clone(); |
61 |
|
reducedMC->Add(rawJZBemMC,-1); |
62 |
|
reducedData->Add(rawJZBemData,-1); |
63 |
|
//this is Kostas' way of doing it - we subtract em to get rid of some of the ttbar contribution (in reality, of flavor-symmetric contribution) |
64 |
< |
MCPeak=find_peak(reducedMC, rawttbarjzbeemmMC, -40, 40, mc, MCPeakError,MCSigma,method); |
65 |
< |
DataPeak=find_peak(reducedData, rawJZBeemmData, -40, 40, data, DataPeakError,DataSigma,method); |
66 |
< |
|
64 |
> |
MCPeak=find_peak(reducedMC, rawttbarjzbeemmMC, -40, 40, mc, MCPeakError,MCSigma,method,saveas); |
65 |
> |
DataPeak=find_peak(reducedData, rawJZBeemmData, -40, 40, data, DataPeakError,DataSigma,method,saveas); |
66 |
> |
delete reducedMC; |
67 |
> |
delete reducedData; |
68 |
|
} |
69 |
|
|
70 |
|
// Revert to original PU reweighting |
72 |
|
|
73 |
|
// MCPeak=find_peak(rawJZBeemmMC, rawttbarjzbeemmMC, -40, 40, mc, MCPeakError,MCSigma,method); |
74 |
|
// DataPeak=find_peak(rawJZBeemmData, rawJZBeemmData, -40, 40, data, DataPeakError,DataSigma,method); |
75 |
< |
dout << "We have found the peak in Data at " << DataPeak << " +/- " << DataPeakError << " with sigma=" << DataSigma << " +/- ?? (not impl.)" << endl; |
76 |
< |
result << "We have found the peak in Data at " << DataPeak << " +/- " << DataPeakError << " with sigma=" << DataSigma << " +/- ?? (not impl.)" << endl; |
77 |
< |
dout << "We have found the peak in MC at " << MCPeak << " +/- " << MCPeakError << " with sigma=" << MCSigma << " +/- ?? (not impl.)" << endl; |
78 |
< |
result << "We have found the peak in MC at " << MCPeak << " +/- " << MCPeakError << " with sigma=" << MCSigma << " +/- ?? (not impl.)" << endl; |
75 |
> |
dout << " We have found the peak in Data at " << DataPeak << " +/- " << DataPeakError << " with sigma=" << DataSigma << " +/- ?? (not impl.)" << endl; |
76 |
> |
result << " We have found the peak in Data at " << DataPeak << " +/- " << DataPeakError << " with sigma=" << DataSigma << " +/- ?? (not impl.)" << endl; |
77 |
> |
dout << " We have found the peak in MC at " << MCPeak << " +/- " << MCPeakError << " with sigma=" << MCSigma << " +/- ?? (not impl.)" << endl; |
78 |
> |
result << " We have found the peak in MC at " << MCPeak << " +/- " << MCPeakError << " with sigma=" << MCSigma << " +/- ?? (not impl.)" << endl; |
79 |
> |
delete rawJZBeemmData; |
80 |
> |
delete rawJZBeemmMC; |
81 |
> |
delete rawJZBemData; |
82 |
> |
delete rawJZBemMC; |
83 |
|
delete tempcan; |
84 |
|
} |
85 |
|
|
86 |
+ |
void find_peaks(float &MCPeak,float &MCPeakError, float &DataPeak, float &DataPeakError, stringstream &result, bool doPUreweighting, stringstream &datajzb, stringstream &mcjzb) |
87 |
+ |
{ |
88 |
+ |
float mcpeak, datapeak; |
89 |
+ |
float mcpeakerr, datapeakerr; |
90 |
+ |
|
91 |
+ |
float mceepeak,mcmmpeak; |
92 |
+ |
float mceepeakerr,mcmmpeakerr; |
93 |
+ |
|
94 |
+ |
float datammpeak,dataeepeak; |
95 |
+ |
float datammpeakerr,dataeepeakerr; |
96 |
+ |
|
97 |
+ |
float mcSigma, dataSigma; |
98 |
+ |
|
99 |
+ |
dout << "Finding global peak : " << endl; |
100 |
+ |
find_one_peak_combination(TCut(""),mcpeak,mcpeakerr, datapeak,datapeakerr,mcSigma, dataSigma,result,doPUreweighting,""); |
101 |
+ |
dout << "Finding peak for electrons : " << endl; |
102 |
+ |
find_one_peak_combination(TCut("id1==0"),mceepeak,mceepeakerr, dataeepeak,dataeepeakerr,mcSigma, dataSigma,result,doPUreweighting,"_ele"); |
103 |
+ |
dout << "Finding peak for muons : " << endl; |
104 |
+ |
find_one_peak_combination(TCut("id1==1"),mcmmpeak,mcmmpeakerr, datammpeak,datammpeakerr,mcSigma, dataSigma,result,doPUreweighting,"_mu"); |
105 |
+ |
|
106 |
+ |
datajzb << "(" << jzbvariabledata; |
107 |
+ |
mcjzb << "(" << jzbvariablemc; |
108 |
+ |
|
109 |
+ |
if(dataeepeak>0) datajzb << "- (id1==id2)*(id1==0)*" << TMath::Abs(dataeepeak) << " "; |
110 |
+ |
else datajzb << "+ (id1==id2)*(id1==0)*" << TMath::Abs(dataeepeak) << " "; |
111 |
+ |
|
112 |
+ |
if(datammpeak>0) datajzb << "- (id1==id2)*(id1==1)*" << TMath::Abs(datammpeak) << " "; |
113 |
+ |
else datajzb << "+ (id1==id2)*(id1==1)*" << TMath::Abs(datammpeak) << " "; |
114 |
+ |
|
115 |
+ |
if(datapeak>0) datajzb << "- (id1!=id2)*" << TMath::Abs(datapeak) << " "; |
116 |
+ |
else datajzb << "+ (id1!=id2)*" << TMath::Abs(datapeak) << " "; |
117 |
+ |
|
118 |
+ |
datajzb << ")"; |
119 |
+ |
|
120 |
+ |
if(mceepeak>0) mcjzb << "- (id1==id2)*(id1==0)*" << TMath::Abs(mceepeak) << " "; |
121 |
+ |
else mcjzb << "+ (id1==id2)*(id1==0)*" << TMath::Abs(mceepeak) << " "; |
122 |
+ |
|
123 |
+ |
if(mcmmpeak>0) mcjzb << "- (id1==id2)*(id1==1)*" << TMath::Abs(mcmmpeak) << " "; |
124 |
+ |
else mcjzb << "+ (id1==id2)*(id1==1)*" << TMath::Abs(mcmmpeak) << " "; |
125 |
+ |
|
126 |
+ |
if(mcpeak>0) mcjzb << "- (id1!=id2)*" << TMath::Abs(mcpeak) << " "; |
127 |
+ |
else mcjzb << "+ (id1!=id2)*" << TMath::Abs(mcpeak) << " "; |
128 |
+ |
|
129 |
+ |
mcjzb << ")"; |
130 |
+ |
} |
131 |
+ |
|
132 |
|
void make_special_obs_pred_mll_plot(string mcjzb, float jzbthreshold) { |
133 |
|
float min=70.0; |
134 |
|
float max=115.0; |
159 |
|
TH1F *RcorrJZBSBeemm; |
160 |
|
TH1F *LcorrJZBSBeemm; |
161 |
|
|
162 |
< |
TH1F *RcorrJZBeemmNoS; |
162 |
> |
// TH1F *RcorrJZBeemmNoS; |
163 |
|
|
164 |
|
if(PlottingSetup::RestrictToMassPeak) { |
165 |
|
RcorrJZBSBem = allsamples.Draw("RcorrJZBSBem", "mll",nbins,min,max, "JZB [GeV]", "events", sidebandcut&&cutOSOF&&cutnJets&&largerzero,data, luminosity); |
246 |
|
TIter next(bgleg->GetListOfPrimitives()); |
247 |
|
TObject* obj; |
248 |
|
// Copy the nice bkgd legend skipping the "data" |
249 |
< |
while ( obj = next() ) |
249 |
> |
while ( (obj = next()) ) |
250 |
|
if ( strcmp(((TLegendEntry*)obj)->GetObject()->GetName(),"datahistoOSSF") ) |
251 |
|
kinleg->GetListOfPrimitives()->Add(obj); |
252 |
|
|
265 |
|
TIter next(hmc.GetHists()); |
266 |
|
TObject* obj; |
267 |
|
TH1* hratio = NULL; |
268 |
< |
while ( obj = next() ) { |
268 |
> |
while ( (obj = next()) ) { |
269 |
|
if ( !hratio ) { |
270 |
|
hratio = (TH1*)obj->Clone(); |
271 |
|
hratio->SetName("hratio"); |
326 |
|
TCut cut; |
327 |
|
cut=cutmass&&cutOSSF&&cutnJets&&ibasiccut; |
328 |
|
if(filename=="nJets") cut=cutmass&&cutOSSF&&ibasiccut; |
329 |
+ |
if(filename=="nJets_osof") cut=cutmass&&cutOSOF&&ibasiccut; |
330 |
|
if(filename=="nJets_nocuts_except_mll_ossf") cut=cutmass&&cutOSSF; |
331 |
|
if(filename=="mll") { |
332 |
|
cut=cutOSSF&&cutnJets&&ibasiccut; |
354 |
|
THStack mcstack = allsamples.DrawStack("mcstack",variable,nbins,min,max, xlabel, "events",cut,mc,luminosity); |
355 |
|
if(variable=="pfJetGoodPt[0]") datahisto->SetMaximum(10*datahisto->GetMaximum()); |
356 |
|
if(variable=="pt") datahisto->SetMaximum(10*datahisto->GetMaximum()); |
357 |
< |
if(filename=="mll_inclusive") datahisto->SetMinimum(1); |
357 |
> |
if(filename=="mll_inclusive"||filename=="mll_inclusive_mm"||filename=="mll_inclusive_ee") datahisto->SetMinimum(1); |
358 |
|
if(filename=="mll_osof") datahisto->SetMaximum(10*datahisto->GetMaximum()); |
359 |
|
if(filename=="mll_osof") datahisto->SetMinimum(9); |
302 |
– |
|
360 |
|
datahisto->SetMaximum(5.3*datahisto->GetMaximum()); |
361 |
|
datahisto->Draw("e1"); |
362 |
|
ckin->Update(); |
487 |
|
bool dolog=true; |
488 |
|
bool nolog=false; |
489 |
|
if(doPF) write_warning(__FUNCTION__,"Please use caution when trying to produce PF plots; not all versions of the JZB trees have these variables!"); |
490 |
< |
float mll_low=40; |
490 |
> |
float mll_low=50; |
491 |
|
float mll_hi=160; |
492 |
|
if(!PlottingSetup::RestrictToMassPeak) { |
493 |
< |
mll_low=10; |
493 |
> |
mll_low=20; |
494 |
|
mll_hi=210; |
495 |
|
} |
496 |
< |
/* |
496 |
> |
|
497 |
> |
make_kin_plot("met[4]","",40,0,200,dolog,"MET [GeV]","met",doPF,true); |
498 |
|
make_kin_plot("mll","",(int)((mll_hi-mll_low)/5),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll",doPF,true); |
499 |
|
make_kin_plot("mll","",(int)((mll_hi-mll_low)/5),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_osof",doPF,true,true); |
500 |
|
make_kin_plot("mll","",(int)((mll_hi-mll_low)/5),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_ee",doPF,true); |
502 |
|
make_kin_plot("mll","",(int)((mll_hi-mll_low)),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_inclusive",doPF,true); |
503 |
|
make_kin_plot("mll","",(int)((mll_hi-mll_low)),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_inclusive_ee",doPF,true); |
504 |
|
make_kin_plot("mll","",(int)((mll_hi-mll_low)),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_inclusive_mm",doPF,true); |
505 |
< |
make_kin_plot("mll","",(int)((mll_hi-mll_low)),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_inclusive_osof",doPF,true); |
506 |
< |
make_kin_plot("mll","",(int)((350-mll_low)),mll_low,350,dolog,"m_{ll} [GeV]","mll_inclusive_highrange",doPF);/* |
505 |
> |
make_kin_plot("mll","",(int)((mll_hi-mll_low))/2,mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_inclusive_osof",doPF,true); |
506 |
> |
make_kin_plot("mll","",(int)((350-mll_low)),mll_low,350,dolog,"m_{ll} [GeV]","mll_inclusive_highrange",doPF); |
507 |
|
make_kin_plot("numVtx","",(int)(30.5-(-0.5)),-0.5,30.5,nolog,"N(Vtx)","numVtx",doPF); |
508 |
< |
make_kin_plot("jetpt[0]","",40,0,200,dolog,"leading jet p_{T} [GeV]","pfJetGoodPt_0",doPF); |
509 |
< |
make_kin_plot("jeteta[0]","",40,-5,5,nolog,"leading jet #eta","pfJetGoodEta_0",doPF); |
508 |
> |
// make_kin_plot("jetpt[0]","",40,0,200,dolog,"leading jet p_{T} [GeV]","pfJetGoodPt_0",doPF); |
509 |
> |
// make_kin_plot("jeteta[0]","",40,-5,5,nolog,"leading jet #eta","pfJetGoodEta_0",doPF); |
510 |
|
make_kin_plot("pt","",50,0,400,dolog,"Z p_{T} [GeV]","Zpt",doPF); |
511 |
|
make_kin_plot("pt1","",50,0,100,nolog,"p_{T} [GeV]","pt1",doPF); |
512 |
|
make_kin_plot("pt2","",50,0,100,nolog,"p_{T} [GeV]","pt2",doPF); |
513 |
|
make_kin_plot("eta1","",40,-5,5,nolog,"#eta_{l}","eta",doPF); |
514 |
|
make_kin_plot("jzb[1]","",100,-150,150,dolog,"JZB [GeV]","jzb_ossf",doPF); |
515 |
< |
make_kin_plot("pfJetGoodNum","",8,0.5,8.5,dolog,"nJets","nJets",doPF); |
516 |
< |
make_kin_plot("pfJetGoodNum","",8,0.5,8.5,dolog,"nJets","nJets_nocuts_except_mll_ossf",doPF);*//* |
515 |
> |
make_kin_plot("pfJetGoodNum40","",9,-0.5,8.5,dolog,"nJets","nJets",doPF); |
516 |
> |
make_kin_plot("pfJetGoodNum40","",9,-0.5,8.5,dolog,"nJets","nJets_osof",doPF); |
517 |
> |
make_kin_plot("pfJetGoodNum40","",9,-0.5,8.5,dolog,"nJets","nJets_nocuts_except_mll_ossf",doPF); |
518 |
|
if(!doPF) make_special_mll_plot((int)((mll_hi-mll_low)/5),mll_low,mll_hi,dolog,"m_{ll} [GeV]"); |
519 |
|
stringstream jzbcut; |
520 |
|
jzbcut << "((is_data&&("<<datajzb<<")>100)||(!is_data&&("<<mcjzb<<")>100))"; |
528 |
|
make_kin_plot("mll",jzbcut2.str(),(int)((mll_hi-mll_low)/5),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_osof_aboveJZB150",doPF,true); |
529 |
|
stringstream jzbcut3; |
530 |
|
jzbcut3 << "((is_data&&("<<datajzb<<")>50)||(!is_data&&("<<mcjzb<<")>50))"; |
531 |
< |
make_kin_plot("mll",jzbcut3.str(),(int)((mll_hi-mll_low)/5),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_aboveJZB50",doPF,true); make_kin_plot("mll",jzbcut3.str(),(int)((mll_hi-mll_low)/5),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_osof_aboveJZB50",doPF,true,true); |
532 |
< |
*/ |
531 |
> |
make_kin_plot("mll",jzbcut3.str(),(int)((mll_hi-mll_low)/5),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_aboveJZB50",doPF,true); |
532 |
> |
make_kin_plot("mll",jzbcut3.str(),(int)((mll_hi-mll_low)/5),mll_low,mll_hi,dolog,"m_{ll} [GeV]","mll_osof_aboveJZB50",doPF,true,true); |
533 |
> |
|
534 |
> |
make_special_obs_pred_mll_plot(mcjzb,0); |
535 |
|
make_special_obs_pred_mll_plot(mcjzb,50); |
536 |
|
make_special_obs_pred_mll_plot(mcjzb,100); |
537 |
|
make_special_obs_pred_mll_plot(mcjzb,150); |
608 |
|
|
609 |
|
|
610 |
|
flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak ---- the arguments changed |
611 |
< |
for(int ijzb=0;ijzb<jzb_cuts.size();ijzb++) { |
611 |
> |
for(int ijzb=0;ijzb<(int)jzb_cuts.size();ijzb++) { |
612 |
|
float jzbcut=jzb_cuts[ijzb]; // Comparison plots are done for this JZB cut |
613 |
|
float mll_low=50;float mll_high=170; |
614 |
|
if(!PlottingSetup::RestrictToMassPeak) { |
915 |
|
return return_functions; |
916 |
|
} |
917 |
|
|
918 |
< |
void do_prediction_plot(string jzb, TCanvas *globalcanvas, float sigma, float high, int use_data, bool overlay_signal = false,string subdir="" ) |
918 |
> |
void do_prediction_plot(string jzb, TCanvas *globalcanvas, float high, int use_data, bool overlay_signal = false,string subdir="" ) |
919 |
|
{ |
859 |
– |
// note: sigma is not used ATM |
920 |
|
switch_overunderflow(true); |
921 |
|
bool is_data=false; |
922 |
|
bool use_signal=false; |
991 |
|
|
992 |
|
TH1F *Bpred = (TH1F*)LcorrJZBeemm->Clone("Bpred"); |
993 |
|
TH1F *JBpred = (TH1F*)JLcorrJZBeemm->Clone("Bpred"); |
994 |
+ |
|
995 |
+ |
TH1F *BpredSys = new TH1F("Bpredsys","Bpredsys",PlottingSetup::global_ratio_binning.size()-1,&PlottingSetup::global_ratio_binning[0]); |
996 |
+ |
ClearHisto(BpredSys); |
997 |
+ |
|
998 |
|
if(PlottingSetup::RestrictToMassPeak) { |
999 |
|
Bpred->Add(RcorrJZBem,1.0/3.); |
1000 |
|
Bpred->Add(LcorrJZBem,-1.0/3.); |
1017 |
|
JBpred->Add(JLcorrJZBSBem,-1.0/3.); |
1018 |
|
JBpred->Add(JRcorrJZBSBeemm,1.0/3.); |
1019 |
|
JBpred->Add(JLcorrJZBSBeemm,-1.0/3.); |
1020 |
+ |
|
1021 |
+ |
//Systematics: |
1022 |
+ |
AddSquared(BpredSys,JLcorrJZBeemm,zjetsestimateuncertONPEAK*zjetsestimateuncertONPEAK); |
1023 |
+ |
AddSquared(BpredSys,JRcorrJZBem,emuncertONPEAK*emuncertONPEAK*(1.0/9)); |
1024 |
+ |
AddSquared(BpredSys,JLcorrJZBem,emuncertONPEAK*emuncertONPEAK*(1.0/9)); |
1025 |
+ |
AddSquared(BpredSys,JRcorrJZBSBem,emsidebanduncertONPEAK*emsidebanduncertONPEAK*(1.0/9)); |
1026 |
+ |
AddSquared(BpredSys,JLcorrJZBSBem,emsidebanduncertONPEAK*emsidebanduncertONPEAK*(1.0/9)); |
1027 |
+ |
AddSquared(BpredSys,JRcorrJZBSBeemm,eemmsidebanduncertONPEAK*eemmsidebanduncertONPEAK*(1.0/9)); |
1028 |
+ |
AddSquared(BpredSys,JLcorrJZBSBeemm,eemmsidebanduncertONPEAK*eemmsidebanduncertONPEAK*(1.0/9)); |
1029 |
|
} else { |
1030 |
|
Bpred->Add(RcorrJZBem,1.0); |
1031 |
|
Bpred->Add(LcorrJZBem,-1.0); |
1035 |
|
//these are for the ratio |
1036 |
|
JBpred->Add(JRcorrJZBem,1.0); |
1037 |
|
JBpred->Add(JLcorrJZBem,-1.0); |
1038 |
+ |
|
1039 |
+ |
//Systematics |
1040 |
+ |
AddSquared(BpredSys,JLcorrJZBeemm,zjetsestimateuncertOFFPEAK*zjetsestimateuncertOFFPEAK); |
1041 |
+ |
AddSquared(BpredSys,JRcorrJZBem,emuncertOFFPEAK*emuncertOFFPEAK); |
1042 |
+ |
AddSquared(BpredSys,JLcorrJZBem,emuncertOFFPEAK*emuncertOFFPEAK); |
1043 |
+ |
|
1044 |
|
} |
1045 |
|
|
1046 |
< |
|
1046 |
> |
SQRT(BpredSys); |
1047 |
> |
BpredSys->Divide(JBpred); |
1048 |
> |
|
1049 |
|
|
1050 |
|
flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak ---- prediction changed |
1051 |
|
TH1F *Tpred = (TH1F*)RcorrJZBem->Clone("Bpred"); |
1102 |
|
lm4RcorrJZBeemm->Draw("hist,same"); |
1103 |
|
legBpred->AddEntry(RcorrJZBeemm,"observed","p"); |
1104 |
|
legBpred->AddEntry(Bpred,"predicted","l"); |
1105 |
< |
legBpred->AddEntry(analytical_function[1],"predicted fit","l"); |
1106 |
< |
legBpred->AddEntry(analytical_function[2],"stat. uncert.","l"); |
1107 |
< |
if ( overlay_signal ) legBpred->AddEntry(lm4RcorrJZBeemm,"LM4","l"); |
1105 |
> |
// legBpred->AddEntry(analytical_function[1],"predicted fit","l"); |
1106 |
> |
// legBpred->AddEntry(analytical_function[2],"stat. uncert.","l"); |
1107 |
> |
legBpred->AddEntry(lm4RcorrJZBeemm,"LM4","l"); |
1108 |
|
legBpred->Draw(); |
1109 |
|
DrawPrelim(); |
1110 |
|
|
1147 |
|
// speciallegBpred->AddEntry(lm4RcorrJZBeemmC,"LM4","l"); |
1148 |
|
speciallegBpred->AddEntry(lm4RcorrJZBeemm,"LM4","l"); |
1149 |
|
speciallegBpred->Draw(); |
1150 |
< |
save_with_ratio(JRcorrJZBeemm,JBpred,predcomppad,subdir+"Bpred_Data_____PredictionComposition",true,use_data!=1,"data/pred"); |
1150 |
> |
save_with_ratio(JRcorrJZBeemm,JBpred,predcomppad,subdir+"Bpred_Data_____PredictionComposition",true,true,"data/pred",BpredSys); |
1151 |
|
|
1152 |
|
TCanvas *specialcanv = new TCanvas("specialcanv","specialcanv"); |
1153 |
|
THStack kostack = allsamples.DrawStack("RcorrJZBeemm",jzb.c_str(),nbins,low,hi, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,!is_data, luminosity,use_signal); |
1173 |
|
Bpred->Draw("hist,same"); |
1174 |
|
RcorrJZBeemm->Draw("e1x0,same");//HAVE IT ON TOP! |
1175 |
|
legBpred->AddEntry(RcorrJZBeemm,"MC true","p"); |
1176 |
< |
legBpred->AddEntry(Bpred,"MC predicted","l"); |
1176 |
> |
legBpred->AddEntry(Bpred,"MC predicted","l"); |
1177 |
|
if(versok) legBpred->AddEntry((TObject*)0,"",""); // Just for alignment // causes seg fault on root v5.18 |
1178 |
|
if(versok) legBpred->AddEntry((TObject*)0,"",""); // causes seg fault on root v5.18 |
1179 |
|
if ( overlay_signal ) legBpred->AddEntry(lm4RcorrJZBeemm,"LM4","l"); |
1225 |
|
string ytitle("ratio"); |
1226 |
|
if ( use_data==1 ) ytitle = "data/pred"; |
1227 |
|
//save_with_ratio(JRcorrJZBeemm,JBpred,kinpad,Bpredsaveas,true,use_data!=1,ytitle); |
1228 |
< |
save_with_ratio(JRcorrJZBeemm,JBpred,kinpad,subdir+Bpredsaveas,true,false,ytitle);//not extending the y range anymore up to 4 |
1228 |
> |
save_with_ratio(JRcorrJZBeemm,JBpred,kinpad,subdir+Bpredsaveas,true,true,ytitle,BpredSys);//not extending the y range anymore up to 4 |
1229 |
|
|
1230 |
|
|
1231 |
|
flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak |
1232 |
< |
/// The part below is meaningless for the offpeak analysis (it's a comparison of the different estimates but there is but one estimate!) |
1232 |
> |
// The part below is meaningless for the offpeak analysis (it's a comparison of the different estimates but there is but one estimate!) |
1233 |
|
if(PlottingSetup::RestrictToMassPeak) { |
1234 |
|
TH1F *Bpredem = (TH1F*)LcorrJZBeemm->Clone("Bpredem"); |
1235 |
|
Bpredem->Add(RcorrJZBem); |
1320 |
|
|
1321 |
|
void do_prediction_plots(string mcjzb, string datajzb, float DataSigma, float MCSigma, bool overlay_signal ) { |
1322 |
|
TCanvas *globalcanvas = new TCanvas("globalcanvas","Prediction Canvas"); |
1323 |
< |
do_prediction_plot(datajzb,globalcanvas,DataSigma,jzbHigh ,data,overlay_signal); |
1323 |
> |
do_prediction_plot(datajzb,globalcanvas,jzbHigh ,data,overlay_signal); |
1324 |
|
if ( !PlottingSetup::Approved ) { |
1325 |
< |
do_prediction_plot(mcjzb,globalcanvas,MCSigma,jzbHigh ,mc,overlay_signal); |
1326 |
< |
do_prediction_plot(mcjzb,globalcanvas,MCSigma,jzbHigh ,mcwithsignal,overlay_signal); |
1325 |
> |
do_prediction_plot(mcjzb,globalcanvas,jzbHigh ,mc,overlay_signal); |
1326 |
> |
do_prediction_plot(mcjzb,globalcanvas,jzbHigh ,mcwithsignal,overlay_signal); |
1327 |
|
} else { |
1328 |
|
write_info(__FUNCTION__,"You set approved to true, therefore not producing prediction/observation plots for MC with and without signal."); |
1329 |
|
} |
1746 |
|
|
1747 |
|
|
1748 |
|
dout << "Going to increase the cross section for diboson samples ... " << endl; |
1749 |
< |
for(int i=0;i<SamplesToBeModified.size();i++) { |
1749 |
> |
for(int i=0;i<(int)SamplesToBeModified.size();i++) { |
1750 |
|
float origxs=(allsamples.collection)[SamplesToBeModified[i]].xs; |
1751 |
|
(allsamples.collection)[SamplesToBeModified[i]].xs=origxs*stretchfactor; |
1752 |
|
dout << " Increased xs for sample " << (allsamples.collection)[SamplesToBeModified[i]].filename << " from " << origxs << " to " << (allsamples.collection)[SamplesToBeModified[i]].xs << " (by a factor of " << stretchfactor << ")" << endl; |
1758 |
|
dout << "Going to produce JZB plots" << endl; |
1759 |
|
produce_stretched_jzb_plots(mcjzb,datajzb,ratio_binning); |
1760 |
|
TCanvas *gloca = new TCanvas("gloca","gloca"); |
1680 |
– |
float sigma=123456; |
1761 |
|
|
1762 |
|
dout << "Going to produce prediction plots" << endl; |
1763 |
< |
do_prediction_plot(mcjzb, gloca, sigma, PlottingSetup::jzbHigh, 0, false,"Dibosons/Bpred/" ); // do only MC plots, no signal |
1764 |
< |
do_prediction_plot(mcjzb, gloca, sigma, PlottingSetup::jzbHigh, 0, false,"Dibosons/Bpred/" ); // do MC plots with signal |
1763 |
> |
do_prediction_plot(mcjzb, gloca, PlottingSetup::jzbHigh, 0, false,"Dibosons/Bpred/" ); // do only MC plots, no signal |
1764 |
> |
do_prediction_plot(mcjzb, gloca, PlottingSetup::jzbHigh, 0, false,"Dibosons/Bpred/" ); // do MC plots with signal |
1765 |
|
delete gloca; |
1766 |
|
|
1767 |
|
dout << "Going to reset the cross section for diboson samples ... " << endl; |
1768 |
< |
for(int i=0;i<SamplesToBeModified.size();i++) { |
1768 |
> |
for(int i=0;i<(int)SamplesToBeModified.size();i++) { |
1769 |
|
float Upxs=(allsamples.collection)[SamplesToBeModified[i]].xs; |
1770 |
|
(allsamples.collection)[SamplesToBeModified[i]].xs=(allsamples.collection)[SamplesToBeModified[i]].xs*(1.0/stretchfactor); |
1771 |
|
string Upname=(allsamples.collection)[SamplesToBeModified[i]].samplename; |
1862 |
|
string posneg[] = {"pos","neg"}; |
1863 |
|
TCut tkRegions[] = {cutOSSF&&cutnJets&&cutmass,cutOSOF&&cutnJets&&cutmass,cutOSSF&&cutnJets&&sidebandcut,cutOSOF&&cutnJets&&sidebandcut}; |
1864 |
|
|
1865 |
< |
for(int ijzb=0;ijzb<jzb_cuts.size();ijzb++) { |
1865 |
> |
for(int ijzb=0;ijzb<(int)jzb_cuts.size();ijzb++) { |
1866 |
|
TCut jzbMC[] = { give_jzb_expression(mcjzb,jzb_cuts[ijzb],"pos"), give_jzb_expression(mcjzb,jzb_cuts[ijzb],"neg") }; |
1867 |
|
dout << "_________________________________________________________" << endl; |
1868 |
|
dout << "Table for JZB> " << jzb_cuts[ijzb] << endl; |
1869 |
< |
for(int isample=0;isample<(allsamples.collection).size();isample++) { |
1869 |
> |
for(int isample=0;isample<(int)(allsamples.collection).size();isample++) { |
1870 |
|
if(!(allsamples.collection)[isample].is_data) dout << (allsamples.collection)[isample].samplename << " & "; |
1871 |
|
else dout << "Sample & "; |
1872 |
|
for(int iregion=0;iregion<nRegions;iregion++) { |
2125 |
|
|
2126 |
|
void draw_ttbar_and_zjets_shape(string mcjzb, string datajzb) { |
2127 |
|
int all_leptons=-1; |
2048 |
– |
int electrons_only=0; |
2049 |
– |
int mu_only=1; |
2050 |
– |
int twojetswith50gev=1; |
2128 |
|
int threejetswith30gev=0; |
2129 |
|
/* |
2130 |
+ |
int twojetswith50gev=1; |
2131 |
+ |
int electrons_only=0; |
2132 |
+ |
int mu_only=1; |
2133 |
+ |
|
2134 |
|
draw_ttbar_and_zjets_shape_for_one_configuration(mcjzb,datajzb,all_leptons,twojetswith50gev); |
2135 |
|
draw_ttbar_and_zjets_shape_for_one_configuration(mcjzb,datajzb,all_leptons,threejetswith30gev); |
2136 |
|
|
2144 |
|
draw_ttbar_and_zjets_shape_for_one_configuration(mcjzb,datajzb,all_leptons,threejetswith30gev,true); |
2145 |
|
} |
2146 |
|
|
2147 |
< |
void find_correction_factors(string &jzbvardata,string &jzbvarmc) { |
2067 |
< |
//first: colorful plots |
2147 |
> |
float find_one_correction_factor(string FindKeyword, TCut SpecialCut, string SaveAs) { |
2148 |
|
TCanvas *cancorr = new TCanvas("cancorr","Canvas for Response Correction"); |
2149 |
|
cancorr->SetLogz(); |
2150 |
|
cancorr->SetRightMargin(0.13); |
2151 |
< |
flag_this_change(__FUNCTION__,__LINE__,true);//PlottingSetup::RestrictToMassPeak |
2152 |
< |
TCut zptforresponsepresentation("pt<600"&&cutmass&&cutOSSF&&"((sumJetPt[1]/pt)<5.0)"); |
2151 |
> |
TCut zptforresponsepresentation("pt<600"&&cutmass&&cutOSSF&&"((sumJetPt[1]/pt)<5.0)"&&SpecialCut); |
2152 |
> |
if(PlottingSetup::DoBTag) zptforresponsepresentation=zptforresponsepresentation&&PlottingSetup::bTagRequirement; |
2153 |
|
TH2F *niceresponseplotd = new TH2F("niceresponseplotd","",100,0,600,100,0,5); |
2154 |
< |
(allsamples.collection)[allsamples.FindSample("Data")[0]].events->Draw("sumJetPt[1]/pt:pt>>niceresponseplotd",zptforresponsepresentation); |
2154 |
> |
vector<int> SampleIndices=allsamples.FindSample(FindKeyword); |
2155 |
> |
for(int iSample=0;iSample<SampleIndices.size();iSample++) { |
2156 |
> |
dout << " Response correction : Using sample " << (allsamples.collection)[SampleIndices[iSample]].filename << endl; |
2157 |
> |
(allsamples.collection)[SampleIndices[iSample]].events->Draw("sumJetPt[1]/pt:pt>>+niceresponseplotd",zptforresponsepresentation*cutWeight); |
2158 |
> |
} |
2159 |
> |
|
2160 |
|
niceresponseplotd->SetStats(0); |
2161 |
|
niceresponseplotd->GetXaxis()->SetTitle("Z p_{T} [GeV]"); |
2162 |
|
niceresponseplotd->GetYaxis()->SetTitle("Response"); |
2165 |
|
niceresponseplotd->Draw("COLZ"); |
2166 |
|
TProfile * profd = (TProfile*)niceresponseplotd->ProfileX(); |
2167 |
|
profd->SetMarkerSize(DataMarkerSize); |
2168 |
< |
profd->Fit("pol0","","same,e1",30,400); |
2168 |
> |
profd->Fit("pol0","","same,e1",100,400); |
2169 |
|
DrawPrelim(); |
2170 |
< |
TText* title = write_text(0.5,0.7,"Data"); |
2170 |
> |
string stitle="Data"; |
2171 |
> |
if(!Contains(FindKeyword,"Data")) stitle="MC"; |
2172 |
> |
TText* title = write_text(0.5,0.7,stitle.c_str()); |
2173 |
|
title->SetTextAlign(12); |
2174 |
|
title->Draw(); |
2175 |
|
TF1 *datapol=(TF1*)profd->GetFunction("pol0"); |
2176 |
< |
float datacorrection=datapol->GetParameter(0); |
2177 |
< |
stringstream dataresstring; |
2178 |
< |
dataresstring<<"Response: "<<std::setprecision(2)<<100*datacorrection<<" %"; |
2179 |
< |
TText* restitle = write_text(0.5,0.65,dataresstring.str()); |
2176 |
> |
float correction=datapol->GetParameter(0); |
2177 |
> |
stringstream resstring; |
2178 |
> |
resstring<<"Response: "<<std::setprecision(2)<<100*correction<<" %"; |
2179 |
> |
TText* restitle = write_text(0.5,0.65,resstring.str()); |
2180 |
|
restitle->SetTextAlign(12); |
2181 |
|
restitle->SetTextSize(0.03); |
2182 |
|
restitle->Draw(); |
2183 |
< |
CompleteSave(cancorr,"ResponseCorrection/Response_Correction_Illustration_Data"); |
2183 |
> |
CompleteSave(cancorr,"ResponseCorrection/Response_Correction_Illustration_New_"+SaveAs); |
2184 |
> |
delete cancorr; |
2185 |
> |
delete niceresponseplotd; |
2186 |
> |
return correction; |
2187 |
> |
} |
2188 |
> |
|
2189 |
> |
void find_correction_factors(string &jzbvardata,string &jzbvarmc) { |
2190 |
|
|
2191 |
< |
TH2F *niceresponseplotm = new TH2F("niceresponseplotm","",100,0,600,100,0,5); |
2192 |
< |
(allsamples.collection)[allsamples.FindSample("DY")[0]].events->Draw("sumJetPt[1]/pt:pt>>niceresponseplotm",zptforresponsepresentation); |
2193 |
< |
niceresponseplotm->SetStats(0); |
2194 |
< |
niceresponseplotm->GetXaxis()->SetTitle("Z p_{T} [GeV]"); |
2195 |
< |
niceresponseplotm->GetYaxis()->SetTitle("Response"); |
2196 |
< |
niceresponseplotm->GetXaxis()->CenterTitle(); |
2197 |
< |
niceresponseplotm->GetYaxis()->CenterTitle(); |
2198 |
< |
niceresponseplotm->Draw("COLZ"); |
2199 |
< |
(allsamples.collection)[allsamples.FindSample("DY")[0]].events->Draw("sumJetPt[1]/pt:pt",zptforresponsepresentation,"PROF,same"); |
2107 |
< |
TProfile * profm = (TProfile*)niceresponseplotm->ProfileX(); |
2108 |
< |
profm->SetMarkerSize(DataMarkerSize); |
2109 |
< |
profm->Fit("pol0","","same,e1",30,400); |
2110 |
< |
DrawMCPrelim(); |
2111 |
< |
title = write_text(0.5,0.7,"MC simulation"); |
2112 |
< |
title->SetTextAlign(12); |
2113 |
< |
title->Draw(); |
2114 |
< |
TF1 *mcpol=(TF1*)profm->GetFunction("pol0"); |
2115 |
< |
float mccorrection=mcpol->GetParameter(0); |
2116 |
< |
stringstream mcresstring; |
2117 |
< |
mcresstring<<"Response: "<<std::setprecision(2)<<100*mccorrection<<" %"; |
2118 |
< |
TText* mcrestitle = write_text(0.5,0.65,mcresstring.str()); |
2119 |
< |
mcrestitle->SetTextAlign(12); |
2120 |
< |
mcrestitle->SetTextSize(0.03); |
2121 |
< |
mcrestitle->Draw(); |
2122 |
< |
CompleteSave(cancorr,"ResponseCorrection/Response_Correction_Illustration_MC"); |
2191 |
> |
//Step 1 : Get results |
2192 |
> |
float datacorrection=find_one_correction_factor("Data","","Data"); |
2193 |
> |
float mccorrection=find_one_correction_factor("DY","","MC"); |
2194 |
> |
|
2195 |
> |
float dataEEcorrection=find_one_correction_factor("Data","id1==0","Data_ee"); |
2196 |
> |
float mcEEcorrection=find_one_correction_factor("DY","id1==0","MC_ee"); |
2197 |
> |
|
2198 |
> |
float dataMMcorrection=find_one_correction_factor("Data","id1==1","Data_mm"); |
2199 |
> |
float mcMMcorrection=find_one_correction_factor("DY","id1==1","MC_mm"); |
2200 |
|
|
2201 |
+ |
cout << "Corrections : " << endl; |
2202 |
+ |
cout << " Data : " << datacorrection << endl; |
2203 |
+ |
cout << " ee (" << dataEEcorrection << ") , mm (" << dataMMcorrection << ")" << endl; |
2204 |
+ |
cout << " MC : " << mccorrection << endl; |
2205 |
+ |
cout << " ee (" << mcEEcorrection << ") , mm (" << mcMMcorrection << ")" << endl; |
2206 |
|
|
2207 |
< |
//Step 2: Getting the result |
2208 |
< |
// TCut zptcutforresponse("pt>30&&pt<300&&TMath::Abs(91.2-mll)<20&&id1==id2&&(ch1*ch2<0)"); |
2207 |
> |
cout << " result from new one for data: " << 0.937097 << endl; |
2208 |
> |
cout << " result from old one for data: " << 0.937097 << endl; |
2209 |
> |
cout << " result from new one for MC: " << 0.96848 << endl; |
2210 |
> |
cout << " result from old one for MC: " << 0.963502 << endl; |
2211 |
> |
|
2212 |
> |
|
2213 |
> |
//Step 2: Processing the result and making it into something useful :-) |
2214 |
|
stringstream jzbvardatas; |
2215 |
< |
if(datacorrection>1) jzbvardatas<<"(jzb[1]-"<<datacorrection-1<<"*pt)"; |
2216 |
< |
if(datacorrection<1) jzbvardatas<<"(jzb[1]+"<<1-datacorrection<<"*pt)"; |
2215 |
> |
jzbvardatas << "("; |
2216 |
> |
|
2217 |
> |
if(dataEEcorrection>=1) jzbvardatas<<"((id1==0&&id1==id2)*(jzb[1]-" << dataEEcorrection-1 << "*pt))"; |
2218 |
> |
if(dataEEcorrection<1) jzbvardatas<<"((id1==0&&id1==id2)*(jzb[1]+" << 1-dataEEcorrection << "*pt))"; |
2219 |
> |
|
2220 |
> |
if(dataMMcorrection>=1) jzbvardatas<<"+((id1==1&&id1==id2)*(jzb[1]-" << dataMMcorrection-1 << "*pt))"; |
2221 |
> |
if(dataMMcorrection<1) jzbvardatas<<"+((id1==1&&id1==id2)*(jzb[1]+" << 1-dataMMcorrection << "*pt))"; |
2222 |
> |
|
2223 |
> |
float averagecorrection=(dataMMcorrection+dataEEcorrection)/2.0; |
2224 |
> |
|
2225 |
> |
if(averagecorrection>1) jzbvardatas<<"+((id1!=id2)*(jzb[1]-" << averagecorrection-1 << "*pt))"; |
2226 |
> |
if(averagecorrection<1) jzbvardatas<<"+((id1!=id2)*(jzb[1]-" << 1-averagecorrection << "*pt))"; |
2227 |
> |
|
2228 |
> |
jzbvardatas << ")"; |
2229 |
|
jzbvardata=jzbvardatas.str(); |
2230 |
+ |
|
2231 |
|
stringstream jzbvarmcs; |
2232 |
< |
if(mccorrection>1) jzbvarmcs<<"(jzb[1]-"<<mccorrection-1<<"*pt)"; |
2233 |
< |
if(mccorrection<1) jzbvarmcs<<"(jzb[1]+"<<1-mccorrection<<"*pt)"; |
2232 |
> |
jzbvarmcs << "("; |
2233 |
> |
|
2234 |
> |
if(mcEEcorrection>=1) jzbvarmcs<<"((id1==0&&id1==id2)*(jzb[1]-" << mcEEcorrection-1 << "*pt))"; |
2235 |
> |
if(mcEEcorrection<1) jzbvarmcs<<"((id1==0&&id1==id2)*(jzb[1]+" << 1-mcEEcorrection << "*pt))"; |
2236 |
> |
|
2237 |
> |
if(mcMMcorrection>=1) jzbvarmcs<<"+((id1==1&&id1==id2)*(jzb[1]-" << mcMMcorrection-1 << "*pt))"; |
2238 |
> |
if(mcMMcorrection<1) jzbvarmcs<<"+((id1==1&&id1==id2)*(jzb[1]+" << 1-mcMMcorrection << "*pt))"; |
2239 |
> |
|
2240 |
> |
float averagemccorrection=(mcMMcorrection+mcEEcorrection)/2.0; |
2241 |
> |
|
2242 |
> |
if(averagemccorrection>1) jzbvarmcs<<"+((id1!=id2)*(jzb[1]-" << averagemccorrection-1 << "*pt))"; |
2243 |
> |
if(averagemccorrection<1) jzbvarmcs<<"+((id1!=id2)*(jzb[1]-" << 1-averagemccorrection << "*pt))"; |
2244 |
> |
|
2245 |
> |
jzbvarmcs << ")"; |
2246 |
|
jzbvarmc=jzbvarmcs.str(); |
2247 |
+ |
|
2248 |
|
dout << "JZB Z pt correction summary : " << endl; |
2249 |
|
dout << " Data: The response is " << datacorrection << " --> jzb variable is now : " << jzbvardata << endl; |
2250 |
|
dout << " MC : The response is " << mccorrection << " --> jzb variable is now : " << jzbvarmc << endl; |
2251 |
+ |
|
2252 |
|
} |
2253 |
|
|
2254 |
|
void pick_up_events(string cut) { |
2256 |
|
allsamples.PickUpEvents(cut); |
2257 |
|
} |
2258 |
|
|
2259 |
< |
void save_template(string mcjzb, string datajzb,vector<float> jzb_cuts,float MCPeakError,float DataPeakError) { |
2259 |
> |
void save_template(string mcjzb, string datajzb,vector<float> jzb_cuts,float MCPeakError,float DataPeakError, vector<float> jzb_shape_limit_bins) { |
2260 |
|
dout << "Saving configuration template!" << endl; |
2261 |
|
ofstream configfile; |
2262 |
|
configfile.open("../DistributedModelCalculations/last_configuration.C"); |
2276 |
|
configfile<<"string datajzb=\"datajzb_ERROR\";\n"; |
2277 |
|
configfile<<"string mcjzb=\"mcjzb_ERROR\";\n"; |
2278 |
|
configfile<<"vector<float>jzb_cuts;\n"; |
2279 |
+ |
configfile<<"vector<float>jzb_shape_limit_bins;\n"; |
2280 |
|
configfile<<"float MCPeakError=-999;\n"; |
2281 |
|
configfile<<"float DataPeakError=-999;\n"; |
2282 |
|
configfile<<"}\n\n"; |
2286 |
|
configfile<<"mcjzb=\""<<mcjzb<<"\";\n\n"; |
2287 |
|
configfile<<"\n\nMCPeakError="<<MCPeakError<<";\n"; |
2288 |
|
configfile<<"DataPeakError="<<DataPeakError<<";\n\n"; |
2289 |
< |
for(int i=0;i<jzb_cuts.size();i++) configfile<<"jzb_cuts.push_back("<<jzb_cuts[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2289 |
> |
for(int i=0;i<(int)jzb_cuts.size();i++) configfile<<"jzb_cuts.push_back("<<jzb_cuts[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2290 |
> |
configfile<<"\n\n"; |
2291 |
> |
for(int i=0;i<(int)Nobs.size();i++) configfile<<"Nobs.push_back("<<Nobs[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2292 |
> |
for(int i=0;i<(int)Npred.size();i++) configfile<<"Npred.push_back("<<Npred[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2293 |
> |
for(int i=0;i<(int)Nprederr.size();i++) configfile<<"Nprederr.push_back("<<Nprederr[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2294 |
|
configfile<<"\n\n"; |
2295 |
< |
for(int i=0;i<Nobs.size();i++) configfile<<"Nobs.push_back("<<Nobs[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2296 |
< |
for(int i=0;i<Npred.size();i++) configfile<<"Npred.push_back("<<Npred[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2297 |
< |
for(int i=0;i<Nprederr.size();i++) configfile<<"Nprederr.push_back("<<Nprederr[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2295 |
> |
for(int i=0;i<(int)flippedNobs.size();i++) configfile<<"flippedNobs.push_back("<<flippedNobs[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2296 |
> |
for(int i=0;i<(int)flippedNpred.size();i++) configfile<<"flippedNpred.push_back("<<flippedNpred[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2297 |
> |
for(int i=0;i<(int)flippedNprederr.size();i++) configfile<<"flippedNprederr.push_back("<<flippedNprederr[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2298 |
> |
for(int i=0;i<(int)jzb_shape_limit_bins.size();i++) configfile<<"jzb_shape_limit_bins.push_back("<<jzb_shape_limit_bins[i]<<"); // JZB shape bin boundary at " << jzb_shape_limit_bins[i] << "\n"; |
2299 |
|
configfile<<"\n\n"; |
2180 |
– |
for(int i=0;i<flippedNobs.size();i++) configfile<<"flippedNobs.push_back("<<flippedNobs[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2181 |
– |
for(int i=0;i<flippedNpred.size();i++) configfile<<"flippedNpred.push_back("<<flippedNpred[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2182 |
– |
for(int i=0;i<flippedNprederr.size();i++) configfile<<"flippedNprederr.push_back("<<flippedNprederr[i]<<"); // JZB cut at " << jzb_cuts[i] << "\n"; |
2300 |
|
configfile<<"\n\n"; |
2301 |
|
configfile<<"luminosity="<<luminosity<<";\n"; |
2302 |
|
configfile<<"RestrictToMassPeak="<<RestrictToMassPeak<<";//defines the type of analysis we're running\n"; |
2322 |
|
float max=histos[0]->GetMaximum(); |
2323 |
|
if(manualmin>=0) min=manualmin; |
2324 |
|
else { |
2325 |
< |
for(int i=1;i<histos.size();i++) { |
2325 |
> |
for(int i=1;i<(int)histos.size();i++) { |
2326 |
|
float curmin=get_nonzero_minimum(histos[i]); |
2327 |
|
float curmax=histos[i]->GetMaximum(); |
2328 |
|
if(curmin<min) min=curmin; |
2333 |
|
histos[0]->Draw(drawoption.c_str()); |
2334 |
|
stringstream drawopt; |
2335 |
|
drawopt << drawoption << ",same"; |
2336 |
< |
for(int i=1;i<histos.size();i++) { |
2336 |
> |
for(int i=1;i<(int)histos.size();i++) { |
2337 |
|
histos[i]->Draw(drawopt.str().c_str()); |
2338 |
|
} |
2339 |
|
} |
2342 |
|
//in the case of the on peak analysis, we compare the 3 control regions to the real value |
2343 |
|
//in the case of the OFF peak analysis, we compare our control region to the real value |
2344 |
|
TCut weightbackup=cutWeight; |
2345 |
< |
cutWeight="1.0"; |
2345 |
> |
// cutWeight="1.0"; |
2346 |
|
float simulatedlumi = luminosity; //in pb please - adjust to your likings |
2347 |
|
|
2348 |
|
|
2400 |
|
leg->Draw("same"); |
2401 |
|
DrawMCPrelim(simulatedlumi); |
2402 |
|
CompleteSave(tcan,"Systematics/"+prestring+"/ttbar_shape_comparison"); |
2286 |
– |
cout << __LINE__ << endl; |
2403 |
|
TH1F *TZemcopy = (TH1F*)TZem->Clone("TZemcopy"); |
2288 |
– |
cout << __LINE__ << endl; |
2404 |
|
TH1F *TZeemmcopy = (TH1F*)TZeemm->Clone("TZeemmcopy"); |
2290 |
– |
cout << __LINE__ << endl; |
2405 |
|
TH1F *TSeemmcopy; |
2406 |
|
TH1F *TSemcopy; |
2407 |
|
if(PlottingSetup::RestrictToMassPeak) { |
2592 |
|
} |
2593 |
|
|
2594 |
|
|
2595 |
< |
void zjets_prediction_comparison() { |
2595 |
> |
void zjets_prediction_comparison(string mcjzbWithPU) { |
2596 |
> |
TCanvas *zcan = new TCanvas("zcan","zcan"); |
2597 |
> |
zcan->SetLogy(1); |
2598 |
> |
TCut weightbackup=cutWeight; |
2599 |
> |
|
2600 |
> |
/* |
2601 |
|
// Do it without PU re-weighting |
2602 |
|
float MCPeakNoPU=0,MCPeakErrorNoPU=0,DataPeakNoPU=0,DataPeakErrorNoPU=0,MCSigma=0,DataSigma=0; |
2603 |
|
stringstream resultsNoPU; |
2604 |
|
|
2605 |
|
stringstream mcjzbnoPU; |
2606 |
< |
find_peaks(MCPeakNoPU,MCPeakErrorNoPU, DataPeakNoPU, DataPeakErrorNoPU,MCSigma,DataSigma,resultsNoPU,false); |
2606 |
> |
find_peaks(MCPeakNoPU,MCPeakErrorNoPU, DataPeakNoPU, DataPeakErrorNoPU,MCSigma,DataSigma,resultsNoPU,true); |
2607 |
|
if(MCPeakNoPU>0) mcjzbnoPU<<"("<<jzbvariablemc<<"-"<<TMath::Abs(MCPeakNoPU)<<")"; |
2608 |
|
else mcjzbnoPU<<"("<<jzbvariablemc<<"+"<<TMath::Abs(MCPeakNoPU)<<")"; |
2609 |
|
|
2610 |
|
string mcjzb = mcjzbnoPU.str(); |
2611 |
|
dout << "The peak corrected JZB expression for MC without pileup is : " << mcjzb << endl; |
2612 |
|
|
2494 |
– |
TCut weightbackup=cutWeight; |
2613 |
|
cutWeight="1.0"; |
2614 |
< |
float sbg_min=0.; |
2615 |
< |
float sbg_max=100.; |
2616 |
< |
int sbg_nbins=5; |
2614 |
> |
*/ |
2615 |
> |
string mcjzb = mcjzbWithPU; // this is with PURW, if you want without it you have to uncomment the part above (and comment out this line) |
2616 |
> |
|
2617 |
> |
vector<float> binning; |
2618 |
> |
binning.push_back(0); |
2619 |
> |
binning.push_back(20); |
2620 |
> |
binning.push_back(40); |
2621 |
> |
binning.push_back(60); |
2622 |
> |
binning.push_back(80); |
2623 |
> |
binning.push_back(100); |
2624 |
> |
// float sbg_min=0.; |
2625 |
> |
// float sbg_max=100.; |
2626 |
> |
// int sbg_nbins=5; |
2627 |
|
float simulatedlumi = luminosity;//in pb please - adjust to your likings |
2628 |
|
|
2629 |
|
TCut kPos((mcjzb+">0").c_str()); |
2630 |
|
TCut kNeg((mcjzb+"<0").c_str()); |
2631 |
|
string var( "abs("+mcjzb+")" ); |
2632 |
|
|
2633 |
< |
TCut kcut(cutmass&&cutOSSF&&"pfJetGoodNum>2"); |
2634 |
< |
TH1F *hJZBpos = systsamples.Draw("hJZBpos",var,sbg_nbins,sbg_min,sbg_max, "JZB [GeV]", "events", |
2635 |
< |
kcut&&kPos,mc,simulatedlumi,systsamples.FindSample("/DY")); |
2508 |
< |
TH1F *hJZBneg = systsamples.Draw("hJZBneg",var,sbg_nbins,sbg_min,sbg_max, "JZB [GeV]", "events", |
2509 |
< |
kcut&&kNeg,mc,simulatedlumi,systsamples.FindSample("/DY")); |
2633 |
> |
TCut kcut(cutmass&&cutOSSF&&cutnJets); |
2634 |
> |
TH1F *hJZBpos = systsamples.Draw("hJZBpos",var,binning, "JZB [GeV]", "events",kcut&&kPos,mc,simulatedlumi,systsamples.FindSample("/DY")); |
2635 |
> |
TH1F *hJZBneg = systsamples.Draw("hJZBneg",var,binning, "JZB [GeV]", "events",kcut&&kNeg,mc,simulatedlumi,systsamples.FindSample("/DY")); |
2636 |
|
hJZBpos->SetLineColor(kBlack); |
2637 |
|
hJZBneg->SetLineColor(kRed); |
2638 |
|
|
2639 |
< |
Int_t nbins = 5; |
2514 |
< |
Float_t xmax = 100.; |
2515 |
< |
|
2516 |
< |
|
2517 |
< |
TCanvas *zcan = new TCanvas("zcan","zcan"); |
2518 |
< |
zcan->SetLogy(1); |
2519 |
< |
|
2639 |
> |
hJZBpos->SetMinimum(1.0); |
2640 |
|
hJZBpos->Draw("e1"); |
2641 |
|
hJZBneg->Draw("same,hist"); |
2642 |
|
hJZBpos->Draw("same,e1"); // So it's on top... |
2656 |
|
hratio->GetYaxis()->SetTitle("Observed/Predicted"); |
2657 |
|
hratio->Draw("e1"); |
2658 |
|
|
2659 |
< |
TLine *top = new TLine(sbg_min,1.25,sbg_max,1.25); |
2660 |
< |
TLine *center = new TLine(sbg_min,1.0,sbg_max,1.0); |
2661 |
< |
TLine *bottom = new TLine(sbg_min,0.75,sbg_max,0.75); |
2659 |
> |
TLine *top = new TLine(binning[0],1.25,binning[binning.size()-1],1.25); |
2660 |
> |
TLine *center = new TLine(binning[0],1.0,binning[binning.size()-1],1.0); |
2661 |
> |
TLine *bottom = new TLine(binning[0],0.75,binning[binning.size()-1],0.75); |
2662 |
|
|
2663 |
|
|
2664 |
|
top->SetLineColor(kBlue);top->SetLineStyle(2); |
2700 |
|
tout << "\\begin{tabular}{l|cc}" << endl; |
2701 |
|
tout << "\\hline" << endl; |
2702 |
|
tout << "& {\\OFZP} & {\\OFSB} \\\\\\hline" << endl; |
2703 |
< |
tout << "\\#(events) & "<<OFSB<<" & "<<OFZP<<"\\\\ \\hline" << endl; |
2703 |
> |
tout << "\\#(events) & "<<OFZP<<" & "<<OFSB<<"\\\\ \\hline" << endl; |
2704 |
|
tout << "\\end{tabular}" << endl; |
2705 |
|
tout << "\\end{center}" << endl; |
2706 |
|
tout << "\\end{table}" << endl; |
2718 |
|
|
2719 |
|
TCanvas *cannie = new TCanvas("cannie","cannie"); |
2720 |
|
|
2721 |
< |
for(int icut=0;icut<jzb_cuts.size();icut++) { |
2721 |
> |
for(int icut=0;icut<(int)jzb_cuts.size();icut++) { |
2722 |
|
float currcut=jzb_cuts[icut]; |
2723 |
|
int nbins=1;float low=currcut; |
2724 |
|
vector<int> mysample; |
2777 |
|
//prediction part |
2778 |
|
if(is_data) cout << "Data prediction & "; |
2779 |
|
if(subselection!="none") cout << subselection << " prediction &"; |
2780 |
< |
for(int ij=0;ij<jzb_cuts.size();ij++) cout << jzbcutprediction[ij] << " vs " << metcutprediction[ij] << " & "; |
2780 |
> |
for(int ij=0;ij<(int)jzb_cuts.size();ij++) cout << jzbcutprediction[ij] << " vs " << metcutprediction[ij] << " & "; |
2781 |
|
|
2782 |
|
cout << endl; |
2783 |
|
//observation part |
2784 |
|
if(is_data) cout << "Data observation & "; |
2785 |
|
if(subselection!="none") cout << subselection << " observation &"; |
2786 |
< |
for(int ij=0;ij<jzb_cuts.size();ij++) cout << jzbcutobservation[ij] << " vs " << metcutobservation[ij] << " & "; |
2786 |
> |
for(int ij=0;ij<(int)jzb_cuts.size();ij++) cout << jzbcutobservation[ij] << " vs " << metcutobservation[ij] << " & "; |
2787 |
|
cout << endl; |
2788 |
|
cout << "_________________________________________________________________" << endl; |
2789 |
|
delete cannie; |
2844 |
|
leg->SetLineColor(kBlack); |
2845 |
|
leg->SetTextFont(62); |
2846 |
|
|
2847 |
< |
for ( int icut=0; icut<jzb_bins.size(); ++icut ) { |
2847 |
> |
for ( int icut=0; icut<(int)jzb_bins.size(); ++icut ) { |
2848 |
|
|
2849 |
|
ostringstream selection; |
2850 |
|
selection << mcjzb << ">" << jzb_bins[icut]; |
2896 |
|
// Calls the above function for each signal sample |
2897 |
|
void plot_jzb_sel_eff(string mcjzb, samplecollection &signalsamples, vector<float> bins ) |
2898 |
|
{ |
2899 |
< |
for (int isignal=0; isignal<signalsamples.collection.size();isignal++) { |
2899 |
> |
for (int isignal=0; isignal<(int)signalsamples.collection.size();isignal++) { |
2900 |
|
dout << "JZB selection efficiency curve: " << std::endl; |
2901 |
|
JZBSelEff(mcjzb,(signalsamples.collection)[isignal].events,(signalsamples.collection)[isignal].samplename,bins); // Only for some selected samples |
2902 |
|
} |
2939 |
|
TH1F *RcorrJZBSBeemm; |
2940 |
|
TH1F *LcorrJZBSBeemm; |
2941 |
|
|
2942 |
< |
TH1F *RcorrJZBeemmNoS; |
2942 |
> |
// TH1F *RcorrJZBeemmNoS; |
2943 |
|
|
2944 |
|
//these are for the ratio |
2945 |
|
TH1F *JRcorrJZBeemm = qcdsamples.Draw("JRcorrJZBeemm",jzb.c_str(),PlottingSetup::global_ratio_binning, "JZB [GeV]", "events", cutmass&&cutOSSF&&cutnJets,is_data, luminosity,use_signal); |
3023 |
|
//3rd last argument: do special bpred ratio, 2nd last argument: extended range!, last: y-axis title |
3024 |
|
string ytitle("ratio"); |
3025 |
|
if ( use_data==1 ) ytitle = "data/pred"; |
3026 |
< |
save_with_ratio(JRcorrJZBeemm,JBpred,kinpad,"QCD/Bpred",true,use_data!=1,ytitle); |
3026 |
> |
save_with_ratio(JRcorrJZBeemm,JBpred,kinpad,"QCD/Bpred",true,false,ytitle); |
3027 |
|
|
3028 |
|
TH1F *allevents = qcdsamples.Draw("allevents","pfJetGoodNum",1,0,100, "internal code", "events", "" ,mc, luminosity); |
3029 |
|
TH1F *ossf = qcdsamples.Draw("ossf","pfJetGoodNum",1,0,100, "internal code", "events", cutOSSF ,mc, luminosity); |
3070 |
|
dout << " Observation increases by : " << RcorrJZBeemm->Integral() << endl; |
3071 |
|
|
3072 |
|
dout << endl; |
3073 |
< |
for(int i=0;i<bins.size();i++) { |
3073 |
> |
for(int i=0;i<(int)bins.size();i++) { |
3074 |
|
dout << " JZB > " << bins[i] << " : " << endl; |
3075 |
|
dout << " Observation increases by : " << RcorrJZBeemm->Integral(RcorrJZBeemm->FindBin(bins[i]),RcorrJZBeemm->GetNbinsX()) << endl; |
3076 |
|
if(PlottingSetup::RestrictToMassPeak) { |
3118 |
|
leg->SetY2(1.0); |
3119 |
|
|
3120 |
|
|
3121 |
< |
for(int isample=0;isample<allsamples.collection.size();isample++) { |
3121 |
> |
for(int isample=0;isample<(int)allsamples.collection.size();isample++) { |
3122 |
|
string nowname=(allsamples.collection)[isample].filename; |
3123 |
|
cout << "Drawing: " << nowname << " (sample " << isample+1 << " / " << allsamples.collection.size() << ")" << endl; |
3124 |
|
individualpt1histos[isample] = allsamples.Draw(GetNumericHistoName(),"pt1",50,0,50, "p_{T,1}", "events",cutOSSF&&cutnJets,mc,luminosity,allsamples.FindSample(nowname)); |
3142 |
|
ptsancan->cd(2); |
3143 |
|
fpt2->Draw(); |
3144 |
|
|
3145 |
< |
for(int isample=0;isample<allsamples.collection.size();isample++) { |
3145 |
> |
for(int isample=0;isample<(int)allsamples.collection.size();isample++) { |
3146 |
|
ptsancan->cd(1); |
3147 |
|
individualpt1histos[isample]->DrawNormalized("same,histo"); |
3148 |
|
ptsancan->cd(2); |
3158 |
|
void do_mlls_plot(string mcjzb) { |
3159 |
|
cout << "At this point we'd plot the mll distribution" << endl; |
3160 |
|
TCanvas *sigcan = new TCanvas("sigcan","sigcan"); |
3161 |
< |
for(int isig=0;isig<(signalsamples.collection).size();isig++) { |
3161 |
> |
for(int isig=0;isig<(int)(signalsamples.collection).size();isig++) { |
3162 |
|
if(!(signalsamples.collection)[isig].events) continue; |
3163 |
|
string nowname=(signalsamples.collection)[isig].filename; |
3164 |
|
TH1F *mll = signalsamples.Draw("mllhisto","mll",150,0,150, "m_{ll}", "events",cutOSSF&&cutnJets,mc,luminosity,signalsamples.FindSample(nowname)); |
3197 |
|
findme.push_back("TTJets"); |
3198 |
|
findme.push_back("LM"); |
3199 |
|
|
3200 |
< |
for(int ifind=0;ifind<findme.size();ifind++) { |
3200 |
> |
for(int ifind=0;ifind<(int)findme.size();ifind++) { |
3201 |
|
vector<int> selsamples = allsamples.FindSample(findme[ifind]); |
3202 |
|
TH2F *metvsjzb = new TH2F("metvsjzb","metvsjzb",200,0,100,400,-100,100); |
3203 |
< |
for(int isel=0;isel<selsamples.size();isel++) { |
3203 |
> |
for(int isel=0;isel<(int)selsamples.size();isel++) { |
3204 |
|
cout << "Producing MET:JZB plot ... working on sample: " << allsamples.collection[selsamples[isel]].filename << endl; |
3205 |
|
allsamples.collection[selsamples[isel]].events->Draw("jzb[1]:met[4]>>+metvsjzb",cutmass&&cutOSSF); |
3206 |
|
} |