Fork me on GitHub

Context Navigation

source: git/examples/Validation.cpp@ 2b1e147

Visit:

ImprovedOutputFile Timing dual_readout llp

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

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

Download in other formats: