ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

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

Comparing UserCode/CmsHi/JetAnalysis/src/HiInclusiveJetAnalyzer.cc (file contents):
Revision 1.1 by yilmaz, Tue Sep 20 19:45:06 2011 UTC vs.
Revision 1.9 by yilmaz, Sun Apr 22 19:12:44 2012 UTC

#	Line 21 | Line 21
21		#include "DataFormats/HeavyIonEvent/interface/Centrality.h"
22
23		#include "DataFormats/CaloTowers/interface/CaloTowerCollection.h"
24	–	#include "DataFormats/PatCandidates/interface/Jet.h"
24		#include "DataFormats/JetReco/interface/CaloJetCollection.h"
25		#include "SimDataFormats/HiGenData/interface/GenHIEvent.h"
26		#include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h"
#	Line 46 | Line 45 | HiInclusiveJetAnalyzer::HiInclusiveJetAn
45
46
47		jetTag_ = iConfig.getParameter<InputTag>("jetTag");
48	<	vtxTag_ = iConfig.getUntrackedParameter<edm::InputTag>("vtxTag",edm::InputTag("hiSelectedVertex"));
48	>	vtxTag_ = iConfig.getUntrackedParameter<edm::InputTag>("vtxTag",edm::InputTag("hiSelectedVertex"));
49
50		isMC_ = iConfig.getUntrackedParameter<bool>("isMC",false);
51	<
51	>
52		if(isMC_){
53		genjetTag_ = iConfig.getParameter<InputTag>("genjetTag");
54		eventInfoTag_ = iConfig.getParameter<InputTag>("eventInfoTag");
#	Line 57 | Line 56 | HiInclusiveJetAnalyzer::HiInclusiveJetAn
56		verbose_ = iConfig.getUntrackedParameter<bool>("verbose",false);
57
58		useCentrality_ = iConfig.getUntrackedParameter<bool>("useCentrality",false);
59	<	useVtx_ = iConfig.getUntrackedParameter<bool>("useVtx",true);
59	>	useVtx_ = iConfig.getUntrackedParameter<bool>("useVtx",false);
60		useJEC_ = iConfig.getUntrackedParameter<bool>("useJEC",true);
61	+	usePat_ = iConfig.getUntrackedParameter<bool>("usePAT",true);
62	+
63	+	doLifeTimeTagging_ = iConfig.getUntrackedParameter<bool>("doLifeTimeTagging",false);
64	+	pfCandidateLabel_ = iConfig.getUntrackedParameter<edm::InputTag>("pfCandidateLabel",edm::InputTag("particleFlowTmp"));
65
66		if(!isMC_){
67		L1gtReadout_ = iConfig.getParameter<edm::InputTag>("L1gtReadout");
#	Line 77 | Line 80 | HiInclusiveJetAnalyzer::HiInclusiveJetAn
80		}
81
82		cout<<" jet collection : "<<jetTag_<<endl;
83	<	if(isMC_)cout<<" genjet collection : "<<genjetTag_<<endl;
84	<
83	>	if(isMC_){
84	>	cout<<" genjet collection : "<<genjetTag_<<endl;
85	>	genPtMin_ = iConfig.getUntrackedParameter<double>("genPtMin",0);
86	>	}
87
88
89		}
#	Line 102 | Line 107 | HiInclusiveJetAnalyzer::beginJob() {
107		string jetTagTitle = jetTag_.label()+" Jet Analysis Tree";
108		t = fs1->make<TTree>("t",jetTagTitle.c_str());
109
105	–
110		//  TTree* t= new TTree("t","Jet Response Analyzer");
111	<	t->Branch("run",&jets_.run,"run/I");
111	>	//t->Branch("run",&jets_.run,"run/I");
112		t->Branch("evt",&jets_.evt,"evt/I");
113	<	t->Branch("lumi",&jets_.lumi,"lumi/I");
113	>	//t->Branch("lumi",&jets_.lumi,"lumi/I");
114		t->Branch("b",&jets_.b,"b/F");
115	<	t->Branch("vx",&jets_.vx,"vx/F");
116	<	t->Branch("vy",&jets_.vy,"vy/F");
117	<	t->Branch("vz",&jets_.vz,"vz/F");
118	<	t->Branch("hf",&jets_.hf,"hf/F");
115	>	if (useVtx_) {
116	>	t->Branch("vx",&jets_.vx,"vx/F");
117	>	t->Branch("vy",&jets_.vy,"vy/F");
118	>	t->Branch("vz",&jets_.vz,"vz/F");
119	>	}
120	>	if (useCentrality_) {
121	>	t->Branch("hf",&jets_.hf,"hf/F");
122	>	t->Branch("bin",&jets_.bin,"bin/I");
123	>	}
124	>
125		t->Branch("nref",&jets_.nref,"nref/I");
116	–	t->Branch("bin",&jets_.bin,"bin/I");
126		t->Branch("rawpt",jets_.rawpt,"rawpt[nref]/F");
127		t->Branch("jtpt",jets_.jtpt,"jtpt[nref]/F");
128		t->Branch("jteta",jets_.jteta,"jteta[nref]/F");
129		t->Branch("jty",jets_.jty,"jty[nref]/F");
130		t->Branch("jtphi",jets_.jtphi,"jtphi[nref]/F");
131	+	t->Branch("jtpu",jets_.jtpu,"jtpu[nref]/F");
132
133	+	// b-jet discriminators
134	+	if (doLifeTimeTagging_) {
135	+
136	+	t->Branch("discr_csvMva",jets_.discr_csvMva,"discr_csvMva[nref]/F");
137	+	t->Branch("discr_csvSimple",jets_.discr_csvSimple,"discr_csvSimple[nref]/F");
138	+	t->Branch("discr_muByIp3",jets_.discr_muByIp3,"discr_muByIp3[nref]/F");
139	+	t->Branch("discr_muByPt",jets_.discr_muByPt,"discr_muByPt[nref]/F");
140	+	t->Branch("discr_prob",jets_.discr_prob,"discr_prob[nref]/F");
141	+	t->Branch("discr_probb",jets_.discr_probb,"discr_probb[nref]/F");
142	+	t->Branch("discr_tcHighEff",jets_.discr_tcHighEff,"discr_tcHighEff[nref]/F");
143	+	t->Branch("discr_tcHighPur",jets_.discr_tcHighPur,"discr_tcHighPur[nref]/F");
144	+
145	+	t->Branch("nsvtx",    jets_.nsvtx,    "nsvtx[nref]/b");
146	+	t->Branch("svtxntrk", jets_.svtxntrk, "svtxntrk[nref]/b");
147	+	t->Branch("svtxdl",   jets_.svtxdl,   "svtxdl[nref]/F");
148	+	t->Branch("svtxdls",  jets_.svtxdls,  "svtxdls[nref]/F");
149	+	t->Branch("svtxm",    jets_.svtxm,    "svtxm[nref]/F");
150	+	t->Branch("svtxpt",   jets_.svtxpt,   "svtxpt[nref]/F");
151	+
152	+	t->Branch("nIPtracks",jets_.nIPtracks,"nIPTracks[nref]/b");
153	+	t->Branch("nselIPtracks",jets_.nselIPtracks,"nselIPTracks[nref]/b");
154	+
155	+	t->Branch("mue",     jets_.mue,     "mue[nref]/F");
156	+	t->Branch("mupt",    jets_.mupt,    "mupt[nref]/F");
157	+	t->Branch("mueta",   jets_.mueta,   "mueta[nref]/F");
158	+	t->Branch("muphi",   jets_.muphi,   "muphi[nref]/F");
159	+	t->Branch("mudr",    jets_.mudr,    "mudr[nref]/F");
160	+	t->Branch("muptrel", jets_.muptrel, "muptrel[nref]/F");
161	+	t->Branch("muchg",   jets_.muchg,   "muchg[nref]/I");
162	+	}
163	+
164		if(isMC_){
165		t->Branch("pthat",&jets_.pthat,"pthat/F");
166
#	Line 132 | Line 173 | HiInclusiveJetAnalyzer::beginJob() {
173		t->Branch("refdrjt",jets_.refdrjt,"refdrjt[nref]/F");
174		// matched parton
175		t->Branch("refparton_pt",jets_.refparton_pt,"refparton_pt[nref]/F");
176	<	t->Branch("refparton_flavor",jets_.refparton_flavor,"refparton_flavor[nref]/F");
176	>	t->Branch("refparton_flavor",jets_.refparton_flavor,"refparton_flavor[nref]/I");
177	>	t->Branch("refparton_flavorForB",jets_.refparton_flavorForB,"refparton_flavorForB[nref]/I");
178
179		// For all gen jets, matched or unmatched
180		t->Branch("ngen",&jets_.ngen,"ngen/I");
#	Line 144 | Line 186 | HiInclusiveJetAnalyzer::beginJob() {
186		t->Branch("gendphijt",jets_.gendphijt,"gendphijt[ngen]/F");
187		t->Branch("gendrjt",jets_.gendrjt,"gendrjt[ngen]/F");
188		}
189	<
189	>	/*
190		if(!isMC_){
191		t->Branch("nL1TBit",&jets_.nL1TBit,"nL1TBit/I");
192		t->Branch("l1TBit",jets_.l1TBit,"l1TBit[nL1TBit]/O");
#	Line 156 | Line 198 | HiInclusiveJetAnalyzer::beginJob() {
198		t->Branch("hltBit",jets_.hltBit,"hltBit[nHLTBit]/O");
199
200		}
201	<
201	>	*/
202		TH1D::SetDefaultSumw2();
203
204
#	Line 177 | Line 219 | HiInclusiveJetAnalyzer::analyze(const Ev
219
220		LogDebug("HiInclusiveJetAnalyzer")<<"START event: "<<event<<" in run "<<run<<endl;
221
180	–
222		int bin = -1;
223		double hf = 0.;
224		double b = 999.;
225
226		if(useCentrality_){
186	–	//if(!isMC_){
227		if(!centrality_) centrality_ = new CentralityProvider(iSetup);
228		centrality_->newEvent(iEvent,iSetup); // make sure you do this first in every event
229		//double c = centrality_->centralityValue();
#	Line 193 | Line 233 | HiInclusiveJetAnalyzer::analyze(const Ev
233
234		bin = centrality_->getBin();
235		b = centrality_->bMean();
196	–	//}
197	–	/*
198	–	else{
199	–
200	–	TFile * centFile = new TFile("/net/hidsk0001/d00/scratch/mnguyen/CMSSW_3_9_1_patch1/src/macros/Hydjet_CentTable.root");
201	–
202	–	edm::Handle<reco::Centrality> cent;
203	–	iEvent.getByLabel(edm::InputTag("hiCentrality"),cent);
204	–	//cout<<" grabbed centrality "<<endl;
205	–	CentralityBins::RunMap cmap = getCentralityFromFile(centFile, "makeCentralityTableTFile", "HFhitsHydjet_2760GeV", 1, 1);
206	–
207	–	// Still getting cent from hits.  When tower tables become available, we need to switch
208	–	hf = cent->EtHFhitSum();
209	–	//cout<<" hf "<<hf<<endl;
210	–	bin = cmap[run]->getBin(hf);
211	–	b = cmap[run]->bMeanOfBin(bin);
212	–	}
213	–	*/
236		}
237
238
#	Line 236 | Line 258 | HiInclusiveJetAnalyzer::analyze(const Ev
258		}
259		}
260
261	<
262	<
263	<	edm::Handle<pat::JetCollection> jets;
261	>	edm::Handle<pat::JetCollection> patjets;
262	>	if(usePat_)iEvent.getByLabel(jetTag_, patjets);
263	>
264	>	edm::Handle<reco::JetView> jets;
265		iEvent.getByLabel(jetTag_, jets);
266
267	<
267	>	edm::Handle<reco::PFCandidateCollection> pfCandidates;
268	>	if(doLifeTimeTagging_){
269	>	iEvent.getByLabel(pfCandidateLabel_,pfCandidates);
270	>	}
271		// FILL JRA TREE
272
273		jets_.b = b;
#	Line 253 | Line 279 | HiInclusiveJetAnalyzer::analyze(const Ev
279		}
280
281		for(unsigned int j = 0 ; j < jets->size(); ++j){
282	<	const pat::Jet& jet = (*jets)[j];
282	>	const reco::Jet& jet = (*jets)[j];
283
284		//cout<<" jet pt "<<jet.pt()<<endl;
285		//if(jet.pt() < jetPtMin) continue;
286	<	if (useJEC_) jets_.rawpt[jets_.nref]=jet.correctedJet("Uncorrected").pt();
286	>	if (useJEC_ && usePat_){
287	>	jets_.rawpt[jets_.nref]=(*patjets)[j].correctedJet("Uncorrected").pt();
288	>	}
289	>
290	>	if(doLifeTimeTagging_){
291	>
292	>	if(jetTag_.label()=="icPu5patJets"){
293	>	jets_.discr_csvMva[jets_.nref]    = (*patjets)[j].bDiscriminator("icPu5CombinedSecondaryVertexMVABJetTags");
294	>	jets_.discr_csvSimple[jets_.nref] = (*patjets)[j].bDiscriminator("icPu5CombinedSecondaryVertexBJetTags");
295	>	jets_.discr_muByIp3[jets_.nref]   = (*patjets)[j].bDiscriminator("icPu5SoftMuonByIP3dBJetTags");
296	>	jets_.discr_muByPt[jets_.nref]    = (*patjets)[j].bDiscriminator("icPu5SoftMuonByPtBJetTags");
297	>	jets_.discr_prob[jets_.nref]      = (*patjets)[j].bDiscriminator("icPu5JetBProbabilityBJetTags");
298	>	jets_.discr_probb[jets_.nref]     = (*patjets)[j].bDiscriminator("icPu5JetProbabilityBJetTags");
299	>	jets_.discr_tcHighEff[jets_.nref] = (*patjets)[j].bDiscriminator("icPu5TrackCountingHighEffBJetTags");
300	>	jets_.discr_tcHighPur[jets_.nref] = (*patjets)[j].bDiscriminator("icPu5TrackCountingHighPurBJetTags");
301	>	}
302	>	else if(jetTag_.label()=="akPu3PFpatJets"){
303	>	jets_.discr_csvMva[jets_.nref]    = (*patjets)[j].bDiscriminator("akPu3PFCombinedSecondaryVertexMVABJetTags");
304	>	jets_.discr_csvSimple[jets_.nref] = (*patjets)[j].bDiscriminator("akPu3PFCombinedSecondaryVertexBJetTags");
305	>	jets_.discr_muByIp3[jets_.nref]   = (*patjets)[j].bDiscriminator("akPu3PFSoftMuonByIP3dBJetTags");
306	>	jets_.discr_muByPt[jets_.nref]    = (*patjets)[j].bDiscriminator("akPu3PFSoftMuonByPtBJetTags");
307	>	jets_.discr_prob[jets_.nref]      = (*patjets)[j].bDiscriminator("akPu3PFJetBProbabilityBJetTags");
308	>	jets_.discr_probb[jets_.nref]     = (*patjets)[j].bDiscriminator("akPu3PFJetProbabilityBJetTags");
309	>	jets_.discr_tcHighEff[jets_.nref] = (*patjets)[j].bDiscriminator("akPu3PFTrackCountingHighEffBJetTags");
310	>	jets_.discr_tcHighPur[jets_.nref] = (*patjets)[j].bDiscriminator("akPu3PFTrackCountingHighPurBJetTags");
311	>	}
312	>	else{
313	>	cout<<" you fail at btagging "<<endl;
314	>	}
315	>
316	>	const reco::SecondaryVertexTagInfo& tagInfoSV=*(*patjets)[j].tagInfoSecondaryVertex();
317	>	if (tagInfoSV.nVertices()>0) {
318	>	Measurement1D m1D = tagInfoSV.flightDistance(0);
319	>	jets_.nsvtx[jets_.nref]     = tagInfoSV.nVertices();
320	>	jets_.svtxntrk[jets_.nref]  = tagInfoSV.nVertexTracks(0);
321	>	jets_.svtxdl[jets_.nref]    = m1D.value();
322	>	jets_.svtxdls[jets_.nref]   = m1D.significance();
323	>
324	>	const reco::Vertex& svtx = tagInfoSV.secondaryVertex(0);
325	>
326	>	jets_.svtxm[jets_.nref]    = svtx.p4().mass();
327	>	jets_.svtxpt[jets_.nref]   = svtx.p4().pt();
328	>
329	>	}
330	>	else {
331	>	jets_.nsvtx[jets_.nref]    =   0;
332	>	jets_.svtxntrk[jets_.nref] =   0;
333	>	jets_.svtxdl[jets_.nref]   = 0.0;
334	>	jets_.svtxdls[jets_.nref]  = 0.0;
335	>	jets_.svtxm[jets_.nref]    = 0.0;
336	>	jets_.svtxpt[jets_.nref]   = 0.0;
337	>	}
338	>	const reco::TrackIPTagInfo& tagInfoIP=*(*patjets)[j].tagInfoTrackIP();
339	>
340	>	jets_.nIPtracks[jets_.nref] = tagInfoIP.tracks().size();
341	>	jets_.nselIPtracks[jets_.nref] = tagInfoIP.selectedTracks().size();
342	>
343	>	// would be nice to save the info for the tracks, but requires a more sophisticated tree
344	>	/*
345	>	//cout << "Jet pt: " << tagInfoIP.jet()->pt() << endl;
346	>	cout << "Tot tracks: " << tagInfoIP.tracks().size() << endl;
347	>	TrackRefVector selTracks=tagInfoIP.selectedTracks();
348	>	int n=selTracks.size();
349	>	cout << "Sel tracks: " << n << endl;
350	>	// false      cout << " Pt  \t d len \t jet dist \t p3d \t p2d\t ip3d \t ip2d " << endl;
351	>	GlobalPoint pv(tagInfoIP.primaryVertex()->position().x(),tagInfoIP.primaryVertex()->position().y(),tagInfoIP.primaryVertex()->position().z());
352	>	cout << pv << " vs " << tagInfoIP.primaryVertex()->position()   << endl;
353	>	for(int j=0;j< n;j++)
354	>	{
355	>	TrackIPTagInfo::TrackIPData data = tagInfoIP.impactParameterData()[j];
356	>	cout << selTracks[j]->pt() << "\t";
357	>	cout << tagInfoIP.probabilities(0)[j] << "\t";
358	>	cout << tagInfoIP.probabilities(1)[j] << "\t";
359	>	cout << data.ip3d.value() << "\t";
360	>	cout << data.ip3d.significance() << "\t";
361	>	cout << data.distanceToJetAxis.value() << "\t";
362	>	cout << data.distanceToJetAxis.significance() << "\t";
363	>	cout << data.distanceToGhostTrack.value() << "\t";
364	>	cout << data.distanceToGhostTrack.significance() << "\t";
365	>	cout << data.closestToJetAxis << "\t";
366	>	cout << (data.closestToJetAxis - pv).mag() << "\t";
367	>	cout << data.closestToGhostTrack << "\t";
368	>	cout << (data.closestToGhostTrack - pv).mag() << "\t";
369	>	cout << data.ip2d.value() << "\t";
370	>	cout << data.ip2d.significance() <<  endl;
371	>	}
372	>	*/
373	>
374	>	const reco::PFCandidateCollection *pfCandidateColl = &(*pfCandidates);
375	>	int pfMuonIndex = getPFJetMuon(jet, pfCandidateColl);
376	>	if(pfMuonIndex >=0){
377	>	const reco::PFCandidate muon = pfCandidateColl->at(pfMuonIndex);
378	>	jets_.mupt[jets_.nref]    =  muon.pt();
379	>	jets_.mueta[jets_.nref]   =  muon.eta();
380	>	jets_.muphi[jets_.nref]   =  muon.phi();
381	>	jets_.mue[jets_.nref]     =  muon.energy();
382	>	jets_.mudr[jets_.nref]    =  reco::deltaR(jet,muon);
383	>	jets_.muptrel[jets_.nref] =  getPtRel(muon, jet);
384	>	jets_.muchg[jets_.nref]   =  muon.charge();
385	>	}
386	>	else{
387	>	jets_.mupt[jets_.nref]    =  0.0;
388	>	jets_.mueta[jets_.nref]   =  0.0;
389	>	jets_.muphi[jets_.nref]   =  0.0;
390	>	jets_.mue[jets_.nref]     =  0.0;
391	>	jets_.mudr[jets_.nref]    =  9.9;
392	>	jets_.muptrel[jets_.nref] =  0.0;
393	>	jets_.muchg[jets_.nref]   = 0;
394	>	}
395	>	}
396	>
397		jets_.jtpt[jets_.nref] = jet.pt();
398		jets_.jteta[jets_.nref] = jet.eta();
399		jets_.jtphi[jets_.nref] = jet.phi();
400		jets_.jty[jets_.nref] = jet.eta();
401	<
401	>	jets_.jtpu[jets_.nref] = jet.pileup();
402
403	<	if(jet.genJet()){
404	<	jets_.refpt[jets_.nref] = jet.genJet()->pt();
405	<	jets_.refeta[jets_.nref] = jet.genJet()->eta();
406	<	jets_.refphi[jets_.nref] = jet.genJet()->phi();
407	<	jets_.refy[jets_.nref] = jet.genJet()->eta();
408	<	jets_.refdphijt[jets_.nref] = reco::deltaPhi(jet.phi(),jet.genJet()->phi());
409	<	jets_.refdrjt[jets_.nref] = reco::deltaR(jet.eta(),jet.phi(),jet.genJet()->eta(),jet.genJet()->phi());
410	<	}
411	<	else{
412	<	jets_.refpt[jets_.nref] = -999.;
413	<	jets_.refeta[jets_.nref] = -999.;
414	<	jets_.refphi[jets_.nref] = -999.;
415	<	jets_.refy[jets_.nref] = -999.;
416	<	jets_.refdphijt[jets_.nref] = -999.;
417	<	jets_.refdrjt[jets_.nref] = -999.;
418	<	}
403	>	if(isMC_ && usePat_){
404	>	jets_.refparton_flavorForB[jets_.nref] = (*patjets)[j].partonFlavour();
405	>	const reco::GenJet* genjet = (*patjets)[j].genJet();
406	>
407	>	if(genjet){
408	>	jets_.refpt[jets_.nref] = genjet->pt();
409	>	jets_.refeta[jets_.nref] = genjet->eta();
410	>	jets_.refphi[jets_.nref] = genjet->phi();
411	>	jets_.refy[jets_.nref] = genjet->eta();
412	>	jets_.refdphijt[jets_.nref] = reco::deltaPhi(jet.phi(), genjet->phi());
413	>	jets_.refdrjt[jets_.nref] = reco::deltaR(jet.eta(),jet.phi(),genjet->eta(),genjet->phi());
414	>	}else{
415	>	jets_.refpt[jets_.nref] = -999.;
416	>	jets_.refeta[jets_.nref] = -999.;
417	>	jets_.refphi[jets_.nref] = -999.;
418	>	jets_.refy[jets_.nref] = -999.;
419	>	jets_.refdphijt[jets_.nref] = -999.;
420	>	jets_.refdrjt[jets_.nref] = -999.;
421	>	}
422
423	<	// matched partons
424	<	if (jet.genParton()) {
425	<	jets_.refparton_pt[jets_.nref] = jet.genParton()->pt();
426	<	jets_.refparton_flavor[jets_.nref] = jet.genParton()->pdgId();
427	<	} else {
423	>	// matched partons
424	>	const reco::GenParticle * parton = (*patjets)[j].genParton();
425	>	if(parton){
426	>	jets_.refparton_pt[jets_.nref] = parton->pt();
427	>	jets_.refparton_flavor[jets_.nref] = parton->pdgId();
428	>	} else {
429		jets_.refparton_pt[jets_.nref] = -999;
430		jets_.refparton_flavor[jets_.nref] = -999;
431	+	}
432		}
292	–
433
434		jets_.nref++;
435
#	Line 316 | Line 456 | HiInclusiveJetAnalyzer::analyze(const Ev
456		float genjet_pt = genjet.pt();
457
458		// threshold to reduce size of output in minbias PbPb
459	<	if(genjet_pt>20.){
459	>	if(genjet_pt>genPtMin_){
460
461		jets_.genpt[jets_.ngen] = genjet_pt;
462		jets_.geneta[jets_.ngen] = genjet.eta();
#	Line 329 | Line 469 | HiInclusiveJetAnalyzer::analyze(const Ev
469		jets_.gendrjt[jets_.ngen] = -1.0;
470		jets_.genmatchindex[jets_.ngen] = -1;
471
332	–
472		for(unsigned int ijet = 0 ; ijet < jets->size(); ++ijet){
473		const pat::Jet& jet = (*jets)[ijet];
474	<
475	<	if(jet.genJet()){
476	<	if(fabs(genjet.pt()-jet.genJet()->pt())<0.0001 &&
477	<	fabs(genjet.eta()-jet.genJet()->eta())<0.0001 &&
478	<	(fabs(genjet.phi()-jet.genJet()->phi())<0.0001 || fabs(genjet.phi()-jet.genJet()->phi() - 2.0*TMath::Pi()) < 0.0001 )){
474	>	const reco::GenJet* matchedGenJet = (*patjets)[ijet].genJet();
475	>	if(matchedGenJet){
476	>	// poor man's matching, someone fix please
477	>	if(fabs(genjet.pt()-matchedGenJet->pt())<0.0001 &&
478	>	fabs(genjet.eta()-matchedGenJet->eta())<0.0001 &&
479	>	(fabs(genjet.phi()-matchedGenJet->phi())<0.0001 || fabs(genjet.phi()-matchedGenJet->phi() - 2.0*TMath::Pi()) < 0.0001 )){
480
481		jets_.genmatchindex[jets_.ngen] = (int)ijet;
482		jets_.gendphijt[jets_.ngen] = reco::deltaPhi(jet.phi(),genjet.phi());
#	Line 452 | Line 592 | inline bool HiInclusiveJetAnalyzer::getP
592		}
593
594
595	+	int
596	+	HiInclusiveJetAnalyzer::getPFJetMuon(const pat::Jet& pfJet, const reco::PFCandidateCollection *pfCandidateColl)
597	+	{
598	+
599	+	int pfMuonIndex = -1;
600	+	float ptMax = 0.;
601
602	+
603	+	for(unsigned icand=0;icand<pfCandidateColl->size(); icand++) {
604	+	const reco::PFCandidate pfCandidate = pfCandidateColl->at(icand);
605	+
606	+	int id = pfCandidate.particleId();
607	+	if(abs(id) != 3) continue;
608	+
609	+	if(reco::deltaR(pfJet,pfCandidate)>0.5) continue;
610	+
611	+	double pt =  pfCandidate.pt();
612	+	if(pt>ptMax){
613	+	ptMax = pt;
614	+	pfMuonIndex = (int) icand;
615	+	}
616	+	}
617	+
618	+	return pfMuonIndex;
619	+
620	+	}
621	+
622	+
623	+	double
624	+	HiInclusiveJetAnalyzer::getPtRel(const reco::PFCandidate lep, const pat::Jet& jet )
625	+	{
626	+
627	+	float lj_x = jet.p4().px();
628	+	float lj_y = jet.p4().py();
629	+	float lj_z = jet.p4().pz();
630	+
631	+	// absolute values squared
632	+	float lj2  = lj_x*lj_x+lj_y*lj_y+lj_z*lj_z;
633	+	float lep2 = lep.px()*lep.px()+lep.py()*lep.py()+lep.pz()*lep.pz();
634	+
635	+	// projection vec(mu) to lepjet axis
636	+	float lepXlj = lep.px()*lj_x+lep.py()*lj_y+lep.pz()*lj_z;
637	+
638	+	// absolute value squared and normalized
639	+	float pLrel2 = lepXlj*lepXlj/lj2;
640	+
641	+	// lep2 = pTrel2 + pLrel2
642	+	float pTrel2 = lep2-pLrel2;
643	+
644	+	return (pTrel2 > 0) ? std::sqrt(pTrel2) : 0.0;
645	+	}
646
647		DEFINE_FWK_MODULE(HiInclusiveJetAnalyzer);

Diff Legend

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