Fork me on GitHub

Context Navigation

source: git/validation/DelphesValidation.cpp@ a6e8029

Visit:

Last change on this file since a6e8029 was e09b9bf, checked in by Pavel Demin <pavel.demin@…>, 4 years ago
add include for TF1.h to DelphesValidation.cpp
Property mode set to `100644`
File size: 107.5 KB

Rev	Line
[734b267]	1	/*
	2	* Delphes: a framework for fast simulation of a generic collider experiment
	3	* Copyright (C) 2012-2014 Universite catholique de Louvain (UCL), Belgium
	4	*
	5	* This program is free software: you can redistribute it and/or modify
	6	* it under the terms of the GNU General Public License as published by
	7	* the Free Software Foundation, either version 3 of the License, or
	8	* (at your option) any later version.
	9	*
	10	* This program is distributed in the hope that it will be useful,
	11	* but WITHOUT ANY WARRANTY; without even the implied warranty of
	12	* MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
	13	* GNU General Public License for more details.
	14	*
	15	* You should have received a copy of the GNU General Public License
	16	* along with this program. If not, see <http://www.gnu.org/licenses/>.
	17	*/
	18
	19	#include <iostream>
[341014c]	20	#include <typeinfo>
[734b267]	21	#include <utility>
	22	#include <vector>
	23
[341014c]	24	#include "TApplication.h"
[734b267]	25	#include "TROOT.h"
	26	#include "TSystem.h"
	27
	28	#include "TString.h"
	29
[341014c]	30	#include "TCanvas.h"
	31	#include "TClonesArray.h"
	32	#include "TGraph.h"
	33	#include "TGraphErrors.h"
[e09b9bf]	34	#include "TF1.h"
[1408174]	35	#include "TH1.h"
[734b267]	36	#include "TH2.h"
	37	#include "THStack.h"
	38	#include "TLegend.h"
	39	#include "TLorentzVector.h"
[341014c]	40	#include "TMath.h"
[1408174]	41	#include "TMultiGraph.h"
[341014c]	42	#include "TPaveText.h"
	43	#include "TStyle.h"
[734b267]	44
	45	#include "classes/DelphesClasses.h"
	46
[341014c]	47	#include "ExRootAnalysis/ExRootResult.h"
	48	#include "ExRootAnalysis/ExRootTreeBranch.h"
[734b267]	49	#include "ExRootAnalysis/ExRootTreeReader.h"
	50	#include "ExRootAnalysis/ExRootTreeWriter.h"
	51	#include "ExRootAnalysis/ExRootUtilities.h"
	52
	53	using namespace std;
	54
	55	//------------------------------------------------------------------------------
	56
[92237ca]	57	static const int Nbins = 50;
[1408174]	58
	59	int objStyle = 1;
	60	int trackStyle = 7;
	61	int towerStyle = 3;
	62
	63	Color_t objColor = kBlack;
	64	Color_t trackColor = kBlack;
	65	Color_t towerColor = kBlack;
	66
	67	double effLegXmin = 0.22;
	68	double effLegXmax = 0.7;
	69	double effLegYmin = 0.22;
	70	double effLegYmax = 0.5;
	71
	72	double resLegXmin = 0.62;
	73	double resLegXmax = 0.9;
	74	double resLegYmin = 0.52;
	75	double resLegYmax = 0.85;
	76
	77	double topLeftLegXmin = 0.22;
	78	double topLeftLegXmax = 0.7;
	79	double topLeftLegYmin = 0.52;
	80	double topLeftLegYmax = 0.85;
	81
[92237ca]	82	unsigned int k;
	83
[1408174]	84	struct resolPlot
	85	{
[341014c]	86	TH1 *resolHist;
[1408174]	87	double ptmin;
	88	double ptmax;
[92237ca]	89	double etamin;
	90	double etamax;
[1408174]	91	double xmin;
	92	double xmax;
	93	TString obj;
	94
	95	resolPlot();
	96	resolPlot(double ptdown, double ptup, TString object);
[92237ca]	97	resolPlot(double etadown, double etaup, double ptdown, double ptup, TString object);
[1408174]	98	void set(double ptdown, double ptup, TString object, double xmin = 0, double xmax = 2);
[92237ca]	99	void set(double etadown, double etaup, double ptdown, double ptup, TString object, double xmin = 0, double xmax = 2);
[341014c]	100	void print() { std::cout << ptmin << std::endl; }
[1408174]	101	};
	102
	103	resolPlot::resolPlot()
	104	{
	105	}
	106
	107	resolPlot::resolPlot(double ptdown, double ptup, TString object)
	108	{
[341014c]	109	this->set(ptdown, ptup, object);
[1408174]	110	}
	111
[92237ca]	112	resolPlot::resolPlot(double etadown, double etaup, double ptdown, double ptup, TString object)
	113	{
	114	this->set(etadown, etaup, ptdown, ptup, object);
	115	}
	116
[1408174]	117	void resolPlot::set(double ptdown, double ptup, TString object, double xmin, double xmax)
	118	{
	119	ptmin = ptdown;
	120	ptmax = ptup;
	121	obj = object;
	122
[341014c]	123	resolHist = new TH1D(obj + "_delta_pt_" + Form("%4.2f", ptmin) + "_" + Form("%4.2f", ptmax), obj + "_delta_pt_" + Form("%4.2f", ptmin) + "_" + Form("%4.2f", ptmax), 1000, xmin, xmax);
[92237ca]	124	}
	125
	126	void resolPlot::set(double etadown, double etaup, double ptdown, double ptup, TString object, double xmin, double xmax)
	127	{
	128	etamin = etadown;
	129	etamax = etaup;
	130	ptmin = ptdown;
	131	ptmax = ptup;
	132	obj = object;
	133
[341014c]	134	resolHist = new TH1D(obj + "_delta_pt_" + Form("%4.2f", ptmin) + "_" + Form("%4.2f", ptmax) + "_" + Form("%4.2f", etamin) + "_" + Form("%4.2f", etamax), obj + "_delta_pt_" + Form("%4.2f", ptmin) + "_" + Form("%4.2f", ptmax) + "_" + Form("%4.2f", etamin) + "_" + Form("%4.2f", etamax), 1000, xmin, xmax);
[1408174]	135	}
	136
	137	void HistogramsCollection(std::vector<resolPlot> *histos, double ptmin, double ptmax, TString obj, double xmin = 0, double xmax = 2)
	138	{
	139	double width;
	140	double ptdown;
	141	double ptup;
	142	resolPlot ptemp;
	143
[341014c]	144	for(int i = 0; i < Nbins; i++)
[1408174]	145	{
	146	width = (ptmax - ptmin) / Nbins;
[341014c]	147	ptdown = TMath::Power(10, ptmin + i * width);
	148	ptup = TMath::Power(10, ptmin + (i + 1) * width);
[1408174]	149	ptemp.set(ptdown, ptup, obj, xmin, xmax);
	150	histos->push_back(ptemp);
	151	}
	152	}
	153
[92237ca]	154	void HistogramsCollectionVsEta(std::vector<resolPlot> *histos, double etamin, double etamax, double ptmin, double ptmax, TString obj, double xmin = 0, double xmax = 2)
	155	{
	156	resolPlot ptemp;
	157	double width;
	158	double etadown;
	159	double etaup;
	160
[341014c]	161	for(int i = 0; i < Nbins; i++)
[92237ca]	162	{
	163	width = (etamax - etamin) / Nbins;
	164	etadown = etamin + i * width;
[341014c]	165	etaup = etamin + (i + 1) * width;
[92237ca]	166
	167	ptemp.set(etadown, etaup, ptmin, ptmax, obj, xmin, xmax);
	168	histos->push_back(ptemp);
	169	}
	170	}
	171
[1408174]	172	//------------------------------------------------------------------------------
	173
	174	class ExRootResult;
	175	class ExRootTreeReader;
	176
	177	//------------------------------------------------------------------------------
	178
[341014c]	179	void BinLogX(TH1 *h)
[1408174]	180	{
	181	TAxis *axis = h->GetXaxis();
	182	int bins = axis->GetNbins();
	183
	184	Axis_t from = axis->GetXmin();
	185	Axis_t to = axis->GetXmax();
	186	Axis_t width = (to - from) / bins;
	187	Axis_t *new_bins = new Axis_t[bins + 1];
	188
[341014c]	189	for(int i = 0; i <= bins; i++)
[1408174]	190	{
	191	new_bins[i] = TMath::Power(10, from + i * width);
	192	}
	193	axis->Set(bins, new_bins);
	194	delete new_bins;
	195	}
	196
	197	//------------------------------------------------------------------------------
	198
[341014c]	199	template <typename T>
	200	TH1D GetEffPt(TClonesArray branchReco, TClonesArray branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader treeReader)
[1408174]	201	{
	202
[341014c]	203	cout << "** Computing Efficiency of reconstructing " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
[1408174]	204
	205	Long64_t allEntries = treeReader->GetEntries();
	206
	207	GenParticle *particle;
	208	T *recoObj;
	209
	210	TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
	211
	212	Float_t deltaR;
	213	Float_t pt, eta;
	214	Long64_t entry;
	215
	216	Int_t i, j;
	217
[341014c]	218	TH1D *histGenPt = new TH1D(name + " gen spectra Pt", name + " gen spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
	219	TH1D *histRecoPt = new TH1D(name + " reco spectra Pt", name + " reco spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
[1408174]	220
[92237ca]	221	histGenPt->SetDirectory(0);
	222	histRecoPt->SetDirectory(0);
[1408174]	223
[92237ca]	224	BinLogX(histGenPt);
	225	BinLogX(histRecoPt);
[1408174]	226
	227	// Loop over all events
	228	for(entry = 0; entry < allEntries; ++entry)
	229	{
	230	// Load selected branches with data from specified event
	231	treeReader->ReadEntry(entry);
	232
	233	// Loop over all generated particle in event
	234	for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
	235	{
	236
[341014c]	237	particle = (GenParticle *)branchParticle->At(i);
[1408174]	238	genMomentum = particle->P4();
	239
	240	deltaR = 999;
	241
[341014c]	242	pt = genMomentum.Pt();
[92237ca]	243	eta = TMath::Abs(genMomentum.Eta());
	244
[341014c]	245	if(eta > etamax \|\| eta < etamin) continue;
[cb159ed]	246
[341014c]	247	if(particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax)
[534d945]	248	//if (TMath::Abs(particle->PID) == pdgID && (particle->Status>20 && particle->Status <30) && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax )
[92237ca]	249	{
[1408174]	250	// Loop over all reco object in event
	251	for(j = 0; j < branchReco->GetEntriesFast(); ++j)
	252	{
[341014c]	253	recoObj = (T *)branchReco->At(j);
[1408174]	254	recoMomentum = recoObj->P4();
[341014c]	255	//if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
[06a83b3]	256
[1408174]	257	// take the closest parton candidate
	258	if(TMath::Abs(pdgID) == 5)
	259	{
	260	Jet jet = (Jet )recoObj;
[341014c]	261	if(!(jet->BTag & (1 << 0))) continue;
[92237ca]	262
	263	//if(jet->BTag != ) continue;
	264	}
	265
[341014c]	266	if(TMath::Abs(pdgID) == 4)
[92237ca]	267	{
	268	Jet jet = (Jet )recoObj;
[341014c]	269	if(!(jet->BTag & (1 << 0))) continue;
[92237ca]	270	}
[cb159ed]	271
[341014c]	272	if(TMath::Abs(pdgID) == 1)
[92237ca]	273	{
	274	Jet jet = (Jet )recoObj;
[341014c]	275	if(!(jet->BTag & (1 << 0))) continue;
[1408174]	276	}
[92237ca]	277
[1408174]	278	if(TMath::Abs(pdgID) == 15)
	279	{
	280	Jet jet = (Jet )recoObj;
	281	if(jet->TauTag != 1) continue;
	282	}
[92237ca]	283
[1408174]	284	if(genMomentum.DeltaR(recoMomentum) < deltaR)
	285	{
	286	deltaR = genMomentum.DeltaR(recoMomentum);
	287	bestRecoMomentum = recoMomentum;
	288	}
	289	}
[92237ca]	290	histGenPt->Fill(pt);
[341014c]	291	if(deltaR < 0.3 && bestRecoMomentum.Pt() > 0.20 * pt)
	292	{
	293	histRecoPt->Fill(pt);
	294	}
[92237ca]	295	}
	296	}
	297	}
	298
	299	histRecoPt->Sumw2();
	300	histGenPt->Sumw2();
	301
	302	histRecoPt->Divide(histGenPt);
	303	histRecoPt->Scale(100.);
	304
	305	return histRecoPt;
	306	}
	307
[341014c]	308	template <typename T>
	309	TH1D GetEffEta(TClonesArray branchReco, TClonesArray branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader treeReader)
[92237ca]	310	{
	311
[341014c]	312	cout << "** Computing Efficiency of reconstructing " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
[92237ca]	313
	314	Long64_t allEntries = treeReader->GetEntries();
	315
	316	GenParticle *particle;
	317	T *recoObj;
	318
	319	TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
	320
	321	Float_t deltaR;
	322	Float_t pt, eta;
	323	Long64_t entry;
	324
	325	Int_t i, j;
	326
[341014c]	327	TH1D *histGenEta = new TH1D(name + " gen spectra Eta", name + " gen spectra", Nbins, etamin, etamax);
	328	TH1D *histRecoEta = new TH1D(name + " reco spectra Eta", name + " reco spectra", Nbins, etamin, etamax);
[92237ca]	329
	330	histGenEta->SetDirectory(0);
	331	histRecoEta->SetDirectory(0);
	332
	333	// Loop over all events
	334	for(entry = 0; entry < allEntries; ++entry)
	335	{
	336	// Load selected branches with data from specified event
	337	treeReader->ReadEntry(entry);
	338
	339	// Loop over all generated particle in event
	340	for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
	341	{
	342
[341014c]	343	particle = (GenParticle *)branchParticle->At(i);
[92237ca]	344	genMomentum = particle->P4();
	345
	346	deltaR = 999;
	347
[341014c]	348	pt = genMomentum.Pt();
[92237ca]	349	eta = genMomentum.Eta();
	350
[341014c]	351	if(pt > ptmax \|\| pt < ptmin) continue;
[92237ca]	352
[341014c]	353	if(particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax)
[534d945]	354	//if (TMath::Abs(particle->PID) == pdgID && (particle->Status>20 && particle->Status <30) && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax )
[92237ca]	355	{
	356	// Loop over all reco object in event
	357	for(j = 0; j < branchReco->GetEntriesFast(); ++j)
[1408174]	358	{
[341014c]	359	recoObj = (T *)branchReco->At(j);
[92237ca]	360	recoMomentum = recoObj->P4();
	361	// this is simply to avoid warnings from initial state particle
	362	// having infite rapidity ...
[341014c]	363	//if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
[1408174]	364
[92237ca]	365	// take the closest parton candidate
	366	if(TMath::Abs(pdgID) == 5)
	367	{
	368	Jet jet = (Jet )recoObj;
[341014c]	369	if(!(jet->BTag & (1 << 0))) continue;
[92237ca]	370	}
	371
	372	if(TMath::Abs(pdgID) == 4)
	373	{
	374	Jet jet = (Jet )recoObj;
[341014c]	375	if(!(jet->BTag & (1 << 0))) continue;
[92237ca]	376	}
	377
[341014c]	378	if(TMath::Abs(pdgID) == 1)
[92237ca]	379	{
	380	Jet jet = (Jet )recoObj;
[341014c]	381	if(!(jet->BTag & (1 << 0))) continue;
[92237ca]	382	}
	383
	384	if(TMath::Abs(pdgID) == 15)
	385	{
	386	Jet jet = (Jet )recoObj;
	387	if(jet->TauTag != 1) continue;
	388	}
	389	if(genMomentum.DeltaR(recoMomentum) < deltaR)
	390	{
	391	deltaR = genMomentum.DeltaR(recoMomentum);
	392	bestRecoMomentum = recoMomentum;
	393	}
[1408174]	394	}
[92237ca]	395
	396	histGenEta->Fill(eta);
[341014c]	397	if(deltaR < 0.3)
	398	{
	399	histRecoEta->Fill(eta);
	400	}
[1408174]	401	}
	402	}
	403	}
	404
[92237ca]	405	histRecoEta->Sumw2();
	406	histGenEta->Sumw2();
[1408174]	407
[92237ca]	408	histRecoEta->Divide(histGenEta);
	409	histRecoEta->Scale(100.);
[1408174]	410
[92237ca]	411	return histRecoEta;
	412	}
[1408174]	413
[534d945]	414	//------------------------------------------------------------------------------
	415
[341014c]	416	template <typename T>
	417	TH1D GetJetEffPt(TClonesArray branchJet, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
[534d945]	418	{
	419
[341014c]	420	cout << "** Computing Efficiency of reconstructing " << branchJet->GetName() << " with PID " << pdgID << endl;
[534d945]	421
	422	Long64_t allEntries = treeReader->GetEntries();
	423
	424	T *recoObj;
	425
	426	TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
	427
	428	Float_t pt, eta;
	429	Long64_t entry;
	430
	431	Int_t j;
	432
[341014c]	433	TH1D *histGenPt = new TH1D(name + " gen spectra Pt", name + " gen spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
	434	TH1D *histRecoPt = new TH1D(name + " reco spectra Pt", name + " reco spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
[534d945]	435
	436	histGenPt->SetDirectory(0);
	437	histRecoPt->SetDirectory(0);
	438
	439	BinLogX(histGenPt);
	440	BinLogX(histRecoPt);
	441
	442	// Loop over all events
	443	for(entry = 0; entry < allEntries; ++entry)
	444	{
	445	// Load selected branches with data from specified event
	446	treeReader->ReadEntry(entry);
	447
	448	// Loop over all reco object in event
	449	for(j = 0; j < branchJet->GetEntriesFast(); ++j)
	450	{
[341014c]	451	recoObj = (T *)branchJet->At(j);
[534d945]	452	recoMomentum = recoObj->P4();
[341014c]	453	pt = recoMomentum.Pt();
[534d945]	454	eta = TMath::Abs(recoMomentum.Eta());
	455	Jet jet = (Jet )recoObj;
[341014c]	456
	457	if(eta > etamax \|\| eta < etamin) continue;
[534d945]	458	if(pt < ptmin \|\| pt > ptmax) continue;
[341014c]	459
[534d945]	460	Int_t flavor = jet->Flavor;
	461	if(flavor == 21) flavor = 0;
[341014c]	462
[534d945]	463	if(TMath::Abs(pdgID) == 1)
	464	{
[341014c]	465
[534d945]	466	if(flavor < 4)
	467	{
	468	histGenPt->Fill(pt);
[341014c]	469	if(jet->BTag & (1 << 0)) histRecoPt->Fill(pt);
[534d945]	470	}
	471	}
	472	if(TMath::Abs(pdgID) == 4)
	473	{
	474	if(flavor == 4)
	475	{
	476	histGenPt->Fill(pt);
[341014c]	477	if(jet->BTag & (1 << 0)) histRecoPt->Fill(pt);
[534d945]	478	}
	479	}
	480	if(TMath::Abs(pdgID) == 5)
	481	{
	482	if(flavor == 5)
	483	{
	484	histGenPt->Fill(pt);
[341014c]	485	if(jet->BTag & (1 << 0)) histRecoPt->Fill(pt);
[534d945]	486	}
	487	}
[341014c]	488	}
[534d945]	489	}
	490
	491	histRecoPt->Sumw2();
	492	histGenPt->Sumw2();
	493
	494	histRecoPt->Divide(histGenPt);
	495	histRecoPt->Scale(100.);
	496
	497	return histRecoPt;
	498	}
	499
	500	// ------------------------------------------------------------------------------------------------------------------------------------------------------
	501
[341014c]	502	template <typename T>
	503	TH1D GetJetEffEta(TClonesArray branchJet, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader *treeReader)
[534d945]	504	{
	505
[341014c]	506	cout << "** Computing Efficiency of reconstructing " << branchJet->GetName() << " with PID " << pdgID << endl;
[534d945]	507
	508	Long64_t allEntries = treeReader->GetEntries();
	509
	510	T *recoObj;
	511
	512	TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
	513
	514	Float_t pt, eta;
	515	Long64_t entry;
	516
	517	Int_t j;
	518
[341014c]	519	TH1D *histGenEta = new TH1D(name + " gen spectra Eta", name + " gen spectra", Nbins, etamin, etamax);
	520	TH1D *histRecoEta = new TH1D(name + " reco spectra Eta", name + " reco spectra", Nbins, etamin, etamax);
[534d945]	521
	522	histGenEta->SetDirectory(0);
	523	histRecoEta->SetDirectory(0);
	524	// Loop over all events
	525	for(entry = 0; entry < allEntries; ++entry)
	526	{
	527	// Load selected branches with data from specified event
	528	treeReader->ReadEntry(entry);
	529
	530	// Loop over all reco object in event
	531	for(j = 0; j < branchJet->GetEntriesFast(); ++j)
	532	{
[341014c]	533	recoObj = (T *)branchJet->At(j);
[534d945]	534	recoMomentum = recoObj->P4();
[341014c]	535	pt = recoMomentum.Pt();
[534d945]	536	eta = recoMomentum.Eta();
	537	Jet jet = (Jet )recoObj;
[341014c]	538
	539	if(eta > etamax \|\| eta < etamin) continue;
[534d945]	540	if(pt < ptmin \|\| pt > ptmax) continue;
[341014c]	541
[534d945]	542	Int_t flavor = jet->Flavor;
	543	if(flavor == 21) flavor = 0;
[341014c]	544
[534d945]	545	if(TMath::Abs(pdgID) == 1)
	546	{
[341014c]	547	if(flavor == 1 \|\| flavor == 21)
[534d945]	548	{
	549	histGenEta->Fill(eta);
[341014c]	550	if(jet->BTag & (1 << 0)) histRecoEta->Fill(eta);
[534d945]	551	}
	552	}
	553	if(TMath::Abs(pdgID) == 4)
	554	{
	555	if(flavor == 4)
	556	{
	557	histGenEta->Fill(eta);
[341014c]	558	if(jet->BTag & (1 << 0)) histRecoEta->Fill(eta);
[534d945]	559	}
	560	}
	561	if(TMath::Abs(pdgID) == 5)
	562	{
	563	if(flavor == 5)
	564	{
	565	histGenEta->Fill(eta);
[341014c]	566	if(jet->BTag & (1 << 0)) histRecoEta->Fill(eta);
[534d945]	567	}
	568	}
[341014c]	569	}
[534d945]	570	}
	571
	572	histRecoEta->Sumw2();
	573	histGenEta->Sumw2();
	574
	575	histRecoEta->Divide(histGenEta);
	576	histRecoEta->Scale(100.);
	577
	578	return histRecoEta;
	579	}
	580
	581	// -----------------------------------------------------------------------------------------------------------------------------------------------------
[1408174]	582
[341014c]	583	template <typename T>
	584	TH1D GetTauEffPt(TClonesArray branchReco, TClonesArray branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader treeReader)
[1408174]	585	{
	586
[341014c]	587	cout << "** Computing Efficiency of reconstructing " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
[92237ca]	588
	589	Long64_t allEntries = treeReader->GetEntries();
[1408174]	590
	591	GenParticle *particle;
[92237ca]	592	T *recoObj;
[1408174]	593
[92237ca]	594	TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
[1408174]	595
	596	Float_t deltaR;
	597	Float_t pt, eta;
	598	Long64_t entry;
	599
[92237ca]	600	Int_t i, j;
	601
[341014c]	602	TH1D *histGenPt = new TH1D(name + " gen spectra Pt", name + " gen spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
	603	TH1D *histRecoPt = new TH1D(name + " reco spectra Pt", name + " reco spectra cen", Nbins, TMath::Log10(ptmin), TMath::Log10(ptmax));
[92237ca]	604
	605	histGenPt->SetDirectory(0);
	606	histRecoPt->SetDirectory(0);
	607
	608	BinLogX(histGenPt);
	609	BinLogX(histRecoPt);
[1408174]	610
	611	// Loop over all events
	612	for(entry = 0; entry < allEntries; ++entry)
	613	{
	614	// Load selected branches with data from specified event
	615	treeReader->ReadEntry(entry);
	616
[92237ca]	617	// Loop over all generated particle in event
	618	for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
[1408174]	619	{
[92237ca]	620
[341014c]	621	particle = (GenParticle *)branchParticle->At(i);
[92237ca]	622	genMomentum = particle->P4();
[1408174]	623
	624	deltaR = 999;
	625
[341014c]	626	pt = genMomentum.Pt();
[92237ca]	627	eta = TMath::Abs(genMomentum.Eta());
[1408174]	628
[341014c]	629	if(eta > etamax \|\| eta < etamin) continue;
[92237ca]	630
[341014c]	631	if(particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax)
[92237ca]	632	{
	633	// Loop over all reco object in event
	634	for(j = 0; j < branchReco->GetEntriesFast(); ++j)
	635	{
[341014c]	636	recoObj = (T *)branchReco->At(j);
[92237ca]	637	recoMomentum = recoObj->P4();
[1408174]	638	// this is simply to avoid warnings from initial state particle
	639	// having infite rapidity ...
[92237ca]	640	//if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
	641
	642	if(TMath::Abs(pdgID) == 1)
	643	{
	644	Jet jet = (Jet )recoObj;
[341014c]	645	if(jet->TauTag != 1) continue;
[92237ca]	646	}
	647
	648	if(TMath::Abs(pdgID) == 15)
	649	{
	650	Jet jet = (Jet )recoObj;
	651	if(jet->TauTag != 1) continue;
	652	}
[1408174]	653
	654	if(genMomentum.DeltaR(recoMomentum) < deltaR)
	655	{
[92237ca]	656	deltaR = genMomentum.DeltaR(recoMomentum);
	657	bestRecoMomentum = recoMomentum;
[1408174]	658	}
	659	}
[92237ca]	660
	661	histGenPt->Fill(pt);
[341014c]	662	if(deltaR < 0.3)
	663	{
	664	histRecoPt->Fill(pt);
	665	}
[1408174]	666	}
[92237ca]	667	}
	668	}
[1408174]	669
[92237ca]	670	histRecoPt->Sumw2();
	671	histGenPt->Sumw2();
[1408174]	672
[92237ca]	673	histRecoPt->Divide(histGenPt);
	674	histRecoPt->Scale(100.);
[341014c]	675	if(TMath::Abs(pdgID) == 15) histRecoPt->Scale(1 / 0.648);
[92237ca]	676
	677	return histRecoPt;
	678	}
	679
[341014c]	680	template <typename T>
	681	TH1D GetTauEffEta(TClonesArray branchReco, TClonesArray branchParticle, TString name, int pdgID, double ptmin, double ptmax, double etamin, double etamax, ExRootTreeReader treeReader)
[92237ca]	682	{
	683
[341014c]	684	cout << "** Computing Efficiency of reconstructing " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
[92237ca]	685
	686	Long64_t allEntries = treeReader->GetEntries();
	687
	688	GenParticle *particle;
	689	T *recoObj;
	690
	691	TLorentzVector recoMomentum, genMomentum, bestRecoMomentum;
	692
	693	Float_t deltaR;
	694	Float_t pt, eta;
	695	Long64_t entry;
	696
	697	Int_t i, j;
	698
[341014c]	699	TH1D *histGenEta = new TH1D(name + " gen spectra Eta", name + " gen spectra", Nbins, etamin, etamax);
	700	TH1D *histRecoEta = new TH1D(name + " reco spectra Eta", name + " reco spectra", Nbins, etamin, etamax);
[92237ca]	701
	702	histGenEta->SetDirectory(0);
	703	histRecoEta->SetDirectory(0);
	704
	705	// Loop over all events
	706	for(entry = 0; entry < allEntries; ++entry)
	707	{
	708	// Load selected branches with data from specified event
	709	treeReader->ReadEntry(entry);
	710
	711	// Loop over all generated particle in event
	712	for(i = 0; i < branchParticle->GetEntriesFast(); ++i)
	713	{
	714
[341014c]	715	particle = (GenParticle *)branchParticle->At(i);
[92237ca]	716	genMomentum = particle->P4();
	717
	718	deltaR = 999;
	719
[341014c]	720	pt = genMomentum.Pt();
[92237ca]	721	eta = genMomentum.Eta();
	722
[341014c]	723	if(pt > ptmax \|\| pt < ptmin) continue;
[92237ca]	724
[341014c]	725	if(particle->PID == pdgID && genMomentum.Pt() > ptmin && genMomentum.Pt() < ptmax)
[92237ca]	726	{
	727	// Loop over all reco object in event
	728	for(j = 0; j < branchReco->GetEntriesFast(); ++j)
[1408174]	729	{
[341014c]	730	recoObj = (T *)branchReco->At(j);
[92237ca]	731	recoMomentum = recoObj->P4();
	732	// this is simply to avoid warnings from initial state particle
	733	// having infite rapidity ...
[341014c]	734	//if(Momentum.Px() == 0 && genMomentum.Py() == 0) continue;
[92237ca]	735
[341014c]	736	if(TMath::Abs(pdgID) == 1)
[92237ca]	737	{
	738	Jet jet = (Jet )recoObj;
[341014c]	739	if(jet->TauTag != 1) continue;
[92237ca]	740	}
	741
	742	if(TMath::Abs(pdgID) == 15)
[1408174]	743	{
[92237ca]	744	Jet jet = (Jet )recoObj;
	745	if(jet->TauTag != 1) continue;
	746	}
	747
	748	if(genMomentum.DeltaR(recoMomentum) < deltaR)
	749	{
	750	deltaR = genMomentum.DeltaR(recoMomentum);
	751	bestRecoMomentum = recoMomentum;
[1408174]	752	}
	753	}
[92237ca]	754
	755	histGenEta->Fill(eta);
[341014c]	756	if(deltaR < 0.3)
	757	{
	758	histRecoEta->Fill(eta);
	759	}
[1408174]	760	}
	761	}
	762	}
[92237ca]	763
	764	histRecoEta->Sumw2();
	765	histGenEta->Sumw2();
	766
	767	histRecoEta->Divide(histGenEta);
	768	histRecoEta->Scale(100.);
[341014c]	769	if(TMath::Abs(pdgID) == 15) histRecoEta->Scale(1 / 0.648);
[92237ca]	770
	771	return histRecoEta;
[1408174]	772	}
[32677fb]	773
[341014c]	774	template <typename T>
[92237ca]	775	void GetPtres(std::vector<resolPlot> histos, TClonesArray branchReco, TClonesArray branchParticle, int pdgID, Double_t etaMin, Double_t etaMax, ExRootTreeReader treeReader)
[32677fb]	776	{
	777	Long64_t allEntries = treeReader->GetEntries();
	778
	779	cout << "** Computing pt resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
	780
	781	GenParticle *particle;
[341014c]	782	T *recoObj;
[32677fb]	783
	784	TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
	785
	786	Float_t deltaR;
	787	Float_t pt, eta;
	788	Long64_t entry;
	789
	790	Int_t i, j, bin;
	791
	792	// Loop over all events
	793	for(entry = 0; entry < allEntries; ++entry)
	794	{
	795	// Load selected branches with data from specified event
	796	treeReader->ReadEntry(entry);
	797
	798	// Loop over all reconstructed jets in event
	799	for(i = 0; i < branchReco->GetEntriesFast(); ++i)
	800	{
[341014c]	801	recoObj = (T *)branchReco->At(i);
[32677fb]	802	recoMomentum = recoObj->P4();
	803
	804	deltaR = 999;
	805
[341014c]	806	// Loop over all hard partons in event
	807	for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
	808	{
	809	particle = (GenParticle *)branchParticle->At(j);
	810	if(particle->PID == pdgID && particle->Status == 1)
[32677fb]	811	{
	812	genMomentum = particle->P4();
	813
	814	// this is simply to avoid warnings from initial state particle
	815	// having infite rapidity ...
	816	if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
	817
	818	// take the closest parton candidate
	819	if(genMomentum.DeltaR(recoMomentum) < deltaR)
	820	{
[341014c]	821	deltaR = genMomentum.DeltaR(recoMomentum);
	822	bestGenMomentum = genMomentum;
[32677fb]	823	}
	824	}
	825	}
	826
	827	if(deltaR < 0.3)
	828	{
[341014c]	829	pt = bestGenMomentum.Pt();
[32677fb]	830	eta = TMath::Abs(bestGenMomentum.Eta());
	831
[341014c]	832	for(bin = 0; bin < Nbins; bin++)
[32677fb]	833	{
[92237ca]	834	if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta > etaMin && eta < etaMax)
[32677fb]	835	{
[341014c]	836	histos->at(bin).resolHist->Fill(recoMomentum.Pt() / bestGenMomentum.Pt());
[32677fb]	837	}
	838	}
	839	}
	840	}
	841	}
	842	}
	843
[341014c]	844	template <typename T>
[92237ca]	845	void GetEres(std::vector<resolPlot> histos, TClonesArray branchReco, TClonesArray branchParticle, int pdgID, Double_t etaMin, Double_t etaMax, ExRootTreeReader treeReader)
[1408174]	846	{
	847	Long64_t allEntries = treeReader->GetEntries();
	848
[92237ca]	849	cout << "** Computing e resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
	850
	851	GenParticle *particle;
[341014c]	852	T *recoObj;
[92237ca]	853
	854	TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
	855
	856	Float_t deltaR;
	857	Float_t e, eta;
	858	Long64_t entry;
	859
	860	Int_t i, j, bin;
	861
	862	// Loop over all events
	863	for(entry = 0; entry < allEntries; ++entry)
	864	{
	865	// Load selected branches with data from specified event
	866	treeReader->ReadEntry(entry);
	867
[341014c]	868	// Loop over all reconstructed jets in event
[92237ca]	869	for(i = 0; i < branchReco->GetEntriesFast(); ++i)
	870	{
[341014c]	871	recoObj = (T *)branchReco->At(i);
[92237ca]	872	recoMomentum = recoObj->P4();
	873
	874	deltaR = 999;
	875
[341014c]	876	// Loop over all hard partons in event
	877	for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
	878	{
	879	particle = (GenParticle *)branchParticle->At(j);
	880	if(particle->PID == pdgID && particle->Status == 1)
[92237ca]	881	{
	882	genMomentum = particle->P4();
	883
	884	// this is simply to avoid warnings from initial state particle
	885	// having infite rapidity ...
	886	if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
	887
	888	// take the closest parton candidate
	889	if(genMomentum.DeltaR(recoMomentum) < deltaR)
	890	{
[341014c]	891	deltaR = genMomentum.DeltaR(recoMomentum);
	892	bestGenMomentum = genMomentum;
[92237ca]	893	}
	894	}
	895	}
	896
	897	if(deltaR < 0.3)
	898	{
[341014c]	899	e = bestGenMomentum.E();
[92237ca]	900	eta = TMath::Abs(bestGenMomentum.Eta());
	901
[341014c]	902	for(bin = 0; bin < Nbins; bin++)
[92237ca]	903	{
	904	if(e > histos->at(bin).ptmin && e < histos->at(bin).ptmax && eta > etaMin && eta < etaMax)
	905	{
[341014c]	906	histos->at(bin).resolHist->Fill(recoMomentum.E() / bestGenMomentum.E());
[92237ca]	907	}
	908	}
	909	}
	910	}
	911	}
	912	}
	913
[341014c]	914	template <typename T>
[92237ca]	915	void GetPtresVsEta(std::vector<resolPlot> histos, TClonesArray branchReco, TClonesArray branchParticle, int pdgID, Double_t ptMin, Double_t ptMax, ExRootTreeReader treeReader)
	916	{
	917	Long64_t allEntries = treeReader->GetEntries();
	918
	919	cout << "** Computing pt resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
	920
	921	GenParticle *particle;
[341014c]	922	T *recoObj;
[92237ca]	923
	924	TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
	925
	926	Float_t deltaR;
	927	Float_t pt, eta;
	928	Long64_t entry;
	929
	930	Int_t i, j, bin;
	931
	932	// Loop over all events
	933	for(entry = 0; entry < allEntries; ++entry)
	934	{
	935	// Load selected branches with data from specified event
	936	treeReader->ReadEntry(entry);
	937
	938	// Loop over all reconstructed jets in event
	939	for(i = 0; i < branchReco->GetEntriesFast(); ++i)
	940	{
[341014c]	941	recoObj = (T *)branchReco->At(i);
[92237ca]	942	recoMomentum = recoObj->P4();
	943
	944	deltaR = 999;
	945
[341014c]	946	// Loop over all hard partons in event
	947	for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
	948	{
	949	particle = (GenParticle *)branchParticle->At(j);
	950	if(particle->PID == pdgID && particle->Status == 1)
[92237ca]	951	{
	952	genMomentum = particle->P4();
	953
	954	// this is simply to avoid warnings from initial state particle
	955	// having infite rapidity ...
	956	if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
	957
	958	// take the closest parton candidate
	959	if(genMomentum.DeltaR(recoMomentum) < deltaR)
	960	{
[341014c]	961	deltaR = genMomentum.DeltaR(recoMomentum);
	962	bestGenMomentum = genMomentum;
[92237ca]	963	}
	964	}
	965	}
	966
	967	if(deltaR < 0.3)
	968	{
[341014c]	969	pt = bestGenMomentum.Pt();
[92237ca]	970	eta = bestGenMomentum.Eta();
	971
[341014c]	972	for(bin = 0; bin < Nbins; bin++)
[92237ca]	973	{
	974	if(eta > histos->at(bin).etamin && eta < histos->at(bin).etamax && pt > ptMin && pt < ptMax)
	975	{
[341014c]	976	histos->at(bin).resolHist->Fill(recoMomentum.Pt() / bestGenMomentum.Pt());
[92237ca]	977	}
	978	}
	979	}
	980	}
	981	}
	982	}
	983
[341014c]	984	template <typename T>
[92237ca]	985	void GetEresVsEta(std::vector<resolPlot> histos, TClonesArray branchReco, TClonesArray branchParticle, int pdgID, Double_t eMin, Double_t eMax, ExRootTreeReader treeReader)
	986	{
	987	Long64_t allEntries = treeReader->GetEntries();
	988
	989	cout << "** Computing E resolution of " << branchReco->GetName() << " induced by " << branchParticle->GetName() << " with PID " << pdgID << endl;
	990
	991	GenParticle *particle;
[341014c]	992	T *recoObj;
[92237ca]	993
	994	TLorentzVector recoMomentum, genMomentum, bestGenMomentum;
	995
	996	Float_t deltaR;
	997	Float_t e, eta;
	998	Long64_t entry;
	999
	1000	Int_t i, j, bin;
	1001
	1002	// Loop over all events
	1003	for(entry = 0; entry < allEntries; ++entry)
	1004	{
	1005	// Load selected branches with data from specified event
	1006	treeReader->ReadEntry(entry);
	1007
	1008	// Loop over all reconstructed jets in event
	1009	for(i = 0; i < branchReco->GetEntriesFast(); ++i)
	1010	{
[341014c]	1011	recoObj = (T *)branchReco->At(i);
[92237ca]	1012	recoMomentum = recoObj->P4();
	1013
	1014	deltaR = 999;
	1015
[341014c]	1016	// Loop over all hard partons in event
	1017	for(j = 0; j < branchParticle->GetEntriesFast(); ++j)
	1018	{
	1019	particle = (GenParticle *)branchParticle->At(j);
	1020	if(particle->PID == pdgID && particle->Status == 1)
[92237ca]	1021	{
	1022	genMomentum = particle->P4();
	1023
	1024	// this is simply to avoid warnings from initial state particle
	1025	// having infite rapidity ...
	1026	if(genMomentum.Px() == 0 && genMomentum.Py() == 0) continue;
	1027
	1028	// take the closest parton candidate
	1029	if(genMomentum.DeltaR(recoMomentum) < deltaR)
	1030	{
[341014c]	1031	deltaR = genMomentum.DeltaR(recoMomentum);
	1032	bestGenMomentum = genMomentum;
[92237ca]	1033	}
	1034	}
	1035	}
	1036
	1037	if(deltaR < 0.3)
	1038	{
[341014c]	1039	e = bestGenMomentum.E();
[92237ca]	1040	eta = bestGenMomentum.Eta();
	1041
[341014c]	1042	for(bin = 0; bin < Nbins; bin++)
[92237ca]	1043	{
	1044	if(eta > histos->at(bin).etamin && eta < histos->at(bin).etamax && e > eMin && e < eMax)
	1045	{
[341014c]	1046	histos->at(bin).resolHist->Fill(recoMomentum.E() / bestGenMomentum.E());
[92237ca]	1047	}
	1048	}
	1049	}
	1050	}
	1051	}
	1052	}
	1053
	1054	void GetJetsEres(std::vector<resolPlot> histos, TClonesArray branchJet, TClonesArray branchGenJet, ExRootTreeReader treeReader, Double_t etaMin, Double_t etaMax)
	1055	{
	1056
	1057	Long64_t allEntries = treeReader->GetEntries();
	1058
	1059	cout << "** Computing resolution of " << branchJet->GetName() << " induced by " << branchGenJet->GetName() << endl;
[1408174]	1060
	1061	Jet jet, genjet;
	1062
[92237ca]	1063	TLorentzVector jetMomentum, genJetMomentum, bestGenJetMomentum;
	1064
	1065	Float_t deltaR;
	1066	Float_t pt, eta;
	1067	Long64_t entry;
	1068
	1069	Int_t i, j, bin;
	1070
	1071	// Loop over all events
	1072	for(entry = 0; entry < allEntries; ++entry)
	1073	{
	1074	// Load selected branches with data from specified event
	1075	treeReader->ReadEntry(entry);
	1076
[341014c]	1077	if(entry % 10000 == 0) cout << "Event number: " << entry << endl;
[92237ca]	1078
	1079	// Loop over all reconstructed jets in event
[341014c]	1080	for(i = 0; i < TMath::Min(2, branchJet->GetEntriesFast()); ++i) //branchJet->GetEntriesFast(); ++i)
[92237ca]	1081	{
	1082
[341014c]	1083	jet = (Jet *)branchJet->At(i);
[92237ca]	1084	jetMomentum = jet->P4();
	1085
	1086	deltaR = 999;
	1087
[341014c]	1088	// Loop over all hard partons in event
	1089	for(j = 0; j < TMath::Min(2, branchGenJet->GetEntriesFast()); ++j)
	1090	{
	1091	genjet = (Jet *)branchGenJet->At(j);
[92237ca]	1092
	1093	genJetMomentum = genjet->P4();
	1094
	1095	// this is simply to avoid warnings from initial state particle
	1096	// having infite rapidity ...
	1097	if(genJetMomentum.Px() == 0 && genJetMomentum.Py() == 0) continue;
	1098
	1099	// take the closest parton candidate
	1100	if(genJetMomentum.DeltaR(jetMomentum) < deltaR)
	1101	{
[341014c]	1102	deltaR = genJetMomentum.DeltaR(jetMomentum);
	1103	bestGenJetMomentum = genJetMomentum;
[92237ca]	1104	}
	1105	}
	1106
	1107	if(deltaR < 0.3)
	1108	{
[341014c]	1109	pt = genJetMomentum.E();
[92237ca]	1110	eta = genJetMomentum.Eta();
	1111
[341014c]	1112	for(bin = 0; bin < Nbins; bin++)
[92237ca]	1113	{
	1114	if(pt > histos->at(bin).ptmin && pt < histos->at(bin).ptmax && eta < etaMax && eta > etaMin)
	1115	{
[341014c]	1116	histos->at(bin).resolHist->Fill(jetMomentum.E() / bestGenJetMomentum.E());
[92237ca]	1117	}
	1118	}
	1119	}
	1120	}
	1121	}
	1122	}
	1123
	1124	void GetJetsEresVsEta(std::vector<resolPlot> histos, TClonesArray branchJet, TClonesArray branchGenJet, ExRootTreeReader treeReader, Double_t eMin, Double_t eMax)
	1125	{
	1126
	1127	Long64_t allEntries = treeReader->GetEntries();
	1128
	1129	cout << "** Computing resolution of " << branchJet->GetName() << " induced by " << branchGenJet->GetName() << endl;
	1130
	1131	Jet jet, genjet;
	1132
	1133	TLorentzVector jetMomentum, genJetMomentum, bestGenJetMomentum;
	1134
	1135	Float_t deltaR;
	1136	Float_t pt, eta;
	1137	Long64_t entry;
	1138
	1139	Int_t i, j, bin;
	1140
	1141	// Loop over all events
	1142	for(entry = 0; entry < allEntries; ++entry)
	1143	{
	1144	// Load selected branches with data from specified event
	1145	treeReader->ReadEntry(entry);
	1146
[341014c]	1147	if(entry % 10000 == 0) cout << "Event number: " << entry << endl;
[92237ca]	1148
	1149	// Loop over all reconstructed jets in event
[341014c]	1150	for(i = 0; i < TMath::Min(2, branchJet->GetEntriesFast()); ++i) //branchJet->GetEntriesFast(); ++i)
[92237ca]	1151	{
	1152
[341014c]	1153	jet = (Jet *)branchJet->At(i);
[92237ca]	1154	jetMomentum = jet->P4();
	1155
	1156	deltaR = 999;
	1157
[341014c]	1158	// Loop over all hard partons in event
	1159	for(j = 0; j < TMath::Min(2, branchGenJet->GetEntriesFast()); ++j)
	1160	{
	1161	genjet = (Jet *)branchGenJet->At(j);
[92237ca]	1162
	1163	genJetMomentum = genjet->P4();
	1164
	1165	// this is simply to avoid warnings from initial state particle
	1166	// having infite rapidity ...
	1167	if(genJetMomentum.Px() == 0 && genJetMomentum.Py() == 0) continue;
	1168
	1169	// take the closest parton candidate
	1170	if(genJetMomentum.DeltaR(jetMomentum) < deltaR)
	1171	{
[341014c]	1172	deltaR = genJetMomentum.DeltaR(jetMomentum);
	1173	bestGenJetMomentum = genJetMomentum;
[92237ca]	1174	}
	1175	}
	1176
	1177	if(deltaR < 0.3)
	1178	{
	1179
[341014c]	1180	pt = genJetMomentum.E();
[92237ca]	1181	eta = genJetMomentum.Eta();
	1182
[341014c]	1183	for(bin = 0; bin < Nbins; bin++)
[92237ca]	1184	{
	1185	if(eta > histos->at(bin).etamin && eta < histos->at(bin).etamax && pt < eMax && pt > eMin)
	1186	{
[341014c]	1187	histos->at(bin).resolHist->Fill(jetMomentum.E() / bestGenJetMomentum.E());
[92237ca]	1188	}
	1189	}
	1190	}
	1191	}
	1192	}
	1193	}
	1194
[341014c]	1195	std::pair<Double_t, Double_t> GausFit(TH1 *hist)
[92237ca]	1196	{
[341014c]	1197	TF1 f1 = new TF1("f1", "gaus", hist->GetMean() - 2 hist->GetRMS(), hist->GetMean() + 2 * hist->GetRMS());
	1198	hist->Fit("f1", "RQ");
[92237ca]	1199
[341014c]	1200	TF1 f2 = new TF1("f2", "gaus", f1->GetParameter(1) - 2 f1->GetParameter(2), f1->GetParameter(1) + 2 * f1->GetParameter(2));
	1201	hist->Fit("f2", "RQ");
[92237ca]	1202
	1203	Double_t sig = f2->GetParameter(2);
	1204	Double_t sigErr = f2->GetParError(2);
	1205
	1206	delete f1;
	1207	delete f2;
[341014c]	1208	return make_pair(sig, sigErr);
[92237ca]	1209	}
	1210
	1211	TGraphErrors EresGraph(std::vector<resolPlot> *histos, bool rms = false)
	1212	{
	1213	Int_t bin;
	1214	Int_t count = 0;
[341014c]	1215	TGraphErrors gr = TGraphErrors(Nbins / 2);
[92237ca]	1216	double val, error;
[341014c]	1217	for(bin = 0; bin < Nbins; bin++)
[92237ca]	1218	{
[341014c]	1219	std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).resolHist);
	1220	if(rms == true)
	1221	{
	1222	gr.SetPoint(count, (histos->at(bin).ptmin + histos->at(bin).ptmax) / 2.0, 100 * histos->at(bin).resolHist->GetRMS());
	1223	//gr.SetPointError(count,0, 100*sigvalues.second); // to correct
	1224	error = 100 * histos->at(bin).resolHist->GetRMSError();
	1225	val = 100 * histos->at(bin).resolHist->GetRMS();
	1226	if(error > 0.2 * val) error = 0.2 * val;
	1227	gr.SetPointError(count, 0, error); // to correct
	1228	}
	1229	else
	1230	{
	1231
	1232	gr.SetPoint(count, (histos->at(bin).ptmin + histos->at(bin).ptmax) / 2.0, 100 * sigvalues.first);
	1233	error = 100 * sigvalues.second;
	1234	val = 100 * sigvalues.first;
	1235	if(error > 0.2 * val) error = 0.2 * val;
	1236	gr.SetPointError(count, 0, error); // to correct
	1237	//gr.SetPointError(count,0, 100*sigvalues.second);
	1238	}
	1239	count++;
[92237ca]	1240	}
	1241
	1242	return gr;
	1243	}
	1244
	1245	TGraphErrors MetResGraph(std::vector<resolPlot> *histos, bool rms = false)
	1246	{
	1247	Int_t bin;
	1248	Int_t count = 0;
[341014c]	1249	TGraphErrors gr = TGraphErrors(Nbins / 2);
[92237ca]	1250	double val, error;
[341014c]	1251	for(bin = 0; bin < Nbins; bin++)
[92237ca]	1252	{
[341014c]	1253	std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).resolHist);
	1254	if(rms == true)
	1255	{
	1256	gr.SetPoint(count, (histos->at(bin).ptmin + histos->at(bin).ptmax) / 2.0, histos->at(bin).resolHist->GetRMS());
	1257	error = histos->at(bin).resolHist->GetRMSError();
	1258	val = histos->at(bin).resolHist->GetRMS();
	1259	if(error > 0.2 * val) error = 0.2 * val;
	1260	gr.SetPointError(count, 0, error); // to correct
	1261	}
	1262	else
	1263	{
	1264
	1265	gr.SetPoint(count, (histos->at(bin).ptmin + histos->at(bin).ptmax) / 2.0, sigvalues.first);
	1266	val = sigvalues.first;
	1267	error = sigvalues.second;
	1268	if(error > 0.2 * val) error = 0.2 * val;
	1269	gr.SetPointError(count, 0, error);
	1270	//gr.SetPointError(count,0, 100*sigvalues.second);
	1271	}
	1272	count++;
[92237ca]	1273	}
	1274
	1275	return gr;
	1276	}
	1277
	1278	TGraphErrors EresGraphVsEta(std::vector<resolPlot> *histos, bool rms = false)
	1279	{
	1280	Int_t bin;
	1281	Int_t count = 0;
[341014c]	1282	TGraphErrors gr = TGraphErrors(Nbins / 2);
[92237ca]	1283	double val, error;
[341014c]	1284	for(bin = 0; bin < Nbins; bin++)
[92237ca]	1285	{
	1286
[341014c]	1287	std::pair<Double_t, Double_t> sigvalues = GausFit(histos->at(bin).resolHist);
	1288	if(rms == true)
	1289	{
	1290	gr.SetPoint(count, (histos->at(bin).etamin + histos->at(bin).etamax) / 2.0, histos->at(bin).resolHist->GetRMS());
	1291	error = 100 * histos->at(bin).resolHist->GetRMSError();
	1292	val = 100 * histos->at(bin).resolHist->GetRMS();
	1293	if(error > 0.2 * val) error = 0.2 * val;
	1294	gr.SetPointError(count, 0, error); // to correct
	1295	}
	1296	else
	1297	{
	1298	gr.SetPoint(count, (histos->at(bin).etamin + histos->at(bin).etamax) / 2.0, 100 * sigvalues.first);
	1299	val = 100 * sigvalues.first;
	1300	error = 100 * sigvalues.second;
	1301	if(error > 0.2 * val) error = 0.2 * val;
	1302	gr.SetPointError(count, 0, error);
	1303	//gr.SetPointError(count,0, 100*sigvalues.second);
	1304	}
	1305	count++;
[92237ca]	1306	}
	1307
	1308	return gr;
	1309	}
	1310
	1311	void GetMetres(std::vector<resolPlot> histos, TClonesArray branchScalarHT, TClonesArray branchMet, TClonesArray branchJet, ExRootTreeReader *treeReader)
	1312	{
	1313
	1314	Long64_t allEntries = treeReader->GetEntries();
	1315
	1316	cout << "** Computing resolution of " << branchMet->GetName() << " vs " << branchScalarHT->GetName() << endl;
	1317
	1318	MissingET *met;
	1319	ScalarHT *scalarHT;
	1320
	1321	Long64_t entry;
	1322
	1323	Int_t bin;
	1324	Double_t ht;
	1325
	1326	Jet *jet;
	1327	TLorentzVector p1, p2;
	1328
	1329	// Loop over all events
	1330	for(entry = 0; entry < allEntries; ++entry)
	1331	{
	1332	// Load selected branches with data from specified event
	1333	treeReader->ReadEntry(entry);
	1334
[341014c]	1335	if(entry % 10000 == 0) cout << "Event number: " << entry << endl;
[92237ca]	1336
[341014c]	1337	if(branchJet->GetEntriesFast() > 1)
[92237ca]	1338	{
	1339
[341014c]	1340	jet = (Jet *)branchJet->At(0);
[92237ca]	1341	p1 = jet->P4();
[341014c]	1342	jet = (Jet *)branchJet->At(1);
[92237ca]	1343	p2 = jet->P4();
	1344
[341014c]	1345	met = (MissingET *)branchMet->At(0);
	1346	scalarHT = (ScalarHT *)branchScalarHT->At(0);
[92237ca]	1347	ht = scalarHT->HT;
	1348
[341014c]	1349	if(p1.Pt() < 0.75 * ht / 2) continue;
	1350	if(p2.Pt() < 0.75 * ht / 2) continue;
[92237ca]	1351
[341014c]	1352	for(bin = 0; bin < Nbins; bin++)
[92237ca]	1353	{
[341014c]	1354	if(ht > histos->at(bin).ptmin && ht < histos->at(bin).ptmax)
[92237ca]	1355	{
	1356	histos->at(bin).resolHist->Fill(met->P4().Px());
	1357	}
	1358	}
	1359	}
	1360	}
	1361	}
	1362
	1363	//------------------------------------------------------------------------------
	1364
	1365	void addResoGraph(TMultiGraph mg, TGraphErrors gr, TLegend *leg, int style, Color_t color, TString text)
	1366	{
	1367
	1368	gr->SetLineWidth(2);
	1369	gr->SetLineColor(color);
	1370	gr->SetMarkerStyle(style);
	1371	gr->SetMarkerColor(color);
	1372	gr->SetMarkerSize(1.5);
	1373
	1374	std::cout << "Adding " << gr->GetName() << std::endl;
	1375	mg->Add(gr);
[341014c]	1376	leg->AddEntry(gr, text, "p");
[92237ca]	1377	}
	1378
	1379	void DrawAxis(TMultiGraph mg, TLegend leg, double xmin, double xmax, double ymin, double ymax, TString tx, TString ty, bool logx = 0, bool logy = 0)
	1380	{
	1381	mg->SetMinimum(ymin);
	1382	mg->SetMaximum(ymax);
[341014c]	1383	mg->GetXaxis()->SetLimits(xmin, xmax);
[92237ca]	1384
	1385	mg->GetXaxis()->SetTitle(tx);
	1386	mg->GetYaxis()->SetTitle(ty);
	1387
	1388	mg->GetYaxis()->SetTitleSize(0.07);
	1389	mg->GetXaxis()->SetTitleSize(0.07);
	1390	mg->GetYaxis()->SetLabelSize(0.06);
	1391	mg->GetXaxis()->SetLabelSize(0.06);
	1392	mg->GetYaxis()->SetLabelOffset(0.03);
	1393	mg->GetYaxis()->SetTitleOffset(1.4);
	1394	mg->GetXaxis()->SetTitleOffset(1.4);
	1395
	1396	mg->GetXaxis()->SetTitleFont(132);
	1397	mg->GetYaxis()->SetTitleFont(132);
	1398	mg->GetXaxis()->SetLabelFont(132);
	1399	mg->GetYaxis()->SetLabelFont(132);
	1400
	1401	mg->GetYaxis()->SetNdivisions(505);
	1402
	1403	leg->SetTextFont(132);
	1404	leg->SetBorderSize(0);
	1405	leg->SetShadowColor(0);
	1406	leg->SetFillColor(0);
	1407	leg->SetFillStyle(0);
	1408
	1409	gStyle->SetOptTitle(0);
	1410
	1411	if(logx) gPad->SetLogx();
	1412	if(logy) gPad->SetLogy();
	1413
	1414	//gPad->SetGridx();
	1415	//gPad->SetGridy();
	1416	gPad->SetBottomMargin(0.2);
	1417	gPad->SetLeftMargin(0.2);
	1418	gPad->Modified();
	1419	gPad->Update();
	1420	}
	1421
[181c061]	1422	void DelphesValidation(
	1423	const char *inputFilePion,
	1424	const char *inputFileElectron,
	1425	const char *inputFileMuon,
	1426	const char *inputFilePhoton,
	1427	const char *inputFileNeutralHadron,
	1428	const char *inputFileJet,
	1429	const char *inputFileBJet,
	1430	const char *inputFileCJet,
	1431	const char *inputFileTauJet,
	1432	const char *outputFile,
	1433	const char *version)
[92237ca]	1434	{
	1435
	1436	TChain *chainPion = new TChain("Delphes");
	1437	chainPion->Add(inputFilePion);
	1438	ExRootTreeReader *treeReaderPion = new ExRootTreeReader(chainPion);
	1439
	1440	TChain *chainElectron = new TChain("Delphes");
	1441	chainElectron->Add(inputFileElectron);
	1442	ExRootTreeReader *treeReaderElectron = new ExRootTreeReader(chainElectron);
	1443
	1444	TChain *chainMuon = new TChain("Delphes");
	1445	chainMuon->Add(inputFileMuon);
	1446	ExRootTreeReader *treeReaderMuon = new ExRootTreeReader(chainMuon);
	1447
	1448	TChain *chainPhoton = new TChain("Delphes");
	1449	chainPhoton->Add(inputFilePhoton);
	1450	ExRootTreeReader *treeReaderPhoton = new ExRootTreeReader(chainPhoton);
	1451
	1452	TChain *chainNeutralHadron = new TChain("Delphes");
	1453	chainNeutralHadron->Add(inputFileNeutralHadron);
	1454	ExRootTreeReader *treeReaderNeutralHadron = new ExRootTreeReader(chainNeutralHadron);
	1455
	1456	TChain *chainJet = new TChain("Delphes");
	1457	chainJet->Add(inputFileJet);
	1458	ExRootTreeReader *treeReaderJet = new ExRootTreeReader(chainJet);
	1459
	1460	TChain *chainBJet = new TChain("Delphes");
	1461	chainBJet->Add(inputFileBJet);
	1462	ExRootTreeReader *treeReaderBJet = new ExRootTreeReader(chainBJet);
	1463
	1464	TChain *chainCJet = new TChain("Delphes");
	1465	chainCJet->Add(inputFileCJet);
	1466	ExRootTreeReader *treeReaderCJet = new ExRootTreeReader(chainCJet);
	1467
	1468	TChain *chainTauJet = new TChain("Delphes");
	1469	chainTauJet->Add(inputFileTauJet);
	1470	ExRootTreeReader *treeReaderTauJet = new ExRootTreeReader(chainTauJet);
	1471
	1472	TClonesArray *branchParticleElectron = treeReaderElectron->UseBranch("Particle");
	1473	TClonesArray *branchTrackElectron = treeReaderElectron->UseBranch("Track");
	1474	TClonesArray *branchElectron = treeReaderElectron->UseBranch("Electron");
[b4ec6ac]	1475	TClonesArray *branchElectronPF = treeReaderElectron->UseBranch("ElectronPF");
[92237ca]	1476
	1477	TClonesArray *branchParticleMuon = treeReaderMuon->UseBranch("Particle");
	1478	TClonesArray *branchTrackMuon = treeReaderMuon->UseBranch("Track");
	1479	TClonesArray *branchMuon = treeReaderMuon->UseBranch("Muon");
	1480
	1481	TClonesArray *branchParticlePion = treeReaderPion->UseBranch("Particle");
	1482	TClonesArray *branchTrackPion = treeReaderPion->UseBranch("Track");
	1483	TClonesArray *branchPion = treeReaderPion->UseBranch("Pion");
	1484
	1485	TClonesArray *branchParticlePhoton = treeReaderPhoton->UseBranch("Particle");
	1486	TClonesArray *branchTowerPhoton = treeReaderPhoton->UseBranch("Tower");
	1487	TClonesArray *branchPhoton = treeReaderPhoton->UseBranch("Photon");
	1488
	1489	TClonesArray *branchParticleNeutralHadron = treeReaderNeutralHadron->UseBranch("Particle");
	1490	TClonesArray *branchTowerNeutralHadron = treeReaderNeutralHadron->UseBranch("Tower");
	1491
	1492	TClonesArray *branchGenJet = treeReaderJet->UseBranch("GenJet");
	1493	TClonesArray *branchParticleJet = treeReaderJet->UseBranch("Particle");
[13cb1b6]	1494	TClonesArray *branchPFJet = treeReaderJet->UseBranch("PFJet");
[92237ca]	1495	TClonesArray *branchCaloJet = treeReaderJet->UseBranch("CaloJet");
[06a83b3]	1496	TClonesArray *branchJet = treeReaderJet->UseBranch("Jet");
[cb159ed]	1497
[92237ca]	1498	TClonesArray *branchParticleBJet = treeReaderBJet->UseBranch("Particle");
	1499	TClonesArray *branchPFBJet = treeReaderBJet->UseBranch("Jet");
	1500
	1501	TClonesArray *branchParticleCJet = treeReaderCJet->UseBranch("Particle");
	1502	TClonesArray *branchPFCJet = treeReaderCJet->UseBranch("Jet");
	1503
	1504	TClonesArray *branchParticleTauJet = treeReaderTauJet->UseBranch("Particle");
	1505	TClonesArray *branchPFTauJet = treeReaderTauJet->UseBranch("Jet");
	1506
	1507	TClonesArray *branchGenScalarHT = treeReaderJet->UseBranch("GenScalarHT");
[13cb1b6]	1508	TClonesArray *branchMet = treeReaderJet->UseBranch("PFMissingET");
[92237ca]	1509	TClonesArray *branchCaloMet = treeReaderJet->UseBranch("CaloMissingET");
	1510
	1511	std::vector<Color_t> colors;
	1512
	1513	colors.push_back(kBlack);
	1514	colors.push_back(kBlue);
	1515	colors.push_back(kRed);
[341014c]	1516	colors.push_back(kGreen + 1);
	1517	colors.push_back(kMagenta + 1);
[92237ca]	1518	colors.push_back(kOrange);
	1519
	1520	std::vector<Int_t> markerStyles;
	1521
	1522	markerStyles.push_back(20);
	1523	markerStyles.push_back(21);
	1524	markerStyles.push_back(22);
	1525	markerStyles.push_back(23);
	1526	markerStyles.push_back(33);
	1527	markerStyles.push_back(34);
	1528
	1529	TString pdfOutput(outputFile);
	1530	pdfOutput.ReplaceAll(".root", ".pdf");
	1531
	1532	TString figPath = inputFilePion;
	1533	figPath.ReplaceAll(".root", "");
	1534	figPath.ReplaceAll("root", "www/fig");
	1535	Int_t lastSlash = figPath.Last('/');
	1536	Int_t sizePath = figPath.Length();
[341014c]	1537	figPath.Remove(lastSlash + 1, sizePath);
[92237ca]	1538
[85ad2b9]	1539	TString header = pdfOutput;
	1540	header.ReplaceAll(".pdf", "");
	1541	header.ReplaceAll("validation_", "");
	1542	lastSlash = header.Last('/');
	1543	sizePath = header.Length();
[341014c]	1544	header.Remove(0, lastSlash + 1);
[85ad2b9]	1545
	1546	TString vrs(version);
	1547
[341014c]	1548	TPaveText *pave = new TPaveText(0.0, 0.89, 0.94, 0.94, "NDC");
[bba3653]	1549	pave->SetTextAlign(30);
[85ad2b9]	1550	pave->SetTextFont(132);
	1551	pave->SetBorderSize(0);
	1552	pave->SetShadowColor(0);
	1553	pave->SetFillColor(0);
	1554	pave->SetFillStyle(0);
[341014c]	1555	pave->AddText("Delphes " + vrs + " - " + header);
[85ad2b9]	1556
[92237ca]	1557	TString s_etaMin, s_etaMax, s_eta, s_pt, s_e;
	1558
	1559	Double_t ptMin = 1.;
	1560	Double_t ptMax = 50000.;
	1561
	1562	Double_t etaMin = -6.;
	1563	Double_t etaMax = 6.;
	1564
	1565	std::vector<Double_t> ptVals;
	1566
	1567	ptVals.push_back(1.);
	1568	ptVals.push_back(10.);
	1569	ptVals.push_back(100.);
	1570	ptVals.push_back(1000.);
	1571	ptVals.push_back(10000.);
	1572
	1573	std::vector<Double_t> etaVals;
	1574
	1575	etaVals.push_back(0.);
	1576	etaVals.push_back(1.5);
	1577	etaVals.push_back(2.5);
	1578	etaVals.push_back(4.0);
	1579	etaVals.push_back(6.0);
	1580
[341014c]	1581	const int n_etabins = etaVals.size() - 1;
[92237ca]	1582	const int n_ptbins = ptVals.size();
	1583
	1584	//////////////////////////
	1585	// Tracking performance //
	1586	//////////////////////////
	1587
	1588	// --------- Pion Tracks --------- //
	1589
[341014c]	1590	TMultiGraph *mg_trkpi_res_pt = new TMultiGraph("", "");
	1591	TMultiGraph *mg_trkpi_eff_pt = new TMultiGraph("", "");
	1592	TMultiGraph *mg_trkpi_res_eta = new TMultiGraph("", "");
	1593	TMultiGraph *mg_trkpi_eff_eta = new TMultiGraph("", "");
[92237ca]	1594
[341014c]	1595	TLegend *leg_trkpi_res_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
	1596	TLegend leg_trkpi_eff_pt = (TLegend )leg_trkpi_res_pt->Clone();
	1597	TLegend leg_trkpi_res_eta = (TLegend )leg_trkpi_res_pt->Clone();
	1598	TLegend leg_trkpi_eff_eta = (TLegend )leg_trkpi_res_eta->Clone();
[92237ca]	1599
[341014c]	1600	TGraphErrors *gr_trkpi_res_pt = new TGraphErrors[n_etabins];
	1601	TGraphErrors *gr_trkpi_eff_pt = new TGraphErrors[n_etabins];
[dd5e213]	1602	TGraphErrors *gr_trkpi_res_eta = new TGraphErrors[n_ptbins];
	1603	TGraphErrors *gr_trkpi_eff_eta = new TGraphErrors[n_ptbins];
[341014c]	1604	TH1D h_trkpi_eff_pt, h_trkpi_eff_eta;
[92237ca]	1605
[341014c]	1606	std::vector<resolPlot> *plots_trkpi_res_pt = new std::vector<resolPlot>[n_etabins];
	1607	std::vector<resolPlot> *plots_trkpi_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]	1608
	1609	// loop over eta bins
[341014c]	1610	for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]	1611	{
[341014c]	1612	HistogramsCollection(&plots_trkpi_res_pt[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkpi");
	1613	GetPtres<Track>(&plots_trkpi_res_pt[k], branchTrackPion, branchParticlePion, 211, etaVals.at(k), etaVals.at(k + 1), treeReaderPion);
	1614	gr_trkpi_res_pt[k] = EresGraph(&plots_trkpi_res_pt[k]);
[92237ca]	1615
[341014c]	1616	h_trkpi_eff_pt = GetEffPt<Track>(branchTrackPion, branchParticlePion, "Pion", 211, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderPion);
	1617	gr_trkpi_eff_pt[k] = TGraphErrors(h_trkpi_eff_pt);
[92237ca]	1618
[341014c]	1619	s_etaMin = Form("%.1f", etaVals.at(k));
	1620	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]	1621
[341014c]	1622	s_eta = "#pi^{ #pm} , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	1623
[341014c]	1624	gr_trkpi_res_pt[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
	1625	gr_trkpi_eff_pt[k].SetName("trkEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]	1626
[341014c]	1627	addResoGraph(mg_trkpi_res_pt, &gr_trkpi_res_pt[k], leg_trkpi_res_pt, markerStyles.at(k), colors.at(k), s_eta);
	1628	addResoGraph(mg_trkpi_eff_pt, &gr_trkpi_eff_pt[k], leg_trkpi_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]	1629	}
	1630
	1631	// loop over pt
[341014c]	1632	for(k = 0; k < ptVals.size(); k++)
[92237ca]	1633	{
[341014c]	1634	HistogramsCollectionVsEta(&plots_trkpi_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkpi", 0.0, 2.0);
	1635	GetPtresVsEta<Track>(&plots_trkpi_res_eta[k], branchTrackPion, branchParticlePion, 211, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderPion);
	1636	gr_trkpi_res_eta[k] = EresGraphVsEta(&plots_trkpi_res_eta[k]);
[92237ca]	1637
[341014c]	1638	h_trkpi_eff_eta = GetEffEta<Track>(branchTrackPion, branchParticlePion, "Pion", 211, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderPion);
	1639	gr_trkpi_eff_eta[k] = TGraphErrors(h_trkpi_eff_eta);
[92237ca]	1640
[341014c]	1641	s_pt = Form("#pi^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
	1642	if(ptVals.at(k) >= 1000.) s_pt = Form("#pi^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	1643
[341014c]	1644	addResoGraph(mg_trkpi_res_eta, &gr_trkpi_res_eta[k], leg_trkpi_res_eta, markerStyles.at(k), colors.at(k), s_pt);
	1645	addResoGraph(mg_trkpi_eff_eta, &gr_trkpi_eff_eta[k], leg_trkpi_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
[92237ca]	1646	}
	1647
[341014c]	1648	TCanvas *c_trkpi_res_pt = new TCanvas("", "", 800, 600);
[92237ca]	1649
	1650	mg_trkpi_res_pt->Draw("APE");
	1651	DrawAxis(mg_trkpi_res_pt, leg_trkpi_res_pt, ptMin, ptMax, 0.01, 100, "p_{T} [GeV]", "(track resolution in p_{T})/p_{T} (%)", true, true);
	1652	leg_trkpi_res_pt->Draw();
[85ad2b9]	1653	pave->Draw();
[92237ca]	1654
[341014c]	1655	c_trkpi_res_pt->Print(pdfOutput + "(", "pdf");
	1656	c_trkpi_res_pt->Print(figPath + "img_trkpi_res_pt.pdf", "pdf");
	1657	c_trkpi_res_pt->Print(figPath + "img_trkpi_res_pt.png", "png");
[92237ca]	1658
[341014c]	1659	TCanvas *c_trkpi_res_eta = new TCanvas("", "", 800, 600);
[92237ca]	1660
	1661	mg_trkpi_res_eta->Draw("APE");
	1662	DrawAxis(mg_trkpi_res_eta, leg_trkpi_res_eta, etaMin, etaMax, 0.01, 100, " #eta ", "(track resolution in p_{T})/p_{T} (%)", false, true);
	1663	leg_trkpi_res_eta->Draw();
[85ad2b9]	1664	pave->Draw();
[92237ca]	1665
[341014c]	1666	c_trkpi_res_eta->Print(pdfOutput, "pdf");
	1667	c_trkpi_res_eta->Print(figPath + "img_trkpi_res_eta.pdf", "pdf");
	1668	c_trkpi_res_eta->Print(figPath + "img_trkpi_res_eta.png", "png");
[92237ca]	1669
[341014c]	1670	TCanvas *c_trkpi_eff_pt = new TCanvas("", "", 800, 600);
[92237ca]	1671
	1672	mg_trkpi_eff_pt->Draw("APE");
	1673	DrawAxis(mg_trkpi_eff_pt, leg_trkpi_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "tracking efficiency (%)", true, false);
	1674	leg_trkpi_eff_pt->Draw();
[85ad2b9]	1675	pave->Draw();
[92237ca]	1676
[341014c]	1677	c_trkpi_eff_pt->Print(pdfOutput, "pdf");
	1678	c_trkpi_eff_pt->Print(figPath + "img_trkpi_eff_pt.pdf", "pdf");
	1679	c_trkpi_eff_pt->Print(figPath + "img_trkpi_eff_pt.png", "png");
[92237ca]	1680
[341014c]	1681	TCanvas *c_trkpi_eff_eta = new TCanvas("", "", 800, 600);
[92237ca]	1682
	1683	mg_trkpi_eff_eta->Draw("APE");
	1684	DrawAxis(mg_trkpi_eff_eta, leg_trkpi_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "tracking efficiency (%)", false, false);
	1685	leg_trkpi_eff_eta->Draw();
[85ad2b9]	1686	pave->Draw();
	1687
[341014c]	1688	c_trkpi_eff_eta->Print(pdfOutput, "pdf");
	1689	c_trkpi_eff_eta->Print(figPath + "img_trkpi_eff_eta.pdf", "pdf");
	1690	c_trkpi_eff_eta->Print(figPath + "img_trkpi_eff_eta.png", "png");
[92237ca]	1691
	1692	// --------- Electron Tracks --------- //
	1693
[341014c]	1694	TMultiGraph *mg_trkele_res_pt = new TMultiGraph("", "");
	1695	TMultiGraph *mg_trkele_eff_pt = new TMultiGraph("", "");
	1696	TMultiGraph *mg_trkele_res_eta = new TMultiGraph("", "");
	1697	TMultiGraph *mg_trkele_eff_eta = new TMultiGraph("", "");
[92237ca]	1698
[341014c]	1699	TLegend *leg_trkele_res_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
	1700	TLegend leg_trkele_eff_pt = (TLegend )leg_trkele_res_pt->Clone();
	1701	TLegend leg_trkele_res_eta = (TLegend )leg_trkele_res_pt->Clone();
	1702	TLegend leg_trkele_eff_eta = (TLegend )leg_trkele_res_eta->Clone();
[92237ca]	1703
[341014c]	1704	TGraphErrors *gr_trkele_res_pt = new TGraphErrors[n_etabins];
	1705	TGraphErrors *gr_trkele_eff_pt = new TGraphErrors[n_etabins];
[dd5e213]	1706	TGraphErrors *gr_trkele_res_eta = new TGraphErrors[n_ptbins];
	1707	TGraphErrors *gr_trkele_eff_eta = new TGraphErrors[n_ptbins];
	1708
[341014c]	1709	TH1D h_trkele_eff_pt, h_trkele_eff_eta;
[92237ca]	1710
[341014c]	1711	std::vector<resolPlot> *plots_trkele_res_pt = new std::vector<resolPlot>[n_etabins];
	1712	std::vector<resolPlot> *plots_trkele_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]	1713
	1714	// loop over eta bins
[341014c]	1715	for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]	1716	{
[341014c]	1717	HistogramsCollection(&plots_trkele_res_pt[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkele");
	1718	GetPtres<Track>(&plots_trkele_res_pt[k], branchTrackElectron, branchParticleElectron, 11, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
	1719	gr_trkele_res_pt[k] = EresGraph(&plots_trkele_res_pt[k]);
[92237ca]	1720
[341014c]	1721	h_trkele_eff_pt = GetEffPt<Track>(branchTrackElectron, branchParticleElectron, "Electron", 11, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
	1722	gr_trkele_eff_pt[k] = TGraphErrors(h_trkele_eff_pt);
[92237ca]	1723
[341014c]	1724	s_etaMin = Form("%.1f", etaVals.at(k));
	1725	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]	1726
[341014c]	1727	s_eta = "e^{ #pm} , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	1728
[341014c]	1729	gr_trkele_res_pt[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
	1730	gr_trkele_eff_pt[k].SetName("trkEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]	1731
[341014c]	1732	addResoGraph(mg_trkele_res_pt, &gr_trkele_res_pt[k], leg_trkele_res_pt, markerStyles.at(k), colors.at(k), s_eta);
	1733	addResoGraph(mg_trkele_eff_pt, &gr_trkele_eff_pt[k], leg_trkele_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]	1734	}
	1735
	1736	// loop over pt
[341014c]	1737	for(k = 0; k < ptVals.size(); k++)
[92237ca]	1738	{
[341014c]	1739	HistogramsCollectionVsEta(&plots_trkele_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkele", 0.0, 2.0);
	1740	GetPtresVsEta<Track>(&plots_trkele_res_eta[k], branchTrackElectron, branchParticleElectron, 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderElectron);
	1741	gr_trkele_res_eta[k] = EresGraphVsEta(&plots_trkele_res_eta[k]);
[92237ca]	1742
[341014c]	1743	h_trkele_eff_eta = GetEffEta<Track>(branchTrackElectron, branchParticleElectron, "Electron", 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderElectron);
	1744	gr_trkele_eff_eta[k] = TGraphErrors(h_trkele_eff_eta);
[92237ca]	1745
[341014c]	1746	s_pt = Form("e^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
	1747	if(ptVals.at(k) >= 1000.) s_pt = Form("e^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	1748
[341014c]	1749	addResoGraph(mg_trkele_res_eta, &gr_trkele_res_eta[k], leg_trkele_res_eta, markerStyles.at(k), colors.at(k), s_pt);
	1750	addResoGraph(mg_trkele_eff_eta, &gr_trkele_eff_eta[k], leg_trkele_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
[92237ca]	1751	}
	1752
[341014c]	1753	TCanvas *c_trkele_res_pt = new TCanvas("", "", 800, 600);
[92237ca]	1754
	1755	mg_trkele_res_pt->Draw("APE");
	1756	DrawAxis(mg_trkele_res_pt, leg_trkele_res_pt, ptMin, ptMax, 0.01, 100, "p_{T} [GeV]", "(track resolution in p_{T})/p_{T} (%)", true, true);
	1757	leg_trkele_res_pt->Draw();
[85ad2b9]	1758	pave->Draw();
[92237ca]	1759
[341014c]	1760	c_trkele_res_pt->Print(pdfOutput, "pdf");
	1761	c_trkele_res_pt->Print(figPath + "img_trkele_res_pt.pdf", "pdf");
	1762	c_trkele_res_pt->Print(figPath + "img_trkele_res_pt.png", "png");
[92237ca]	1763
[341014c]	1764	TCanvas *c_trkele_res_eta = new TCanvas("", "", 800, 600);
[92237ca]	1765
	1766	mg_trkele_res_eta->Draw("APE");
	1767	DrawAxis(mg_trkele_res_eta, leg_trkele_res_eta, etaMin, etaMax, 0.01, 100, " #eta ", "(track resolution in p_{T})/p_{T} (%)", false, true);
	1768	leg_trkele_res_eta->Draw();
[85ad2b9]	1769	pave->Draw();
[92237ca]	1770
[341014c]	1771	c_trkele_res_eta->Print(pdfOutput, "pdf");
	1772	c_trkele_res_eta->Print(figPath + "img_trkele_res_eta.pdf", "pdf");
	1773	c_trkele_res_eta->Print(figPath + "img_trkele_res_eta.png", "png");
[92237ca]	1774
[341014c]	1775	TCanvas *c_trkele_eff_pt = new TCanvas("", "", 800, 600);
[92237ca]	1776
	1777	mg_trkele_eff_pt->Draw("APE");
	1778	DrawAxis(mg_trkele_eff_pt, leg_trkele_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "tracking efficiency (%)", true, false);
	1779	leg_trkele_eff_pt->Draw();
[85ad2b9]	1780	pave->Draw();
[92237ca]	1781
[341014c]	1782	c_trkele_eff_pt->Print(pdfOutput, "pdf");
	1783	c_trkele_eff_pt->Print(figPath + "img_trkele_eff_pt.pdf", "pdf");
	1784	c_trkele_eff_pt->Print(figPath + "img_trkele_eff_pt.png", "png");
[92237ca]	1785
[341014c]	1786	TCanvas *c_trkele_eff_eta = new TCanvas("", "", 800, 600);
[92237ca]	1787
	1788	mg_trkele_eff_eta->Draw("APE");
	1789	DrawAxis(mg_trkele_eff_eta, leg_trkele_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "tracking efficiency (%)", false, false);
	1790	leg_trkele_eff_eta->Draw();
[85ad2b9]	1791	pave->Draw();
[92237ca]	1792
[341014c]	1793	c_trkele_eff_eta->Print(pdfOutput, "pdf");
	1794	c_trkele_eff_eta->Print(figPath + "img_trkele_eff_eta.pdf", "pdf");
	1795	c_trkele_eff_eta->Print(figPath + "img_trkele_eff_eta.png", "png");
[92237ca]	1796
	1797	// --------- Muon Tracks --------- //
	1798
[341014c]	1799	TMultiGraph *mg_trkmu_res_pt = new TMultiGraph("", "");
	1800	TMultiGraph *mg_trkmu_eff_pt = new TMultiGraph("", "");
	1801	TMultiGraph *mg_trkmu_res_eta = new TMultiGraph("", "");
	1802	TMultiGraph *mg_trkmu_eff_eta = new TMultiGraph("", "");
[92237ca]	1803
[341014c]	1804	TLegend *leg_trkmu_res_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
	1805	TLegend leg_trkmu_eff_pt = (TLegend )leg_trkmu_res_pt->Clone();
[92237ca]	1806
[341014c]	1807	TLegend leg_trkmu_res_eta = (TLegend )leg_trkmu_res_pt->Clone();
	1808	TLegend leg_trkmu_eff_eta = (TLegend )leg_trkmu_res_eta->Clone();
[92237ca]	1809
[341014c]	1810	TGraphErrors *gr_trkmu_res_pt = new TGraphErrors[n_etabins];
	1811	TGraphErrors *gr_trkmu_eff_pt = new TGraphErrors[n_etabins];
[dd5e213]	1812	TGraphErrors *gr_trkmu_res_eta = new TGraphErrors[n_ptbins];
	1813	TGraphErrors *gr_trkmu_eff_eta = new TGraphErrors[n_ptbins];
	1814
[341014c]	1815	TH1D h_trkmu_eff_pt, h_trkmu_eff_eta;
[92237ca]	1816
[341014c]	1817	std::vector<resolPlot> *plots_trkmu_res_pt = new std::vector<resolPlot>[n_etabins];
	1818	std::vector<resolPlot> *plots_trkmu_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]	1819
	1820	// loop over eta bins
[341014c]	1821	for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]	1822	{
[341014c]	1823	HistogramsCollection(&plots_trkmu_res_pt[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkmu");
	1824	GetPtres<Track>(&plots_trkmu_res_pt[k], branchTrackMuon, branchParticleMuon, 13, etaVals.at(k), etaVals.at(k + 1), treeReaderMuon);
	1825	gr_trkmu_res_pt[k] = EresGraph(&plots_trkmu_res_pt[k]);
[92237ca]	1826
[341014c]	1827	h_trkmu_eff_pt = GetEffPt<Track>(branchTrackMuon, branchParticleMuon, "Muon", 13, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderMuon);
	1828	gr_trkmu_eff_pt[k] = TGraphErrors(h_trkmu_eff_pt);
[92237ca]	1829
[341014c]	1830	s_etaMin = Form("%.1f", etaVals.at(k));
	1831	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]	1832
[341014c]	1833	s_eta = "#mu^{ #pm} , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	1834
[341014c]	1835	gr_trkmu_res_pt[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
	1836	gr_trkmu_eff_pt[k].SetName("trkEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]	1837
[341014c]	1838	addResoGraph(mg_trkmu_res_pt, &gr_trkmu_res_pt[k], leg_trkmu_res_pt, markerStyles.at(k), colors.at(k), s_eta);
	1839	addResoGraph(mg_trkmu_eff_pt, &gr_trkmu_eff_pt[k], leg_trkmu_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]	1840	}
	1841
	1842	// loop over pt
[341014c]	1843	for(k = 0; k < ptVals.size(); k++)
[92237ca]	1844	{
[341014c]	1845	HistogramsCollectionVsEta(&plots_trkmu_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkmu", 0.0, 2.0);
	1846	GetPtresVsEta<Track>(&plots_trkmu_res_eta[k], branchTrackMuon, branchParticleMuon, 13, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderMuon);
	1847	gr_trkmu_res_eta[k] = EresGraphVsEta(&plots_trkmu_res_eta[k]);
[92237ca]	1848
[341014c]	1849	h_trkmu_eff_eta = GetEffEta<Track>(branchTrackMuon, branchParticleMuon, "Muon", 13, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderMuon);
	1850	gr_trkmu_eff_eta[k] = TGraphErrors(h_trkmu_eff_eta);
[92237ca]	1851
[341014c]	1852	s_pt = Form("#mu^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
	1853	if(ptVals.at(k) >= 1000.) s_pt = Form("#mu^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	1854
[341014c]	1855	addResoGraph(mg_trkmu_res_eta, &gr_trkmu_res_eta[k], leg_trkmu_res_eta, markerStyles.at(k), colors.at(k), s_pt);
	1856	addResoGraph(mg_trkmu_eff_eta, &gr_trkmu_eff_eta[k], leg_trkmu_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
[92237ca]	1857	}
	1858
[341014c]	1859	TCanvas *c_trkmu_res_pt = new TCanvas("", "", 800, 600);
[92237ca]	1860
	1861	mg_trkmu_res_pt->Draw("APE");
	1862	DrawAxis(mg_trkmu_res_pt, leg_trkmu_res_pt, ptMin, ptMax, 0.01, 100, "p_{T} [GeV]", "(track resolution in p_{T})/p_{T} (%)", true, true);
	1863	leg_trkmu_res_pt->Draw();
[85ad2b9]	1864	pave->Draw();
[92237ca]	1865
[341014c]	1866	c_trkmu_res_pt->Print(pdfOutput, "pdf");
	1867	c_trkmu_res_pt->Print(figPath + "img_trkmu_res_pt.pdf", "pdf");
	1868	c_trkmu_res_pt->Print(figPath + "img_trkmu_res_pt.png", "png");
[92237ca]	1869
[341014c]	1870	TCanvas *c_trkmu_res_eta = new TCanvas("", "", 800, 600);
[92237ca]	1871
	1872	mg_trkmu_res_eta->Draw("APE");
	1873	DrawAxis(mg_trkmu_res_eta, leg_trkmu_res_eta, etaMin, etaMax, 0.01, 100, " #eta ", "(track resolution in p_{T})/p_{T} (%)", false, true);
	1874	leg_trkmu_res_eta->Draw();
[85ad2b9]	1875	pave->Draw();
[92237ca]	1876
[341014c]	1877	c_trkmu_res_eta->Print(pdfOutput, "pdf");
	1878	c_trkmu_res_eta->Print(figPath + "img_trkmu_res_eta.pdf", "pdf");
	1879	c_trkmu_res_eta->Print(figPath + "img_trkmu_res_eta.png", "png");
[92237ca]	1880
[341014c]	1881	TCanvas *c_trkmu_eff_pt = new TCanvas("", "", 800, 600);
[92237ca]	1882
	1883	mg_trkmu_eff_pt->Draw("APE");
	1884	DrawAxis(mg_trkmu_eff_pt, leg_trkmu_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "tracking efficiency (%)", true, false);
	1885	leg_trkmu_eff_pt->Draw();
[85ad2b9]	1886	pave->Draw();
[92237ca]	1887
[341014c]	1888	c_trkmu_eff_pt->Print(pdfOutput, "pdf");
	1889	c_trkmu_eff_pt->Print(figPath + "img_trkmu_eff_pt.pdf", "pdf");
	1890	c_trkmu_eff_pt->Print(figPath + "img_trkmu_eff_pt.png", "png");
[92237ca]	1891
[341014c]	1892	TCanvas *c_trkmu_eff_eta = new TCanvas("", "", 800, 600);
[92237ca]	1893
	1894	mg_trkmu_eff_eta->Draw("APE");
	1895	DrawAxis(mg_trkmu_eff_eta, leg_trkmu_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "tracking efficiency (%)", false, false);
	1896	leg_trkmu_eff_eta->Draw();
[85ad2b9]	1897	pave->Draw();
[92237ca]	1898
[341014c]	1899	c_trkmu_eff_eta->Print(pdfOutput, "pdf");
	1900	c_trkmu_eff_eta->Print(figPath + "img_trkmu_eff_eta.pdf", "pdf");
	1901	c_trkmu_eff_eta->Print(figPath + "img_trkmu_eff_eta.png", "png");
[92237ca]	1902
	1903	//////////////////////
	1904	// Ecal performance //
	1905	//////////////////////
	1906
[341014c]	1907	TMultiGraph *mg_ecal_res_e = new TMultiGraph("", "");
	1908	TMultiGraph *mg_ecal_res_eta = new TMultiGraph("", "");
[92237ca]	1909
[341014c]	1910	TLegend *leg_ecal_res_e = new TLegend(0.55, 0.64, 0.90, 0.90);
	1911	TLegend *leg_ecal_res_eta = new TLegend(0.60, 0.59, 0.95, 0.90);
[92237ca]	1912
[341014c]	1913	TGraphErrors *gr_ecal_res_e = new TGraphErrors[n_etabins];
	1914	TGraphErrors *gr_ecal_res_eta = new TGraphErrors[n_ptbins];
[92237ca]	1915
[341014c]	1916	std::vector<resolPlot> *plots_ecal_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]	1917	std::vector<resolPlot> *plots_ecal_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]	1918
	1919	// loop over eta bins
[341014c]	1920	for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]	1921	{
[341014c]	1922	HistogramsCollection(&plots_ecal_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "ecal");
	1923	GetEres<Tower>(&plots_ecal_res_e[k], branchTowerPhoton, branchParticlePhoton, 22, etaVals.at(k), etaVals.at(k + 1), treeReaderPhoton);
	1924	gr_ecal_res_e[k] = EresGraph(&plots_ecal_res_e[k]);
[92237ca]	1925
[341014c]	1926	s_etaMin = Form("%.1f", etaVals.at(k));
	1927	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]	1928
[341014c]	1929	s_eta = "#gamma , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	1930
[341014c]	1931	gr_ecal_res_e[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
[92237ca]	1932
[341014c]	1933	addResoGraph(mg_ecal_res_e, &gr_ecal_res_e[k], leg_ecal_res_e, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]	1934	}
	1935
	1936	// loop over pt
[341014c]	1937	for(k = 0; k < ptVals.size(); k++)
[92237ca]	1938	{
[341014c]	1939	HistogramsCollectionVsEta(&plots_ecal_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "ecal", 0.0, 2.0);
	1940	GetEresVsEta<Tower>(&plots_ecal_res_eta[k], branchTowerPhoton, branchParticlePhoton, 22, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderPhoton);
	1941	gr_ecal_res_eta[k] = EresGraphVsEta(&plots_ecal_res_eta[k]);
[92237ca]	1942
[341014c]	1943	s_e = Form("#gamma , E = %.0f GeV", ptVals.at(k));
	1944	if(ptVals.at(k) >= 1000.) s_e = Form("#gamma , E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	1945
[341014c]	1946	addResoGraph(mg_ecal_res_eta, &gr_ecal_res_eta[k], leg_ecal_res_eta, markerStyles.at(k), colors.at(k), s_e);
[92237ca]	1947	}
	1948
[341014c]	1949	TCanvas *c_ecal_res_e = new TCanvas("", "", 800, 600);
[92237ca]	1950
	1951	mg_ecal_res_e->Draw("APE");
[341014c]	1952	// DrawAxis(mg_ecal_res_e, leg_ecal_res_e, ptMin, ptMax, 0.5, 100, "E [GeV]", "(ECAL resolution in E)/E (%)", true, true);
[92237ca]	1953	DrawAxis(mg_ecal_res_e, leg_ecal_res_e, ptMin, ptMax, 0.0, 20, "E [GeV]", "(ECAL resolution in E)/E (%)", true, false);
	1954	leg_ecal_res_e->Draw();
[85ad2b9]	1955	pave->Draw();
[92237ca]	1956
[341014c]	1957	c_ecal_res_e->Print(pdfOutput, "pdf");
	1958	c_ecal_res_e->Print(figPath + "img_ecal_res_e.pdf", "pdf");
	1959	c_ecal_res_e->Print(figPath + "img_ecal_res_e.png", "png");
[92237ca]	1960
[341014c]	1961	TCanvas *c_ecal_res_eta = new TCanvas("", "", 800, 600);
[92237ca]	1962
	1963	mg_ecal_res_eta->Draw("APE");
	1964	//DrawAxis(mg_ecal_res_eta, leg_ecal_res_eta, etaMin, etaMax, 0.5, 100, " #eta ", "(ECAL resolution in E)/E (%)", false, true);
	1965	DrawAxis(mg_ecal_res_eta, leg_ecal_res_eta, etaMin, etaMax, 0.0, 20, " #eta ", "(ECAL resolution in E)/E (%)", false, false);
	1966	leg_ecal_res_eta->Draw();
[85ad2b9]	1967	pave->Draw();
[92237ca]	1968
[341014c]	1969	c_ecal_res_eta->Print(pdfOutput, "pdf");
	1970	c_ecal_res_eta->Print(figPath + "img_ecal_res_eta.pdf", "pdf");
	1971	c_ecal_res_eta->Print(figPath + "img_ecal_res_eta.png", "png");
[92237ca]	1972
	1973	//////////////////////
	1974	// Hcal performance //
	1975	//////////////////////
	1976
[341014c]	1977	TMultiGraph *mg_hcal_res_e = new TMultiGraph("", "");
	1978	TMultiGraph *mg_hcal_res_eta = new TMultiGraph("", "");
[92237ca]	1979
[341014c]	1980	TLegend *leg_hcal_res_e = new TLegend(0.55, 0.64, 0.90, 0.90);
	1981	TLegend *leg_hcal_res_eta = new TLegend(0.60, 0.59, 0.95, 0.90);
[92237ca]	1982
[341014c]	1983	TGraphErrors *gr_hcal_res_e = new TGraphErrors[n_etabins];
[dd5e213]	1984	TGraphErrors *gr_hcal_res_eta = new TGraphErrors[n_ptbins];
[92237ca]	1985
[341014c]	1986	std::vector<resolPlot> *plots_hcal_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]	1987	std::vector<resolPlot> *plots_hcal_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]	1988
	1989	// loop over eta bins
[341014c]	1990	for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]	1991	{
[341014c]	1992	HistogramsCollection(&plots_hcal_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "hcal");
	1993	GetEres<Tower>(&plots_hcal_res_e[k], branchTowerNeutralHadron, branchParticleNeutralHadron, 2112, etaVals.at(k), etaVals.at(k + 1), treeReaderNeutralHadron);
[92237ca]	1994
[341014c]	1995	gr_hcal_res_e[k] = EresGraph(&plots_hcal_res_e[k]);
[92237ca]	1996
[341014c]	1997	s_etaMin = Form("%.1f", etaVals.at(k));
	1998	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]	1999
[341014c]	2000	s_eta = "n , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	2001
[341014c]	2002	gr_hcal_res_e[k].SetName("trkRes_" + s_etaMin + "_" + s_etaMax);
[92237ca]	2003
[341014c]	2004	addResoGraph(mg_hcal_res_e, &gr_hcal_res_e[k], leg_hcal_res_e, markerStyles.at(k), colors.at(k), s_eta);
[92237ca]	2005	}
	2006
	2007	// loop over pt
[341014c]	2008	for(k = 0; k < ptVals.size(); k++)
[92237ca]	2009	{
[341014c]	2010	HistogramsCollectionVsEta(&plots_hcal_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "hcal", 0.0, 2.0);
	2011	GetEresVsEta<Tower>(&plots_hcal_res_eta[k], branchTowerNeutralHadron, branchParticleNeutralHadron, 2112, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderNeutralHadron);
	2012	gr_hcal_res_eta[k] = EresGraphVsEta(&plots_hcal_res_eta[k]);
[92237ca]	2013
[341014c]	2014	s_e = Form("n , E = %.0f GeV", ptVals.at(k));
	2015	if(ptVals.at(k) >= 1000.) s_e = Form("n , E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	2016
[341014c]	2017	addResoGraph(mg_hcal_res_eta, &gr_hcal_res_eta[k], leg_hcal_res_eta, markerStyles.at(k), colors.at(k), s_e);
[92237ca]	2018	}
	2019
[341014c]	2020	TCanvas *c_hcal_res_e = new TCanvas("", "", 800, 600);
[92237ca]	2021
	2022	mg_hcal_res_e->Draw("APE");
	2023	//DrawAxis(mg_hcal_res_e, leg_hcal_res_e, ptMin, ptMax, 1, 100, "E [GeV]", "(HCAL resolution in E)/E (%)", true, true);
	2024	DrawAxis(mg_hcal_res_e, leg_hcal_res_e, ptMin, ptMax, 0.0, 50, "E [GeV]", "(HCAL resolution in E)/E (%)", true, false);
	2025	leg_hcal_res_e->Draw();
[85ad2b9]	2026	pave->Draw();
[92237ca]	2027
[341014c]	2028	c_hcal_res_e->Print(pdfOutput, "pdf");
	2029	c_hcal_res_e->Print(figPath + "img_hcal_res_e.pdf", "pdf");
	2030	c_hcal_res_e->Print(figPath + "img_hcal_res_e.png", "png");
[92237ca]	2031
[341014c]	2032	TCanvas *c_hcal_res_eta = new TCanvas("", "", 800, 600);
[92237ca]	2033
	2034	mg_hcal_res_eta->Draw("APE");
	2035	//DrawAxis(mg_hcal_res_eta, leg_hcal_res_eta, etaMin, etaMax, 1, 100, " #eta ", "(HCAL resolution in E)/E (%)", false, true);
	2036	DrawAxis(mg_hcal_res_eta, leg_hcal_res_eta, etaMin, etaMax, 0.0, 50, " #eta ", "(HCAL resolution in E)/E (%)", false, false);
	2037	leg_hcal_res_eta->Draw();
[85ad2b9]	2038	pave->Draw();
[92237ca]	2039
[341014c]	2040	c_hcal_res_eta->Print(pdfOutput, "pdf");
	2041	c_hcal_res_eta->Print(figPath + "img_hcal_res_eta.pdf", "pdf");
	2042	c_hcal_res_eta->Print(figPath + "img_hcal_res_eta.png", "png");
[92237ca]	2043
	2044	////////////////////
	2045	// PF - electrons //
	2046	////////////////////
	2047
	2048	TMultiGraph *mg_pfele_res_e[n_etabins];
	2049	TMultiGraph *mg_pfele_res_eta[n_ptbins];
	2050
	2051	TLegend *leg_pfele_res_e[n_etabins];
	2052	TLegend *leg_pfele_res_eta[n_ptbins];
	2053
[dd5e213]	2054	TGraphErrors *gr_pfele_res_e = new TGraphErrors[n_etabins];
	2055	TGraphErrors *gr_pfele_res_eta = new TGraphErrors[n_ptbins];
	2056	TGraphErrors *gr_trkele_res_e = new TGraphErrors[n_etabins];
	2057	TGraphErrors *gr_trkele_res_eeta = new TGraphErrors[n_ptbins];
[cb159ed]	2058
[341014c]	2059	std::vector<resolPlot> *plots_pfele_res_e = new std::vector<resolPlot>[n_etabins];
	2060	std::vector<resolPlot> *plots_pfele_res_eta = new std::vector<resolPlot>[n_ptbins];
	2061	std::vector<resolPlot> *plots_trkele_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]	2062	std::vector<resolPlot> *plots_trkele_res_eeta = new std::vector<resolPlot>[n_ptbins];
[92237ca]	2063
	2064	TCanvas *c_pfele_res_e[n_etabins];
	2065	TCanvas *c_pfele_res_eta[n_ptbins];
	2066
	2067	// loop over eta bins
[341014c]	2068	for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]	2069	{
[341014c]	2070	mg_pfele_res_e[k] = new TMultiGraph("", "");
	2071	leg_pfele_res_e[k] = new TLegend(0.40, 0.60, 0.75, 0.90);
[92237ca]	2072
[341014c]	2073	HistogramsCollection(&plots_pfele_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "pfele");
	2074	GetEres<Electron>(&plots_pfele_res_e[k], branchElectronPF, branchParticleElectron, 11, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
	2075	gr_pfele_res_e[k] = EresGraph(&plots_pfele_res_e[k]);
[92237ca]	2076
[341014c]	2077	HistogramsCollection(&plots_trkele_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkele");
	2078	GetEres<Track>(&plots_trkele_res_e[k], branchTrackElectron, branchParticleElectron, 11, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
	2079	gr_trkele_res_e[k] = EresGraph(&plots_trkele_res_e[k]);
[cb159ed]	2080
[341014c]	2081	s_etaMin = Form("%.1f", etaVals.at(k));
	2082	s_etaMax = Form("%.1f", etaVals.at(k + 1));
	2083	s_eta = "e^{ #pm}, " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	2084
[341014c]	2085	leg_pfele_res_e[k]->SetTextFont(132);
	2086	leg_pfele_res_e[k]->SetHeader(s_eta);
[92237ca]	2087
[341014c]	2088	addResoGraph(mg_pfele_res_e[k], &gr_ecal_res_e[k], leg_pfele_res_e[k], markerStyles.at(0), colors.at(0), "ECAL");
	2089	addResoGraph(mg_pfele_res_e[k], &gr_trkele_res_e[k], leg_pfele_res_e[k], markerStyles.at(1), colors.at(1), "Track");
	2090	addResoGraph(mg_pfele_res_e[k], &gr_pfele_res_e[k], leg_pfele_res_e[k], markerStyles.at(2), colors.at(2), "Particle-flow");
[92237ca]	2091
[341014c]	2092	c_pfele_res_e[k] = new TCanvas("", "", 800, 600);
[92237ca]	2093
[341014c]	2094	mg_pfele_res_e[k]->Draw("APE");
	2095	//DrawAxis(mg_pfele_res_e[k], leg_pfele_res_e[k], ptMin, ptMax, 0.1, 100, "E [GeV]", "(resolution in E)/E (%)", true, true);
	2096	DrawAxis(mg_pfele_res_e[k], leg_pfele_res_e[k], ptMin, ptMax, 0.0, 20, "E [GeV]", "(resolution in E)/E (%)", true, false);
	2097	leg_pfele_res_e[k]->Draw();
	2098	pave->Draw();
[92237ca]	2099
[341014c]	2100	TString s_etarange = "eta_" + s_etaMin + "_" + s_etaMax + "_";
[92237ca]	2101
[341014c]	2102	c_pfele_res_e[k]->Print(pdfOutput, "pdf");
	2103	c_pfele_res_e[k]->Print(figPath + "img_pfele_res_" + s_etarange + "e.pdf", "pdf");
	2104	c_pfele_res_e[k]->Print(figPath + "img_pfele_res_" + s_etarange + "e.png", "png");
[92237ca]	2105	}
	2106
[341014c]	2107	// loop over eta bins
	2108	for(k = 0; k < ptVals.size(); k++)
[92237ca]	2109	{
	2110
[341014c]	2111	mg_pfele_res_eta[k] = new TMultiGraph("", "");
	2112	leg_pfele_res_eta[k] = new TLegend(0.40, 0.60, 0.75, 0.90);
[cb159ed]	2113
[341014c]	2114	HistogramsCollectionVsEta(&plots_pfele_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "pfele", 0.0, 2.0);
	2115	GetEresVsEta<Electron>(&plots_pfele_res_eta[k], branchElectronPF, branchParticleElectron, 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderElectron);
	2116	gr_pfele_res_eta[k] = EresGraphVsEta(&plots_pfele_res_eta[k]);
[92237ca]	2117
[341014c]	2118	HistogramsCollectionVsEta(&plots_trkele_res_eeta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkele", 0.0, 2.0);
	2119	GetEresVsEta<Track>(&plots_trkele_res_eeta[k], branchTrackElectron, branchParticleElectron, 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderElectron);
	2120	gr_trkele_res_eeta[k] = EresGraphVsEta(&plots_trkele_res_eeta[k]);
[92237ca]	2121
[341014c]	2122	s_e = Form("e^{ #pm}, E = %.0f GeV", ptVals.at(k));
	2123	if(ptVals.at(k) >= 1000.) s_e = Form("e^{ #pm}, E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	2124
[341014c]	2125	leg_pfele_res_eta[k]->SetTextFont(132);
	2126	leg_pfele_res_eta[k]->SetHeader(s_e);
[92237ca]	2127
[341014c]	2128	addResoGraph(mg_pfele_res_eta[k], &gr_ecal_res_eta[k], leg_pfele_res_eta[k], markerStyles.at(0), colors.at(0), "ECAL");
	2129	addResoGraph(mg_pfele_res_eta[k], &gr_trkele_res_eeta[k], leg_pfele_res_eta[k], markerStyles.at(1), colors.at(1), "Track");
	2130	addResoGraph(mg_pfele_res_eta[k], &gr_pfele_res_eta[k], leg_pfele_res_eta[k], markerStyles.at(2), colors.at(2), "Particle-flow");
[92237ca]	2131
[341014c]	2132	c_pfele_res_eta[k] = new TCanvas("", "", 800, 600);
[92237ca]	2133
[341014c]	2134	mg_pfele_res_eta[k]->Draw("APE");
	2135	//DrawAxis(mg_pfele_res_eta[k], leg_pfele_res_eta[k], etaMin, etaMax, 0.1, 1000, "#eta", "(resolution in E)/E (%)", false, true);
	2136	DrawAxis(mg_pfele_res_eta[k], leg_pfele_res_eta[k], etaMin, etaMax, 0.0, 50, "#eta", "(resolution in E)/E (%)", false, false);
	2137	leg_pfele_res_eta[k]->Draw();
	2138	pave->Draw();
[92237ca]	2139
[341014c]	2140	TString s_ptrange = Form("pt_%.0f_", ptVals.at(k));
[92237ca]	2141
[341014c]	2142	c_pfele_res_eta[k]->Print(pdfOutput, "pdf");
	2143	c_pfele_res_eta[k]->Print(figPath + "img_pfele_res_" + s_ptrange + "eta.pdf", "pdf");
	2144	c_pfele_res_eta[k]->Print(figPath + "img_pfele_res_" + s_ptrange + "eta.png", "png");
[92237ca]	2145	}
	2146
	2147	/////////////////
	2148	// PF - Pions //
	2149	/////////////////
	2150
	2151	TMultiGraph *mg_pfpi_res_e[n_etabins];
	2152	TMultiGraph *mg_pfpi_res_eta[n_ptbins];
	2153
	2154	TLegend *leg_pfpi_res_e[n_etabins];
	2155	TLegend *leg_pfpi_res_eta[n_ptbins];
	2156
[341014c]	2157	TGraphErrors *gr_pfpi_res_e = new TGraphErrors[n_etabins];
	2158	TGraphErrors *gr_pfpi_res_eta = new TGraphErrors[n_ptbins];
[92237ca]	2159
[341014c]	2160	TGraphErrors *gr_trkpi_res_e = new TGraphErrors[n_etabins];
[dd5e213]	2161	TGraphErrors *gr_trkpi_res_eeta = new TGraphErrors[n_ptbins];
[cb159ed]	2162
[341014c]	2163	std::vector<resolPlot> *plots_pfpi_res_e = new std::vector<resolPlot>[n_etabins];
	2164	std::vector<resolPlot> *plots_pfpi_res_eta = new std::vector<resolPlot>[n_ptbins];
	2165	std::vector<resolPlot> *plots_trkpi_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]	2166	std::vector<resolPlot> *plots_trkpi_res_eeta = new std::vector<resolPlot>[n_ptbins];
[92237ca]	2167
	2168	TCanvas *c_pfpi_res_e[n_etabins];
	2169	TCanvas *c_pfpi_res_eta[n_ptbins];
	2170
	2171	// loop over eta bins
[341014c]	2172	for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]	2173	{
[341014c]	2174	mg_pfpi_res_e[k] = new TMultiGraph("", "");
	2175	leg_pfpi_res_e[k] = new TLegend(0.40, 0.60, 0.75, 0.90);
[cb159ed]	2176
[341014c]	2177	HistogramsCollection(&plots_pfpi_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "pfpi");
	2178	GetEres<Track>(&plots_pfpi_res_e[k], branchPion, branchParticlePion, 211, etaVals.at(k), etaVals.at(k + 1), treeReaderPion);
	2179	gr_pfpi_res_e[k] = EresGraph(&plots_pfpi_res_e[k]);
[cb159ed]	2180
[341014c]	2181	HistogramsCollection(&plots_trkpi_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "trkpi");
	2182	GetEres<Track>(&plots_trkpi_res_e[k], branchTrackPion, branchParticlePion, 211, etaVals.at(k), etaVals.at(k + 1), treeReaderPion);
	2183	gr_trkpi_res_e[k] = EresGraph(&plots_trkpi_res_e[k]);
[92237ca]	2184
[341014c]	2185	s_etaMin = Form("%.1f", etaVals.at(k));
	2186	s_etaMax = Form("%.1f", etaVals.at(k + 1));
	2187	s_eta = "#pi^{ #pm}, " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	2188
[341014c]	2189	leg_pfpi_res_e[k]->SetTextFont(132);
	2190	leg_pfpi_res_e[k]->SetHeader(s_eta);
[92237ca]	2191
[341014c]	2192	addResoGraph(mg_pfpi_res_e[k], &gr_hcal_res_e[k], leg_pfpi_res_e[k], markerStyles.at(0), colors.at(0), "HCAL");
	2193	addResoGraph(mg_pfpi_res_e[k], &gr_trkpi_res_e[k], leg_pfpi_res_e[k], markerStyles.at(1), colors.at(1), "Track");
	2194	addResoGraph(mg_pfpi_res_e[k], &gr_pfpi_res_e[k], leg_pfpi_res_e[k], markerStyles.at(2), colors.at(2), "Particle-flow");
[92237ca]	2195
[341014c]	2196	c_pfpi_res_e[k] = new TCanvas("", "", 800, 600);
[92237ca]	2197
[341014c]	2198	mg_pfpi_res_e[k]->Draw("APE");
	2199	//DrawAxis(mg_pfpi_res_e[k], leg_pfpi_res_e[k], ptMin, ptMax, 0.1, 100, "E [GeV]", "(resolution in E)/E (%)", true, true);
	2200	DrawAxis(mg_pfpi_res_e[k], leg_pfpi_res_e[k], ptMin, ptMax, 0.1, 50, "E [GeV]", "(resolution in E)/E (%)", true, false);
	2201	leg_pfpi_res_e[k]->Draw();
	2202	pave->Draw();
[92237ca]	2203
[341014c]	2204	TString s_etarange = "eta_" + s_etaMin + "_" + s_etaMax + "_";
[92237ca]	2205
[341014c]	2206	c_pfpi_res_e[k]->Print(pdfOutput, "pdf");
	2207	c_pfpi_res_e[k]->Print(figPath + "img_pfpi_res_" + s_etarange + "e.pdf", "pdf");
	2208	c_pfpi_res_e[k]->Print(figPath + "img_pfpi_res_" + s_etarange + "e.png", "png");
[92237ca]	2209	}
	2210
[341014c]	2211	// loop over eta bins
	2212	for(k = 0; k < ptVals.size(); k++)
[92237ca]	2213	{
	2214
[341014c]	2215	mg_pfpi_res_eta[k] = new TMultiGraph("", "");
	2216	leg_pfpi_res_eta[k] = new TLegend(0.40, 0.60, 0.75, 0.90);
[cb159ed]	2217
[341014c]	2218	HistogramsCollectionVsEta(&plots_pfpi_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "pfpi", 0.0, 2.0);
	2219	GetEresVsEta<Track>(&plots_pfpi_res_eta[k], branchPion, branchParticlePion, 211, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderPion);
	2220	gr_pfpi_res_eta[k] = EresGraphVsEta(&plots_pfpi_res_eta[k]);
[cb159ed]	2221
[341014c]	2222	HistogramsCollectionVsEta(&plots_trkpi_res_eeta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "trkpi", 0.0, 2.0);
	2223	GetEresVsEta<Track>(&plots_trkpi_res_eeta[k], branchPion, branchParticlePion, 211, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), treeReaderPion);
	2224	gr_trkpi_res_eeta[k] = EresGraphVsEta(&plots_trkpi_res_eeta[k]);
[92237ca]	2225
[341014c]	2226	s_e = Form("#pi^{ #pm}, E = %.0f GeV", ptVals.at(k));
	2227	if(ptVals.at(k) >= 1000.) s_e = Form("#pi^{ #pm}, E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	2228
[341014c]	2229	leg_pfpi_res_eta[k]->SetTextFont(132);
	2230	leg_pfpi_res_eta[k]->SetHeader(s_e);
[92237ca]	2231
[341014c]	2232	addResoGraph(mg_pfpi_res_eta[k], &gr_hcal_res_eta[k], leg_pfpi_res_eta[k], markerStyles.at(0), colors.at(0), "HCAL");
	2233	addResoGraph(mg_pfpi_res_eta[k], &gr_trkpi_res_eeta[k], leg_pfpi_res_eta[k], markerStyles.at(1), colors.at(1), "Track");
	2234	addResoGraph(mg_pfpi_res_eta[k], &gr_pfpi_res_eta[k], leg_pfpi_res_eta[k], markerStyles.at(2), colors.at(2), "Particle-flow");
[92237ca]	2235
[341014c]	2236	c_pfpi_res_eta[k] = new TCanvas("", "", 800, 600);
[92237ca]	2237
[341014c]	2238	mg_pfpi_res_eta[k]->Draw("APE");
	2239	//DrawAxis(mg_pfpi_res_eta[k], leg_pfpi_res_eta[k], etaMin, etaMax, 0.1, 1000, "#eta", "(resolution in E)/E (%)", false, true);
	2240	DrawAxis(mg_pfpi_res_eta[k], leg_pfpi_res_eta[k], etaMin, etaMax, 0.0, 50, "#eta", "(resolution in E)/E (%)", false, false);
	2241	leg_pfpi_res_eta[k]->Draw();
	2242	pave->Draw();
[92237ca]	2243
[341014c]	2244	TString s_ptrange = Form("pt_%.0f_", ptVals.at(k));
[92237ca]	2245
[341014c]	2246	c_pfpi_res_eta[k]->Print(pdfOutput, "pdf");
	2247	c_pfpi_res_eta[k]->Print(figPath + "img_pfpi_res_" + s_ptrange + "eta.pdf", "pdf");
	2248	c_pfpi_res_eta[k]->Print(figPath + "img_pfpi_res_" + s_ptrange + "eta.png", "png");
[92237ca]	2249	}
	2250
	2251	/////////////////
	2252	// PF - jets //
	2253	/////////////////
	2254
	2255	TMultiGraph *mg_pfjet_res_e[n_etabins];
	2256	TMultiGraph *mg_pfjet_res_eta[n_ptbins];
	2257
	2258	TLegend *leg_pfjet_res_e[n_etabins];
	2259	TLegend *leg_pfjet_res_eta[n_ptbins];
	2260
[341014c]	2261	TGraphErrors *gr_pfjet_res_e = new TGraphErrors[n_etabins];
[dd5e213]	2262	TGraphErrors *gr_pfjet_res_eta = new TGraphErrors[n_ptbins];
[92237ca]	2263
[341014c]	2264	TGraphErrors *gr_cajet_res_e = new TGraphErrors[n_etabins];
[dd5e213]	2265	TGraphErrors *gr_cajet_res_eta = new TGraphErrors[n_ptbins];
[92237ca]	2266
[341014c]	2267	std::vector<resolPlot> *plots_pfjet_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]	2268	std::vector<resolPlot> *plots_pfjet_res_eta = new std::vector<resolPlot>[n_ptbins];
[341014c]	2269	std::vector<resolPlot> *plots_cajet_res_e = new std::vector<resolPlot>[n_etabins];
[dd5e213]	2270	std::vector<resolPlot> *plots_cajet_res_eta = new std::vector<resolPlot>[n_ptbins];
[92237ca]	2271
	2272	TCanvas *c_pfjet_res_e[n_etabins];
	2273	TCanvas *c_pfjet_res_eta[n_ptbins];
	2274
	2275	// loop over eta bins
[341014c]	2276	for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]	2277	{
	2278
[341014c]	2279	mg_pfjet_res_e[k] = new TMultiGraph("", "");
	2280	leg_pfjet_res_e[k] = new TLegend(0.40, 0.70, 0.90, 0.90);
[92237ca]	2281
[341014c]	2282	HistogramsCollection(&plots_pfjet_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "pfjet");
	2283	HistogramsCollection(&plots_cajet_res_e[k], TMath::Log10(ptMin), TMath::Log10(ptMax), "cajet");
[92237ca]	2284
[341014c]	2285	GetJetsEres(&plots_pfjet_res_e[k], branchPFJet, branchGenJet, treeReaderJet, etaVals.at(k), etaVals.at(k + 1));
	2286	GetJetsEres(&plots_cajet_res_e[k], branchCaloJet, branchGenJet, treeReaderJet, etaVals.at(k), etaVals.at(k + 1));
[92237ca]	2287
[341014c]	2288	gr_pfjet_res_e[k] = EresGraph(&plots_pfjet_res_e[k]);
	2289	gr_cajet_res_e[k] = EresGraph(&plots_cajet_res_e[k]);
[92237ca]	2290
[341014c]	2291	s_etaMin = Form("%.1f", etaVals.at(k));
	2292	s_etaMax = Form("%.1f", etaVals.at(k + 1));
	2293	s_eta = "anti-k_{T}, R = 0.4, " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	2294
[341014c]	2295	leg_pfjet_res_e[k]->SetTextFont(132);
	2296	leg_pfjet_res_e[k]->SetHeader(s_eta);
[92237ca]	2297
[341014c]	2298	addResoGraph(mg_pfjet_res_e[k], &gr_cajet_res_e[k], leg_pfjet_res_e[k], markerStyles.at(0), colors.at(0), "Calorimeter Jets");
	2299	addResoGraph(mg_pfjet_res_e[k], &gr_pfjet_res_e[k], leg_pfjet_res_e[k], markerStyles.at(1), colors.at(1), "Particle-flow Jets");
[92237ca]	2300
[341014c]	2301	c_pfjet_res_e[k] = new TCanvas("", "", 800, 600);
[92237ca]	2302
[341014c]	2303	mg_pfjet_res_e[k]->Draw("APE");
	2304	//DrawAxis(mg_pfjet_res_e[k], leg_pfjet_res_e[k], 10, ptMax, 0.5, 100, "E [GeV]", "(resolution in E)/E (%)", true, true);
	2305	DrawAxis(mg_pfjet_res_e[k], leg_pfjet_res_e[k], 10, ptMax, 0.0, 30, "E [GeV]", "(resolution in E)/E (%)", true, false);
	2306	leg_pfjet_res_e[k]->Draw();
	2307	pave->Draw();
[92237ca]	2308
[341014c]	2309	TString s_etarange = "eta_" + s_etaMin + "_" + s_etaMax + "_";
[92237ca]	2310
[341014c]	2311	c_pfjet_res_e[k]->Print(pdfOutput, "pdf");
	2312	c_pfjet_res_e[k]->Print(figPath + "img_pfjet_res_" + s_etarange + "e.pdf", "pdf");
	2313	c_pfjet_res_e[k]->Print(figPath + "img_pfjet_res_" + s_etarange + "e.png", "png");
[92237ca]	2314	}
	2315
[341014c]	2316	// loop over eta bins
	2317	for(k = 0; k < ptVals.size(); k++)
[92237ca]	2318	{
	2319
[341014c]	2320	mg_pfjet_res_eta[k] = new TMultiGraph("", "");
	2321	leg_pfjet_res_eta[k] = new TLegend(0.30, 0.70, 0.85, 0.90);
[92237ca]	2322
[341014c]	2323	HistogramsCollectionVsEta(&plots_pfjet_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "pfjet", 0.0, 2.0);
	2324	HistogramsCollectionVsEta(&plots_cajet_res_eta[k], etaMin, etaMax, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), "cajet", 0.0, 2.0);
[92237ca]	2325
[341014c]	2326	GetJetsEresVsEta(&plots_pfjet_res_eta[k], branchPFJet, branchGenJet, treeReaderJet, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k));
	2327	GetJetsEresVsEta(&plots_cajet_res_eta[k], branchCaloJet, branchGenJet, treeReaderJet, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k));
[92237ca]	2328
[341014c]	2329	gr_pfjet_res_eta[k] = EresGraphVsEta(&plots_pfjet_res_eta[k]);
	2330	gr_cajet_res_eta[k] = EresGraphVsEta(&plots_cajet_res_eta[k]);
[92237ca]	2331
[341014c]	2332	s_e = Form("anti-k_{T}, R = 0.4, jets, E = %.0f GeV", ptVals.at(k));
	2333	if(ptVals.at(k) >= 1000.) s_e = Form("anti-k_{T}, R = 0.4, E = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	2334
[341014c]	2335	leg_pfjet_res_eta[k]->SetTextFont(132);
	2336	leg_pfjet_res_eta[k]->SetHeader(s_e);
[92237ca]	2337
[341014c]	2338	addResoGraph(mg_pfjet_res_eta[k], &gr_cajet_res_eta[k], leg_pfjet_res_eta[k], markerStyles.at(0), colors.at(0), "Calorimeter Jets");
	2339	addResoGraph(mg_pfjet_res_eta[k], &gr_pfjet_res_eta[k], leg_pfjet_res_eta[k], markerStyles.at(1), colors.at(1), "Particle-flow Jets");
[92237ca]	2340
[341014c]	2341	c_pfjet_res_eta[k] = new TCanvas("", "", 800, 600);
[92237ca]	2342
[341014c]	2343	mg_pfjet_res_eta[k]->Draw("APE");
	2344	//DrawAxis(mg_pfjet_res_eta[k], leg_pfjet_res_eta[k], etaMin, etaMax, 0.1, 1000, "#eta", "(resolution in E)/E (%)", false, true);
	2345	DrawAxis(mg_pfjet_res_eta[k], leg_pfjet_res_eta[k], etaMin, etaMax, 0.0, 50, "#eta", "(resolution in E)/E (%)", false, false);
	2346	leg_pfjet_res_eta[k]->Draw();
	2347	pave->Draw();
[92237ca]	2348
[341014c]	2349	TString s_ptrange = Form("pt_%.0f_", ptVals.at(k));
[92237ca]	2350
[341014c]	2351	c_pfjet_res_eta[k]->Print(pdfOutput, "pdf");
	2352	c_pfjet_res_eta[k]->Print(figPath + "img_pfjet_res_" + s_ptrange + "eta.pdf", "pdf");
	2353	c_pfjet_res_eta[k]->Print(figPath + "img_pfjet_res_" + s_ptrange + "eta.png", "png");
[92237ca]	2354	}
	2355
[341014c]	2356	/////////////////////
	2357	// PF Missing ET ///
	2358	/////////////////////
[92237ca]	2359
[341014c]	2360	TMultiGraph *mg_met_res_ht = new TMultiGraph("", "");
	2361	TLegend *leg_met_res_ht = new TLegend(0.60, 0.22, 0.90, 0.42);
[92237ca]	2362
[341014c]	2363	std::vector<resolPlot> plots_pfmet, plots_camet;
[92237ca]	2364
[341014c]	2365	HistogramsCollection(&plots_pfmet, TMath::Log10(ptMin), TMath::Log10(ptMax), "pfMET", -500, 500);
	2366	HistogramsCollection(&plots_camet, TMath::Log10(ptMin), TMath::Log10(ptMax), "caMET", -500, 500);
[92237ca]	2367
[341014c]	2368	GetMetres(&plots_pfmet, branchGenScalarHT, branchMet, branchPFJet, treeReaderJet);
	2369	GetMetres(&plots_camet, branchGenScalarHT, branchCaloMet, branchCaloJet, treeReaderJet);
[92237ca]	2370
[341014c]	2371	TGraphErrors gr_pfmet_res_ht = MetResGraph(&plots_pfmet, true);
	2372	TGraphErrors gr_camet_res_ht = MetResGraph(&plots_camet, true);
[92237ca]	2373
[341014c]	2374	addResoGraph(mg_met_res_ht, &gr_camet_res_ht, leg_met_res_ht, markerStyles.at(0), colors.at(0), "Calorimeter E_{T}^{miss}");
	2375	addResoGraph(mg_met_res_ht, &gr_pfmet_res_ht, leg_met_res_ht, markerStyles.at(1), colors.at(1), "Particle-flow E_{T}^{miss}");
[92237ca]	2376
[341014c]	2377	TCanvas *c_met_res_ht = new TCanvas("", "", 800, 600);
[92237ca]	2378
[341014c]	2379	mg_met_res_ht->Draw("APE");
	2380	DrawAxis(mg_met_res_ht, leg_met_res_ht, 1000, 100000, 0.1, 1000, " #sum p_{T} [GeV]", "resolution in E_{x,y}^{miss} [GeV]", true, true);
[92237ca]	2381
[341014c]	2382	leg_met_res_ht->Draw();
	2383	pave->Draw();
	2384	c_met_res_ht->Print(pdfOutput, "pdf");
	2385	c_met_res_ht->Print(figPath + "img_met_res_ht.pdf", "pdf");
	2386	c_met_res_ht->Print(figPath + "img_met_res_ht.png", "png");
[92237ca]	2387
[341014c]	2388	/////////////////////////////////////////
	2389	// Electron Reconstruction Efficiency ///
	2390	/////////////////////////////////////////
[92237ca]	2391
[341014c]	2392	TMultiGraph *mg_recele_eff_pt = new TMultiGraph("", "");
	2393	TMultiGraph *mg_recele_eff_eta = new TMultiGraph("", "");
[92237ca]	2394
[341014c]	2395	TLegend *leg_recele_eff_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
	2396	TLegend *leg_recele_eff_eta = new TLegend(0.55, 0.22, 0.90, 0.48);
[92237ca]	2397
[341014c]	2398	TGraphErrors *gr_recele_eff_pt = new TGraphErrors[n_etabins];
	2399	TGraphErrors *gr_recele_eff_eta = new TGraphErrors[n_ptbins];
	2400	TH1D h_recele_eff_pt, h_recele_eff_eta;
[92237ca]	2401
[341014c]	2402	// loop over eta bins
	2403	for(k = 0; k < etaVals.size() - 1; k++)
	2404	{
[92237ca]	2405
[341014c]	2406	h_recele_eff_pt = GetEffPt<Electron>(branchElectron, branchParticleElectron, "Electron", 11, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderElectron);
	2407	gr_recele_eff_pt[k] = TGraphErrors(h_recele_eff_pt);
[92237ca]	2408
[341014c]	2409	s_etaMin = Form("%.1f", etaVals.at(k));
	2410	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]	2411
[341014c]	2412	s_eta = "e^{ #pm} , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	2413
[341014c]	2414	gr_recele_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]	2415
[341014c]	2416	addResoGraph(mg_recele_eff_pt, &gr_recele_eff_pt[k], leg_recele_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
	2417	}
[92237ca]	2418
[341014c]	2419	// loop over pt
	2420	for(k = 0; k < ptVals.size(); k++)
	2421	{
	2422	h_recele_eff_eta = GetEffEta<Electron>(branchElectron, branchParticleElectron, "Electron", 11, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderElectron);
	2423	gr_recele_eff_eta[k] = TGraphErrors(h_recele_eff_eta);
[92237ca]	2424
[341014c]	2425	s_pt = Form("e^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
	2426	if(ptVals.at(k) >= 1000.) s_pt = Form("e^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	2427
[341014c]	2428	addResoGraph(mg_recele_eff_eta, &gr_recele_eff_eta[k], leg_recele_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
	2429	}
[92237ca]	2430
[341014c]	2431	TCanvas *c_recele_eff_pt = new TCanvas("", "", 800, 600);
[1408174]	2432
[341014c]	2433	mg_recele_eff_pt->Draw("APE");
	2434	DrawAxis(mg_recele_eff_pt, leg_recele_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "reconstruction efficiency (%)", true, false);
	2435	leg_recele_eff_pt->Draw();
	2436	pave->Draw();
[1408174]	2437
[341014c]	2438	c_recele_eff_pt->Print(pdfOutput, "pdf");
	2439	c_recele_eff_pt->Print(figPath + "img_recele_eff_pt.pdf", "pdf");
	2440	c_recele_eff_pt->Print(figPath + "img_recele_eff_pt.png", "png");
[1408174]	2441
[341014c]	2442	TCanvas *c_recele_eff_eta = new TCanvas("", "", 800, 600);
[1408174]	2443
[341014c]	2444	mg_recele_eff_eta->Draw("APE");
	2445	DrawAxis(mg_recele_eff_eta, leg_recele_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "reconstruction efficiency (%)", false, false);
	2446	leg_recele_eff_eta->Draw();
	2447	pave->Draw();
[1408174]	2448
[341014c]	2449	c_recele_eff_eta->Print(pdfOutput, "pdf");
	2450	c_recele_eff_eta->Print(figPath + "img_recele_eff_eta.pdf", "pdf");
	2451	c_recele_eff_eta->Print(figPath + "img_recele_eff_eta.png", "png");
[92237ca]	2452
[341014c]	2453	/////////////////////////////////////////
	2454	// Muon Reconstruction Efficiency ///
	2455	/////////////////////////////////////////
[92237ca]	2456
[341014c]	2457	TMultiGraph *mg_recmu_eff_pt = new TMultiGraph("", "");
	2458	TMultiGraph *mg_recmu_eff_eta = new TMultiGraph("", "");
[92237ca]	2459
[341014c]	2460	TLegend *leg_recmu_eff_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
	2461	TLegend *leg_recmu_eff_eta = new TLegend(0.55, 0.22, 0.90, 0.48);
[1408174]	2462
[341014c]	2463	TGraphErrors *gr_recmu_eff_pt = new TGraphErrors[n_etabins];
	2464	TGraphErrors *gr_recmu_eff_eta = new TGraphErrors[n_ptbins];
	2465	TH1D h_recmu_eff_pt, h_recmu_eff_eta;
[1408174]	2466
[341014c]	2467	// loop over eta bins
	2468	for(k = 0; k < etaVals.size() - 1; k++)
	2469	{
[1408174]	2470
[341014c]	2471	h_recmu_eff_pt = GetEffPt<Muon>(branchMuon, branchParticleMuon, "muon", 13, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderMuon);
	2472	gr_recmu_eff_pt[k] = TGraphErrors(h_recmu_eff_pt);
[1408174]	2473
[341014c]	2474	s_etaMin = Form("%.1f", etaVals.at(k));
	2475	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]	2476
[341014c]	2477	s_eta = "#mu^{ #pm} , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[1408174]	2478
[341014c]	2479	gr_recmu_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[1408174]	2480
[341014c]	2481	addResoGraph(mg_recmu_eff_pt, &gr_recmu_eff_pt[k], leg_recmu_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
	2482	}
[1408174]	2483
[341014c]	2484	// loop over pt
	2485	for(k = 0; k < ptVals.size(); k++)
	2486	{
	2487	h_recmu_eff_eta = GetEffEta<Muon>(branchMuon, branchParticleMuon, "muon", 13, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderMuon);
	2488	gr_recmu_eff_eta[k] = TGraphErrors(h_recmu_eff_eta);
[1408174]	2489
[341014c]	2490	s_pt = Form("#mu^{ #pm} , p_{T} = %.0f GeV", ptVals.at(k));
	2491	if(ptVals.at(k) >= 1000.) s_pt = Form("#mu^{ #pm} , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1408174]	2492
[341014c]	2493	addResoGraph(mg_recmu_eff_eta, &gr_recmu_eff_eta[k], leg_recmu_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
	2494	}
[1408174]	2495
[341014c]	2496	TCanvas *c_recmu_eff_pt = new TCanvas("", "", 800, 600);
[1408174]	2497
[341014c]	2498	mg_recmu_eff_pt->Draw("APE");
	2499	DrawAxis(mg_recmu_eff_pt, leg_recmu_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "reconstruction efficiency (%)", true, false);
	2500	leg_recmu_eff_pt->Draw();
	2501	pave->Draw();
[1408174]	2502
[341014c]	2503	c_recmu_eff_pt->Print(pdfOutput, "pdf");
	2504	c_recmu_eff_pt->Print(figPath + "img_recmu_eff_pt.pdf", "pdf");
	2505	c_recmu_eff_pt->Print(figPath + "img_recmu_eff_pt.png", "png");
[1408174]	2506
[341014c]	2507	TCanvas *c_recmu_eff_eta = new TCanvas("", "", 800, 600);
[1408174]	2508
[341014c]	2509	mg_recmu_eff_eta->Draw("APE");
	2510	DrawAxis(mg_recmu_eff_eta, leg_recmu_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "reconstruction efficiency (%)", false, false);
	2511	leg_recmu_eff_eta->Draw();
	2512	pave->Draw();
[1408174]	2513
[341014c]	2514	c_recmu_eff_eta->Print(pdfOutput, "pdf");
	2515	c_recmu_eff_eta->Print(figPath + "img_recmu_eff_eta.pdf", "pdf");
	2516	c_recmu_eff_eta->Print(figPath + "img_recmu_eff_eta.png", "png");
[1408174]	2517
[341014c]	2518	/////////////////////////////////////////
	2519	// Photon Reconstruction Efficiency ///
	2520	/////////////////////////////////////////
[1408174]	2521
[341014c]	2522	TMultiGraph *mg_recpho_eff_pt = new TMultiGraph("", "");
	2523	TMultiGraph *mg_recpho_eff_eta = new TMultiGraph("", "");
[92237ca]	2524
[341014c]	2525	TLegend *leg_recpho_eff_pt = new TLegend(0.55, 0.22, 0.90, 0.48);
	2526	TLegend *leg_recpho_eff_eta = new TLegend(0.55, 0.22, 0.90, 0.48);
[92237ca]	2527
[341014c]	2528	TGraphErrors *gr_recpho_eff_pt = new TGraphErrors[n_etabins];
	2529	TGraphErrors *gr_recpho_eff_eta = new TGraphErrors[n_ptbins];
	2530	TH1D h_recpho_eff_pt, h_recpho_eff_eta;
[1408174]	2531
[341014c]	2532	// loop over eta bins
	2533	for(k = 0; k < etaVals.size() - 1; k++)
	2534	{
[1408174]	2535
[341014c]	2536	h_recpho_eff_pt = GetEffPt<Photon>(branchPhoton, branchParticlePhoton, "Photon", 22, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderPhoton);
	2537	gr_recpho_eff_pt[k] = TGraphErrors(h_recpho_eff_pt);
[1408174]	2538
[341014c]	2539	s_etaMin = Form("%.1f", etaVals.at(k));
	2540	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]	2541
[341014c]	2542	s_eta = "#gamma , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	2543
[341014c]	2544	gr_recpho_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[1408174]	2545
[341014c]	2546	addResoGraph(mg_recpho_eff_pt, &gr_recpho_eff_pt[k], leg_recpho_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
	2547	}
[1160e4f]	2548
[341014c]	2549	// loop over pt
	2550	for(k = 0; k < ptVals.size(); k++)
	2551	{
	2552	h_recpho_eff_eta = GetEffEta<Photon>(branchPhoton, branchParticlePhoton, "Photon", 22, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderPhoton);
	2553	gr_recpho_eff_eta[k] = TGraphErrors(h_recpho_eff_eta);
[1160e4f]	2554
[341014c]	2555	s_pt = Form("#gamma , p_{T} = %.0f GeV", ptVals.at(k));
	2556	if(ptVals.at(k) >= 1000.) s_pt = Form("#gamma , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1160e4f]	2557
[341014c]	2558	addResoGraph(mg_recpho_eff_eta, &gr_recpho_eff_eta[k], leg_recpho_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
	2559	}
[1408174]	2560
[341014c]	2561	TCanvas *c_recpho_eff_pt = new TCanvas("", "", 800, 600);
[1408174]	2562
[341014c]	2563	mg_recpho_eff_pt->Draw("APE");
	2564	DrawAxis(mg_recpho_eff_pt, leg_recpho_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "reconstruction efficiency (%)", true, false);
	2565	leg_recpho_eff_pt->Draw();
	2566	pave->Draw();
[92237ca]	2567
[341014c]	2568	c_recpho_eff_pt->Print(pdfOutput, "pdf");
	2569	c_recpho_eff_pt->Print(figPath + "img_recpho_eff_pt.pdf", "pdf");
	2570	c_recpho_eff_pt->Print(figPath + "img_recpho_eff_pt.png", "png");
[92237ca]	2571
[341014c]	2572	TCanvas *c_recpho_eff_eta = new TCanvas("", "", 800, 600);
[92237ca]	2573
[341014c]	2574	mg_recpho_eff_eta->Draw("APE");
	2575	DrawAxis(mg_recpho_eff_eta, leg_recpho_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "reconstruction efficiency (%)", false, false);
	2576	leg_recpho_eff_eta->Draw();
	2577	pave->Draw();
[92237ca]	2578
[341014c]	2579	c_recpho_eff_eta->Print(pdfOutput, "pdf");
	2580	c_recpho_eff_eta->Print(figPath + "img_recpho_eff_eta.pdf", "pdf");
	2581	c_recpho_eff_eta->Print(figPath + "img_recpho_eff_eta.png", "png");
[92237ca]	2582
[341014c]	2583	/////////////////////////////////////////
	2584	// B-jets Efficiency/ mistag rates ///
	2585	/////////////////////////////////////////
[92237ca]	2586
[341014c]	2587	TMultiGraph *mg_recbjet_eff_pt = new TMultiGraph("", "");
	2588	TMultiGraph *mg_recbjet_eff_eta = new TMultiGraph("", "");
[92237ca]	2589
[341014c]	2590	TLegend *leg_recbjet_eff_pt = new TLegend(0.50, 0.22, 0.90, 0.48);
	2591	TLegend *leg_recbjet_eff_eta = new TLegend(0.50, 0.22, 0.90, 0.48);
[92237ca]	2592
[341014c]	2593	TGraphErrors *gr_recbjet_eff_pt = new TGraphErrors[n_etabins];
	2594	TGraphErrors *gr_recbjet_eff_eta = new TGraphErrors[n_ptbins];
	2595	TH1D h_recbjet_eff_pt, h_recbjet_eff_eta;
[92237ca]	2596
[341014c]	2597	// loop over eta bins
	2598	for(k = 0; k < etaVals.size() - 1; k++)
	2599	{
[92237ca]	2600
[341014c]	2601	h_recbjet_eff_pt = GetJetEffPt<Jet>(branchPFBJet, "BJet", 5, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderBJet);
	2602	//h_recbjet_eff_pt = GetEffPt<Jet>(branchPFBJet, branchParticleBJet, "BJet", 5, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderBJet);
	2603	gr_recbjet_eff_pt[k] = TGraphErrors(h_recbjet_eff_pt);
[92237ca]	2604
[341014c]	2605	s_etaMin = Form("%.1f", etaVals.at(k));
	2606	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]	2607
[341014c]	2608	s_eta = "b-jet , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	2609
[341014c]	2610	gr_recbjet_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]	2611
[341014c]	2612	addResoGraph(mg_recbjet_eff_pt, &gr_recbjet_eff_pt[k], leg_recbjet_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
	2613	}
[1408174]	2614
[341014c]	2615	// loop over pt
	2616	for(k = 0; k < ptVals.size(); k++)
	2617	{
	2618	h_recbjet_eff_eta = GetJetEffEta<Jet>(branchPFBJet, "BJet", 5, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderBJet);
	2619	//h_recbjet_eff_eta = GetEffEta<Jet>(branchPFBJet, branchParticleBJet, "BJet", 5, 0.5ptVals.at(k), 2.0ptVals.at(k) ,etaMin, etaMax , treeReaderBJet);
	2620	gr_recbjet_eff_eta[k] = TGraphErrors(h_recbjet_eff_eta);
[734b267]	2621
[341014c]	2622	s_pt = Form("b-jet , p_{T} = %.0f GeV", ptVals.at(k));
	2623	if(ptVals.at(k) >= 1000.) s_pt = Form("b-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1408174]	2624
[341014c]	2625	addResoGraph(mg_recbjet_eff_eta, &gr_recbjet_eff_eta[k], leg_recbjet_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
	2626	}
[1408174]	2627
[341014c]	2628	TCanvas *c_recbjet_eff_pt = new TCanvas("", "", 800, 600);
[1408174]	2629
[341014c]	2630	mg_recbjet_eff_pt->Draw("APE");
	2631	DrawAxis(mg_recbjet_eff_pt, leg_recbjet_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "b - tag efficiency (%)", true, false);
	2632	leg_recbjet_eff_pt->Draw();
	2633	pave->Draw();
[1408174]	2634
[341014c]	2635	c_recbjet_eff_pt->Print(pdfOutput, "pdf");
	2636	c_recbjet_eff_pt->Print(figPath + "img_recbjet_eff_pt.pdf", "pdf");
	2637	c_recbjet_eff_pt->Print(figPath + "img_recbjet_eff_pt.png", "png");
[1408174]	2638
[341014c]	2639	TCanvas *c_recbjet_eff_eta = new TCanvas("", "", 800, 600);
[1408174]	2640
[341014c]	2641	mg_recbjet_eff_eta->Draw("APE");
	2642	DrawAxis(mg_recbjet_eff_eta, leg_recbjet_eff_eta, etaMin, etaMax, 0.0, 100, " #eta ", "b - tag efficiency (%)", false, false);
	2643	leg_recbjet_eff_eta->Draw();
	2644	pave->Draw();
[1408174]	2645
[341014c]	2646	c_recbjet_eff_eta->Print(pdfOutput, "pdf");
	2647	c_recbjet_eff_eta->Print(figPath + "img_recbjet_eff_eta.pdf", "pdf");
	2648	c_recbjet_eff_eta->Print(figPath + "img_recbjet_eff_eta.png", "png");
[1408174]	2649
[341014c]	2650	// ------ c - mistag ------
[1408174]	2651
[341014c]	2652	TMultiGraph *mg_recbjet_cmis_pt = new TMultiGraph("", "");
	2653	TMultiGraph *mg_recbjet_cmis_eta = new TMultiGraph("", "");
[1408174]	2654
[341014c]	2655	TLegend *leg_recbjet_cmis_pt = new TLegend(0.50, 0.64, 0.90, 0.90);
	2656	TLegend *leg_recbjet_cmis_eta = new TLegend(0.50, 0.64, 0.90, 0.90);
[1408174]	2657
[341014c]	2658	TGraphErrors *gr_recbjet_cmis_pt = new TGraphErrors[n_etabins];
	2659	TGraphErrors *gr_recbjet_cmis_eta = new TGraphErrors[n_ptbins];
	2660	TH1D h_recbjet_cmis_pt, h_recbjet_cmis_eta;
[1408174]	2661
[341014c]	2662	// loop over eta bins
	2663	for(k = 0; k < etaVals.size() - 1; k++)
	2664	{
[1408174]	2665
[341014c]	2666	h_recbjet_cmis_pt = GetJetEffPt<Jet>(branchPFCJet, "CJet", 4, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderCJet);
	2667	//h_recbjet_cmis_pt = GetEffPt<Jet>(branchPFCJet, branchParticleCJet, "CJet", 4, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderCJet);
	2668	gr_recbjet_cmis_pt[k] = TGraphErrors(h_recbjet_cmis_pt);
[1408174]	2669
[341014c]	2670	s_etaMin = Form("%.1f", etaVals.at(k));
	2671	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]	2672
[341014c]	2673	s_eta = "c-jet , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[1408174]	2674
[341014c]	2675	gr_recbjet_cmis_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[1408174]	2676
[341014c]	2677	addResoGraph(mg_recbjet_cmis_pt, &gr_recbjet_cmis_pt[k], leg_recbjet_cmis_pt, markerStyles.at(k), colors.at(k), s_eta);
	2678	}
[1408174]	2679
[341014c]	2680	// loop over pt
	2681	for(k = 0; k < ptVals.size(); k++)
	2682	{
	2683	h_recbjet_cmis_eta = GetJetEffEta<Jet>(branchPFCJet, "CJet", 4, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderCJet);
	2684	//h_recbjet_cmis_eta = GetEffEta<Jet>(branchPFCJet, branchParticleCJet, "CJet", 4, 0.5ptVals.at(k), 2.0ptVals.at(k) ,etaMin, etaMax , treeReaderCJet);
	2685	gr_recbjet_cmis_eta[k] = TGraphErrors(h_recbjet_cmis_eta);
[1408174]	2686
[341014c]	2687	s_pt = Form("c-jet , p_{T} = %.0f GeV", ptVals.at(k));
	2688	if(ptVals.at(k) >= 1000.) s_pt = Form("c-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1408174]	2689
[341014c]	2690	addResoGraph(mg_recbjet_cmis_eta, &gr_recbjet_cmis_eta[k], leg_recbjet_cmis_eta, markerStyles.at(k), colors.at(k), s_pt);
	2691	}
[1408174]	2692
[341014c]	2693	TCanvas *c_recbjet_cmis_pt = new TCanvas("", "", 800, 600);
[1408174]	2694
[341014c]	2695	mg_recbjet_cmis_pt->Draw("APE");
	2696	DrawAxis(mg_recbjet_cmis_pt, leg_recbjet_cmis_pt, ptMin, ptMax, 0.0, 20, "p_{T} [GeV]", "c - mistag rate (%)", true, false);
	2697	leg_recbjet_cmis_pt->Draw();
	2698	pave->Draw();
[1408174]	2699
[341014c]	2700	c_recbjet_cmis_pt->Print(pdfOutput, "pdf");
	2701	c_recbjet_cmis_pt->Print(figPath + "img_recbjet_cmis_pt.pdf", "pdf");
	2702	c_recbjet_cmis_pt->Print(figPath + "img_recbjet_cmis_pt.png", "png");
[1408174]	2703
[341014c]	2704	TCanvas *c_recbjet_cmis_eta = new TCanvas("", "", 800, 600);
[1408174]	2705
[341014c]	2706	mg_recbjet_cmis_eta->Draw("APE");
	2707	DrawAxis(mg_recbjet_cmis_eta, leg_recbjet_cmis_eta, etaMin, etaMax, 0.0, 20, " #eta ", "c - mistag rate (%)", false, false);
	2708	leg_recbjet_cmis_eta->Draw();
	2709	pave->Draw();
[1408174]	2710
[341014c]	2711	c_recbjet_cmis_eta->Print(pdfOutput, "pdf");
	2712	c_recbjet_cmis_eta->Print(figPath + "img_recbjet_cmis_eta.pdf", "pdf");
	2713	c_recbjet_cmis_eta->Print(figPath + "img_recbjet_cmis_eta.png", "png");
[1408174]	2714
[341014c]	2715	// ------ light - mistag ------
[1408174]	2716
[341014c]	2717	TMultiGraph *mg_recbjet_lmis_pt = new TMultiGraph("", "");
	2718	TMultiGraph *mg_recbjet_lmis_eta = new TMultiGraph("", "");
[cb159ed]	2719
[341014c]	2720	TLegend *leg_recbjet_lmis_pt = new TLegend(0.50, 0.64, 0.90, 0.90);
	2721	TLegend *leg_recbjet_lmis_eta = new TLegend(0.50, 0.64, 0.90, 0.90);
[1408174]	2722
[341014c]	2723	TGraphErrors *gr_recbjet_lmis_pt = new TGraphErrors[n_etabins];
	2724	TGraphErrors *gr_recbjet_lmis_eta = new TGraphErrors[n_ptbins];
	2725	TH1D h_recbjet_lmis_pt, h_recbjet_lmis_eta;
[1408174]	2726
[341014c]	2727	// loop over eta bins
	2728	for(k = 0; k < etaVals.size() - 1; k++)
	2729	{
[1408174]	2730
[341014c]	2731	h_recbjet_lmis_pt = GetJetEffPt<Jet>(branchJet, "Jet", 1, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderJet);
	2732	//h_recbjet_lmis_pt = GetEffPt<Jet>(branchJet, branchParticleJet, "Jet", 1, ptMin, ptMax, etaVals.at(k), etaVals.at(k+1), treeReaderJet);
	2733	gr_recbjet_lmis_pt[k] = TGraphErrors(h_recbjet_lmis_pt);
[1408174]	2734
[341014c]	2735	s_etaMin = Form("%.1f", etaVals.at(k));
	2736	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]	2737
[341014c]	2738	s_eta = "uds-jet , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	2739
[341014c]	2740	gr_recbjet_lmis_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]	2741
[341014c]	2742	addResoGraph(mg_recbjet_lmis_pt, &gr_recbjet_lmis_pt[k], leg_recbjet_lmis_pt, markerStyles.at(k), colors.at(k), s_eta);
	2743	}
[92237ca]	2744
[341014c]	2745	// loop over pt
	2746	for(k = 0; k < ptVals.size(); k++)
	2747	{
	2748	h_recbjet_lmis_eta = GetJetEffEta<Jet>(branchJet, "Jet", 1, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderJet);
	2749	//h_recbjet_lmis_eta = GetEffEta<Jet>(branchJet, branchParticleJet, "Jet", 1, 0.5ptVals.at(k), 2.0ptVals.at(k) ,etaMin, etaMax , treeReaderJet);
	2750	gr_recbjet_lmis_eta[k] = TGraphErrors(h_recbjet_lmis_eta);
[92237ca]	2751
[341014c]	2752	s_pt = Form("uds-jet , p_{T} = %.0f GeV", ptVals.at(k));
	2753	if(ptVals.at(k) >= 1000.) s_pt = Form("uds-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	2754
[341014c]	2755	addResoGraph(mg_recbjet_lmis_eta, &gr_recbjet_lmis_eta[k], leg_recbjet_lmis_eta, markerStyles.at(k), colors.at(k), s_pt);
	2756	}
[92237ca]	2757
[341014c]	2758	TCanvas *c_recbjet_lmis_pt = new TCanvas("", "", 800, 600);
[92237ca]	2759
[341014c]	2760	mg_recbjet_lmis_pt->Draw("APE");
[92237ca]	2761
[341014c]	2762	DrawAxis(mg_recbjet_lmis_pt, leg_recbjet_lmis_pt, ptMin, ptMax, 0.0, 1.0, "p_{T} [GeV]", "light - mistag rate (%)", true, false);
[92237ca]	2763
[341014c]	2764	leg_recbjet_lmis_pt->Draw();
	2765	pave->Draw();
[92237ca]	2766
[341014c]	2767	c_recbjet_lmis_pt->Print(pdfOutput, "pdf");
	2768	c_recbjet_lmis_pt->Print(figPath + "img_recbjet_lmis_pt.pdf", "pdf");
	2769	c_recbjet_lmis_pt->Print(figPath + "img_recbjet_lmis_pt.png", "png");
[92237ca]	2770
[341014c]	2771	TCanvas *c_recbjet_lmis_eta = new TCanvas("", "", 800, 600);
[fa068d3]	2772
[341014c]	2773	mg_recbjet_lmis_eta->Draw("APE");
	2774	DrawAxis(mg_recbjet_lmis_eta, leg_recbjet_lmis_eta, etaMin, etaMax, 0.0, 1.0, " #eta ", "light - mistag rate (%)", false, false);
	2775	leg_recbjet_lmis_eta->Draw();
	2776	pave->Draw();
[1408174]	2777
[341014c]	2778	c_recbjet_lmis_eta->Print(pdfOutput, "pdf");
	2779	c_recbjet_lmis_eta->Print(figPath + "img_recbjet_lmis_eta.pdf", "pdf");
	2780	c_recbjet_lmis_eta->Print(figPath + "img_recbjet_lmis_eta.png", "png");
[1408174]	2781
[341014c]	2782	///////////////////////////////////////////
	2783	// tau-jets Efficiency/ mistag rates ///
	2784	///////////////////////////////////////////
[92237ca]	2785
[341014c]	2786	TMultiGraph *mg_rectaujet_eff_pt = new TMultiGraph("", "");
	2787	TMultiGraph *mg_rectaujet_eff_eta = new TMultiGraph("", "");
[92237ca]	2788
[341014c]	2789	TLegend *leg_rectaujet_eff_pt = new TLegend(0.50, 0.22, 0.90, 0.48);
	2790	TLegend *leg_rectaujet_eff_eta = new TLegend(0.50, 0.22, 0.90, 0.48);
[92237ca]	2791
[341014c]	2792	TGraphErrors *gr_rectaujet_eff_pt = new TGraphErrors[n_etabins];
	2793	TGraphErrors *gr_rectaujet_eff_eta = new TGraphErrors[n_ptbins];
	2794	TH1D h_rectaujet_eff_pt, h_rectaujet_eff_eta;
[92237ca]	2795
[341014c]	2796	// loop over eta bins
	2797	for(k = 0; k < etaVals.size() - 1; k++)
	2798	{
[92237ca]	2799
[341014c]	2800	h_rectaujet_eff_pt = GetTauEffPt<Jet>(branchPFTauJet, branchParticleTauJet, "TauJet", 15, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderTauJet);
	2801	gr_rectaujet_eff_pt[k] = TGraphErrors(h_rectaujet_eff_pt);
[92237ca]	2802
[341014c]	2803	s_etaMin = Form("%.1f", etaVals.at(k));
	2804	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[92237ca]	2805
[341014c]	2806	s_eta = "#tau-jet , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[92237ca]	2807
[341014c]	2808	gr_rectaujet_eff_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[92237ca]	2809
[341014c]	2810	addResoGraph(mg_rectaujet_eff_pt, &gr_rectaujet_eff_pt[k], leg_rectaujet_eff_pt, markerStyles.at(k), colors.at(k), s_eta);
	2811	}
[92237ca]	2812
[341014c]	2813	// loop over pt
	2814	for(k = 0; k < ptVals.size(); k++)
	2815	{
	2816	h_rectaujet_eff_eta = GetTauEffEta<Jet>(branchPFTauJet, branchParticleTauJet, "TauJet", 15, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderTauJet);
	2817	gr_rectaujet_eff_eta[k] = TGraphErrors(h_rectaujet_eff_eta);
[92237ca]	2818
[341014c]	2819	s_pt = Form("#tau-jet , p_{T} = %.0f GeV", ptVals.at(k));
	2820	if(ptVals.at(k) >= 1000.) s_pt = Form("#tau-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[92237ca]	2821
[341014c]	2822	addResoGraph(mg_rectaujet_eff_eta, &gr_rectaujet_eff_eta[k], leg_rectaujet_eff_eta, markerStyles.at(k), colors.at(k), s_pt);
	2823	}
[1408174]	2824
[341014c]	2825	TCanvas *c_rectaujet_eff_pt = new TCanvas("", "", 800, 600);
[1408174]	2826
[341014c]	2827	mg_rectaujet_eff_pt->Draw("APE");
	2828	DrawAxis(mg_rectaujet_eff_pt, leg_rectaujet_eff_pt, ptMin, ptMax, 0.0, 100, "p_{T} [GeV]", "#tau - tag efficiency (%)", true, false);
	2829	leg_rectaujet_eff_pt->Draw();
	2830	pave->Draw();
[1408174]	2831
[341014c]	2832	c_rectaujet_eff_pt->Print(pdfOutput, "pdf");
	2833	c_rectaujet_eff_pt->Print(figPath + "img_rectaujet_eff_pt.pdf", "pdf");
	2834	c_rectaujet_eff_pt->Print(figPath + "img_rectaujet_eff_pt.png", "png");
[1408174]	2835
[341014c]	2836	TCanvas *c_rectaujet_eff_eta = new TCanvas("", "", 800, 600);
[1408174]	2837
[341014c]	2838	mg_rectaujet_eff_eta->Draw("APE");
	2839	DrawAxis(mg_rectaujet_eff_eta, leg_rectaujet_eff_eta, etaMin, etaMax, 0.0, 100., " #eta ", "#tau - tag efficiency (%)", false, false);
	2840	leg_rectaujet_eff_eta->Draw();
	2841	pave->Draw();
[1408174]	2842
[341014c]	2843	c_rectaujet_eff_eta->Print(pdfOutput, "pdf");
	2844	c_rectaujet_eff_eta->Print(figPath + "img_rectaujet_eff_eta.pdf", "pdf");
	2845	c_rectaujet_eff_eta->Print(figPath + "img_rectaujet_eff_eta.png", "png");
[1408174]	2846
[341014c]	2847	//--------------- tau mistag rate ----------
[1408174]	2848
[341014c]	2849	TMultiGraph *mg_rectaujet_mis_pt = new TMultiGraph("", "");
	2850	TMultiGraph *mg_rectaujet_mis_eta = new TMultiGraph("", "");
[1408174]	2851
[341014c]	2852	TLegend *leg_rectaujet_mis_pt = new TLegend(0.50, 0.64, 0.90, 0.90);
	2853	TLegend *leg_rectaujet_mis_eta = new TLegend(0.50, 0.64, 0.90, 0.90);
[1408174]	2854
[341014c]	2855	TGraphErrors *gr_rectaujet_mis_pt = new TGraphErrors[n_etabins];
	2856	TGraphErrors *gr_rectaujet_mis_eta = new TGraphErrors[n_ptbins];
	2857	TH1D h_rectaujet_mis_pt, h_rectaujet_mis_eta;
[1408174]	2858
[341014c]	2859	// loop over eta bins
	2860	for(k = 0; k < etaVals.size() - 1; k++)
	2861	{
[1408174]	2862
[341014c]	2863	h_rectaujet_mis_pt = GetTauEffPt<Jet>(branchJet, branchParticleJet, "TauJet", 1, ptMin, ptMax, etaVals.at(k), etaVals.at(k + 1), treeReaderJet);
	2864	gr_rectaujet_mis_pt[k] = TGraphErrors(h_rectaujet_mis_pt);
[1408174]	2865
[341014c]	2866	s_etaMin = Form("%.1f", etaVals.at(k));
	2867	s_etaMax = Form("%.1f", etaVals.at(k + 1));
[1408174]	2868
[341014c]	2869	s_eta = "uds-jet , " + s_etaMin + " < \| #eta \| < " + s_etaMax;
[1408174]	2870
[341014c]	2871	gr_rectaujet_mis_pt[k].SetName("recEff_" + s_etaMin + "_" + s_etaMax);
[1408174]	2872
[341014c]	2873	addResoGraph(mg_rectaujet_mis_pt, &gr_rectaujet_mis_pt[k], leg_rectaujet_mis_pt, markerStyles.at(k), colors.at(k), s_eta);
	2874	}
[1408174]	2875
[341014c]	2876	// loop over pt
	2877	for(k = 0; k < ptVals.size(); k++)
	2878	{
	2879	h_rectaujet_mis_eta = GetTauEffEta<Jet>(branchJet, branchParticleJet, "TauJet", 1, 0.5 * ptVals.at(k), 2.0 * ptVals.at(k), etaMin, etaMax, treeReaderJet);
	2880	gr_rectaujet_mis_eta[k] = TGraphErrors(h_rectaujet_mis_eta);
[1408174]	2881
[341014c]	2882	s_pt = Form("uds-jet , p_{T} = %.0f GeV", ptVals.at(k));
	2883	if(ptVals.at(k) >= 1000.) s_pt = Form("uds-jet , p_{T} = %.0f TeV", ptVals.at(k) / 1000.);
[1408174]	2884
[341014c]	2885	addResoGraph(mg_rectaujet_mis_eta, &gr_rectaujet_mis_eta[k], leg_rectaujet_mis_eta, markerStyles.at(k), colors.at(k), s_pt);
	2886	}
[1408174]	2887
[341014c]	2888	TCanvas *c_rectaujet_mis_pt = new TCanvas("", "", 800, 600);
[1408174]	2889
[341014c]	2890	mg_rectaujet_mis_pt->Draw("APE");
	2891	DrawAxis(mg_rectaujet_mis_pt, leg_rectaujet_mis_pt, ptMin, ptMax, 0.0, 5., "p_{T} [GeV]", "#tau - mistag(%)", true, false);
	2892	leg_rectaujet_mis_pt->Draw();
	2893	pave->Draw();
[1408174]	2894
[341014c]	2895	c_rectaujet_mis_pt->Print(pdfOutput, "pdf");
	2896	c_rectaujet_mis_pt->Print(figPath + "img_rectaujet_mis_pt.pdf", "pdf");
	2897	c_rectaujet_mis_pt->Print(figPath + "img_rectaujet_mis_pt.png", "png");
[1408174]	2898
[341014c]	2899	TCanvas *c_rectaujet_mis_eta = new TCanvas("", "", 800, 600);
[1408174]	2900
[341014c]	2901	mg_rectaujet_mis_eta->Draw("APE");
	2902	DrawAxis(mg_rectaujet_mis_eta, leg_rectaujet_mis_eta, etaMin, etaMax, 0.0, 5., " #eta ", "#tau - mistag (%)", false, false);
	2903	leg_rectaujet_mis_eta->Draw();
	2904	pave->Draw();
[1408174]	2905
[341014c]	2906	c_rectaujet_mis_eta->Print(pdfOutput + ")", "pdf");
	2907	c_rectaujet_mis_eta->Print(figPath + "img_rectaujet_mis_eta.pdf", "pdf");
	2908	c_rectaujet_mis_eta->Print(figPath + "img_rectaujet_mis_eta.png", "png");
[1408174]	2909
[341014c]	2910	// ---- store resolution histograms in the output (for leave efficiencies out) ---
[1408174]	2911
[341014c]	2912	TFile *fout = new TFile(outputFile, "recreate");
[1408174]	2913
[341014c]	2914	for(int bin = 0; bin < Nbins; bin++)
[1408174]	2915	{
	2916
[341014c]	2917	for(k = 0; k < etaVals.size() - 1; k++)
[92237ca]	2918	{
[341014c]	2919	plots_trkpi_res_pt[k].at(bin).resolHist->Write();
	2920	plots_trkele_res_pt[k].at(bin).resolHist->Write();
	2921	plots_trkmu_res_pt[k].at(bin).resolHist->Write();
	2922	plots_ecal_res_e[k].at(bin).resolHist->Write();
	2923	plots_hcal_res_e[k].at(bin).resolHist->Write();
	2924	plots_pfele_res_e[k].at(bin).resolHist->Write();
	2925	plots_pfpi_res_e[k].at(bin).resolHist->Write();
	2926	plots_pfjet_res_e[k].at(bin).resolHist->Write();
	2927	plots_cajet_res_e[k].at(bin).resolHist->Write();
[92237ca]	2928	}
	2929
[341014c]	2930	for(k = 0; k < ptVals.size(); k++)
[92237ca]	2931	{
[341014c]	2932	plots_trkpi_res_eta[k].at(bin).resolHist->Write();
	2933	plots_trkele_res_eta[k].at(bin).resolHist->Write();
	2934	plots_trkmu_res_eta[k].at(bin).resolHist->Write();
	2935	plots_ecal_res_eta[k].at(bin).resolHist->Write();
	2936	plots_hcal_res_eta[k].at(bin).resolHist->Write();
	2937	plots_pfele_res_eta[k].at(bin).resolHist->Write();
	2938	plots_pfpi_res_eta[k].at(bin).resolHist->Write();
	2939	plots_pfjet_res_eta[k].at(bin).resolHist->Write();
	2940	plots_cajet_res_eta[k].at(bin).resolHist->Write();
[92237ca]	2941	}
[1408174]	2942
[92237ca]	2943	plots_pfmet.at(bin).resolHist->Write();
	2944	plots_camet.at(bin).resolHist->Write();
[1408174]	2945	}
	2946
	2947	fout->Write();
	2948
	2949	cout << "** Exiting..." << endl;
	2950
	2951	delete treeReaderElectron;
	2952	delete treeReaderMuon;
	2953	delete treeReaderPhoton;
	2954	delete treeReaderJet;
	2955	delete treeReaderBJet;
	2956	delete treeReaderTauJet;
	2957	delete chainElectron;
	2958	delete chainMuon;
	2959	delete chainPhoton;
	2960	delete chainJet;
	2961	delete chainBJet;
	2962	delete chainTauJet;
	2963	}
[734b267]	2964
	2965	//------------------------------------------------------------------------------
	2966
	2967	int main(int argc, char *argv[])
	2968	{
[181c061]	2969	char *appName = "DelphesValidation";
[734b267]	2970
[85ad2b9]	2971	if(argc != 12)
[734b267]	2972	{
[181c061]	2973	cout << " Usage: " << appName << " input_file_electron input_file_muon input_file_photon input_file_jet input_file_bjet input_file_taujet output_file version" << endl;
	2974	cout << " input_file_pion - input file in ROOT format ('Delphes' tree)," << endl;
	2975	cout << " input_file_electron - input file in ROOT format ('Delphes' tree)," << endl;
[1408174]	2976	cout << " input_file_muon - input file in ROOT format ('Delphes' tree)," << endl;
	2977	cout << " input_file_photon - input file in ROOT format ('Delphes' tree)," << endl;
[92237ca]	2978	cout << " input_file_neutralhadron - input file in ROOT format ('Delphes' tree)," << endl;
[1408174]	2979	cout << " input_file_jet - input file in ROOT format ('Delphes' tree)," << endl;
	2980	cout << " input_file_bjet - input file in ROOT format ('Delphes' tree)," << endl;
[92237ca]	2981	cout << " input_file_cjet - input file in ROOT format ('Delphes' tree)," << endl;
[1408174]	2982	cout << " input_file_taujet - input file in ROOT format ('Delphes' tree)," << endl;
[181c061]	2983	cout << " output_file - output file in ROOT format," << endl;
	2984	cout << " version - Delphes version" << endl;
[85ad2b9]	2985
[734b267]	2986	return 1;
	2987	}
	2988
	2989	gROOT->SetBatch();
	2990
	2991	int appargc = 1;
	2992	char *appargv[] = {appName};
	2993	TApplication app(appName, &appargc, appargv);
	2994
[181c061]	2995	DelphesValidation(argv[1], argv[2], argv[3], argv[4], argv[5], argv[6], argv[7], argv[8], argv[9], argv[10], argv[11]);
[734b267]	2996	}

Note: See TracBrowser for help on using the repository browser.

Download in other formats: