ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/cbrown/Development/Plotting/Modules/ZbTools.C

Comparing UserCode/cbrown/Development/Plotting/Modules/ZbTools.C (file contents):
Revision 1.2 by buchmann, Tue Jul 10 16:00:42 2012 UTC vs.
Revision 1.13 by buchmann, Fri Nov 23 12:10:58 2012 UTC

#	Line 23 | Line 23 | using namespace PlottingSetup;
23		namespace ZbData {
24		vector<float> data_over_mc;
25
26	<	TCut ZplusBsel("pt1>20&&pt2>20&&mll>60&&mll<120");
27	<	TCut LeadingB("bTagProbCSVBP[0]>0.679");
28	<	TCut EtaB("abs(pfJetGoodEta[0])<1.3");
29	<	TCut PhiZcut("pfJetDphiZ[0]>2.7");
26	>	TCut ZplusBsel("pt1>20&&pt2>20&&mll>60&&mll<120&&id1==id2");
27	>	TCut LeadingB("Zb3010_bTagProbCSVBP[0]>0.898");
28	>	TCut EtaB("abs(Zb3010_pfJetGoodEta[0])<1.3");
29	>	TCut PhiZcut("abs(Zb3010_pfJetDphiZ[0])>2.7");
30
31		}
32
33		using namespace ZbData;
34
35	+
36	+	//*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*  DELETE EVERYTHING BETWEEN THIS LINE /*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*//
37	+	void SpecialCutFlow(string identifier);
38	+
39	+	//*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*            AND THIS LINE             /*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*//
40	+
41	+
42		void print_yield(TCut cut) {
43		TH1F *data  = allsamples.Draw("data",    "mll",1,0,10000000, "m_{ll} [GeV]", "events", cut,data,luminosity);
44	<	cout << "  Yields for " << (const char*) cut << endl;
45	<	cout << "     data: " << data->Integral() << endl;
44	>	dout << "  Yields for " << (const char*) cut << endl;
45	>	dout << "     data: " << data->Integral() << endl;
46		THStack mcm = allsamples.DrawStack("mc", "mll",1,0,10000000, "m_{ll} [GeV]", "events", cut,mc,luminosity);
47		int nhists = mcm.GetHists()->GetEntries();
48	<	cout << "Going to loop over " << nhists << " histograms " << endl;
48	>	dout << "Going to loop over " << nhists << " histograms " << endl;
49		TFile *test = new TFile("test.root","RECREATE");
50		mcm.Write();
51		test->Close();
#	Line 46 | Line 53 | void print_yield(TCut cut) {
53		TH1* hist = static_cast<TH1*>(mcm.GetHists()->At(i));
54		cout << "     " << hist->GetName() << ": " << hist->Integral() << endl;
55		}
56	+
57		cout << endl << endl;
58	+
59		delete data;
60		}
61
#	Line 100 | Line 109 | Value getfrom2Dmap(TH2F *map, int ixbin,
109		return sum;
110		}
111
103	–
104	–	void data_mc_agreement_2d() {
105	–
106	–	vector<int> emptyvector;
107	–	float binX[5]={10,50,100,200,1000};
108	–	float binY[5]={0.0,0.1,0.15,0.2,0.3};
109	–	vector<float>binningx;
110	–	vector<float>binningy;
111	–	binningx.assign(binX,binX+5);
112	–	binningy.assign(binY,binY+5);
113	–
114	–	TCut cut(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
115	–
116	–	TCanvas *agrcan = new TCanvas("agrcan","agrcan");
117	–
118	–	cout << "Drawing 2d for data " << endl;
119	–	TH2F *alphazptplane_data = allsamples.Draw("alphazptplane_data","alpha:pt",binningx,binningy,"Z p_{T}","#alpha",cut,data,luminosity,emptyvector,false);
120	–	cout << "Drawing 2d for MC " << endl;
121	–	TH2F *alphazptplane_mc   = allsamples.Draw("alphazptplane_mc",  "alpha:pt",binningx,binningy,"Z p_{T}","#alpha",cut,mc,  luminosity,emptyvector,false);
122	–
123	–
124	–	TFile *f = new TFile("test.root","RECREATE");
125	–	alphazptplane_data->Write();
126	–	f->Close();
127	–
128	–	for(int ix=1;ix<binningx.size();ix++) {
129	–	cout << "pt in " << binningx[ix-1] << " , " << binningx[ix] << endl;
130	–	TGraphErrors gr;
131	–	for(int iy=1;iy<binningy.size();iy++) {
132	–	Value mc   = Value((getfrom2Dmap(alphazptplane_mc,ix,iy)).getValue(),0.0);// not taking MC statistics into account
133	–	Value data = getfrom2Dmap(alphazptplane_data,ix,iy);
134	–	cout << " For an alpha cut at " << binningy[iy] << " we get : mc (" << mc << ") , data (" << data << ") --> data/mc = " << data/mc << endl;
135	–	gr.SetPoint(iy-1,binningy[iy],(data/mc).getValue());
136	–	gr.SetPointError(iy-1,0,(data/mc).getError());
137	–	}
138	–
139	–	gr.Draw("AP*");
140	–	gr.GetYaxis()->SetTitle("data / mc");
141	–	gr.GetYaxis()->CenterTitle();
142	–	gr.GetXaxis()->SetTitle("#alpha");
143	–	gr.GetXaxis()->CenterTitle();
144	–	gr.GetXaxis()->SetRangeUser(0,0.4);
145	–	gr.GetYaxis()->SetRangeUser(0.5,1.5);
146	–	gr.Draw("AP*");
147	–	DrawPrelim();
148	–	stringstream ptbin;
149	–	ptbin << binningx[ix-1] << " < p_{T}^{Z} < " << binningx[ix] << " GeV";
150	–	TText *ptt = write_title(ptbin.str());
151	–	ptt->SetX(0.7);
152	–	ptt->SetY(0.7);
153	–	ptt->SetTextSize(0.03);
154	–	ptt->Draw();
155	–	TLine *l = new TLine(0.0,1.0,0.4,1.0);
156	–	l->SetLineStyle(2);
157	–	l->SetLineColor(kBlue);
158	–	l->Draw();
159	–
160	–	stringstream saveas;
161	–	saveas << "Zplusb_data_over_mc___" << binningx[ix-1] << "to" << binningx[ix];
162	–	CompleteSave(agrcan,saveas.str());
163	–	//    agrcan->SaveAs("tester.png");
164	–	}
165	–
166	–	delete agrcan;
167	–
168	–	}
169	–
112		void print_all_b_yields() {
113		cout << "Basic selection with a b jet" << endl;
114	<	print_yield(ZplusBsel&&"pfJetGoodNumBtag>0");
114	>	print_yield(ZplusBsel&&"Zb3010_pfJetGoodNumBtag>0");
115		cout << "Leading jet is a b-jet " << endl;
116		print_yield(ZplusBsel&&LeadingB);
117		cout << "|#eta_{b}|<1.3 " << endl;
#	Line 178 | Line 120 | void print_all_b_yields() {
120		print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
121		cout << "10<ptZ<1000 GeV" << endl;
122		print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000");
123	+	cout << "#alpha < 0.35" << endl;
124	+	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.35");
125		cout << "#alpha < 0.3" << endl;
126	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.3");
126	>	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.3");
127		cout << "#alpha < 0.2" << endl;
128	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.2");
128	>	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.2");
129		cout << "#alpha < 0.15" << endl;
130	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.15");
130	>	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.15");
131		cout << "#alpha < 0.1" << endl;
132	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.1");
133	<	cout << "#alpha < 0.05" << endl;
134	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.05");
132	>	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.1");
133	>	//  cout << "#alpha < 0.05" << endl;
134	>	//  print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.05");
135		}
136
137	+	void DrawEvilCutFlow() {
138	+	stringstream MegaCut;
139	+
140	+	MegaCut << "   1*(" << (const char*) (ZplusBsel&&TCut("Zb3010_pfJetGoodNumBtag>0")) << ")";
141	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB) << ")";
142	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB) << ")";
143	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut) << ")";
144	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000")) << ")";
145	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.3")) << ")";
146	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.2")) << ")";
147	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.15")) << ")";
148	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.1")) << ")";
149	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.05")) << ")";
150	+	MegaCut << " ";
151	+
152	+	TCanvas *can = new TCanvas("can","can");
153	+	can->SetLogy(1);
154	+	TCut basecut("mll>6||mll<10");
155	+
156	+
157	+
158	+	TLegend *leg = allsamples.allbglegend();
159	+
160	+
161	+	TH1F *data  = allsamples.Draw("data",    MegaCut.str(),11,-0.5,10.5, "", "events", basecut,data,luminosity);
162	+	THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
163	+
164	+	float runningsum=0;
165	+	for(int i=data->GetNbinsX();i>0;i--) {
166	+	runningsum+=data->GetBinContent(i);
167	+	data->SetBinContent(i,runningsum);
168	+	}
169	+
170	+	float minimum=data->GetMaximum();
171	+
172	+	TH1F* h;
173	+	TIter nextOF(mcm.GetHists());
174	+	while ( h = (TH1F*)nextOF() ) {
175	+	float runningsum=0;
176	+	minimum=data->GetMaximum();//done deliberately at each step so get the minimum of the last (!) contribution
177	+	for(int i=h->GetNbinsX();i>0;i--) {
178	+	runningsum+=h->GetBinContent(i);
179	+	h->SetBinContent(i,runningsum);
180	+	if(runningsum<minimum&&runningsum>0) minimum=runningsum;
181	+	}
182	+	}
183	+
184	+
185	+
186	+	data->SetMinimum(0.05*minimum);
187	+	can->cd(1)->SetBottomMargin(0.25);
188	+	data->GetXaxis()->SetBinLabel(1,"Pre-selection");
189	+	data->GetXaxis()->SetBinLabel(2,"Contains b jet");
190	+	data->GetXaxis()->SetBinLabel(3,"Leading jet is b-jet");
191	+	data->GetXaxis()->SetBinLabel(4,"|#eta_{b}|<1.3 ");
192	+	data->GetXaxis()->SetBinLabel(5,"|#delta#phi(b,Z)|>2.7");
193	+	data->GetXaxis()->SetBinLabel(6,"10<p_{T}^{Z}<1000 GeV");
194	+	data->GetXaxis()->SetBinLabel(7,"#alpha < 0.3");
195	+	data->GetXaxis()->SetBinLabel(8,"#alpha < 0.2");
196	+	data->GetXaxis()->SetBinLabel(9,"#alpha < 0.15");
197	+	data->GetXaxis()->SetBinLabel(10,"#alpha < 0.1");
198	+	data->GetXaxis()->SetBinLabel(11,"#alpha < 0.05");
199	+	data->GetXaxis()->LabelsOption("v");
200	+
201	+	data->Draw("e1");
202	+	mcm.Draw("same");
203	+	data->Draw("same,e1,axis");
204	+	data->Draw("same,e1");
205	+	//  data->Draw("same,TEXT");
206	+
207	+	leg->Draw();
208	+
209	+	DrawPrelim();
210	+
211	+	CompleteSave(can,"CutFlow");
212	+	}
213	+
214		void draw_Zb_kin_vars() {
215	<	draw_kin_variable("pfBJetDphiZ",ZplusBsel&&LeadingB,100,0,3.2,"#delta#phi (Z,b lead)","DeltaPhiZBlead",1);
216	<	draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",50,0,200,"Z p_{T}","Zpt",1);
217	<	draw_kin_variable("pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",50,0,200,"Sub-Leading Jet Pt","Jet2Pt",1);
218	<	draw_kin_variable("pfJetGoodPtlBtag[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",50,0,200,"Leading Jet Pt (B)","LeadingJetPt",1);
219	<	draw_kin_variable("alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,2,"#alpha","alpha",1);
215	>	draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Z p_{T}","Official/Zpt",1);
216	>	draw_kin_variable("Zb3010_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt",1);
217	>	draw_kin_variable("Zb3010_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt",1);
218	>	draw_kin_variable("Zb3010_pfJetGoodPt[1]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,2,"p_{T}^{J2} / p_{T}^{Z}","Official/SubLeadingJetPt_Over_ZPt",1);
219	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<50",   40,0,2,"#alpha","Official/alpha_pt_30_to_50",1);
220	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>50&&pt<75",  40,0,2,"#alpha","Official/alpha_50_to_75",1);
221	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>75&&pt<125", 40,0,2,"#alpha","Official/alpha_pt_75_to_125",1);
222	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>125&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_125_to_1000",1);
223	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_30_to_1000",1);
224	>
225	>	draw_kin_variable("Zb3010_pfBJetDphiZ[0]",ZplusBsel&&LeadingB&&"pt>10&&pt<1000&&pfBJetDphiZ[0]>0",30,0,3.2,"#delta#phi (Z,b lead)","DeltaPhiZBlead",0);
226	>
227	>	/*
228	>	draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,200,"p_{T}^{l1}","Lep/pt1",1);
229	>	draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==0",20,0,200,"p_{T}^{e1}","Lep/pt1_e",1);
230	>	draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==1",20,0,200,"p_{T}^{#mu1}","Lep/pt1_m",1);
231	>	draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,200,"p_{T}^{l2}","Lep/pt2",1);
232	>	draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==0",20,0,200,"p_{T}^{e2}","Lep/pt2_e",1);
233	>	draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==1",20,0,200,"p_{T}^{#mu2}","Lep/pt2_m",1);
234	>	draw_kin_variable("eta1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,-3.1,3.1,"#eta^{l1}","Lep/eta1",1);
235	>	draw_kin_variable("eta1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==0",20,-3.1,3.1,"#eta^{e1}","Lep/eta1_e",1);
236	>	draw_kin_variable("eta1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==1",20,-3.1,3.1,"#eta^{m1}","Lep/eta1_m",1);
237	>	draw_kin_variable("eta2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,-3.1,3.1,"#eta^{l1}","Lep/eta2",1);
238	>	draw_kin_variable("eta2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==0",20,-3.1,3.1,"#eta^{e2}","Lep/eta2_e",1);
239	>	draw_kin_variable("eta2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==1",20,-3.1,3.1,"#eta^{#mu2}","Lep/eta2_m",1);
240	>	draw_kin_variable("numVtx",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",30,0,30,"N(Vtx)","Lep/nVtx",1);
241	>	*/
242		}
243
244		void draw_mpf_vars() {
245		ZbData::data_over_mc.clear();
246	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.3",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p3",0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
247	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.2",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p2",0,"#alpha<0.2, 10 GeV < p_{T}^{Z} < 1000 GeV");
248	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.15",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p15",0,"#alpha<0.15, 10 GeV < p_{T}^{Z} < 1000 GeV");
249	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.1",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p1",0,"#alpha<0.1, 10 GeV < p_{T}^{Z} < 1000 GeV");
250	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.05",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
251	<
252	<	draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.3",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p3",0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
253	<	draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.2",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p2",0,"#alpha<0.2, 10 GeV < p_{T}^{Z} < 1000 GeV");
254	<	draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.15",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p15",0,"#alpha<0.15, 10 GeV < p_{T}^{Z} < 1000 GeV");
255	<	draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.1",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p1",0,"#alpha<0.1, 10 GeV < p_{T}^{Z} < 1000 GeV");
256	<	draw_kin_variable("pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"alpha<0.05",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz__mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
246	>	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.3",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p3",0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
247	>	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.2",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p2",0,"#alpha<0.2, 10 GeV < p_{T}^{Z} < 1000 GeV");
248	>	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.15",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p15",0,"#alpha<0.15, 10 GeV < p_{T}^{Z} < 1000 GeV");
249	>	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.1",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p1",0,"#alpha<0.1, 10 GeV < p_{T}^{Z} < 1000 GeV");
250	>	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.05",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
251	>
252	>	draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.3",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p3",0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
253	>	draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.2",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p2",0,"#alpha<0.2, 10 GeV < p_{T}^{Z} < 1000 GeV");
254	>	draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.15",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p15",0,"#alpha<0.15, 10 GeV < p_{T}^{Z} < 1000 GeV");
255	>	draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.1",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz___alpha_smaller_0p1",0,"#alpha<0.1, 10 GeV < p_{T}^{Z} < 1000 GeV");
256	>	draw_kin_variable("Zb3010_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb3010_alpha<0.05",20,0,2,"p_{T} b-jet / p_{T} Z","ptb_over_ptz__mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
257		}
258
259	<	Value get_Zb_data_over_mc(string variable, TCut cut) {
259	>	/*double GetMedian(TH1F *histo) {
260	>	int numBins = histo->GetXaxis()->GetNbins();
261	>	Double_t x[numBins];
262	>	Double_t y[numBins];
263	>	for (int i = 1; i <= numBins; i++) {
264	>	x[i] = histo->GetBinCenter(i);
265	>	y[i] = histo->GetBinContent(i);
266	>	}
267	>	return TMath::Median(numBins, &x, &y);
268	>	}*/
269	>
270	>	Value get_Zb_data_over_mc(string variable, string variable2, TCut cut, string saveas) {
271		//  write_warning(__FUNCTION__,"Debugging this function, therefore always returning 3.1415+/-1.010101"); return Value(3.1415,1.010101);
272	<	TH1F *hdata  = allsamples.Draw("hdata",variable,1000,0,30, "#alpha [GeV]", "events", cut,data,luminosity);
273	<	TH1F *hmc    = allsamples.Draw("hmc" , variable,1000,0,30, "#alpha [GeV]", "events", cut,mc,  luminosity);
272	>	TCanvas *cn = new TCanvas("cn","cn");
273	>	string varname="MPF";
274	>	if(Contains(variable,"pfJetGood")) varname="R_{abs}";
275	>	TH1F *hdata  = allsamples.Draw("hdata",variable,1200,0,30, varname, "events",  cut,data,luminosity);
276	>	TH1F *hmc    = allsamples.Draw("hmc" , variable,1200,0,30, varname, "events", cut, mc,  luminosity);
277
278		float a=hdata->GetMean();
279		float b=hmc->GetMean();
#	Line 224 | Line 281 | Value get_Zb_data_over_mc(string variabl
281		float db=hmc->GetMeanError();
282		float factor = a / b;
283		float error = TMath::Sqrt(  (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
284	+
285	+	cn->cd();
286	+	hdata->GetYaxis()->SetTitle("events / 0.1");
287	+	hdata->Rebin(4);
288	+	hmc->Rebin(4);
289	+	hdata->GetXaxis()->SetRangeUser(0,2);
290	+	hdata->Draw("e1");
291	+	hmc->Draw("histo,same");
292	+	hdata->Draw("e1,same");
293	+	stringstream summary;
294	+	summary << "#splitline{" << varname << "}{#splitline{data: " << std::setprecision(4) << a << " #pm " << std::setprecision(4) << da << "}{#splitline{mc:   " << std::setprecision(4) << b << " #pm " << std::setprecision(4) << db << "}{ratio: " << factor;
295	+	summary << " #pm " << error << "}}}";
296	+	TText *infobox = write_title(summary.str());
297	+	infobox->SetX(0.75);
298	+	infobox->SetTextSize(0.03);
299	+	infobox->SetY(0.75);
300	+	infobox->Draw();
301	+
302	+
303	+	DrawPrelim();
304	+	CompleteSave(cn,"ResponsePlots/"+saveas);
305	+
306	+	cout << "Data;" << a << ";"<<da<<";MC"<<b<<";"<<db<<endl;
307	+	delete cn;
308		delete hdata;
309		delete hmc;
310	+
311	+	//    cout << a << "+/-" << da << "        " << b << "=/-" << db << endl;
312	+
313		return Value(factor,error);
314		}
315
316
317
318		void new_data_mc_agreement_2d() {
319	+
320	+	write_error(__FUNCTION__,"need to make alpha<0.3 inclusive again!!!! IMPORTANT SINCE THAT'S THE FINAL RESULT!!!");
321	+
322	+	gStyle->SetOptFit(0);
323		TCut basecut(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
324	+	const int nalphacuts=6;
325	+	float alphacuts[nalphacuts] = {0.1,0.15,0.2,0.25,0.3,0.35};
326	+	int PointThree=-1;
327	+	for(int i=0;i<nalphacuts;i++) {
328	+	if((alphacuts[i]-0.3)<0.01) PointThree=i-1;
329	+	}
330	+	if(PointThree<0) {
331	+	write_warning(__FUNCTION__,"alpha=0.3 MUST BE ONE OF THE ALPHA VALUES!");
332	+	assert(PointThree<0);
333	+	}
334
237	–	const int nalphacuts=5;
238	–	float alphacuts[nalphacuts] = {0.1,0.15,0.2,0.25,0.3};
335		const int nptcuts=5;
336	<	float ptcuts[nptcuts]={10,50,100,200,1000};
336	>	float ptcuts[nptcuts]={30,50,75,125,1000};
337
338		float MPF_Results[nalphacuts][nptcuts];
339		float MPF_Errors[nalphacuts][nptcuts];
#	Line 248 | Line 344 | void new_data_mc_agreement_2d() {
344		for(int ia=0;ia<nalphacuts;ia++) {
345		for(int ipt=0;ipt<nptcuts-1;ipt++) {
346		stringstream specialcut;
347	<	specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (alpha<" << alphacuts[ia] << "))";
348	<	Value MPF_data_over_mc = get_Zb_data_over_mc("mpf",TCut(basecut && specialcut.str().c_str()));
349	<	Value RABS_data_over_mc = get_Zb_data_over_mc("pfJetGoodPt[0]/pt",TCut(basecut && specialcut.str().c_str()));
347	>	if(ia>0) specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<" << alphacuts[ia] << " && Zb3010_alpha>" << alphacuts[ia-1] << "))";
348	>	else     specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<" << alphacuts[ia] << "))";
349	>	cout << specialcut.str() << endl;
350	>	Value MPF_data_over_mc = get_Zb_data_over_mc("mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(alphacuts[ia]));
351	>	Value RABS_data_over_mc = get_Zb_data_over_mc("Zb3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__alpha_"+any2string(alphacuts[ia]));
352
353		MPF_Results[ia][ipt]=MPF_data_over_mc.getValue();
354		MPF_Errors[ia][ipt]=MPF_data_over_mc.getError();
#	Line 258 | Line 356 | void new_data_mc_agreement_2d() {
356		RABS_Results[ia][ipt]=RABS_data_over_mc.getValue();
357		RABS_Errors[ia][ipt]=RABS_data_over_mc.getError();
358
359	<	cout << "alpha cut at " << alphacuts[ia] << " and pt range " << ptcuts[ipt] << " , " << ptcuts[ipt+1] << " \t : \t" << specialcut.str() << " \t \t --> " << MPF_data_over_mc << " (MPF) " << RABS_data_over_mc << " (RABS)" << endl;
359	>	cout << "alpha cut at " << alphacuts[ia+1] << " and pt range " << ptcuts[ipt] << " , " << ptcuts[ipt+1] << " \t : \t" << specialcut.str() << " \t \t --> " << MPF_data_over_mc << " (MPF) " << RABS_data_over_mc << " (RABS)" << endl;
360		}
361		}
362
#	Line 269 | Line 367 | void new_data_mc_agreement_2d() {
367
368		TGraphErrors *MPF_FinalGraph = new TGraphErrors();
369		TGraphErrors *RABS_FinalGraph = new TGraphErrors();
370	+	TGraphErrors *MPF_FaceValueAtPoint3 = new TGraphErrors();
371	+	TGraphErrors *RABS_FaceValueAtPoint3 = new TGraphErrors();
372
373		for(int ipt=0;ipt<nptcuts-1;ipt++) {
374
#	Line 277 | Line 377 | void new_data_mc_agreement_2d() {
377		cout << "Going to create histo for " << filename.str() << endl;
378		TGraphErrors *mpf_gr =new TGraphErrors();
379		TGraphErrors *rabs_gr =new TGraphErrors();
380	+
381	+	int rabs_counter=0;
382	+	int mpf_counter=0;
383	+
384		for(int ia=0;ia<nalphacuts;ia++) {
385	<	mpf_gr->SetPoint(ia,alphacuts[ia],MPF_Results[ia][ipt]);
386	<	mpf_gr->SetPointError(ia,0,MPF_Errors[ia][ipt]);
387	<	rabs_gr->SetPoint(ia,alphacuts[ia],RABS_Results[ia][ipt]);
388	<	rabs_gr->SetPointError(ia,0,RABS_Errors[ia][ipt]);
385	>	cout << "Setting a point for an alpha cut at " << alphacuts[ia] << " with value " << MPF_Results[ia][ipt] << " (MPF) and " << RABS_Results[ia][ipt] << " (RABS)" << endl;
386	>	if(MPF_Results[ia][ipt]==MPF_Results[ia][ipt] && MPF_Errors[ia][ipt]==MPF_Errors[ia][ipt]) { // remove nan's
387	>	mpf_gr->SetPoint(mpf_counter,alphacuts[ia],MPF_Results[ia][ipt]);
388	>	mpf_gr->SetPointError(mpf_counter,0,MPF_Errors[ia][ipt]);
389	>	mpf_counter++;
390	>	} else {
391	>	dout << "Detected BS for MPF method (bin ignored) : " << endl;
392	>	dout << "       alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
393	>	cout << "       value (MPF) : " << MPF_Results[ia][ipt] << " +/- " << MPF_Errors[ia][ipt] << endl;
394	>	}
395	>	if(RABS_Results[ia][ipt]==RABS_Results[ia][ipt] && RABS_Errors[ia][ipt]==RABS_Errors[ia][ipt]) {
396	>	rabs_gr->SetPoint(rabs_counter,alphacuts[ia],RABS_Results[ia][ipt]);
397	>	rabs_gr->SetPointError(rabs_counter,0,RABS_Errors[ia][ipt]);
398	>	rabs_counter++;
399	>	} else {
400	>	dout << "Detected BS for RABS method (bin ignored) : " << endl;
401	>	dout << "       alpha: " << alphacuts[ia] << " , pt range: " << ptcuts[ipt] << " < pt < " << ptcuts[ipt+1] << endl;
402	>	cout << "       value (RABS) : " << RABS_Results[ia][ipt] << " +/- " << RABS_Errors[ia][ipt] << endl;
403	>	}
404		}
405	+	cout << "About to add a point in pt (" << (ptcuts[ipt]+ptcuts[ipt+1])/2 << " ) to face value : " << MPF_Results[PointThree][ipt] << endl;
406	+
407	+	stringstream specialcut;
408	+	specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<0.3))";
409	+	Value MPF_data_over_mc = get_Zb_data_over_mc("mpf","",TCut(basecut && specialcut.str().c_str()),"MPF___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_0.3");
410	+	Value RABS_data_over_mc = get_Zb_data_over_mc("Zb3010_pfJetGoodPt[0]/pt","",TCut(basecut && specialcut.str().c_str()),"Rabs___pt_"+any2string(ptcuts[ipt])+"_"+any2string(ptcuts[ipt+1])+"__INCLUSIVEalpha_0.3");
411	+	MPF_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,MPF_data_over_mc.getValue());
412	+	MPF_FaceValueAtPoint3->SetPointError(ipt,0,MPF_data_over_mc.getError());
413	+	RABS_FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1])/2,RABS_data_over_mc.getValue());
414	+	RABS_FaceValueAtPoint3->SetPointError(ipt,0,RABS_data_over_mc.getError());
415
416	+	can->cd();
417		mpf_gr->SetMarkerStyle(21);
418		mpf_gr->SetMarkerColor(kBlue);
419		mpf_gr->Draw("AP");
420	<	mpf_gr->Fit("pol1");
420	>	mpf_gr->Fit("pol1","");
421	>	mpf_gr->GetXaxis()->SetTitle("#alpha");
422	>	mpf_gr->GetXaxis()->CenterTitle();
423	>	mpf_gr->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
424	>	mpf_gr->GetYaxis()->CenterTitle();
425		TF1 *mpf_pol=(TF1*)mpf_gr->GetFunction("pol1");
426	+	TF1 *mpf_pol_complete = new TF1("mpf_pol_complete","[0]+[1]*x");
427	+	mpf_pol_complete->SetParameters(mpf_pol->GetParameters());
428	+	mpf_pol_complete->SetLineWidth(1);
429	+	mpf_pol_complete->SetLineColor(kBlue);
430	+	float min,max;
431	+	FindMinMax(mpf_gr,min,max);
432	+	TH1F *histo = new TH1F("histo","histo",1,0,0.4);
433	+	histo->SetStats(0);
434	+	histo->GetXaxis()->SetRangeUser(0,0.4);
435	+	histo->GetYaxis()->SetRangeUser(min,max);
436	+	mpf_gr->GetYaxis()->SetRangeUser(min,max);
437	+
438	+	TPolyLine *mpf_fit_uncert = GetFitUncertaintyShape(mpf_gr, "pol1", min, max,0.0,0.4);
439	+	mpf_pol->SetLineColor(TColor::GetColor("#0101DF"));
440	+	mpf_fit_uncert->SetFillColor(TColor::GetColor("#0101DF"));
441	+	histo->Draw();
442	+	mpf_fit_uncert->Draw("same,f");
443	+	mpf_gr->Draw("P");
444	+	mpf_pol_complete->Draw("same");
445	+
446		float mpf_a=mpf_pol->GetParameter(0);
447		float mpf_b=mpf_pol->GetParameter(1);
448		MPF_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),mpf_a);
#	Line 297 | Line 451 | void new_data_mc_agreement_2d() {
451		MPF_ExtraPolatedResults[ipt]=mpf_a;
452
453		stringstream mpf_info;
454	<	mpf_info << "#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (MPF) }}{Extrapolated value: " << std::setprecision(3) << mpf_a << "}";
454	>	mpf_info <<  "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (MPF) }}{Extrapolated value: " << std::setprecision(3) << mpf_a << "}}{#chi^{2}/NDF : " << mpf_pol->GetChisquare() << " / " << mpf_pol->GetNDF() << "}";
455	>
456		TText *mpf_mark = write_title(mpf_info.str());
457		mpf_mark->SetX(0.75);
458		mpf_mark->SetTextSize(0.03);
#	Line 306 | Line 461 | void new_data_mc_agreement_2d() {
461
462		string filenamebkp=filename.str();
463		filename << "__MPF";
464	+	DrawPrelim();
465		CompleteSave(can,filename.str());
466		cout << "MPF : " << mpf_a << " + " << mpf_b << " * alpha " << endl;
467
468		filename.str("");
469
470		rabs_gr->SetMarkerStyle(21);
471	<	rabs_gr->SetMarkerColor(kBlue);
471	>	rabs_gr->SetMarkerColor(kRed);
472		rabs_gr->Draw("AP");
473	<	rabs_gr->Fit("pol1");
473	>	rabs_gr->Fit("pol1","");
474	>	rabs_gr->GetXaxis()->SetTitle("#alpha");
475	>	rabs_gr->GetXaxis()->CenterTitle();
476	>	rabs_gr->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
477	>	rabs_gr->GetYaxis()->CenterTitle();
478		TF1 *rabs_pol=(TF1*)rabs_gr->GetFunction("pol1");
479	+	TF1 *rabs_pol_complete = new TF1("rabs_pol_complete","[0]+[1]*x");
480	+	rabs_pol_complete->SetParameters(rabs_pol->GetParameters());
481	+	rabs_pol_complete->SetLineColor(kRed);
482	+	rabs_pol_complete->SetLineWidth(1);
483	+	FindMinMax(rabs_gr,min,max);
484	+	rabs_gr->GetYaxis()->SetRangeUser(min,max);
485	+	histo->GetYaxis()->SetRangeUser(min,max);
486	+	TPolyLine *rabs_fit_uncert = GetFitUncertaintyShape(rabs_gr, "pol1", min, max,0.0,0.4);
487	+	rabs_pol->SetLineColor(TColor::GetColor("#0101DF"));
488	+	rabs_pol->SetFillColor(TColor::GetColor("#0101DF"));
489	+	histo->Draw();
490	+	rabs_fit_uncert->Draw("same,f");
491	+	rabs_gr->Draw("P");
492	+	rabs_pol_complete->Draw("same");
493	+
494	+
495		float rabs_a=rabs_pol->GetParameter(0);
496		float rabs_b=rabs_pol->GetParameter(1);
497		RABS_FinalGraph->SetPoint(ipt,0.5*(ptcuts[ipt]+ptcuts[ipt+1]),rabs_a);
498		RABS_FinalGraph->SetPointError(ipt,0,rabs_pol->GetParError(0));
323	–	cout << "!+!+!+!+!+!+!+!+!+ Just added a point to the final plot : " << rabs_a << " +/- " << rabs_pol->GetParError(0) << endl;
499
500	+	cout << "!+!+!+!+!+!+!+!+!+ Just added a point to the final plot : " << rabs_a << " +/- " << rabs_pol->GetParError(0) << endl;
501
502		stringstream rabs_info;
503	<	rabs_info << "#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (R_{abs}) }}{Extrapolated value: " << std::setprecision(3) << rabs_a << "}";
503	>	rabs_info << "#splitline{#splitline{#splitline{" << ptcuts[ipt] << " GeV < p_{T}^{Z} < " << ptcuts[ipt+1] << " GeV }{ (R_{abs}) }}{Extrapolated value: " << std::setprecision(3) << rabs_a << "}}{#chi^{2}/NDF : " << rabs_pol->GetChisquare() << " / " << rabs_pol->GetNDF() << "}";
504		TText *rabs_mark = write_title(rabs_info.str());
505		rabs_mark->SetX(0.75);
506		rabs_mark->SetTextSize(0.03);
507		rabs_mark->SetY(0.75);
508		rabs_mark->Draw();
509	<
509	>
510	>
511		RABS_ExtraPolatedResults[ipt]=rabs_a;
512		filename << filenamebkp << "__RABS";
513	+	DrawPrelim();
514		CompleteSave(can,filename.str());
515		cout << "RABS : " << rabs_a << " + " << rabs_b << " * alpha " << endl;
516		}
517
518		MPF_FinalGraph->SetMarkerStyle(21);
519	<	MPF_FinalGraph->SetMarkerColor(kRed);
520	<	MPF_FinalGraph->Fit("pol0");
519	>	MPF_FinalGraph->SetMarkerColor(kBlue);
520	>	MPF_FinalGraph->Fit("pol0",""); // quiet, use minos
521	>
522		TF1 *mpf_pol0=(TF1*)MPF_FinalGraph->GetFunction("pol0");
523	+	TF1 *mpf_pol0_complete = new TF1("mpf_pol0_complete","[0]+[1]*x");
524	+	mpf_pol0_complete->SetLineWidth(1);
525	+	mpf_pol0_complete->SetLineColor(kBlue);
526	+	mpf_pol0->SetLineColor(kBlue);
527	+	MPF_FinalGraph->Draw("AP");
528	+	MPF_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
529	+	MPF_FinalGraph->GetXaxis()->CenterTitle();
530	+	MPF_FinalGraph->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
531	+	MPF_FinalGraph->GetYaxis()->CenterTitle();
532		MPF_FinalGraph->Draw("AP");
533	+	mpf_pol0_complete->Draw("same");
534	+
535		stringstream mpf_result;
536	<	mpf_result << "C_{abs}= " << mpf_pol0->GetParameter(0) << " #pm " << mpf_pol0->GetParError(0);
536	>	mpf_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/mpf_pol0->GetParameter(0) << " #pm " << std::setprecision(2) << mpf_pol0->GetParError(0)/(mpf_pol0->GetParameter(0)*mpf_pol0->GetParameter(0)) << "}{ (MPF)}";
537		TText *rmpf_final_mark = write_title(mpf_result.str());
538		rmpf_final_mark->SetX(0.75);
539		rmpf_final_mark->SetY(0.75);
540		rmpf_final_mark->SetTextSize(0.03);
541		rmpf_final_mark->Draw();
542	+	DrawPrelim();
543	+	float MPFResult=1/mpf_pol0->GetParameter(0);
544	+	float MPFResultError=mpf_pol0->GetParError(0)/(mpf_pol0->GetParameter(0)*mpf_pol0->GetParameter(0));
545
546		stringstream filename2;
547		filename2 << "Extrapolation/FINAL_RESULT_MPF";
#	Line 358 | Line 551 | void new_data_mc_agreement_2d() {
551		RABS_FinalGraph->SetMarkerStyle(21);
552		RABS_FinalGraph->SetMarkerStyle(21);
553		RABS_FinalGraph->SetMarkerColor(kRed);
554	<	RABS_FinalGraph->Fit("pol0");
554	>	RABS_FinalGraph->Fit("pol0","QE"); // quiet, use minos
555	>	RABS_FinalGraph->Draw("AP");
556	>	RABS_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
557	>	RABS_FinalGraph->GetXaxis()->CenterTitle();
558	>	RABS_FinalGraph->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
559	>	RABS_FinalGraph->GetYaxis()->CenterTitle();
560		RABS_FinalGraph->Draw("AP");
561		TF1 *rabs_pol0=(TF1*)RABS_FinalGraph->GetFunction("pol0");
562		stringstream rabs_result;
563	<	rabs_result << "C_{abs}= " << rabs_pol0->GetParameter(0) << " #pm " << rabs_pol0->GetParError(0);
563	>	rabs_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/rabs_pol0->GetParameter(0) << " #pm " << std::setprecision(2) << rabs_pol0->GetParError(0)/(rabs_pol0->GetParameter(0)*rabs_pol0->GetParameter(0)) << "}{ (R_{abs})}";;
564		TText *rabs_final_mark = write_title(rabs_result.str());
565		rabs_final_mark->SetX(0.75);
566		rabs_final_mark->SetY(0.75);
567		rabs_final_mark->SetTextSize(0.03);
568		rabs_final_mark->Draw();
569	+	DrawPrelim();
570	+
571	+	float RabsResult=1/rabs_pol0->GetParameter(0);
572	+	float RabsResultError=rabs_pol0->GetParError(0)/(rabs_pol0->GetParameter(0)*rabs_pol0->GetParameter(0));
573
574		filename2.str("");
575		filename2 << "Extrapolation/FINAL_RESULT_RABS";
576		CompleteSave(can,filename2.str());
577
578	+	can->cd();
579	+	MPF_FaceValueAtPoint3->SetMarkerStyle(21);
580	+	MPF_FaceValueAtPoint3->SetMarkerColor(kBlue);
581	+	MPF_FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
582	+	MPF_FaceValueAtPoint3->Draw("AP");
583	+	MPF_FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
584	+	MPF_FaceValueAtPoint3->GetXaxis()->CenterTitle();
585	+	MPF_FaceValueAtPoint3->GetYaxis()->SetTitle("<MPF>__{data}/<MPF>_{mc} for #alpha<0.3");
586	+	MPF_FaceValueAtPoint3->GetYaxis()->CenterTitle();
587	+	MPF_FaceValueAtPoint3->Draw("AP");
588	+	TF1 *faceval_pol0=(TF1*)MPF_FaceValueAtPoint3->GetFunction("pol0");
589	+	faceval_pol0->SetLineColor(kBlue);
590	+	faceval_pol0->SetLineWidth(1);
591	+	faceval_pol0->Draw("same");
592	+	TF1 *faceval_pol0_complete = new TF1("faceval_pol0_complete","[0]+[1]*x");
593	+	faceval_pol0_complete->SetParameters(faceval_pol0->GetParameters());
594	+	stringstream faceval_result;
595	+	faceval_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/faceval_pol0->GetParameter(0) << " #pm " << std::setprecision(2) << faceval_pol0->GetParError(0)/(faceval_pol0->GetParameter(0)*faceval_pol0->GetParameter(0)) << "}{ (MPF at #alpha<0.3)}";;
596	+	TText *faceval_final_mark = write_title(faceval_result.str());
597	+	faceval_final_mark->SetX(0.75);
598	+	faceval_final_mark->SetY(0.75);
599	+	faceval_final_mark->SetTextSize(0.03);
600	+	faceval_final_mark->Draw();
601	+	DrawPrelim();
602	+
603	+	float FaceValResult=1/faceval_pol0->GetParameter(0);
604	+	float FaceValResultError=faceval_pol0->GetParError(0)/(faceval_pol0->GetParameter(0)*faceval_pol0->GetParameter(0));
605	+
606	+	filename2.str("");
607	+	filename2 << "Extrapolation/FINAL_RESULT_MPF_FaceValp3";
608	+	CompleteSave(can,filename2.str());
609	+
610	+	can->cd();
611	+	RABS_FaceValueAtPoint3->SetMarkerStyle(21);
612	+	RABS_FaceValueAtPoint3->SetMarkerColor(kRed);
613	+	RABS_FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
614	+	RABS_FaceValueAtPoint3->Draw("AP");
615	+	RABS_FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
616	+	RABS_FaceValueAtPoint3->GetXaxis()->CenterTitle();
617	+	RABS_FaceValueAtPoint3->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc} for #alpha<0.3");
618	+	RABS_FaceValueAtPoint3->GetYaxis()->CenterTitle();
619	+	RABS_FaceValueAtPoint3->Draw("AP");
620	+	TF1 *faceval_pol0RABS=(TF1*)RABS_FaceValueAtPoint3->GetFunction("pol0");
621	+	TF1 *faceval_pol0RABS_complete = new TF1("faceval_pol0_complete","[0]+[1]*x");
622	+	faceval_pol0RABS_complete->SetParameters(faceval_pol0RABS->GetParameters());
623	+	faceval_pol0RABS->SetLineWidth(1);
624	+	faceval_pol0RABS->SetLineColor(kRed);
625	+	faceval_pol0RABS_complete->Draw("same");
626	+	stringstream RABSfaceval_result;
627	+	RABSfaceval_result << "#splitline{C_{abs}= " << std::setprecision(5) << 1/faceval_pol0RABS->GetParameter(0) << " #pm " << std::setprecision(2) << faceval_pol0RABS->GetParError(0)/(faceval_pol0RABS->GetParameter(0)*faceval_pol0RABS->GetParameter(0)) << "}{ (R_{abs} at #alpha<0.3)}";;
628	+	TText *RABSfaceval_final_mark = write_title(RABSfaceval_result.str());
629	+	RABSfaceval_final_mark->SetX(0.75);
630	+	RABSfaceval_final_mark->SetY(0.75);
631	+	RABSfaceval_final_mark->SetTextSize(0.03);
632	+	RABSfaceval_final_mark->Draw();
633	+	DrawPrelim();
634	+
635	+	float RABSFaceValResult=1/faceval_pol0RABS->GetParameter(0);
636	+	float RABSFaceValResultError=faceval_pol0RABS->GetParError(0)/(faceval_pol0RABS->GetParameter(0)*faceval_pol0RABS->GetParameter(0));
637	+
638	+	filename2.str("");
639	+	filename2 << "Extrapolation/FINAL_RESULT_RABS_FaceValp3";
640	+	CompleteSave(can,filename2.str());
641	+
642	+	cout << "FINAL RESULTS: " << endl;
643	+	cout << "MPF : " << MPFResult  << " +/- " << MPFResultError  << endl;
644	+	cout << "Rabs: " << RabsResult << " +/- " << RabsResultError << endl;
645	+	cout << "Face value at 0.3 :" << FaceValResult << " +/- " << FaceValResultError << "  (MPF) " << endl;
646	+	cout << "Face value at 0.3 :" << RABSFaceValResult << " +/- " << RABSFaceValResultError << "  (RABS) " << endl;
647	+
648	+	delete can;
649	+
650	+
651	+	}
652	+
653	+
654	+	void DoPUStudy(string identifier) {
655	+
656	+	float numVtxcuts[7]={0,10,15,20};
657	+
658	+	TCanvas *can = new TCanvas("can","can");
659	+	TH1F *malpha[8];
660	+	TH1F *dalpha[8];
661	+	cout << identifier << ";";
662	+	cout << endl;
663	+	for(int i=1;i<5;i++) {
664	+	stringstream sSpecialCut;
665	+	if(i<4) {
666	+	sSpecialCut << "numVtx>" << numVtxcuts[i-1] << " && numVtx<" << numVtxcuts[i];
667	+	cout << numVtxcuts[i-1] << "  < nVtx < " << numVtxcuts[i] << ";";
668	+	} else {
669	+	sSpecialCut << "numVtx>" << numVtxcuts[i-1];
670	+	cout << numVtxcuts[i-1] << ";";
671	+	}
672	+
673	+	TCut SpecialCut(sSpecialCut.str().c_str());
674	+
675	+
676	+	malpha[i] = allsamples.Draw("malpha"+any2string(i),"mpf",1000,0,10,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,mc,luminosity);
677	+	dalpha[i] = allsamples.Draw("dalpha"+any2string(i),"mpf",1000,0,10,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,data,luminosity);
678	+
679	+
680	+	float a=dalpha[i]->GetMean();
681	+	float b=malpha[i]->GetMean();
682	+	float da=dalpha[i]->GetMeanError();
683	+	float db=malpha[i]->GetMeanError();
684	+	float factor = a / b;
685	+	float error = TMath::Sqrt(  (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
686	+
687	+	cout << malpha[i]->GetMean() << ";" << malpha[i]->GetMeanError() << ";" << dalpha[i]->GetMean() << ";" << dalpha[i]->GetMeanError() << ";" << malpha[i]->Integral()<<";"<<dalpha[i]->Integral() << ";"<<factor << ";" << error <<  endl;
688	+	malpha[i]->SetLineColor(i+1);
689	+	if(i==0) malpha[i]->Draw("histo");
690	+	else malpha[i]->Draw("histo,same");
691	+	}
692	+
693	+	CompleteSave(can,"AlphaOverview_"+identifier);
694	+	delete can;
695	+
696	+	/*  for(int i=1;i<8;i++) {
697	+	delete malpha[i];
698	+	delete dalpha[i];
699	+	}
700	+	*/
701	+
702	+	}
703	+
704	+	void ScenarioComparison() {
705	+	//very dumb way to do this.
706	+	TGraphAsymmErrors *gr[5];
707	+	TF1 *fit[5];
708	+	bool dofit=false;
709	+
710	+	TCanvas *can = new TCanvas("can","can",500,500);
711	+
712	+	TLegend *leg = new TLegend(0.2,0.2,0.5,0.45);
713	+	leg->SetBorderSize(0);
714	+	leg->SetFillColor(kWhite);
715	+
716	+	for(int i=0;i<5;i++) {
717	+	gr[i] = new TGraphAsymmErrors();
718	+	double x[4]   = {5,12.5,17.5,20};
719	+	double dxl[4] = {5,2.5,2.5,2.5};
720	+	double dxh[4] = {5,2.5,2.5,10};
721	+
722	+	string name="";
723	+	string color="";
724	+
725	+	if(i==0) {
726	+	double y[4]   = {1.008,1.013,1.003,1.007};
727	+	double dyl[4] = {0.011,0.009,0.013,0.020};
728	+	double dyh[4] = {0.011,0.009,0.013,0.020};
729	+	name="30/30";
730	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
731	+	color="#a01d99";
732	+	}
733	+	if(i==1) {
734	+	double y[4]   = {1.008+0.001,1.013+0.001,1.003+0.001,1.007+0.001};
735	+	double dyl[4] = {0.011,0.009,0.013,0.020};
736	+	double dyh[4] = {0.011,0.009,0.013,0.020};
737	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
738	+	name="30/15";
739	+	color="#2464bc";
740	+	}
741	+	if(i==2) {
742	+	double y[4]   = {1.014,1.026,1.014,1.059};
743	+	double dyl[4] = {0.011,0.010,0.014,0.023};
744	+	double dyh[4] = {0.011,0.010,0.014,0.023};
745	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
746	+	name="15/15";
747	+	color="#f49230";
748	+	}
749	+	if(i==3) {
750	+	double y[4]   = {1.015,1.020,1.008,1.050};
751	+	double dyl[4] = {0.011,0.010,0.014,0.022};
752	+	double dyh[4] = {0.011,0.010,0.014,0.022};
753	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
754	+	name="20/15";
755	+	color="#c72f3c";
756	+	}
757	+	if(i==4) {
758	+	double y[4]   = {1.015+0.001,1.020+0.001,1.008+0.001,1.050+0.001};
759	+	double dyl[4] = {0.011,0.010,0.014,0.022};
760	+	double dyh[4] = {0.011,0.010,0.014,0.022};
761	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
762	+	name="20/20";
763	+	color="#7cb433";
764	+	}
765	+
766	+
767	+	gr[i]->SetTitle();
768	+	gr[i]->GetXaxis()->SetTitle("N(Vtx)");
769	+	gr[i]->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
770	+	gr[i]->GetXaxis()->CenterTitle();
771	+	gr[i]->GetYaxis()->CenterTitle();
772	+	gr[i]->SetName(name.c_str());
773	+	gr[i]->SetMarkerSize(0.00001);
774	+	gr[i]->GetYaxis()->SetRangeUser(0.90,1.10);
775	+	gr[i]->SetLineColor(TColor::GetColor(color.c_str()));
776	+	gr[i]->SetMarkerColor(TColor::GetColor(color.c_str()));
777	+	if(dofit) {
778	+	gr[i]->Fit("pol1");
779	+	fit[i] = (TF1*)gr[i]->GetFunction("pol1");
780	+	fit[i]->SetLineColor(TColor::GetColor(color.c_str()));
781	+	}
782	+	leg->AddEntry(gr[i],name.c_str(),"l");
783	+	if(i==0) gr[i]->Draw("AP*");
784	+	else gr[i]->Draw("P");
785	+	}
786	+	leg->Draw();
787	+	DrawPrelim();
788	+
789	+	CompleteSave(can,"ScenarioComparison");
790	+	delete can;
791	+
792	+	}
793	+
794	+	void ScenarioComparisonInclusive() {
795	+	//dumbest way ever to do this. but ok we only need to do it once.
796	+	TGraphAsymmErrors *gr[5];
797	+	TF1 *fit[5];
798	+	bool dofit=true;
799	+
800	+	TCanvas *can = new TCanvas("can","can",500,500);
801	+
802	+	TLegend *leg = new TLegend(0.2,0.2,0.5,0.45);
803	+	leg->SetBorderSize(0);
804	+	leg->SetFillColor(kWhite);
805	+
806	+	for(int i=0;i<5;i++) {
807	+	gr[i] = new TGraphAsymmErrors();
808	+	double x[4]   = {5,12.5,17.5,20};
809	+	double dxl[4] = {5,2.5,2.5,2.5};
810	+	double dxh[4] = {5,2.5,2.5,10};
811	+
812	+	string name="";
813	+	string color="";
814	+
815	+	if(i==0) {
816	+	double y[4]   = {1.01103,1.01496,1.01486,1.02785};
817	+	double dyl[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
818	+	double dyh[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
819	+	name="30/30";
820	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
821	+	color="#a01d99";
822	+	}
823	+	if(i==1) {
824	+	double y[4]   = {1.01103+0.001,1.01496+0.001,1.01486+0.001,1.02785+0.001};
825	+	double dyl[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
826	+	double dyh[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
827	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
828	+	name="30/15";
829	+	color="#2464bc";
830	+	}
831	+	if(i==2) {
832	+	double y[4]   = {1.01575,1.02879,1.0325,1.05456};
833	+	double dyl[4] = {0.00288281,0.00278536,0.00367451,0.00590421};
834	+	double dyh[4] = {0.00288281,0.00278536,0.00367451,0.00590421};
835	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
836	+	name="15/15";
837	+	color="#f49230";
838	+	}
839	+	if(i==3) {
840	+	double y[4]   = {1.01221,1.01922,1.02396,1.04677};
841	+	double dyl[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
842	+	double dyh[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
843	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
844	+	name="20/15";
845	+	color="#c72f3c";
846	+	}
847	+	if(i==4) {
848	+	double y[4]   = {1.01221+0.001,1.01922+0.001,1.02396+0.001,1.04677+0.001};
849	+	double dyl[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
850	+	double dyh[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
851	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
852	+	name="20/20";
853	+	color="#7cb433";
854	+	}
855	+
856	+
857	+	gr[i]->SetTitle();
858	+	gr[i]->GetXaxis()->SetTitle("N(Vtx)");
859	+	gr[i]->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
860	+	gr[i]->GetXaxis()->CenterTitle();
861	+	gr[i]->GetYaxis()->CenterTitle();
862	+	gr[i]->SetName(name.c_str());
863	+	gr[i]->SetMarkerSize(0.00001);
864	+	gr[i]->GetYaxis()->SetRangeUser(0.90,1.10);
865	+	gr[i]->SetLineColor(TColor::GetColor(color.c_str()));
866	+	gr[i]->SetMarkerColor(TColor::GetColor(color.c_str()));
867	+	if(dofit) {
868	+	gr[i]->Fit("pol1");
869	+	fit[i] = (TF1*)gr[i]->GetFunction("pol1");
870	+	fit[i]->SetLineColor(TColor::GetColor(color.c_str()));
871	+	}
872	+	leg->AddEntry(gr[i],name.c_str(),"l");
873	+	if(i==0) gr[i]->Draw("AP*");
874	+	else gr[i]->Draw("P");
875	+	}
876	+	leg->Draw();
877	+	DrawPrelim();
878
879	+	CompleteSave(can,"ScenarioComparisonInclusive");
880		delete can;
881
882	+	}
883	+
884	+	void compare_selection(string identifier) {
885	+	bool recognized_scenario=false;
886	+
887	+	cout << "Running with identifier " << identifier << endl;
888	+	if(identifier=="Zb3010") {
889	+	LeadingB=TCut ("Zb3010_bTagProbCSVBP[0]>0.898");
890	+	EtaB=TCut("abs(Zb3010_pfJetGoodEta[0])<1.3");
891	+	PhiZcut=TCut("abs(Zb3010_pfJetDphiZ[0])>2.7");
892	+	recognized_scenario=true;
893	+	}
894	+
895	+	if(identifier=="Zb3010") {
896	+	LeadingB=TCut ("Zb3010_bTagProbCSVBP[0]>0.898");
897	+	EtaB=TCut("abs(Zb3010_pfJetGoodEta[0])<1.3");
898	+	PhiZcut=TCut("abs(Zb3010_pfJetDphiZ[0])>2.7");
899	+	recognized_scenario=true;
900	+	}
901	+
902	+	if(identifier=="Zb40") {
903	+	LeadingB=TCut ("Zb40_bTagProbCSVBP[0]>0.898");
904	+	EtaB=TCut("abs(Zb40_pfJetGoodEta[0])<1.3");
905	+	PhiZcut=TCut("abs(Zb40_pfJetDphiZ[0])>2.7");
906	+	recognized_scenario=true;
907	+	}
908	+
909	+	if(identifier=="Zb301030") {
910	+	LeadingB=TCut ("Zb301030_bTagProbCSVBP[0]>0.898");
911	+	EtaB=TCut("abs(Zb301030_pfJetGoodEta[0])<1.3");
912	+	PhiZcut=TCut("abs(Zb301030_pfJetDphiZ[0])>2.7");
913	+	recognized_scenario=true;
914	+	}
915	+
916	+	if(identifier=="Zb1530") {
917	+	LeadingB=TCut ("Zb1530_bTagProbCSVBP[0]>0.898");
918	+	EtaB=TCut("abs(Zb1530_pfJetGoodEta[0])<1.3");
919	+	PhiZcut=TCut("abs(Zb1530_pfJetDphiZ[0])>2.7");
920	+	recognized_scenario=true;
921	+	}
922	+
923	+	if(identifier=="Zb301010") {
924	+	LeadingB=TCut ("Zb301010_bTagProbCSVBP[0]>0.898");
925	+	EtaB=TCut("abs(Zb301010_pfJetGoodEta[0])<1.3");
926	+	PhiZcut=TCut("abs(Zb301010_pfJetDphiZ[0])>2.7");
927	+	recognized_scenario=true;
928	+	}
929	+
930	+	if(identifier=="Zb1510") {
931	+	LeadingB=TCut ("Zb1510_bTagProbCSVBP[0]>0.898");
932	+	EtaB=TCut("abs(Zb1510_pfJetGoodEta[0])<1.3");
933	+	PhiZcut=TCut("abs(Zb1510_pfJetDphiZ[0])>2.7");
934	+	recognized_scenario=true;
935	+	}
936	+
937	+	if(identifier=="Zb301010") {
938	+	LeadingB=TCut ("Zb301010_bTagProbCSVBP[0]>0.898");
939	+	EtaB=TCut("abs(Zb301010_pfJetGoodEta[0])<1.3");
940	+	PhiZcut=TCut("abs(Zb301010_pfJetDphiZ[0])>2.7");
941	+	recognized_scenario=true;
942	+	}
943	+
944	+	if(identifier=="ZbMikko") {
945	+	LeadingB=TCut ("ZbMikko_bTagProbCSVBP[0]>0.898 && ZbMikko__IsClean");
946	+	EtaB=TCut("abs(ZbMikko_pfJetGoodEta[0])<1.3 && ZbMikko__IsClean");
947	+	PhiZcut=TCut("abs(ZbMikko_pfJetDphiZ[0])>2.7 && ZbMikko__IsClean");
948	+	recognized_scenario=true;
949	+	}
950	+
951	+	if(identifier=="Zb3010_SecEta3") {
952	+	LeadingB=TCut ("Zb3010_SecEta3_bTagProbCSVBP[0]>0.898");
953	+	EtaB=TCut("abs(Zb3010_SecEta3_pfJetGoodEta[0])<1.3");
954	+	PhiZcut=TCut("abs(Zb3010_SecEta3_pfJetDphiZ[0])>2.7");
955	+	recognized_scenario=true;
956	+	}
957	+
958	+	if(identifier=="Zb3010_SecEta5") {
959	+	LeadingB=TCut ("Zb3010_SecEta5_bTagProbCSVBP[0]>0.898");
960	+	EtaB=TCut("abs(Zb3010_SecEta5_pfJetGoodEta[0])<1.3");
961	+	PhiZcut=TCut("abs(Zb3010_SecEta5_pfJetDphiZ[0])>2.7");
962	+	recognized_scenario=true;
963	+	}
964	+
965	+	if(identifier=="Zb3010_p5Clean") {
966	+	LeadingB=TCut ("Zb3010_p5Clean_bTagProbCSVBP[0]>0.898 && Zb3010_p5Clean_IsClean");
967	+	EtaB=TCut("abs(Zb3010_p5Clean_pfJetGoodEta[0])<1.3 && Zb3010_p5Clean_IsClean");
968	+	PhiZcut=TCut("abs(Zb3010_p5Clean_pfJetDphiZ[0])>2.7 && Zb3010_p5Clean_IsClean");
969	+	recognized_scenario=true;
970	+	}
971	+
972	+	if(identifier=="Zb3010CHS") {
973	+	LeadingB=TCut ("Zb3010CHS_bTagProbCSVBP[0]>0.898");
974	+	EtaB=TCut("abs(Zb3010CHS_pfJetGoodEta[0])<1.3");
975	+	PhiZcut=TCut("abs(Zb3010CHS_pfJetDphiZ[0])>2.7");
976	+	recognized_scenario=true;
977	+	}
978	+
979	+	assert(recognized_scenario);
980	+
981	+	cout << "Cuts: " << endl;
982	+	cout << "   " << (const char*) LeadingB << endl;
983	+	cout << "   " << (const char*) EtaB << endl;
984	+	cout << "   " << (const char*) PhiZcut << endl;
985	+
986	+	cout << endl << endl;
987	+
988	+	//    ScenarioComparison();
989	+	//    ScenarioComparisonInclusive();
990	+
991	+
992	+	//     DoPUStudy(identifier);
993	+
994	+	//   SpecialCutFlow(identifier);
995	+	//    draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",17000,30,200,"Z p_{T}","Official/"+identifier+"/Zpt__"+identifier,1);
996	+	//   draw_kin_variable(identifier+"_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",40,0,200,"Leading Jet Pt (B)","Official/"+identifier+"/LeadingJetPt__"+identifier,1);
997	+	//   draw_kin_variable(identifier+"_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",40,0,200,"Sub-Leading Jet Pt [GeV]","Official/"+identifier+"/Jet2Pt__"+identifier,1);
998	+	//   draw_kin_variable(identifier+"_pfJetGoodPt[1]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,2,"p_{T}^{J2} / p_{T}^{Z}","Official/"+identifier+"/SubLeadingJetPt_Over_ZPt__"+identifier,1);
999	+	draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>30&&pt<50"),20,0,2,"#alpha","Official/"+identifier+"/alpha_pt_30_to_50__"+identifier,1);
1000	+	draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>50&&pt<75"),20,0,2,"#alpha","Official/"+identifier+"/alpha_50_to_75__"+identifier,1);
1001	+	draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>75&&pt<125"),20,0,2,"#alpha","Official/"+identifier+"/alpha_pt_75_to_125__"+identifier,1);
1002	+	draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>125&&pt<1000"),20,0,2,"#alpha","Official/"+identifier+"/alpha_pt_125_to_1000__"+identifier,1);
1003	+	draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut,20,0,2,"#alpha","Official/"+identifier+"/alpha_full__"+identifier,1);
1004	+	//   draw_kin_variable(identifier+"_pfBJetDphiZ[0]",ZplusBsel&&LeadingB&&"pt>10&&pt<1000&&pfBJetDphiZ[0]>0",30,0,3.2,"#delta#phi (Z,b lead)","Official/"+identifier+"/DeltaPhiZBlead__"+identifier,0);
1005	+
1006	+	//   draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.3").c_str()),20,0,2,"Z+b MPF","Official/"+identifier+"/mpf_alpha_smaller_0p3___"+identifier,0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
1007	+	//   draw_kin_variable(identifier+"_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.3").c_str()),20,0,2,"p_{T} b-jet / p_{T} Z","Official/"+identifier+"/ptb_over_ptz___alpha_smaller_0p3______"+identifier,0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
1008	+
1009	+
1010	+	//   */
1011
1012		}
1013
1014	+	void GetNumberEventsInsideOutsideAlphaWindow() {
1015	+
1016	+	TCut cut =  ZplusBsel&&LeadingB&&EtaB&&PhiZcut;
1017	+	TCut in  =  TCut("Zb3010_alpha<=0.3");
1018	+	TCut out =  TCut("Zb3010_alpha>0.3");
1019
1020	+	TH1F *data_IN  = allsamples.Draw("data_IN",    "mll",1,0,150, "nothing [GeV]", "events", cut&&in, data,luminosity);
1021	+	TH1F *data_OUT = allsamples.Draw("data_OUT",   "mll",1,0,150, "nothing [GeV]", "events", cut&&out,data,luminosity);
1022	+	TH1F *mc_IN    = allsamples.Draw("mc_IN",      "mll",1,0,150, "nothing [GeV]", "events", cut&&in, mc  ,luminosity);
1023	+	TH1F *mc_OUT   = allsamples.Draw("mc_OUT",     "mll",1,0,150, "nothing [GeV]", "events", cut&&out,mc  ,luminosity);
1024	+
1025	+	dout << "  Data: In " << data_IN->Integral() << " vs out " << data_OUT->Integral() << endl;
1026	+	dout << "  MC  : In " << mc_IN->Integral()   << " vs out " << mc_OUT->Integral()   << endl;
1027	+	}
1028
1029	<	void do_basic_ZB_analysis() {
1029	>
1030	>	void do_basic_ZB_analysis(int do_inclusive) {
1031	>
1032	>	//https://twiki.cern.ch/twiki/bin/viewauth/CMS/BtagPOG#2011_Data_and_MC
1033	>	//https://twiki.cern.ch/twiki/pub/CMS/BtagPOG/SFb-ttbar_payload.txt
1034	>
1035	>	if(do_inclusive) {
1036	>	LeadingB=TCut("Zb3010_bTagProbCSVBP[0]>=-10000");
1037	>	cutWeight=TCut("(weight*(weight<1000)*(is_data+(!is_data)*((id1==id2&&id1==0)*0.95+(id1==id2&&id1==1)*0.88+(id1!=id2)*0.92)))");
1038	>	write_info(__FUNCTION__,"Inclusive selection - no btag requirement, no btag efficiency/mistag correction (weight)");
1039	>	} else {
1040	>	LeadingB=TCut("Zb3010_bTagProbCSVBP[0]>0.898");
1041	>	cutWeight=TCut("(weight*(weight<1000)*(is_data+(!is_data)*((id1==id2&&id1==0)*0.95+(id1==id2&&id1==1)*0.88+(id1!=id2)*0.92)*Zb3010_BTagWgt))");
1042	>	write_info(__FUNCTION__,"Exclusive selection - switched on btag requirement and btag efficiency/mistag correction (weight)");
1043	>	}
1044	>
1045	>	bool doquick=true;
1046	>
1047		TCanvas *zbcanvas = new TCanvas("zbcanvas","zbcanvas");
386	–	write_warning(__FUNCTION__,"TESTING NEW DATA MC AGREEMENT");
1048
1049	<	/*
1050	<	write_warning(__FUNCTION__,"Not doing b yields (already implemented)");
1051	<	//  print_all_b_yields();
1052	<	write_warning(__FUNCTION__,"Not drawing mpf variables (already implemented)");
1053	<	//  draw_mpf_vars();
1054	<	write_warning(__FUNCTION__,"Not drawing mpf variables (already implemented)");
1055	<	draw_Zb_kin_vars();
1049	>	//     compare_selection("Zb3010_p5Clean");
1050	>	//     compare_selection("Zb3010CHS");
1051	>	//     compare_selection("Zb3010");
1052	>	//     compare_selection("Zb301030");
1053	>	//     compare_selection("Zb1530");
1054	>	//     compare_selection("Zb301010");
1055	>	//     compare_selection("Zb3010_SecEta3");
1056	>	//     compare_selection("Zb3010_SecEta5");
1057	>	//     compare_selection("Zb1510");
1058	>	//     compare_selection("Zb301010");
1059	>	//     compare_selection("ZbMikko");
1060	>	//     compare_selection("Zb3010");
1061	>	//     compare_selection("Zb40");
1062	>
1063	>	if(!doquick) {
1064	>	print_all_b_yields();
1065	>	draw_mpf_vars();
1066	>	//  DrawEvilCutFlow();
1067	>
1068	>	draw_Zb_kin_vars();
1069	>
1070	>	}
1071
1072	<	data_mc_agreement_2d();
1073	<	*/
398	<	new_data_mc_agreement_2d();
1072	>	GetNumberEventsInsideOutsideAlphaWindow();
1073	>	//  new_data_mc_agreement_2d();
1074
1075	+	write_error(__FUNCTION__,"Need to implement switch for inclusive result; inclusive result needs to have btag efficiency reweighting deactivated!");
1076		delete zbcanvas;
1077		}
1078	+
1079	+
1080	+
1081	+
1082	+	//*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*  DELETE EVERYTHING BELOW THIS LINE /*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*//
1083	+
1084	+
1085	+	const char* concatenate(string bla1, string bla2) {
1086	+	stringstream bla;
1087	+	bla << bla1 << bla2;
1088	+	return bla.str().c_str();
1089	+	}
1090	+
1091	+	void SpecialCutFlow(string identifier) {
1092	+	stringstream MegaCut;
1093	+
1094	+
1095	+
1096	+	MegaCut << "   1*(" << (const char*) (ZplusBsel&&TCut(concatenate(identifier,"_pfJetGoodNumBtag>0"))) << ")";
1097	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB) << ")";
1098	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB) << ")";
1099	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut) << ")";
1100	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000").c_str())) << ")";
1101	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.3").c_str())) << ")";
1102	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.2").c_str())) << ")";
1103	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.15").c_str())) << ")";
1104	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.1").c_str())) << ")";
1105	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.05").c_str())) << ")";
1106	+	MegaCut << " ";
1107	+
1108	+	TCanvas *can = new TCanvas("can","can");
1109	+	can->SetLogy(1);
1110	+	TCut basecut(ZplusBsel);
1111	+
1112	+
1113	+
1114	+	TLegend *leg = allsamples.allbglegend();
1115	+
1116	+
1117	+	TH1F *data  = allsamples.Draw("data",    MegaCut.str(),11,-0.5,10.5, "", "events", basecut,data,luminosity);
1118	+	THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
1119	+
1120	+	float runningsum=0;
1121	+	for(int i=data->GetNbinsX();i>0;i--) {
1122	+	runningsum+=data->GetBinContent(i);
1123	+	data->SetBinContent(i,runningsum);
1124	+	}
1125	+
1126	+	float minimum=data->GetMaximum();
1127	+
1128	+	TH1F* h;
1129	+	TIter nextOF(mcm.GetHists());
1130	+
1131	+	float mcsum=0;
1132	+	float zb=0;
1133	+	while ( h = (TH1F*)nextOF() ) {
1134	+	float runningsum=0;
1135	+	minimum=data->GetMaximum();//done deliberately at each step so get the minimum of the last (!) contribution
1136	+	for(int i=h->GetNbinsX();i>0;i--) {
1137	+	runningsum+=h->GetBinContent(i);
1138	+	h->SetBinContent(i,runningsum);
1139	+	if(runningsum<minimum&&runningsum>0) minimum=runningsum;
1140	+	}
1141	+	if(Contains(h->GetName(),"Z_b_")) {
1142	+	zb=h->GetBinContent(h->GetNbinsX());
1143	+	}
1144	+	mcsum+=h->GetBinContent(h->GetNbinsX());
1145	+	}
1146	+
1147	+
1148	+
1149	+	data->SetMinimum(0.05*minimum);
1150	+	can->cd(1)->SetBottomMargin(0.25);
1151	+	data->GetXaxis()->SetBinLabel(1,"Pre-selection");
1152	+	data->GetXaxis()->SetBinLabel(2,"Contains b jet");
1153	+	data->GetXaxis()->SetBinLabel(3,"Leading jet is b-jet");
1154	+	data->GetXaxis()->SetBinLabel(4,"|#eta_{b}|<1.3 ");
1155	+	data->GetXaxis()->SetBinLabel(5,"|#delta#phi(b,Z)|>2.7");
1156	+	data->GetXaxis()->SetBinLabel(6,"10<p_{T}^{Z}<1000 GeV");
1157	+	data->GetXaxis()->SetBinLabel(7,"#alpha < 0.3");
1158	+	data->GetXaxis()->SetBinLabel(8,"#alpha < 0.2");
1159	+	data->GetXaxis()->SetBinLabel(9,"#alpha < 0.15");
1160	+	data->GetXaxis()->SetBinLabel(10,"#alpha < 0.1");
1161	+	data->GetXaxis()->SetBinLabel(11,"#alpha < 0.05");
1162	+	data->GetXaxis()->LabelsOption("v");
1163	+
1164	+	stringstream purityinfo;
1165	+	purityinfo << "Purity: " << std::setprecision(3) << 100*zb/mcsum << " %";
1166	+	stringstream neventsinfo;
1167	+	neventsinfo << "Nevents: " << data->GetBinContent(data->GetNbinsX());
1168	+	TH1F *crap = new TH1F("crap","",1,0,1);
1169	+	crap->SetLineColor(kWhite);
1170	+	leg->AddEntry(crap,purityinfo.str().c_str(),"l");
1171	+	leg->AddEntry(crap,neventsinfo.str().c_str(),"l");
1172	+
1173	+	data->Draw("e1");
1174	+	mcm.Draw("same");
1175	+	data->Draw("same,e1,axis");
1176	+	data->Draw("same,e1");
1177	+	//  data->Draw("same,TEXT");
1178	+
1179	+	leg->Draw();
1180	+
1181	+	DrawPrelim();
1182	+
1183	+	CompleteSave(can,"CutFlow___"+identifier);
1184	+
1185	+	delete crap;
1186	+	delete data;
1187	+
1188	+
1189	+	}

Diff Legend

–	Removed lines
+	Added lines
<	Changed lines
>	Changed lines