[ViewVC] Diff of: cvsroot/UserCode/SchieferD/jetcalib/mcfakerate

Comparing UserCode/SchieferD/jetcalib/mcfakerate_x.cpp (file contents):
Revision 1.1 by schiefer, Wed Aug 15 17:20:06 2007 UTC vs.
Revision 1.4 by schiefer, Tue Sep 4 13:31:20 2007 UTC

#	Line 17 \| Line 17
17		#include <TH1F.h>
18		#include <TCanvas.h>
19		#include <TLegend.h>
20	+	#include <TMath.h>
21
22		#include <iostream>
23		#include <iomanip>
#	Line 33 \| Line 34 \| using namespace std;
34		////////////////////////////////////////////////////////////////////////////////
35
36		//______________________________________________________________________________
37	<	void setBinomialErrors(TH1F* hEff,const TH1F* hFake, const TH1F* hAll);
37	>	void setBinomialErrors(TH1F* hEff,const TH1F* hFake, const TH1F* hAll);
38	>	string getLegLabel(const string& jetalg);
39
40		////////////////////////////////////////////////////////////////////////////////
41		// main
#	Line 47 \| Line 49 \| int main(int argc,charargv)**
49
50		vector<string> input = cl.get_vector<string>("input");
51		string treename = cl.get_value<string> ("treename", "t");
52	<	double drmax = cl.get_value<double> ("drmax", 0.3);
53	<	int nptbins = cl.get_value<int> ("nptbins", 50);
52	>	float drmax = cl.get_value<float> ("drmax", 0.3);
53	>	int netbins = cl.get_value<int> ("netbins", 50);
54		int netabins = cl.get_value<int> ("netabins", 25);
55		int nphibins = cl.get_value<int> ("nphibins", 25);
56	<	double ptmin = cl.get_value<double> ("ptmin", 0.0);
57	<	double ptmax = cl.get_value<double> ("ptmax", 200.0);
58	<	double etamin = cl.get_value<double> ("etamin", -5.0);
59	<	double etamax = cl.get_value<double> ("etamax", 5.0);
60	<	double phimin = cl.get_value<double> ("phimin", -3.2);
61	<	double phimax = cl.get_value<double> ("phimax", 3.2);
56	>	int nemfbins = cl.get_value<int> ("nemfbins", 15);
57	>	float etmin = cl.get_value<float> ("etmin", 0.0);
58	>	float etmax = cl.get_value<float> ("etmax", 200.0);
59	>	float etamin = cl.get_value<float> ("etamin", -5.0);
60	>	float etamax = cl.get_value<float> ("etamax", 5.0);
61	>	float phimin = cl.get_value<float> ("phimin", -3.2);
62	>	float phimax = cl.get_value<float> ("phimax", 3.2);
63	>	float emfmin = cl.get_value<float> ("emfmin", 0.0);
64	>	float emfmax = cl.get_value<float> ("emfmax", 1.0);
65	>	short applyjes = cl.get_value<short> ("applyjes", 0);
66
67		if (!cl.check()) return 0;
68		cl.print();
#	Line 73 \| Line 79 \| int main(int argc,charargv)**
79		}
80
81		// loop over all input samples
82	<	vector<TH1F*> ptHistos;
82	>	vector<TH1F*> etHistos;
83		vector<TH1F*> etaHistos;
84		vector<TH1F*> phiHistos;
85	+	vector<TH1F*> emfHistos;
86		for (unsigned int i=0;i<input.size();++i) {
87		TFile* f = new TFile(input[i].c_str(),"READ");
88		TTree* t = (TTree*)f->Get(treename.c_str());
89		cout<<jetalgs[i]<<" sample has "<<t->GetEntries()<<" entries."<<endl;
90
91
92	<	string namePtAll = "PtAll_"+jetalgs[i];
93	<	string namePtFake = "PtFake_"+jetalgs[i];
94	<	string namePtEff = "PtEff_"+jetalgs[i];
95	<
96	<	TH1F* hPtAll = new TH1F(namePtAll.c_str(),namePtAll.c_str(),
97	<	nptbins,ptmin,ptmax);
98	<	TH1F* hPtFake = new TH1F(namePtFake.c_str(),namePtFake.c_str(),
99	<	nptbins,ptmin,ptmax);
100	<	TH1F* hPtEff = new TH1F(namePtEff.c_str(),namePtEff.c_str(),
101	<	nptbins,ptmin,ptmax);
102	<
103	<	hPtAll->Sumw2();
104	<	hPtFake->Sumw2();
105	<	hPtEff->Sumw2();
92	>	string nameEtAll = "EtAll_"+jetalgs[i];
93	>	string nameEtFake = "EtFake_"+jetalgs[i];
94	>	string nameEtEff = "EtEff_"+jetalgs[i];
95	>
96	>	TH1F* hEtAll = new TH1F(nameEtAll.c_str(),nameEtAll.c_str(),
97	>	netbins,etmin,etmax);
98	>	TH1F* hEtFake = new TH1F(nameEtFake.c_str(),nameEtFake.c_str(),
99	>	netbins,etmin,etmax);
100	>	TH1F* hEtEff = new TH1F(nameEtEff.c_str(),nameEtEff.c_str(),
101	>	netbins,etmin,etmax);
102	>
103	>	hEtAll->Sumw2();
104	>	hEtFake->Sumw2();
105	>	hEtEff->Sumw2();
106
107		string nameEtaAll = "EtaAll_"+jetalgs[i];
108		string nameEtaFake = "EtaFake_"+jetalgs[i];
#	Line 127 \| Line 134 \| int main(int argc,charargv)**
134		hPhiFake->Sumw2();
135		hPhiEff->Sumw2();
136
137	+	string nameEmfAll = "EmfAll_"+jetalgs[i];
138	+	string nameEmfFake = "EmfFake_"+jetalgs[i];
139	+	string nameEmfEff = "EmfEff_"+jetalgs[i];
140	+
141	+	TH1F* hEmfAll = new TH1F(nameEmfAll.c_str(),nameEmfAll.c_str(),
142	+	nemfbins,emfmin,emfmax);
143	+	TH1F* hEmfFake = new TH1F(nameEmfFake.c_str(),nameEmfFake.c_str(),
144	+	nemfbins,emfmin,emfmax);
145	+	TH1F* hEmfEff = new TH1F(nameEmfEff.c_str(),nameEmfEff.c_str(),
146	+	nemfbins,emfmin,emfmax);
147	+
148	+	hEmfAll->Sumw2();
149	+	hEmfFake->Sumw2();
150	+	hEmfEff->Sumw2();
151	+
152		float weight; t->SetBranchAddress("weight",&weight);
153		char njt; t->SetBranchAddress("njt", &njt);
154	<	float jtpt[100]; t->SetBranchAddress("jtpt", jtpt);
154	>	float jtet[100]; t->SetBranchAddress("jtet", jtet);
155		float jteta[100]; t->SetBranchAddress("jteta", jteta);
156		float jtphi[100]; t->SetBranchAddress("jtphi", jtphi);
157	+	float jtemf[100]; t->SetBranchAddress("jtemf", jtemf);
158		float jtgendr[100]; t->SetBranchAddress("jtgendr",jtgendr);
159	<
159	>	float jtjes[100][3];
160	>	if (applyjes>0&&t->FindBranch("jtjes")) t->SetBranchAddress("jtjes",jtjes);
161	>
162		int nevts=t->GetEntries();
163
164		for (int ievt=0;ievt<nevts;ievt++) {
165	+
166		t->GetEntry(ievt);
167	+
168		for (int ijt=0;ijt<njt;ijt++) {
142	–	hPtAll ->Fill(jtpt[ijt], weight);
143	–	hEtaAll->Fill(jteta[ijt],weight);
144	–	hPhiAll->Fill(jtphi[ijt],weight);
169
170	<	if (jtgendr[ijt]<0.0\|\|jtgendr[ijt]>drmax) {
171	<	hPtFake ->Fill(jtpt[ijt], weight);
172	<	hEtaFake->Fill(jteta[ijt],weight);
173	<	hPhiFake->Fill(jtphi[ijt],weight);
170	>	float dr = jtgendr[ijt];
171	>	float et = (applyjes==0) ? jtet[ijt]:jtet[ijt]*jtjes[ijt][applyjes-1];
172	>	float eta = jteta[ijt];
173	>	float phi = jtphi[ijt];
174	>	float emf = jtemf[ijt]; if (TMath::IsNaN(emf)) emf = 0.0;
175	>
176	>	if (et <etmin \|\|et >etmax) continue;
177	>	if (eta<etamin\|\|eta>etamax) continue;
178	>	if (phi<phimin\|\|phi>phimax) continue;
179	>	if (emf<emfmin\|\|emf>emfmax) continue;
180	>
181	>	hEtAll ->Fill(et, weight);
182	>	hEtaAll->Fill(eta,weight);
183	>	hPhiAll->Fill(phi,weight);
184	>	hEmfAll->Fill(emf,weight);
185	>
186	>	if (dr<0.0\|\|dr>drmax) {
187	>	hEtFake ->Fill(et, weight);
188	>	hEtaFake->Fill(eta,weight);
189	>	hPhiFake->Fill(phi,weight);
190	>	hEmfFake->Fill(emf,weight);
191		}
192		}
193		}
194
195	<	hPtEff->Divide(hPtFake,hPtAll,1.0,1.0,"B");
196	<	setBinomialErrors(hPtEff,hPtFake,hPtAll);
197	<	ptHistos.push_back(hPtEff);
195	>	hEtEff->Divide(hEtFake,hEtAll,1.0,1.0,"B");
196	>	setBinomialErrors(hEtEff,hEtFake,hEtAll);
197	>	etHistos.push_back(hEtEff);
198
199		hEtaEff->Divide(hEtaFake,hEtaAll,1.0,1.0,"B");
200		setBinomialErrors(hEtaEff,hEtaFake,hEtaAll);
#	Line 162 \| Line 203 \| int main(int argc,charargv)**
203		hPhiEff->Divide(hPhiFake,hPhiAll,1.0,1.0,"B");
204		setBinomialErrors(hPhiEff,hPhiFake,hPhiAll);
205		phiHistos.push_back(hPhiEff);
206	+
207	+	hEmfEff->Divide(hEmfFake,hEmfAll,1.0,1.0,"B");
208	+	setBinomialErrors(hEmfEff,hEmfFake,hEmfAll);
209	+	emfHistos.push_back(hEmfEff);
210		}
211
212
#	Line 173 \| Line 218 \| int main(int argc,charargv)**
218
219		Color_t color;
220
221	+	// emf canvas
222	+	TCanvas* cEmfEff = new TCanvas("FakeVsEmf", "FakeVsEmf",10,420,400,400);
223	+	TLegend* lEmfEff = new TLegend(0.65,0.85-emfHistos.size()*0.05,0.86,0.85);
224	+	cEmfEff->cd();
225	+	lEmfEff->SetFillColor(10);
226	+
227	+	color = kBlack;
228	+	for (unsigned int i=0;i<emfHistos.size();i++) {
229	+	TH1F* h = emfHistos[i];
230	+	lEmfEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
231	+	h->SetLineWidth(2);
232	+	h->SetMarkerColor(color);
233	+	h->SetLineColor(color++);
234	+	h->SetMarkerStyle(kFullCircle);
235	+	h->SetMarkerSize(0.7);
236	+	if (i==0) {
237	+	h->SetTitle("#epsilon_{FAKE} vs emf");
238	+	h->SetXTitle("emf");
239	+	h->SetYTitle("#epsilon_{FAKE}");
240	+	h->Draw("E");
241	+	h->SetMinimum(0.0);
242	+	h->SetMaximum(1.05);
243	+	}
244	+	else h->Draw("ESAME");
245	+	}
246	+	lEmfEff->Draw();
247	+
248		// phi canvas
249		TCanvas* cPhiEff = new TCanvas("FakeVsPhi", "FakeVsPhi",830,10,400,400);
250		TLegend* lPhiEff = new TLegend(0.65,0.85-phiHistos.size()*0.05,0.86,0.85);
#	Line 182 \| Line 254 \| int main(int argc,charargv)**
254		color = kBlack;
255		for (unsigned int i=0;i<phiHistos.size();i++) {
256		TH1F* h = phiHistos[i];
257	<	lPhiEff->AddEntry(h,jetalgs[i].c_str(),"l");
257	>	lPhiEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
258	>	h->SetLineWidth(2);
259		h->SetMarkerColor(color);
260		h->SetLineColor(color++);
261		h->SetMarkerStyle(kFullCircle);
262	<	h->SetMarkerSize(0.5);
262	>	h->SetMarkerSize(0.7);
263		if (i==0) {
264	<	h->SetTitle("#epsilon_{fake} vs #phi");
265	<	h->SetXTitle("jet #phi");
266	<	h->SetYTitle("#epsilon_{fake}");
264	>	h->SetTitle("#epsilon_{FAKE} vs #phi");
265	>	h->SetXTitle("#phi");
266	>	h->SetYTitle("#epsilon_{FAKE}");
267		h->Draw("E");
268		h->SetMinimum(0.0);
269		h->SetMaximum(1.05);
#	Line 200 \| Line 273 \| int main(int argc,charargv)**
273		lPhiEff->Draw();
274
275		// eta canvas
276	<	TCanvas* cEtaEff = new TCanvas("EffVsEta", "EffVsEta",420,10,400,400);
276	>	TCanvas* cEtaEff = new TCanvas("FakeVsEta", "FakeVsEta",420,10,400,400);
277		TLegend* lEtaEff = new TLegend(0.65,0.85-etaHistos.size()*0.05,0.86,0.85);
278		cEtaEff->cd();
279		lEtaEff->SetFillColor(10);
#	Line 208 \| Line 281 \| int main(int argc,charargv)**
281		color = kBlack;
282		for (unsigned int i=0;i<etaHistos.size();i++) {
283		TH1F* h = etaHistos[i];
284	<	lEtaEff->AddEntry(h,jetalgs[i].c_str(),"l");
284	>	lEtaEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
285	>	h->SetLineWidth(2);
286		h->SetMarkerColor(color);
287		h->SetLineColor(color++);
288		h->SetMarkerStyle(kFullCircle);
289	<	h->SetMarkerSize(0.5);
289	>	h->SetMarkerSize(0.7);
290		if (i==0) {
291	<	h->SetTitle("#epsilon_{fake} vs #eta");
292	<	h->SetXTitle("jet #eta");
293	<	h->SetYTitle("#epsilon_{fake}");
291	>	h->SetTitle("#epsilon_{FAKE} vs #eta");
292	>	h->SetXTitle("#eta");
293	>	h->SetYTitle("#epsilon_{FAKE}");
294		h->Draw("E");
295		h->SetMinimum(0.0);
296		h->SetMaximum(1.05);
#	Line 226 \| Line 300 \| int main(int argc,charargv)**
300		}
301		lEtaEff->Draw();
302
303	<	// pt canvas
304	<	TCanvas* cPtEff = new TCanvas("EffVsPt", "EffVsPt",10,10,400,400);
305	<	TLegend* lPtEff = new TLegend(0.65,0.85-ptHistos.size()*0.05,0.86,0.85);
306	<	cPtEff->cd();
307	<	lPtEff->SetFillColor(10);
303	>	// et canvas
304	>	TCanvas* cEtEff = new TCanvas("FakeVsEt", "FakeVsEt",10,10,400,400);
305	>	TLegend* lEtEff = new TLegend(0.65,0.85-etHistos.size()*0.05,0.86,0.85);
306	>	cEtEff->cd();
307	>	lEtEff->SetFillColor(10);
308
309		color = kBlack;
310	<	for (unsigned int i=0;i<ptHistos.size();i++) {
311	<	TH1F* h = ptHistos[i];
312	<	lPtEff->AddEntry(h,jetalgs[i].c_str(),"l");
310	>	for (unsigned int i=0;i<etHistos.size();i++) {
311	>	TH1F* h = etHistos[i];
312	>	lEtEff->AddEntry(h,getLegLabel(jetalgs[i]).c_str(),"lp");
313	>	h->SetLineWidth(2);
314		h->SetMarkerColor(color);
315		h->SetLineColor(color++);
316		h->SetMarkerStyle(kFullCircle);
317	<	h->SetMarkerSize(0.5);
317	>	h->SetMarkerSize(0.7);
318		if (i==0) {
319	<	h->SetTitle("#epsilon_{fake} vs p_{T}");
320	<	h->SetXTitle("jet p_{T} [GeV]");
321	<	h->SetYTitle("#epsilon_{fake}");
319	>	h->SetTitle("#epsilon_{FAKE} vs E_{T}");
320	>	h->SetXTitle("E_{T} [GeV]");
321	>	h->SetYTitle("#epsilon_{FAKE}");
322		h->Draw("E");
323		h->SetMinimum(0.0);
324		h->SetMaximum(1.05);
325		}
326		else h->Draw("ESAME");
327		}
328	<	lPtEff->Draw();
328	>	lEtEff->Draw();
329
330		app->Run();
331
#	Line 266 \| Line 341 \| int main(int argc,charargv)**
341		void setBinomialErrors(TH1F* hEff,const TH1F* hFake, const TH1F* hAll)
342		{
343		for (int i=1;i<hEff->GetNbinsX();i++) {
344	<	double nfake = hFake->GetBinContent(i);
345	<	double nall = hAll ->GetBinContent(i);
346	<	double eeff = (nall>0.0) ? std::sqrt(nfake/(nallnall)(1-nfake/nall)) : 0.0;
344	>	float nfake = hFake->GetBinContent(i);
345	>	float nall = hAll ->GetBinContent(i);
346	>	float eeff = (nall>0.0) ? std::sqrt(nfake/(nallnall)(1-nfake/nall)):0.0;
347		hEff->SetBinError(i,eeff);
348		}
349		}
350	+
351	+
352	+	//______________________________________________________________________________
353	+	string getLegLabel(const string& jetalg)
354	+	{
355	+	string result;
356	+	string tmp(jetalg);
357	+
358	+	if (tmp.find("fastjet")==0) {
359	+	result = "fast k_{T}, D=";
360	+	tmp=tmp.substr(7);
361	+	}
362	+	else if (tmp.find("icone")==0) {
363	+	result = "icone, R=";
364	+	tmp = tmp.substr(5);
365	+	}
366	+	else if (tmp.find("mcone")==0) {
367	+	result = "mcone, R=";
368	+	tmp = tmp.substr(5);
369	+	}
370	+	else return result;
371	+
372	+	stringstream ssparam;
373	+	ssparam<<tmp;
374	+	float param;
375	+	ssparam>>param;
376	+	param = 0.1*param;
377	+	stringstream ssparam2;
378	+	ssparam2<<param;
379	+	result += ssparam2.str();
380	+
381	+	return result;
382	+	}

Diff Legend

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

Comparing UserCode/SchieferD/jetcalib/mcfakerate_x.cpp (file contents): Revision 1.1 by schiefer, Wed Aug 15 17:20:06 2007 UTC vs. Revision 1.4 by schiefer, Tue Sep 4 13:31:20 2007 UTC

Diff Legend

Comparing UserCode/SchieferD/jetcalib/mcfakerate_x.cpp (file contents):
Revision 1.1 by schiefer, Wed Aug 15 17:20:06 2007 UTC vs.
Revision 1.4 by schiefer, Tue Sep 4 13:31:20 2007 UTC