ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/CmsHi/JetAnalysis/src/RecHitComparison.cc

Comparing UserCode/CmsHi/JetAnalysis/src/RecHitComparison.cc (file contents):
Revision 1.1 by yilmaz, Thu Sep 23 13:28:28 2010 UTC vs.
Revision 1.2 by yilmaz, Mon Oct 18 16:13:37 2010 UTC

#	Line 52 | Line 52
52		#include "SimDataFormats/HiGenData/interface/GenHIEvent.h"
53		#include "DataFormats/PatCandidates/interface/Jet.h"
54		#include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
55	+	#include "DataFormats/HcalRecHit/interface/HcalRecHitCollections.h"
56	+	#include "DataFormats/HcalDetId/interface/HcalDetId.h"
57
58		#include "TNtuple.h"
59
#	Line 79 | Line 81 | class RecHitComparison : public edm::EDA
81		edm::Handle<vector<double> > ktRhos;
82		edm::Handle<vector<double> > akRhos;
83
84	<	edm::Handle<edm::SortedCollection<EcalRecHit,edm::StrictWeakOrdering<EcalRecHit> > > jets1;
85	<	edm::Handle<edm::SortedCollection<EcalRecHit,edm::StrictWeakOrdering<EcalRecHit> > > jets2;
84	>	edm::Handle<edm::SortedCollection<EcalRecHit,edm::StrictWeakOrdering<EcalRecHit> > > ebHits1;
85	>	edm::Handle<edm::SortedCollection<EcalRecHit,edm::StrictWeakOrdering<EcalRecHit> > > ebHits2;
86	>	edm::Handle<edm::SortedCollection<EcalRecHit,edm::StrictWeakOrdering<EcalRecHit> > > eeHits1;
87	>	edm::Handle<edm::SortedCollection<EcalRecHit,edm::StrictWeakOrdering<EcalRecHit> > > eeHits2;
88	>
89	>	edm::Handle<HFRecHitCollection> hfHits1;
90	>	edm::Handle<HFRecHitCollection> hfHits2;
91	>	edm::Handle<HBHERecHitCollection> hbheHits1;
92	>	edm::Handle<HBHERecHitCollection> hbheHits2;
93
85	–	//   typedef edm::Handle<edm::SortedCollection<EcalRecHit,edm::StrictWeakOrdering<EcalRecHit> > >::const_iterator EcalIterator;
94		typedef vector<EcalRecHit>::const_iterator EcalIterator;
95	+	typedef vector<HFRecHit>::const_iterator HFIterator;
96	+	typedef vector<HBHERecHit>::const_iterator HBHEIterator;
97
98		edm::Handle<reco::CaloJetCollection> signalJets;
99
100	<	TNtuple* nthit;
100	>	edm::InputTag HcalRecHitHFSrc1_;
101	>	edm::InputTag HcalRecHitHFSrc2_;
102	>	edm::InputTag HcalRecHitHBHESrc1_;
103	>	edm::InputTag HcalRecHitHBHESrc2_;
104	>	edm::InputTag EBSrc1_;
105	>	edm::InputTag EBSrc2_;
106	>	edm::InputTag EESrc1_;
107	>	edm::InputTag EESrc2_;
108	>
109	>	TNtuple* ntEB;
110	>	TNtuple* ntEE;
111	>	TNtuple* ntHBHE;
112	>	TNtuple* ntHF;
113		TNtuple* ntjet;
114
115		double cone;
#	Line 114 | Line 136 | RecHitComparison::RecHitComparison(const
136		cone(0.5)
137		{
138		//now do what ever initialization is needed
139	+	HcalRecHitHFSrc1_ = iConfig.getUntrackedParameter<edm::InputTag>("hcalHFRecHitSrc1",edm::InputTag("hfreco"));
140	+	HcalRecHitHFSrc2_ = iConfig.getUntrackedParameter<edm::InputTag>("hcalHFRecHitSrc2",edm::InputTag("hfreco"));
141	+	HcalRecHitHBHESrc1_ = iConfig.getUntrackedParameter<edm::InputTag>("hcalHBHERecHitSrc1",edm::InputTag("hbhereco"));
142	+	HcalRecHitHBHESrc2_ = iConfig.getUntrackedParameter<edm::InputTag>("hcalHBHERecHitSrc2",edm::InputTag("hbhereco"));
143	+	EBSrc1_ = iConfig.getUntrackedParameter<edm::InputTag>("EBRecHitSrc1",edm::InputTag("ecalRecHit","EcalRecHitsEB","RECO"));
144	+	EBSrc2_ = iConfig.getUntrackedParameter<edm::InputTag>("EBRecHitSrc2",edm::InputTag("ecalRecHit","EcalRecHitsEB","SIGNALONLY"));
145	+	EESrc1_ = iConfig.getUntrackedParameter<edm::InputTag>("EERecHitSrc1",edm::InputTag("ecalRecHit","EcalRecHitsEE","RECO"));
146	+	EESrc2_ = iConfig.getUntrackedParameter<edm::InputTag>("EERecHitSrc2",edm::InputTag("ecalRecHit","EcalRecHitsEE","SIGNALONLY"));
147
148		}
149
#	Line 143 | Line 173 | RecHitComparison::analyze(const edm::Eve
173		geo = pGeo.product();
174		}
175
146	–	edm::InputTag jetTag1("ecalRecHit","EcalRecHitsEB","SIGNALONLY");
147	–	edm::InputTag jetTag2("ecalRecHit","EcalRecHitsEB","RECO");
148	–
176		edm::InputTag signalTag("iterativeConePu5CaloJets","","SIGNALONLY");
177		edm::InputTag centTag("hiCentrality","","RECO");
178
179		using namespace edm;
180
181		ev.getByLabel(centTag,cent);
182	<	ev.getByLabel(jetTag1,jets1);
183	<	ev.getByLabel(jetTag2,jets2);
182	>	ev.getByLabel(EBSrc1_,ebHits1);
183	>	ev.getByLabel(EBSrc1_,ebHits2);
184		ev.getByLabel(signalTag,signalJets);
185
186	+	ev.getByLabel(HcalRecHitHFSrc1_,hfHits1);
187	+	ev.getByLabel(HcalRecHitHFSrc2_,hfHits2);
188	+	ev.getByLabel(HcalRecHitHBHESrc1_,hbheHits1);
189	+	ev.getByLabel(HcalRecHitHBHESrc2_,hbheHits2);
190	+
191	+	ev.getByLabel(EESrc1_,eeHits1);
192	+	ev.getByLabel(EESrc2_,eeHits2);
193	+
194		double hf = cent->EtHFhitSum();
195		double sumet = cent->EtMidRapiditySum();
196		int run = ev.id().run();
#	Line 190 | Line 225 | RecHitComparison::analyze(const edm::Eve
225		f3.push_back(0);
226		}
227
228	<	for(unsigned int j1 = 0 ; j1 < jets1->size(); ++j1){
228	>	for(unsigned int j1 = 0 ; j1 < ebHits1->size(); ++j1){
229
230	<	const EcalRecHit& jet1 = (*jets1)[j1];
230	>	const EcalRecHit& jet1 = (*ebHits1)[j1];
231		double e1 = jet1.energy();
232
233		const GlobalPoint& pos1=geo->getPosition(jet1.id());
#	Line 238 | Line 273 | RecHitComparison::analyze(const edm::Eve
273
274		GlobalPoint pos2;
275		double e2 = -1;
276	<	int siz = jets2->size();
277	<
243	<	EcalIterator hitit = jets2->find(jet1.id());
244	<	if(hitit != jets2->end()){
276	>	EcalIterator hitit = ebHits2->find(jet1.id());
277	>	if(hitit != ebHits2->end()){
278		e2 = hitit->energy();
279		pos2=geo->getPosition(hitit->id());
280		}else{
#	Line 252 | Line 285 | RecHitComparison::analyze(const edm::Eve
285		double eta2 = pos2.eta();
286		double phi2 = pos2.eta();
287		double et2 = e2*sin(pos2.theta());
288	<
289	<	if(isjet) nthit->Fill(e1,et1,e2,et2,eta2,phi2,sumet,hf,bin,jetpt,drjet);
290	<
288	>	if(isjet) ntEB->Fill(e1,et1,e2,et2,eta2,phi2,sumet,hf,bin,jetpt,drjet);
289	>	}
290	>
291	>	for(unsigned int i = 0; i < eeHits1->size(); ++i){
292	>	const EcalRecHit & jet1= (*eeHits1)[i];
293	>	double e1 = jet1.energy();
294	>	const GlobalPoint& pos1=geo->getPosition(jet1.id());
295	>	double eta1 = pos1.eta();
296	>	double phi1 = pos1.phi();
297	>	double et1 = e1*sin(pos1.theta());
298	>	double drjet = -1;
299	>	double jetpt = -1;
300	>	bool isjet = false;
301	>	int matchedJet = -1;
302	>	for(unsigned int j = 0 ; j < signalJets->size(); ++j){
303	>	const reco::CaloJet & jet = (*signalJets)[j];
304	>	double dr = reco::deltaR(eta1,phi1,jet.eta(),jet.phi());
305	>	if(dr < cone){
306	>	jetpt = jet.pt();
307	>	drjet = dr;
308	>	isjet = true;
309	>	matchedJet = j;
310	>	}
311	>	}
312	>	GlobalPoint pos2;
313	>	double e2 = -1;
314	>	EcalIterator hitit = eeHits2->find(jet1.id());
315	>	if(hitit != eeHits2->end()){
316	>	e2 = hitit->energy();
317	>	pos2=geo->getPosition(hitit->id());
318	>	}else{
319	>	e2 = 0;
320	>	pos2 = pos1;
321	>	}
322	>	double eta2 = pos2.eta();
323	>	double phi2 = pos2.eta();
324	>	double et2 = e2*sin(pos2.theta());
325	>	if(isjet) ntEE->Fill(e1,et1,e2,et2,eta2,phi2,sumet,hf,bin,jetpt,drjet);
326	>	}
327	>
328	>	for(unsigned int i = 0; i < hbheHits1->size(); ++i){
329	>	const HBHERecHit & jet1= (*hbheHits1)[i];
330	>	double e1 = jet1.energy();
331	>	const GlobalPoint& pos1=geo->getPosition(jet1.id());
332	>	double eta1 = pos1.eta();
333	>	double phi1 = pos1.phi();
334	>	double et1 = e1*sin(pos1.theta());
335	>	double drjet = -1;
336	>	double jetpt = -1;
337	>	bool isjet = false;
338	>	int matchedJet = -1;
339	>	for(unsigned int j = 0 ; j < signalJets->size(); ++j){
340	>	const reco::CaloJet & jet = (*signalJets)[j];
341	>	double dr = reco::deltaR(eta1,phi1,jet.eta(),jet.phi());
342	>	if(dr < cone){
343	>	jetpt = jet.pt();
344	>	drjet = dr;
345	>	isjet = true;
346	>	matchedJet = j;
347	>	}
348	>	}
349	>	GlobalPoint pos2;
350	>	double e2 = -1;
351	>	HBHEIterator hitit = hbheHits2->find(jet1.id());
352	>	if(hitit != hbheHits2->end()){
353	>	e2 = hitit->energy();
354	>	pos2=geo->getPosition(hitit->id());
355	>	}else{
356	>	e2 = 0;
357	>	pos2 = pos1;
358	>	}
359	>	double eta2 = pos2.eta();
360	>	double phi2 = pos2.eta();
361	>	double et2 = e2*sin(pos2.theta());
362	>	if(isjet) ntHBHE->Fill(e1,et1,e2,et2,eta2,phi2,sumet,hf,bin,jetpt,drjet);
363	>	}
364	>
365	>	for(unsigned int i = 0; i < hfHits1->size(); ++i){
366	>	const HFRecHit & jet1= (*hfHits1)[i];
367	>	double e1 = jet1.energy();
368	>	const GlobalPoint& pos1=geo->getPosition(jet1.id());
369	>	double eta1 = pos1.eta();
370	>	double phi1 = pos1.phi();
371	>	double et1 = e1*sin(pos1.theta());
372	>	double drjet = -1;
373	>	double jetpt = -1;
374	>	bool isjet = false;
375	>	int matchedJet = -1;
376	>	for(unsigned int j = 0 ; j < signalJets->size(); ++j){
377	>	const reco::CaloJet & jet = (*signalJets)[j];
378	>	double dr = reco::deltaR(eta1,phi1,jet.eta(),jet.phi());
379	>	if(dr < cone){
380	>	jetpt = jet.pt();
381	>	drjet = dr;
382	>	isjet = true;
383	>	matchedJet = j;
384	>	}
385	>	}
386	>	GlobalPoint pos2;
387	>	double e2 = -1;
388	>	HFIterator hitit = hfHits2->find(jet1.id());
389	>	if(hitit != hfHits2->end()){
390	>	e2 = hitit->energy();
391	>	pos2=geo->getPosition(hitit->id());
392	>	}else{
393	>	e2 = 0;
394	>	pos2 = pos1;
395	>	}
396	>	double eta2 = pos2.eta();
397	>	double phi2 = pos2.eta();
398	>	double et2 = e2*sin(pos2.theta());
399	>	if(isjet) ntHF->Fill(e1,et1,e2,et2,eta2,phi2,sumet,hf,bin,jetpt,drjet);
400		}
401
402		for(unsigned int j1 = 0 ; j1 < signalJets->size(); ++j1){
#	Line 273 | Line 415 | RecHitComparison::analyze(const edm::Eve
415		void
416		RecHitComparison::beginJob()
417		{
418	<	nthit = fs->make<TNtuple>("nthit","","e1:et1:e2:et2:eta:phi:sumet:hf:bin:ptjet:drjet");
418	>	ntEB = fs->make<TNtuple>("ntEB","","e1:et1:e2:et2:eta:phi:sumet:hf:bin:ptjet:drjet");
419	>	ntEE = fs->make<TNtuple>("ntEE","","e1:et1:e2:et2:eta:phi:sumet:hf:bin:ptjet:drjet");
420	>	ntHBHE = fs->make<TNtuple>("ntHBHE","","e1:et1:e2:et2:eta:phi:sumet:hf:bin:ptjet:drjet");
421	>	ntHF = fs->make<TNtuple>("ntHF","","e1:et1:e2:et2:eta:phi:sumet:hf:bin:ptjet:drjet");
422		ntjet = fs->make<TNtuple>("ntjet","","bin:pt:eta:ethit:f05:f1:f15:f2:f25:f3:em:emf");
423
424		}

Diff Legend

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