Fork me on GitHub

Context Navigation

source: git/validation/DelphesValidation.cpp@ 861ad5a

Visit:

ImprovedOutputFile Timing

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

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

Download in other formats: