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.4 by buchmann, Fri Aug 10 11:47:53 2012 UTC vs.
Revision 1.9 by buchmann, Wed Nov 14 10:46:45 2012 UTC

#	Line 24 | Line 24 | namespace ZbData {
24		vector<float> data_over_mc;
25
26		TCut ZplusBsel("pt1>20&&pt2>20&&mll>60&&mll<120");
27	<	TCut LeadingB("Zb30_bTagProbCSVBP[0]>0.679");
28	<	TCut EtaB("abs(Zb30_pfJetGoodEta[0])<1.3");
29	<	TCut PhiZcut("Zb30_pfJetDphiZ[0]>2.7");
27	>	TCut LeadingB("Zb3010_bTagProbCSVBP[0]>0.898");
28	>	TCut EtaB("abs(Zb3010_pfJetGoodEta[0])<1.3");
29	>	TCut PhiZcut("abs(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	–	DrawPrelim();
163	–	CompleteSave(agrcan,saveas.str());
164	–	//    agrcan->SaveAs("tester.png");
165	–	}
166	–
167	–	delete agrcan;
168	–
169	–	}
170	–
112		void print_all_b_yields() {
113		cout << "Basic selection with a b jet" << endl;
114	<	print_yield(ZplusBsel&&"Zb30_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 180 | Line 121 | void print_all_b_yields() {
121		cout << "10<ptZ<1000 GeV" << endl;
122		print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000");
123		cout << "#alpha < 0.3" << endl;
124	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.3");
124	>	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.3");
125		cout << "#alpha < 0.2" << endl;
126	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.2");
126	>	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.2");
127		cout << "#alpha < 0.15" << endl;
128	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.15");
128	>	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.15");
129		cout << "#alpha < 0.1" << endl;
130	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.1");
131	<	cout << "#alpha < 0.05" << endl;
132	<	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.05");
130	>	print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.1");
131	>	//  cout << "#alpha < 0.05" << endl;
132	>	//  print_yield(ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"(Zb3010_alpha)<0.05");
133		}
134
135	+	void DrawEvilCutFlow() {
136	+	stringstream MegaCut;
137	+
138	+	MegaCut << "   1*(" << (const char*) (ZplusBsel&&TCut("Zb3010_pfJetGoodNumBtag>0")) << ")";
139	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB) << ")";
140	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB) << ")";
141	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut) << ")";
142	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000")) << ")";
143	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.3")) << ")";
144	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.2")) << ")";
145	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.15")) << ")";
146	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.1")) << ")";
147	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut("pt>10&&pt<1000&&Zb3010_alpha<0.05")) << ")";
148	+	MegaCut << " ";
149	+
150	+	TCanvas *can = new TCanvas("can","can");
151	+	can->SetLogy(1);
152	+	TCut basecut("mll>6||mll<10");
153	+
154	+
155	+
156	+	TLegend *leg = allsamples.allbglegend();
157	+
158	+
159	+	TH1F *data  = allsamples.Draw("data",    MegaCut.str(),11,-0.5,10.5, "", "events", basecut,data,luminosity);
160	+	THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
161	+
162	+	float runningsum=0;
163	+	for(int i=data->GetNbinsX();i>0;i--) {
164	+	runningsum+=data->GetBinContent(i);
165	+	data->SetBinContent(i,runningsum);
166	+	}
167	+
168	+	float minimum=data->GetMaximum();
169	+
170	+	TH1F* h;
171	+	TIter nextOF(mcm.GetHists());
172	+	while ( h = (TH1F*)nextOF() ) {
173	+	float runningsum=0;
174	+	minimum=data->GetMaximum();//done deliberately at each step so get the minimum of the last (!) contribution
175	+	for(int i=h->GetNbinsX();i>0;i--) {
176	+	runningsum+=h->GetBinContent(i);
177	+	h->SetBinContent(i,runningsum);
178	+	if(runningsum<minimum&&runningsum>0) minimum=runningsum;
179	+	}
180	+	}
181	+
182	+
183	+
184	+	data->SetMinimum(0.05*minimum);
185	+	can->cd(1)->SetBottomMargin(0.25);
186	+	data->GetXaxis()->SetBinLabel(1,"Pre-selection");
187	+	data->GetXaxis()->SetBinLabel(2,"Contains b jet");
188	+	data->GetXaxis()->SetBinLabel(3,"Leading jet is b-jet");
189	+	data->GetXaxis()->SetBinLabel(4,"|#eta_{b}|<1.3 ");
190	+	data->GetXaxis()->SetBinLabel(5,"|#delta#phi(b,Z)|>2.7");
191	+	data->GetXaxis()->SetBinLabel(6,"10<p_{T}^{Z}<1000 GeV");
192	+	data->GetXaxis()->SetBinLabel(7,"#alpha < 0.3");
193	+	data->GetXaxis()->SetBinLabel(8,"#alpha < 0.2");
194	+	data->GetXaxis()->SetBinLabel(9,"#alpha < 0.15");
195	+	data->GetXaxis()->SetBinLabel(10,"#alpha < 0.1");
196	+	data->GetXaxis()->SetBinLabel(11,"#alpha < 0.05");
197	+	data->GetXaxis()->LabelsOption("v");
198	+
199	+	data->Draw("e1");
200	+	mcm.Draw("same");
201	+	data->Draw("same,e1,axis");
202	+	data->Draw("same,e1");
203	+	//  data->Draw("same,TEXT");
204	+
205	+	leg->Draw();
206	+
207	+	DrawPrelim();
208	+
209	+	CompleteSave(can,"CutFlow");
210	+	}
211	+
212		void draw_Zb_kin_vars() {
213	<	draw_kin_variable("Zb30_pfBJetDphiZ",ZplusBsel&&LeadingB,100,0,3.2,"#delta#phi (Z,b lead)","DeltaPhiZBlead",1);
214	<	draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",50,0,200,"Z p_{T}","Zpt",1);
215	<	draw_kin_variable("Zb30_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",50,0,200,"Sub-Leading Jet Pt","Jet2Pt",1);
216	<	draw_kin_variable("Zb30_pfJetGoodPtlBtag[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",50,0,200,"Leading Jet Pt (B)","LeadingJetPt",1);
217	<	draw_kin_variable("Zb30_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,2,"#alpha","alpha",1);
213	>	draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Z p_{T}","Official/Zpt",1);
214	>	draw_kin_variable("Zb3010_pfJetGoodPt[1]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Sub-Leading Jet Pt","Official/Jet2Pt",1);
215	>	draw_kin_variable("Zb3010_pfJetGoodPt[0]",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",25,0,200,"Leading Jet Pt (B)","Official/LeadingJetPt",1);
216	>	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);
217	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<50",   40,0,2,"#alpha","Official/alpha_pt_30_to_50",1);
218	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>50&&pt<75",  40,0,2,"#alpha","Official/alpha_50_to_75",1);
219	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>75&&pt<125", 40,0,2,"#alpha","Official/alpha_pt_75_to_125",1);
220	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>125&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_125_to_1000",1);
221	>	draw_kin_variable("Zb3010_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>30&&pt<1000",40,0,2,"#alpha","Official/alpha_pt_30_to_1000",1);
222	>
223	>	draw_kin_variable("pfBJetDphiZ[0]",ZplusBsel&&LeadingB&&"pt>10&&pt<1000&&pfBJetDphiZ[0]>0",30,0,3.2,"#delta#phi (Z,b lead)","DeltaPhiZBlead",0);
224	>
225	>	/*
226	>	draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,200,"p_{T}^{l1}","Lep/pt1",1);
227	>	draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==0",20,0,200,"p_{T}^{e1}","Lep/pt1_e",1);
228	>	draw_kin_variable("pt1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id1==1",20,0,200,"p_{T}^{#mu1}","Lep/pt1_m",1);
229	>	draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,0,200,"p_{T}^{l2}","Lep/pt2",1);
230	>	draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==0",20,0,200,"p_{T}^{e2}","Lep/pt2_e",1);
231	>	draw_kin_variable("pt2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000&&id2==1",20,0,200,"p_{T}^{#mu2}","Lep/pt2_m",1);
232	>	draw_kin_variable("eta1",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,-3.1,3.1,"#eta^{l1}","Lep/eta1",1);
233	>	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);
234	>	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);
235	>	draw_kin_variable("eta2",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",20,-3.1,3.1,"#eta^{l1}","Lep/eta2",1);
236	>	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);
237	>	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);
238	>	draw_kin_variable("numVtx",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",30,0,30,"N(Vtx)","Lep/nVtx",1);
239	>	*/
240		}
241
242		void draw_mpf_vars() {
243		ZbData::data_over_mc.clear();
244	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.3",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p3",0,"#alpha<0.3, 10 GeV < p_{T}^{Z} < 1000 GeV");
245	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.2",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p2",0,"#alpha<0.2, 10 GeV < p_{T}^{Z} < 1000 GeV");
246	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.15",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p15",0,"#alpha<0.15, 10 GeV < p_{T}^{Z} < 1000 GeV");
247	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.1",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p1",0,"#alpha<0.1, 10 GeV < p_{T}^{Z} < 1000 GeV");
248	<	draw_kin_variable("mpf",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_alpha<0.05",20,0,2,"Z+b MPF","mpf_alpha_smaller_0p05",0,"#alpha<0.05, 10 GeV < p_{T}^{Z} < 1000 GeV");
249	<
250	<	draw_kin_variable("Zb30_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_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");
251	<	draw_kin_variable("Zb30_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_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");
252	<	draw_kin_variable("Zb30_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_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");
253	<	draw_kin_variable("Zb30_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_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");
254	<	draw_kin_variable("Zb30_pfJetGoodPt[0]/pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000"&&"Zb30_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");
244	>	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");
245	>	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");
246	>	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");
247	>	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");
248	>	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");
249	>
250	>	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");
251	>	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");
252	>	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");
253	>	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");
254	>	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");
255		}
256
257		Value get_Zb_data_over_mc(string variable, TCut cut) {
#	Line 219 | Line 259 | Value get_Zb_data_over_mc(string variabl
259		TH1F *hdata  = allsamples.Draw("hdata",variable,1000,0,30, "#alpha [GeV]", "events", cut,data,luminosity);
260		TH1F *hmc    = allsamples.Draw("hmc" , variable,1000,0,30, "#alpha [GeV]", "events", cut,mc,  luminosity);
261
262	<	float a=hdata->GetMean();
262	>	hdata->Fit("gaus","QE");
263	>	hmc->Fit("gaus","QE");
264	>	TF1 *fdata = (TF1*) hdata->GetFunction("gaus");
265	>	TF1 *fmc   = (TF1*) hmc->GetFunction("gaus");
266	>
267	>	/*  float a=hdata->GetMean();
268		float b=hmc->GetMean();
269		float da=hdata->GetMeanError();
270	<	float db=hmc->GetMeanError();
270	>	float db=hmc->GetMeanError();*/
271	>	float a = fdata->GetParameter(1);
272	>	float da = fdata->GetParError(1);
273	>	float b = fmc->GetParameter(1);
274	>	float db = fmc->GetParError(1);
275		float factor = a / b;
276		float error = TMath::Sqrt(  (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
277		delete hdata;
278		delete hmc;
279	+
280	+	//    cout << a << "+/-" << da << "        " << b << "=/-" << db << endl;
281	+
282		return Value(factor,error);
283		}
284
#	Line 235 | Line 287 | Value get_Zb_data_over_mc(string variabl
287		void new_data_mc_agreement_2d() {
288		gStyle->SetOptFit(0);
289		TCut basecut(ZplusBsel&&LeadingB&&EtaB&&PhiZcut);
238	–
290		const int nalphacuts=5;
291		float alphacuts[nalphacuts] = {0.1,0.15,0.2,0.25,0.3};
292	+	int PointThree=-1;
293	+	for(int i=0;i<nalphacuts;i++) {
294	+	if((alphacuts[i]-0.3)<0.01) PointThree=i-1;
295	+	}
296	+	if(PointThree<0) {
297	+	write_warning(__FUNCTION__,"alpha=0.3 MUST BE ONE OF THE ALPHA VALUES!");
298	+	assert(PointThree<0);
299	+	}
300	+
301		const int nptcuts=5;
302	<	float ptcuts[nptcuts]={10,50,100,200,1000};
302	>	float ptcuts[nptcuts]={30,50,75,125,1000};
303
304		float MPF_Results[nalphacuts][nptcuts];
305		float MPF_Errors[nalphacuts][nptcuts];
#	Line 250 | Line 310 | void new_data_mc_agreement_2d() {
310		for(int ia=0;ia<nalphacuts;ia++) {
311		for(int ipt=0;ipt<nptcuts-1;ipt++) {
312		stringstream specialcut;
313	<	specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (alpha<" << alphacuts[ia] << "))";
313	>	specialcut << "((pt>" << ptcuts[ipt] << " && pt< " << ptcuts[ipt+1] << ") && (Zb3010_alpha<" << alphacuts[ia] << "))";
314	>	cout << specialcut.str() << endl;
315		Value MPF_data_over_mc = get_Zb_data_over_mc("mpf",TCut(basecut && specialcut.str().c_str()));
316	<	Value RABS_data_over_mc = get_Zb_data_over_mc("Zb30_pfJetGoodPt[0]/pt",TCut(basecut && specialcut.str().c_str()));
316	>	Value RABS_data_over_mc = get_Zb_data_over_mc("Zb3010_pfJetGoodPt[0]/pt",TCut(basecut && specialcut.str().c_str()));
317
318		MPF_Results[ia][ipt]=MPF_data_over_mc.getValue();
319		MPF_Errors[ia][ipt]=MPF_data_over_mc.getError();
#	Line 260 | Line 321 | void new_data_mc_agreement_2d() {
321		RABS_Results[ia][ipt]=RABS_data_over_mc.getValue();
322		RABS_Errors[ia][ipt]=RABS_data_over_mc.getError();
323
324	<	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;
324	>	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;
325		}
326		}
327
#	Line 271 | Line 332 | void new_data_mc_agreement_2d() {
332
333		TGraphErrors *MPF_FinalGraph = new TGraphErrors();
334		TGraphErrors *RABS_FinalGraph = new TGraphErrors();
335	+	TGraphErrors *FaceValueAtPoint3 = new TGraphErrors();
336
337		for(int ipt=0;ipt<nptcuts-1;ipt++) {
338
#	Line 280 | Line 342 | void new_data_mc_agreement_2d() {
342		TGraphErrors *mpf_gr =new TGraphErrors();
343		TGraphErrors *rabs_gr =new TGraphErrors();
344		for(int ia=0;ia<nalphacuts;ia++) {
345	<	mpf_gr->SetPoint(ia,alphacuts[ia],MPF_Results[ia][ipt]);
346	<	mpf_gr->SetPointError(ia,0,MPF_Errors[ia][ipt]);
347	<	rabs_gr->SetPoint(ia,alphacuts[ia],RABS_Results[ia][ipt]);
348	<	rabs_gr->SetPointError(ia,0,RABS_Errors[ia][ipt]);
345	>	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;
346	>	mpf_gr->SetPoint(ia,alphacuts[ia],MPF_Results[ia][ipt]);     // avoiding empty bins (nan's)
347	>	mpf_gr->SetPointError(ia,0,MPF_Errors[ia][ipt]);             // avoiding empty bins (nan's)
348	>	rabs_gr->SetPoint(ia,alphacuts[ia],RABS_Results[ia][ipt]); // avoiding empty bins (nan's)
349	>	rabs_gr->SetPointError(ia,0,RABS_Errors[ia][ipt]);         // avoiding empty bins (nan's)
350		}
351	+	cout << "About to add a point in pt (" << ptcuts[ipt]+ptcuts[ipt+1] << " ) to face value : " << MPF_Results[PointThree][ipt] << endl;
352	+	FaceValueAtPoint3->SetPoint(ipt,(ptcuts[ipt]+ptcuts[ipt+1]),MPF_Results[PointThree][ipt]);
353	+	FaceValueAtPoint3->SetPointError(ipt,0,MPF_Errors[PointThree][ipt]);
354
355		mpf_gr->SetMarkerStyle(21);
356		mpf_gr->SetMarkerColor(kBlue);
#	Line 292 | Line 358 | void new_data_mc_agreement_2d() {
358		mpf_gr->Fit("pol1","QE"); // quiet, use minos
359		mpf_gr->GetXaxis()->SetTitle("#alpha");
360		mpf_gr->GetXaxis()->CenterTitle();
361	<	mpf_gr->GetYaxis()->SetTitle("data/mc");
361	>	mpf_gr->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
362		mpf_gr->GetYaxis()->CenterTitle();
363		TF1 *mpf_pol=(TF1*)mpf_gr->GetFunction("pol1");
364		float mpf_a=mpf_pol->GetParameter(0);
#	Line 325 | Line 391 | void new_data_mc_agreement_2d() {
391		rabs_gr->Draw("AP");
392		rabs_gr->GetXaxis()->SetTitle("#alpha");
393		rabs_gr->GetXaxis()->CenterTitle();
394	<	rabs_gr->GetYaxis()->SetTitle("data/mc");
394	>	rabs_gr->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
395		rabs_gr->GetYaxis()->CenterTitle();
396
397		TF1 *rabs_pol=(TF1*)rabs_gr->GetFunction("pol1");
#	Line 362 | Line 428 | void new_data_mc_agreement_2d() {
428		MPF_FinalGraph->Draw("AP");
429		MPF_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
430		MPF_FinalGraph->GetXaxis()->CenterTitle();
431	<	MPF_FinalGraph->GetYaxis()->SetTitle("data/mc");
431	>	MPF_FinalGraph->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
432		MPF_FinalGraph->GetYaxis()->CenterTitle();
433		MPF_FinalGraph->Draw("AP");
434
#	Line 374 | Line 440 | void new_data_mc_agreement_2d() {
440		rmpf_final_mark->SetTextSize(0.03);
441		rmpf_final_mark->Draw();
442		DrawPrelim();
443	+	float MPFResult=1/mpf_pol0->GetParameter(0);
444	+	float MPFResultError=mpf_pol0->GetParError(0)/(mpf_pol0->GetParameter(0)*mpf_pol0->GetParameter(0));
445
446		stringstream filename2;
447		filename2 << "Extrapolation/FINAL_RESULT_MPF";
#	Line 387 | Line 455 | void new_data_mc_agreement_2d() {
455		RABS_FinalGraph->Draw("AP");
456		RABS_FinalGraph->GetXaxis()->SetTitle("p_{T}^{Z}");
457		RABS_FinalGraph->GetXaxis()->CenterTitle();
458	<	RABS_FinalGraph->GetYaxis()->SetTitle("data/mc");
458	>	RABS_FinalGraph->GetYaxis()->SetTitle("<R_{abs}>_{data}/<R_{abs}>_{mc}");
459		RABS_FinalGraph->GetYaxis()->CenterTitle();
460		RABS_FinalGraph->Draw("AP");
461		TF1 *rabs_pol0=(TF1*)RABS_FinalGraph->GetFunction("pol0");
#	Line 400 | Line 468 | void new_data_mc_agreement_2d() {
468		rabs_final_mark->Draw();
469		DrawPrelim();
470
471	+	float RabsResult=1/rabs_pol0->GetParameter(0);
472	+	float RabsResultError=rabs_pol0->GetParError(0)/(rabs_pol0->GetParameter(0)*rabs_pol0->GetParameter(0));
473	+
474		filename2.str("");
475		filename2 << "Extrapolation/FINAL_RESULT_RABS";
476		CompleteSave(can,filename2.str());
477
478	+	can->cd();
479	+	TCanvas *FaceValueCanvas = new TCanvas("FaceValueCanvas","FaceValueCanvas");
480	+	FaceValueAtPoint3->SetMarkerStyle(21);
481	+	FaceValueAtPoint3->SetMarkerStyle(21);
482	+	FaceValueAtPoint3->Fit("pol0","QE"); // quiet, use minos
483	+	FaceValueAtPoint3->Draw("AP");
484	+	FaceValueAtPoint3->GetXaxis()->SetTitle("p_{T}^{Z}");
485	+	FaceValueAtPoint3->GetXaxis()->CenterTitle();
486	+	FaceValueAtPoint3->GetYaxis()->SetTitle("data/mc for #alpha<0.3");
487	+	FaceValueAtPoint3->GetYaxis()->CenterTitle();
488	+	FaceValueAtPoint3->Draw("AP");
489	+	TF1 *faceval_pol0=(TF1*)FaceValueAtPoint3->GetFunction("pol0");
490	+	stringstream faceval_result;
491	+	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)}";;
492	+	TText *faceval_final_mark = write_title(faceval_result.str());
493	+	faceval_final_mark->SetX(0.75);
494	+	faceval_final_mark->SetY(0.75);
495	+	faceval_final_mark->SetTextSize(0.03);
496	+	faceval_final_mark->Draw();
497	+	DrawPrelim();
498	+
499	+	float FaceValResult=1/faceval_pol0->GetParameter(0);
500	+	float FaceValResultError=faceval_pol0->GetParError(0)/(faceval_pol0->GetParameter(0)*faceval_pol0->GetParameter(0));
501	+
502	+	filename2.str("");
503	+	filename2 << "Extrapolation/FINAL_RESULT_MPF_FaceValp3";
504	+	CompleteSave(can,filename2.str());
505	+
506	+	cout << "FINAL RESULTS: " << endl;
507	+	cout << "MPF : " << MPFResult  << " +/- " << MPFResultError  << endl;
508	+	cout << "Rabs: " << RabsResult << " +/- " << RabsResultError << endl;
509	+	cout << "Face value at 0.3 :" << FaceValResult << " +/- " << FaceValResultError << endl;
510	+
511	+	delete can;
512	+
513	+
514	+	}
515	+
516	+
517	+	void DoPUStudy(string identifier) {
518	+
519	+	float numVtxcuts[7]={0,10,15,20};
520	+
521	+	TCanvas *can = new TCanvas("can","can");
522	+	TH1F *malpha[8];
523	+	TH1F *dalpha[8];
524	+	cout << identifier << ";";
525	+	cout << endl;
526	+	for(int i=1;i<5;i++) {
527	+	stringstream sSpecialCut;
528	+	if(i<4) {
529	+	sSpecialCut << "numVtx>" << numVtxcuts[i-1] << " && numVtx<" << numVtxcuts[i];
530	+	cout << numVtxcuts[i-1] << "  < nVtx < " << numVtxcuts[i] << ";";
531	+	} else {
532	+	sSpecialCut << "numVtx>" << numVtxcuts[i-1];
533	+	cout << numVtxcuts[i-1] << ";";
534	+	}
535	+
536	+	TCut SpecialCut(sSpecialCut.str().c_str());
537	+
538	+
539	+	malpha[i] = allsamples.Draw("malpha"+any2string(i),"mpf",1000,0,10,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,mc,luminosity);
540	+	dalpha[i] = allsamples.Draw("dalpha"+any2string(i),"mpf",1000,0,10,"MPF","events",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&SpecialCut,data,luminosity);
541	+
542	+
543	+	float a=dalpha[i]->GetMean();
544	+	float b=malpha[i]->GetMean();
545	+	float da=dalpha[i]->GetMeanError();
546	+	float db=malpha[i]->GetMeanError();
547	+	float factor = a / b;
548	+	float error = TMath::Sqrt(  (1/(b*b)) * (da*da) + ((a*a)/(b*b))*db*db);
549	+
550	+	cout << malpha[i]->GetMean() << ";" << malpha[i]->GetMeanError() << ";" << dalpha[i]->GetMean() << ";" << dalpha[i]->GetMeanError() << ";" << malpha[i]->Integral()<<";"<<dalpha[i]->Integral() << ";"<<factor << ";" << error <<  endl;
551	+	malpha[i]->SetLineColor(i+1);
552	+	if(i==0) malpha[i]->Draw("histo");
553	+	else malpha[i]->Draw("histo,same");
554	+	}
555
556	+	CompleteSave(can,"AlphaOverview_"+identifier);
557		delete can;
558
559	+	/*  for(int i=1;i<8;i++) {
560	+	delete malpha[i];
561	+	delete dalpha[i];
562	+	}
563	+	*/
564	+
565	+	}
566	+
567	+	void ScenarioComparison() {
568	+
569	+	TGraphAsymmErrors *gr[5];
570	+	TF1 *fit[5];
571	+	bool dofit=false;
572	+
573	+	TCanvas *can = new TCanvas("can","can",500,500);
574	+
575	+	TLegend *leg = new TLegend(0.2,0.2,0.5,0.45);
576	+	leg->SetBorderSize(0);
577	+	leg->SetFillColor(kWhite);
578	+
579	+	for(int i=0;i<5;i++) {
580	+	gr[i] = new TGraphAsymmErrors();
581	+	double x[4]   = {5,12.5,17.5,20};
582	+	double dxl[4] = {5,2.5,2.5,2.5};
583	+	double dxh[4] = {5,2.5,2.5,10};
584	+
585	+	string name="";
586	+	string color="";
587	+
588	+	if(i==0) {
589	+	double y[4]   = {1.008,1.013,1.003,1.007};
590	+	double dyl[4] = {0.011,0.009,0.013,0.020};
591	+	double dyh[4] = {0.011,0.009,0.013,0.020};
592	+	name="30/30";
593	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
594	+	color="#a01d99";
595	+	}
596	+	if(i==1) {
597	+	double y[4]   = {1.008+0.001,1.013+0.001,1.003+0.001,1.007+0.001};
598	+	double dyl[4] = {0.011,0.009,0.013,0.020};
599	+	double dyh[4] = {0.011,0.009,0.013,0.020};
600	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
601	+	name="30/10";
602	+	color="#2464bc";
603	+	}
604	+	if(i==2) {
605	+	double y[4]   = {1.014,1.026,1.014,1.059};
606	+	double dyl[4] = {0.011,0.010,0.014,0.023};
607	+	double dyh[4] = {0.011,0.010,0.014,0.023};
608	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
609	+	name="15/15";
610	+	color="#f49230";
611	+	}
612	+	if(i==3) {
613	+	double y[4]   = {1.015,1.020,1.008,1.050};
614	+	double dyl[4] = {0.011,0.010,0.014,0.022};
615	+	double dyh[4] = {0.011,0.010,0.014,0.022};
616	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
617	+	name="20/15";
618	+	color="#c72f3c";
619	+	}
620	+	if(i==4) {
621	+	double y[4]   = {1.015+0.001,1.020+0.001,1.008+0.001,1.050+0.001};
622	+	double dyl[4] = {0.011,0.010,0.014,0.022};
623	+	double dyh[4] = {0.011,0.010,0.014,0.022};
624	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
625	+	name="20/20";
626	+	color="#7cb433";
627	+	}
628	+
629	+
630	+	gr[i]->SetTitle();
631	+	gr[i]->GetXaxis()->SetTitle("N(Vtx)");
632	+	gr[i]->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
633	+	gr[i]->GetXaxis()->CenterTitle();
634	+	gr[i]->GetYaxis()->CenterTitle();
635	+	gr[i]->SetName(name.c_str());
636	+	gr[i]->SetMarkerSize(0.00001);
637	+	gr[i]->GetYaxis()->SetRangeUser(0.90,1.10);
638	+	gr[i]->SetLineColor(TColor::GetColor(color.c_str()));
639	+	gr[i]->SetMarkerColor(TColor::GetColor(color.c_str()));
640	+	if(dofit) {
641	+	gr[i]->Fit("pol1");
642	+	fit[i] = (TF1*)gr[i]->GetFunction("pol1");
643	+	fit[i]->SetLineColor(TColor::GetColor(color.c_str()));
644	+	}
645	+	leg->AddEntry(gr[i],name.c_str(),"l");
646	+	if(i==0) gr[i]->Draw("AP*");
647	+	else gr[i]->Draw("P");
648	+	}
649	+	leg->Draw();
650	+	DrawPrelim();
651	+
652	+	CompleteSave(can,"ScenarioComparison");
653	+	delete can;
654	+
655	+	}
656	+
657	+	void ScenarioComparisonInclusive() {
658	+
659	+	TGraphAsymmErrors *gr[5];
660	+	TF1 *fit[5];
661	+	bool dofit=true;
662	+
663	+	TCanvas *can = new TCanvas("can","can",500,500);
664	+
665	+	TLegend *leg = new TLegend(0.2,0.2,0.5,0.45);
666	+	leg->SetBorderSize(0);
667	+	leg->SetFillColor(kWhite);
668	+
669	+	for(int i=0;i<5;i++) {
670	+	gr[i] = new TGraphAsymmErrors();
671	+	double x[4]   = {5,12.5,17.5,20};
672	+	double dxl[4] = {5,2.5,2.5,2.5};
673	+	double dxh[4] = {5,2.5,2.5,10};
674	+
675	+	string name="";
676	+	string color="";
677	+
678	+	if(i==0) {
679	+	double y[4]   = {1.01103,1.01496,1.01486,1.02785};
680	+	double dyl[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
681	+	double dyh[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
682	+	name="30/30";
683	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
684	+	color="#a01d99";
685	+	}
686	+	if(i==1) {
687	+	double y[4]   = {1.01103+0.001,1.01496+0.001,1.01486+0.001,1.02785+0.001};
688	+	double dyl[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
689	+	double dyh[4] = {0.00235044,0.00225264,0.00313574,0.00541047};
690	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
691	+	name="30/10";
692	+	color="#2464bc";
693	+	}
694	+	if(i==2) {
695	+	double y[4]   = {1.01575,1.02879,1.0325,1.05456};
696	+	double dyl[4] = {0.00288281,0.00278536,0.00367451,0.00590421};
697	+	double dyh[4] = {0.00288281,0.00278536,0.00367451,0.00590421};
698	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
699	+	name="15/15";
700	+	color="#f49230";
701	+	}
702	+	if(i==3) {
703	+	double y[4]   = {1.01221,1.01922,1.02396,1.04677};
704	+	double dyl[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
705	+	double dyh[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
706	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
707	+	name="20/15";
708	+	color="#c72f3c";
709	+	}
710	+	if(i==4) {
711	+	double y[4]   = {1.01221+0.001,1.01922+0.001,1.02396+0.001,1.04677+0.001};
712	+	double dyl[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
713	+	double dyh[4] = {0.00263821,0.00260853,0.00353673,0.00584728};
714	+	gr[i] = new TGraphAsymmErrors(4,x,y,dxl,dxh,dyl,dyh);
715	+	name="20/20";
716	+	color="#7cb433";
717	+	}
718	+
719	+
720	+	gr[i]->SetTitle();
721	+	gr[i]->GetXaxis()->SetTitle("N(Vtx)");
722	+	gr[i]->GetYaxis()->SetTitle("<MPF>_{data}/<MPF>_{mc}");
723	+	gr[i]->GetXaxis()->CenterTitle();
724	+	gr[i]->GetYaxis()->CenterTitle();
725	+	gr[i]->SetName(name.c_str());
726	+	gr[i]->SetMarkerSize(0.00001);
727	+	gr[i]->GetYaxis()->SetRangeUser(0.90,1.10);
728	+	gr[i]->SetLineColor(TColor::GetColor(color.c_str()));
729	+	gr[i]->SetMarkerColor(TColor::GetColor(color.c_str()));
730	+	if(dofit) {
731	+	gr[i]->Fit("pol1");
732	+	fit[i] = (TF1*)gr[i]->GetFunction("pol1");
733	+	fit[i]->SetLineColor(TColor::GetColor(color.c_str()));
734	+	}
735	+	leg->AddEntry(gr[i],name.c_str(),"l");
736	+	if(i==0) gr[i]->Draw("AP*");
737	+	else gr[i]->Draw("P");
738	+	}
739	+	leg->Draw();
740	+	DrawPrelim();
741	+
742	+	CompleteSave(can,"ScenarioComparisonInclusive");
743	+	delete can;
744
745		}
746
747	+	void compare_selection(string identifier) {
748	+	bool recognized_scenario=false;
749	+
750	+	cout << "Running with identifier " << identifier << endl;
751	+	if(identifier=="Zb3010") {
752	+	LeadingB=TCut ("Zb3010_bTagProbCSVBP[0]>0.898");
753	+	EtaB=TCut("abs(Zb3010_pfJetGoodEta[0])<1.3");
754	+	PhiZcut=TCut("abs(Zb3010_pfJetDphiZ[0])>2.7");
755	+	recognized_scenario=true;
756	+	}
757	+
758	+	if(identifier=="Zb3010") {
759	+	LeadingB=TCut ("Zb3010_bTagProbCSVBP[0]>0.898");
760	+	EtaB=TCut("abs(Zb3010_pfJetGoodEta[0])<1.3");
761	+	PhiZcut=TCut("abs(Zb3010_pfJetDphiZ[0])>2.7");
762	+	recognized_scenario=true;
763	+	}
764	+
765	+	if(identifier=="Zb40") {
766	+	LeadingB=TCut ("Zb40_bTagProbCSVBP[0]>0.898");
767	+	EtaB=TCut("abs(Zb40_pfJetGoodEta[0])<1.3");
768	+	PhiZcut=TCut("abs(Zb40_pfJetDphiZ[0])>2.7");
769	+	recognized_scenario=true;
770	+	}
771	+
772	+	if(identifier=="Zb301030") {
773	+	LeadingB=TCut ("Zb301030_bTagProbCSVBP[0]>0.898");
774	+	EtaB=TCut("abs(Zb301030_pfJetGoodEta[0])<1.3");
775	+	PhiZcut=TCut("abs(Zb301030_pfJetDphiZ[0])>2.7");
776	+	recognized_scenario=true;
777	+	}
778	+
779	+	if(identifier=="Zb1530") {
780	+	LeadingB=TCut ("Zb1530_bTagProbCSVBP[0]>0.898");
781	+	EtaB=TCut("abs(Zb1530_pfJetGoodEta[0])<1.3");
782	+	PhiZcut=TCut("abs(Zb1530_pfJetDphiZ[0])>2.7");
783	+	recognized_scenario=true;
784	+	}
785	+
786	+	if(identifier=="Zb301010") {
787	+	LeadingB=TCut ("Zb301010_bTagProbCSVBP[0]>0.898");
788	+	EtaB=TCut("abs(Zb301010_pfJetGoodEta[0])<1.3");
789	+	PhiZcut=TCut("abs(Zb301010_pfJetDphiZ[0])>2.7");
790	+	recognized_scenario=true;
791	+	}
792	+
793	+	if(identifier=="Zb1510") {
794	+	LeadingB=TCut ("Zb1510_bTagProbCSVBP[0]>0.898");
795	+	EtaB=TCut("abs(Zb1510_pfJetGoodEta[0])<1.3");
796	+	PhiZcut=TCut("abs(Zb1510_pfJetDphiZ[0])>2.7");
797	+	recognized_scenario=true;
798	+	}
799	+
800	+	if(identifier=="Zb301010") {
801	+	LeadingB=TCut ("Zb301010_bTagProbCSVBP[0]>0.898");
802	+	EtaB=TCut("abs(Zb301010_pfJetGoodEta[0])<1.3");
803	+	PhiZcut=TCut("abs(Zb301010_pfJetDphiZ[0])>2.7");
804	+	recognized_scenario=true;
805	+	}
806	+
807	+	if(identifier=="ZbMikko") {
808	+	LeadingB=TCut ("ZbMikko_bTagProbCSVBP[0]>0.898 && ZbMikko__IsClean");
809	+	EtaB=TCut("abs(ZbMikko_pfJetGoodEta[0])<1.3 && ZbMikko__IsClean");
810	+	PhiZcut=TCut("abs(ZbMikko_pfJetDphiZ[0])>2.7 && ZbMikko__IsClean");
811	+	recognized_scenario=true;
812	+	}
813	+
814	+	if(identifier=="Zb3010_SecEta3") {
815	+	LeadingB=TCut ("Zb3010_SecEta3_bTagProbCSVBP[0]>0.898");
816	+	EtaB=TCut("abs(Zb3010_SecEta3_pfJetGoodEta[0])<1.3");
817	+	PhiZcut=TCut("abs(Zb3010_SecEta3_pfJetDphiZ[0])>2.7");
818	+	recognized_scenario=true;
819	+	}
820	+
821	+	if(identifier=="Zb3010_SecEta5") {
822	+	LeadingB=TCut ("Zb3010_SecEta5_bTagProbCSVBP[0]>0.898");
823	+	EtaB=TCut("abs(Zb3010_SecEta5_pfJetGoodEta[0])<1.3");
824	+	PhiZcut=TCut("abs(Zb3010_SecEta5_pfJetDphiZ[0])>2.7");
825	+	recognized_scenario=true;
826	+	}
827	+
828	+	if(identifier=="Zb3010_p5Clean") {
829	+	LeadingB=TCut ("Zb3010_p5Clean_bTagProbCSVBP[0]>0.898 && Zb3010_p5Clean_IsClean");
830	+	EtaB=TCut("abs(Zb3010_p5Clean_pfJetGoodEta[0])<1.3 && Zb3010_p5Clean_IsClean");
831	+	PhiZcut=TCut("abs(Zb3010_p5Clean_pfJetDphiZ[0])>2.7 && Zb3010_p5Clean_IsClean");
832	+	recognized_scenario=true;
833	+	}
834	+
835	+	if(identifier=="Zb3010CHS") {
836	+	LeadingB=TCut ("Zb3010CHS_bTagProbCSVBP[0]>0.898");
837	+	EtaB=TCut("abs(Zb3010CHS_pfJetGoodEta[0])<1.3");
838	+	PhiZcut=TCut("abs(Zb3010CHS_pfJetDphiZ[0])>2.7");
839	+	recognized_scenario=true;
840	+	}
841	+
842	+	assert(recognized_scenario);
843	+
844	+	cout << "Cuts: " << endl;
845	+	cout << "   " << (const char*) LeadingB << endl;
846	+	cout << "   " << (const char*) EtaB << endl;
847	+	cout << "   " << (const char*) PhiZcut << endl;
848	+
849	+	cout << endl << endl;
850	+
851	+	//    ScenarioComparison();
852	+	//    ScenarioComparisonInclusive();
853	+
854	+
855	+	//     DoPUStudy(identifier);
856	+
857	+	//   SpecialCutFlow(identifier);
858	+	//    draw_kin_variable("pt",ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&"pt>10&&pt<1000",17000,30,200,"Z p_{T}","Official/"+identifier+"/Zpt__"+identifier,1);
859	+	//   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);
860	+	//   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);
861	+	//   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);
862	+	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);
863	+	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);
864	+	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);
865	+	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);
866	+	draw_kin_variable(identifier+"_alpha",ZplusBsel&&LeadingB&&EtaB&&PhiZcut,20,0,2,"#alpha","Official/"+identifier+"/alpha_full__"+identifier,1);
867	+	//   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);
868	+
869	+	//   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");
870	+	//   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");
871
872	+
873	+	//   */
874	+
875	+	}
876
877		void do_basic_ZB_analysis() {
878	+
879	+
880		TCanvas *zbcanvas = new TCanvas("zbcanvas","zbcanvas");
417	–	write_warning(__FUNCTION__,"TESTING NEW DATA MC AGREEMENT");
881
882	<	/*
882	>	//     compare_selection("Zb3010_p5Clean");
883	>	//     compare_selection("Zb3010CHS");
884	>	//     compare_selection("Zb3010");
885	>	//     compare_selection("Zb301030");
886	>	//     compare_selection("Zb1530");
887	>	//     compare_selection("Zb301010");
888	>	//     compare_selection("Zb3010_SecEta3");
889	>	//     compare_selection("Zb3010_SecEta5");
890	>	//     compare_selection("Zb1510");
891	>	//     compare_selection("Zb301010");
892	>	//     compare_selection("ZbMikko");
893	>	//     compare_selection("Zb3010");
894	>	//     compare_selection("Zb40");
895	>
896	>
897		//  print_all_b_yields();
898		//  draw_mpf_vars();
899	<	draw_Zb_kin_vars();
899	>	//  DrawEvilCutFlow();
900	>
901	>	//  draw_Zb_kin_vars();
902
424	–	data_mc_agreement_2d();
425	–	*/
903		new_data_mc_agreement_2d();
904
905	+
906		delete zbcanvas;
907		}
908	+
909	+
910	+
911	+
912	+	//*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*  DELETE EVERYTHING BELOW THIS LINE /*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*/*//
913	+
914	+
915	+	const char* concatenate(string bla1, string bla2) {
916	+	stringstream bla;
917	+	bla << bla1 << bla2;
918	+	return bla.str().c_str();
919	+	}
920	+
921	+	void SpecialCutFlow(string identifier) {
922	+	stringstream MegaCut;
923	+
924	+
925	+
926	+	MegaCut << "   1*(" << (const char*) (ZplusBsel&&TCut(concatenate(identifier,"_pfJetGoodNumBtag>0"))) << ")";
927	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB) << ")";
928	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB) << ")";
929	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut) << ")";
930	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000").c_str())) << ")";
931	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.3").c_str())) << ")";
932	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.2").c_str())) << ")";
933	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.15").c_str())) << ")";
934	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.1").c_str())) << ")";
935	+	MegaCut << " + 1*(" << (const char*) (ZplusBsel&&LeadingB&&EtaB&&PhiZcut&&TCut(((string)"pt>10&&pt<1000&&"+identifier+"_alpha<0.05").c_str())) << ")";
936	+	MegaCut << " ";
937	+
938	+	TCanvas *can = new TCanvas("can","can");
939	+	can->SetLogy(1);
940	+	TCut basecut(ZplusBsel);
941	+
942	+
943	+
944	+	TLegend *leg = allsamples.allbglegend();
945	+
946	+
947	+	TH1F *data  = allsamples.Draw("data",    MegaCut.str(),11,-0.5,10.5, "", "events", basecut,data,luminosity);
948	+	THStack mcm = allsamples.DrawStack("mc", MegaCut.str(),11,-0.5,10.5, "", "events", basecut,mc,luminosity);
949	+
950	+	float runningsum=0;
951	+	for(int i=data->GetNbinsX();i>0;i--) {
952	+	runningsum+=data->GetBinContent(i);
953	+	data->SetBinContent(i,runningsum);
954	+	}
955	+
956	+	float minimum=data->GetMaximum();
957	+
958	+	TH1F* h;
959	+	TIter nextOF(mcm.GetHists());
960	+
961	+	float mcsum=0;
962	+	float zb=0;
963	+	while ( h = (TH1F*)nextOF() ) {
964	+	float runningsum=0;
965	+	minimum=data->GetMaximum();//done deliberately at each step so get the minimum of the last (!) contribution
966	+	for(int i=h->GetNbinsX();i>0;i--) {
967	+	runningsum+=h->GetBinContent(i);
968	+	h->SetBinContent(i,runningsum);
969	+	if(runningsum<minimum&&runningsum>0) minimum=runningsum;
970	+	}
971	+	if(Contains(h->GetName(),"Z_b_")) {
972	+	zb=h->GetBinContent(h->GetNbinsX());
973	+	}
974	+	mcsum+=h->GetBinContent(h->GetNbinsX());
975	+	}
976	+
977	+
978	+
979	+	data->SetMinimum(0.05*minimum);
980	+	can->cd(1)->SetBottomMargin(0.25);
981	+	data->GetXaxis()->SetBinLabel(1,"Pre-selection");
982	+	data->GetXaxis()->SetBinLabel(2,"Contains b jet");
983	+	data->GetXaxis()->SetBinLabel(3,"Leading jet is b-jet");
984	+	data->GetXaxis()->SetBinLabel(4,"|#eta_{b}|<1.3 ");
985	+	data->GetXaxis()->SetBinLabel(5,"|#delta#phi(b,Z)|>2.7");
986	+	data->GetXaxis()->SetBinLabel(6,"10<p_{T}^{Z}<1000 GeV");
987	+	data->GetXaxis()->SetBinLabel(7,"#alpha < 0.3");
988	+	data->GetXaxis()->SetBinLabel(8,"#alpha < 0.2");
989	+	data->GetXaxis()->SetBinLabel(9,"#alpha < 0.15");
990	+	data->GetXaxis()->SetBinLabel(10,"#alpha < 0.1");
991	+	data->GetXaxis()->SetBinLabel(11,"#alpha < 0.05");
992	+	data->GetXaxis()->LabelsOption("v");
993	+
994	+	stringstream purityinfo;
995	+	purityinfo << "Purity: " << std::setprecision(3) << 100*zb/mcsum << " %";
996	+	stringstream neventsinfo;
997	+	neventsinfo << "Nevents: " << data->GetBinContent(data->GetNbinsX());
998	+	TH1F *crap = new TH1F("crap","",1,0,1);
999	+	crap->SetLineColor(kWhite);
1000	+	leg->AddEntry(crap,purityinfo.str().c_str(),"l");
1001	+	leg->AddEntry(crap,neventsinfo.str().c_str(),"l");
1002	+
1003	+	data->Draw("e1");
1004	+	mcm.Draw("same");
1005	+	data->Draw("same,e1,axis");
1006	+	data->Draw("same,e1");
1007	+	//  data->Draw("same,TEXT");
1008	+
1009	+	leg->Draw();
1010	+
1011	+	DrawPrelim();
1012	+
1013	+	CompleteSave(can,"CutFlow___"+identifier);
1014	+
1015	+	delete crap;
1016	+	delete data;
1017	+
1018	+
1019	+	}

Diff Legend

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