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.2 by yilmaz, Tue Sep 20 20:06:06 2011 UTC vs.
Revision 1.10 by mnguyen, Sat Apr 28 12:15:17 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"
27	+	#include "SimDataFormats/GeneratorProducts/interface/HepMCProduct.h"
28		#include "DataFormats/JetReco/interface/GenJetCollection.h"
29		#include "DataFormats/VertexReco/interface/Vertex.h"
30		#include "DataFormats/Math/interface/deltaPhi.h"
#	Line 46 | Line 46 | HiInclusiveJetAnalyzer::HiInclusiveJetAn
46
47
48		jetTag_ = iConfig.getParameter<InputTag>("jetTag");
49	<	vtxTag_ = iConfig.getUntrackedParameter<edm::InputTag>("vtxTag",edm::InputTag("hiSelectedVertex"));
49	>	vtxTag_ = iConfig.getUntrackedParameter<edm::InputTag>("vtxTag",edm::InputTag("hiSelectedVertex"));
50
51		isMC_ = iConfig.getUntrackedParameter<bool>("isMC",false);
52	<
52	>
53		if(isMC_){
54		genjetTag_ = iConfig.getParameter<InputTag>("genjetTag");
55		eventInfoTag_ = iConfig.getParameter<InputTag>("eventInfoTag");
#	Line 57 | Line 57 | HiInclusiveJetAnalyzer::HiInclusiveJetAn
57		verbose_ = iConfig.getUntrackedParameter<bool>("verbose",false);
58
59		useCentrality_ = iConfig.getUntrackedParameter<bool>("useCentrality",false);
60	<	useVtx_ = iConfig.getUntrackedParameter<bool>("useVtx",true);
60	>	useVtx_ = iConfig.getUntrackedParameter<bool>("useVtx",false);
61		useJEC_ = iConfig.getUntrackedParameter<bool>("useJEC",true);
62		usePat_ = iConfig.getUntrackedParameter<bool>("usePAT",true);
63
64	+	doLifeTimeTagging_ = iConfig.getUntrackedParameter<bool>("doLifeTimeTagging",false);
65	+	doLifeTimeTaggingExtras_ = iConfig.getUntrackedParameter<bool>("doLifeTimeTaggingExtras",true);
66	+	pfCandidateLabel_ = iConfig.getUntrackedParameter<edm::InputTag>("pfCandidateLabel",edm::InputTag("particleFlowTmp"));
67	+
68		if(!isMC_){
69		L1gtReadout_ = iConfig.getParameter<edm::InputTag>("L1gtReadout");
70		hltResName_ = iConfig.getUntrackedParameter<string>("hltTrgResults","TriggerResults::HLT");
#	Line 78 | Line 82 | HiInclusiveJetAnalyzer::HiInclusiveJetAn
82		}
83
84		cout<<" jet collection : "<<jetTag_<<endl;
85	<	if(isMC_)cout<<" genjet collection : "<<genjetTag_<<endl;
86	<
85	>	if(isMC_){
86	>	cout<<" genjet collection : "<<genjetTag_<<endl;
87	>	genPtMin_ = iConfig.getUntrackedParameter<double>("genPtMin",0);
88	>	}
89
90
91		}
#	Line 104 | Line 110 | HiInclusiveJetAnalyzer::beginJob() {
110		t = fs1->make<TTree>("t",jetTagTitle.c_str());
111
112		//  TTree* t= new TTree("t","Jet Response Analyzer");
113	<	t->Branch("run",&jets_.run,"run/I");
113	>	//t->Branch("run",&jets_.run,"run/I");
114		t->Branch("evt",&jets_.evt,"evt/I");
115	<	t->Branch("lumi",&jets_.lumi,"lumi/I");
115	>	//t->Branch("lumi",&jets_.lumi,"lumi/I");
116		t->Branch("b",&jets_.b,"b/F");
117	<	t->Branch("vx",&jets_.vx,"vx/F");
118	<	t->Branch("vy",&jets_.vy,"vy/F");
119	<	t->Branch("vz",&jets_.vz,"vz/F");
120	<	t->Branch("hf",&jets_.hf,"hf/F");
117	>	if (useVtx_) {
118	>	t->Branch("vx",&jets_.vx,"vx/F");
119	>	t->Branch("vy",&jets_.vy,"vy/F");
120	>	t->Branch("vz",&jets_.vz,"vz/F");
121	>	}
122	>	if (useCentrality_) {
123	>	t->Branch("hf",&jets_.hf,"hf/F");
124	>	t->Branch("bin",&jets_.bin,"bin/I");
125	>	}
126	>
127		t->Branch("nref",&jets_.nref,"nref/I");
116	–	t->Branch("bin",&jets_.bin,"bin/I");
128		t->Branch("rawpt",jets_.rawpt,"rawpt[nref]/F");
129		t->Branch("jtpt",jets_.jtpt,"jtpt[nref]/F");
130		t->Branch("jteta",jets_.jteta,"jteta[nref]/F");
#	Line 121 | Line 132 | HiInclusiveJetAnalyzer::beginJob() {
132		t->Branch("jtphi",jets_.jtphi,"jtphi[nref]/F");
133		t->Branch("jtpu",jets_.jtpu,"jtpu[nref]/F");
134
135	+	// b-jet discriminators
136	+	if (doLifeTimeTagging_) {
137	+
138	+	t->Branch("discr_csvMva",jets_.discr_csvMva,"discr_csvMva[nref]/F");
139	+	t->Branch("discr_csvSimple",jets_.discr_csvSimple,"discr_csvSimple[nref]/F");
140	+	t->Branch("discr_muByIp3",jets_.discr_muByIp3,"discr_muByIp3[nref]/F");
141	+	t->Branch("discr_muByPt",jets_.discr_muByPt,"discr_muByPt[nref]/F");
142	+	t->Branch("discr_prob",jets_.discr_prob,"discr_prob[nref]/F");
143	+	t->Branch("discr_probb",jets_.discr_probb,"discr_probb[nref]/F");
144	+	t->Branch("discr_tcHighEff",jets_.discr_tcHighEff,"discr_tcHighEff[nref]/F");
145	+	t->Branch("discr_tcHighPur",jets_.discr_tcHighPur,"discr_tcHighPur[nref]/F");
146	+
147	+	t->Branch("nsvtx",    jets_.nsvtx,    "nsvtx[nref]/I");
148	+	t->Branch("svtxntrk", jets_.svtxntrk, "svtxntrk[nref]/I");
149	+	t->Branch("svtxdl",   jets_.svtxdl,   "svtxdl[nref]/F");
150	+	t->Branch("svtxdls",  jets_.svtxdls,  "svtxdls[nref]/F");
151	+	t->Branch("svtxm",    jets_.svtxm,    "svtxm[nref]/F");
152	+	t->Branch("svtxpt",   jets_.svtxpt,   "svtxpt[nref]/F");
153	+
154	+	t->Branch("nIPtrk",jets_.nIPtrk,"nIPtrk[nref]/I");
155	+	t->Branch("nselIPtrk",jets_.nselIPtrk,"nselIPtrk[nref]/I");
156	+
157	+	if (doLifeTimeTaggingExtras_) {
158	+	t->Branch("ipJetIndex",jets_.ipJetIndex,"ipJetIndex[nselIPtrk]/I");
159	+	t->Branch("ipPt",jets_.ipPt,"ipPt[nselIPtrk]/F");
160	+	t->Branch("ipProb0",jets_.ipProb0,"ipProb0[nselIPtrk]/F");
161	+	t->Branch("ipProb1",jets_.ipProb1,"ipProb1[nselIPtrk]/F");
162	+	t->Branch("ip2d",jets_.ip2d,"ip2d[nselIPtrk]/F");
163	+	t->Branch("ip2dSig",jets_.ip2dSig,"ip2dSig[nselIPtrk]/F");
164	+	t->Branch("ip3d",jets_.ip3d,"ip3d[nselIPtrk]/F");
165	+	t->Branch("ip3dSig",jets_.ip3dSig,"ip3dSig[nselIPtrk]/F");
166	+	t->Branch("ipDist2Jet",jets_.ipDist2Jet,"ipDist2Jet[nselIPtrk]/F");
167	+	t->Branch("ipDist2JetSig",jets_.ipDist2JetSig,"ipDist2JetSig[nselIPtrk]/F");
168	+	t->Branch("ipClosest2Jet",jets_.ipClosest2Jet,"ipClosest2Jet[nselIPtrk]/F");
169	+	}
170	+
171	+	t->Branch("mue",     jets_.mue,     "mue[nref]/F");
172	+	t->Branch("mupt",    jets_.mupt,    "mupt[nref]/F");
173	+	t->Branch("mueta",   jets_.mueta,   "mueta[nref]/F");
174	+	t->Branch("muphi",   jets_.muphi,   "muphi[nref]/F");
175	+	t->Branch("mudr",    jets_.mudr,    "mudr[nref]/F");
176	+	t->Branch("muptrel", jets_.muptrel, "muptrel[nref]/F");
177	+	t->Branch("muchg",   jets_.muchg,   "muchg[nref]/I");
178	+	}
179	+
180		if(isMC_){
181	+	t->Branch("beamId1",&jets_.beamId1,"beamId1/I");
182	+	t->Branch("beamId2",&jets_.beamId2,"beamId2/I");
183	+
184		t->Branch("pthat",&jets_.pthat,"pthat/F");
185
186		// Only matched gen jets
#	Line 133 | Line 192 | HiInclusiveJetAnalyzer::beginJob() {
192		t->Branch("refdrjt",jets_.refdrjt,"refdrjt[nref]/F");
193		// matched parton
194		t->Branch("refparton_pt",jets_.refparton_pt,"refparton_pt[nref]/F");
195	<	t->Branch("refparton_flavor",jets_.refparton_flavor,"refparton_flavor[nref]/F");
195	>	t->Branch("refparton_flavor",jets_.refparton_flavor,"refparton_flavor[nref]/I");
196	>	t->Branch("refparton_flavorForB",jets_.refparton_flavorForB,"refparton_flavorForB[nref]/I");
197
198		// For all gen jets, matched or unmatched
199		t->Branch("ngen",&jets_.ngen,"ngen/I");
#	Line 145 | Line 205 | HiInclusiveJetAnalyzer::beginJob() {
205		t->Branch("gendphijt",jets_.gendphijt,"gendphijt[ngen]/F");
206		t->Branch("gendrjt",jets_.gendrjt,"gendrjt[ngen]/F");
207		}
208	<
208	>	/*
209		if(!isMC_){
210		t->Branch("nL1TBit",&jets_.nL1TBit,"nL1TBit/I");
211		t->Branch("l1TBit",jets_.l1TBit,"l1TBit[nL1TBit]/O");
#	Line 157 | Line 217 | HiInclusiveJetAnalyzer::beginJob() {
217		t->Branch("hltBit",jets_.hltBit,"hltBit[nHLTBit]/O");
218
219		}
220	<
220	>	*/
221		TH1D::SetDefaultSumw2();
222
223
#	Line 178 | Line 238 | HiInclusiveJetAnalyzer::analyze(const Ev
238
239		LogDebug("HiInclusiveJetAnalyzer")<<"START event: "<<event<<" in run "<<run<<endl;
240
181	–
241		int bin = -1;
242		double hf = 0.;
243		double b = 999.;
244
245		if(useCentrality_){
187	–	//if(!isMC_){
246		if(!centrality_) centrality_ = new CentralityProvider(iSetup);
247		centrality_->newEvent(iEvent,iSetup); // make sure you do this first in every event
248		//double c = centrality_->centralityValue();
#	Line 194 | Line 252 | HiInclusiveJetAnalyzer::analyze(const Ev
252
253		bin = centrality_->getBin();
254		b = centrality_->bMean();
197	–	//}
198	–	/*
199	–	else{
200	–
201	–	TFile * centFile = new TFile("/net/hidsk0001/d00/scratch/mnguyen/CMSSW_3_9_1_patch1/src/macros/Hydjet_CentTable.root");
202	–
203	–	edm::Handle<reco::Centrality> cent;
204	–	iEvent.getByLabel(edm::InputTag("hiCentrality"),cent);
205	–	//cout<<" grabbed centrality "<<endl;
206	–	CentralityBins::RunMap cmap = getCentralityFromFile(centFile, "makeCentralityTableTFile", "HFhitsHydjet_2760GeV", 1, 1);
207	–
208	–	// Still getting cent from hits.  When tower tables become available, we need to switch
209	–	hf = cent->EtHFhitSum();
210	–	//cout<<" hf "<<hf<<endl;
211	–	bin = cmap[run]->getBin(hf);
212	–	b = cmap[run]->bMeanOfBin(bin);
213	–	}
214	–	*/
255		}
256
257
#	Line 243 | Line 283 | HiInclusiveJetAnalyzer::analyze(const Ev
283		edm::Handle<reco::JetView> jets;
284		iEvent.getByLabel(jetTag_, jets);
285
286	+	edm::Handle<reco::PFCandidateCollection> pfCandidates;
287	+	if(doLifeTimeTagging_){
288	+	iEvent.getByLabel(pfCandidateLabel_,pfCandidates);
289	+	}
290		// FILL JRA TREE
291
292		jets_.b = b;
#	Line 261 | Line 305 | HiInclusiveJetAnalyzer::analyze(const Ev
305		if (useJEC_ && usePat_){
306		jets_.rawpt[jets_.nref]=(*patjets)[j].correctedJet("Uncorrected").pt();
307		}
308	+
309	+	if(doLifeTimeTagging_){
310	+
311	+	if(jetTag_.label()=="icPu5patJets"){
312	+	jets_.discr_csvMva[jets_.nref]    = (*patjets)[j].bDiscriminator("icPu5CombinedSecondaryVertexMVABJetTags");
313	+	jets_.discr_csvSimple[jets_.nref] = (*patjets)[j].bDiscriminator("icPu5CombinedSecondaryVertexBJetTags");
314	+	jets_.discr_muByIp3[jets_.nref]   = (*patjets)[j].bDiscriminator("icPu5SoftMuonByIP3dBJetTags");
315	+	jets_.discr_muByPt[jets_.nref]    = (*patjets)[j].bDiscriminator("icPu5SoftMuonByPtBJetTags");
316	+	jets_.discr_prob[jets_.nref]      = (*patjets)[j].bDiscriminator("icPu5JetBProbabilityBJetTags");
317	+	jets_.discr_probb[jets_.nref]     = (*patjets)[j].bDiscriminator("icPu5JetProbabilityBJetTags");
318	+	jets_.discr_tcHighEff[jets_.nref] = (*patjets)[j].bDiscriminator("icPu5TrackCountingHighEffBJetTags");
319	+	jets_.discr_tcHighPur[jets_.nref] = (*patjets)[j].bDiscriminator("icPu5TrackCountingHighPurBJetTags");
320	+	}
321	+	else if(jetTag_.label()=="akPu3PFpatJets"){
322	+	jets_.discr_csvMva[jets_.nref]    = (*patjets)[j].bDiscriminator("akPu3PFCombinedSecondaryVertexMVABJetTags");
323	+	jets_.discr_csvSimple[jets_.nref] = (*patjets)[j].bDiscriminator("akPu3PFCombinedSecondaryVertexBJetTags");
324	+	jets_.discr_muByIp3[jets_.nref]   = (*patjets)[j].bDiscriminator("akPu3PFSoftMuonByIP3dBJetTags");
325	+	jets_.discr_muByPt[jets_.nref]    = (*patjets)[j].bDiscriminator("akPu3PFSoftMuonByPtBJetTags");
326	+	jets_.discr_prob[jets_.nref]      = (*patjets)[j].bDiscriminator("akPu3PFJetBProbabilityBJetTags");
327	+	jets_.discr_probb[jets_.nref]     = (*patjets)[j].bDiscriminator("akPu3PFJetProbabilityBJetTags");
328	+	jets_.discr_tcHighEff[jets_.nref] = (*patjets)[j].bDiscriminator("akPu3PFTrackCountingHighEffBJetTags");
329	+	jets_.discr_tcHighPur[jets_.nref] = (*patjets)[j].bDiscriminator("akPu3PFTrackCountingHighPurBJetTags");
330	+	}
331	+	else{
332	+	cout<<" b-tagging variables not filled for this collection, turn of doLifeTimeTagging "<<endl;
333	+	}
334	+
335	+	const reco::SecondaryVertexTagInfo& tagInfoSV=*(*patjets)[j].tagInfoSecondaryVertex();
336	+
337	+	jets_.nsvtx[jets_.nref]     = tagInfoSV.nVertices();
338	+
339	+	if (tagInfoSV.nVertices()>0) {
340	+	jets_.svtxntrk[jets_.nref]  = tagInfoSV.nVertexTracks(0);
341	+	// this is the 3d flight distance, for 2-D use (0,true)
342	+	Measurement1D m1D = tagInfoSV.flightDistance(0);
343	+	jets_.svtxdl[jets_.nref]    = m1D.value();
344	+	jets_.svtxdls[jets_.nref]   = m1D.significance();
345	+
346	+	const reco::Vertex& svtx = tagInfoSV.secondaryVertex(0);
347	+	jets_.svtxm[jets_.nref]    = svtx.p4().mass();
348	+	jets_.svtxpt[jets_.nref]   = svtx.p4().pt();
349	+	}
350	+
351	+	const reco::TrackIPTagInfo& tagInfoIP=*(*patjets)[j].tagInfoTrackIP();
352	+
353	+	jets_.nIPtrk[jets_.nref] = tagInfoIP.tracks().size();
354	+	jets_.nselIPtrk[jets_.nref] = tagInfoIP.selectedTracks().size();
355	+
356	+	if (doLifeTimeTaggingExtras_) {
357	+
358	+	TrackRefVector selTracks=tagInfoIP.selectedTracks();
359	+
360	+	GlobalPoint pv(tagInfoIP.primaryVertex()->position().x(),tagInfoIP.primaryVertex()->position().y(),tagInfoIP.primaryVertex()->position().z());
361	+
362	+	for(int it=0;it<jets_.nselIPtrk[jets_.nref] ;it++)
363	+	{
364	+	jets_.ipJetIndex[it]= jets_.nref;
365	+	TrackIPTagInfo::TrackIPData data = tagInfoIP.impactParameterData()[it];
366	+	jets_.ipPt[it] = selTracks[it]->pt();
367	+	jets_.ipProb0[it] = tagInfoIP.probabilities(0)[it];
368	+	jets_.ipProb1[it] = tagInfoIP.probabilities(1)[it];
369	+	jets_.ip2d[it] = data.ip2d.value();
370	+	jets_.ip2dSig[it] = data.ip2d.significance();
371	+	jets_.ip3d[it] = data.ip3d.value();
372	+	jets_.ip3dSig[it] = data.ip3d.significance();
373	+	jets_.ipDist2Jet[it] = data.distanceToJetAxis.value();
374	+	jets_.ipDist2JetSig[it] = data.distanceToJetAxis.significance();
375	+	jets_.ipClosest2Jet[it] = (data.closestToJetAxis - pv).mag();      //decay length
376	+	}
377	+	}
378	+
379	+	const reco::PFCandidateCollection *pfCandidateColl = &(*pfCandidates);
380	+	int pfMuonIndex = getPFJetMuon(jet, pfCandidateColl);
381	+	if(pfMuonIndex >=0){
382	+	const reco::PFCandidate muon = pfCandidateColl->at(pfMuonIndex);
383	+	jets_.mupt[jets_.nref]    =  muon.pt();
384	+	jets_.mueta[jets_.nref]   =  muon.eta();
385	+	jets_.muphi[jets_.nref]   =  muon.phi();
386	+	jets_.mue[jets_.nref]     =  muon.energy();
387	+	jets_.mudr[jets_.nref]    =  reco::deltaR(jet,muon);
388	+	jets_.muptrel[jets_.nref] =  getPtRel(muon, jet);
389	+	jets_.muchg[jets_.nref]   =  muon.charge();
390	+	}
391	+	}
392	+
393		jets_.jtpt[jets_.nref] = jet.pt();
394		jets_.jteta[jets_.nref] = jet.eta();
395		jets_.jtphi[jets_.nref] = jet.phi();
#	Line 268 | Line 397 | HiInclusiveJetAnalyzer::analyze(const Ev
397		jets_.jtpu[jets_.nref] = jet.pileup();
398
399		if(isMC_ && usePat_){
400	+	jets_.refparton_flavorForB[jets_.nref] = (*patjets)[j].partonFlavour();
401		const reco::GenJet* genjet = (*patjets)[j].genJet();
402	+
403		if(genjet){
404		jets_.refpt[jets_.nref] = genjet->pt();
405		jets_.refeta[jets_.nref] = genjet->eta();
#	Line 288 | Line 419 | HiInclusiveJetAnalyzer::analyze(const Ev
419		// matched partons
420		const reco::GenParticle * parton = (*patjets)[j].genParton();
421		if(parton){
422	<	jets_.refparton_pt[jets_.nref] = parton->pt();
423	<	jets_.refparton_flavor[jets_.nref] = parton->pdgId();
424	<	} else {
422	>	jets_.refparton_pt[jets_.nref] = parton->pt();
423	>	jets_.refparton_flavor[jets_.nref] = parton->pdgId();
424	>	} else {
425		jets_.refparton_pt[jets_.nref] = -999;
426		jets_.refparton_flavor[jets_.nref] = -999;
427	<	}
427	>	}
428		}
429
430		jets_.nref++;
#	Line 304 | Line 435 | HiInclusiveJetAnalyzer::analyze(const Ev
435
436		if(isMC_){
437
438	+	edm::Handle<HepMCProduct> hepMCProduct;
439	+	iEvent.getByLabel(eventInfoTag_,hepMCProduct);
440	+	const HepMC::GenEvent* MCEvt = hepMCProduct->GetEvent();
441	+	std::pair<HepMC::GenParticle*,HepMC::GenParticle*> beamParticles = MCEvt->beam_particles();
442	+	jets_.beamId1 = beamParticles.first->pdg_id();
443	+	jets_.beamId2 = beamParticles.second->pdg_id();
444	+
445		edm::Handle<GenEventInfoProduct> hEventInfo;
446		iEvent.getByLabel(eventInfoTag_,hEventInfo);
447		//jets_.pthat = hEventInfo->binningValues()[0];
#	Line 321 | Line 459 | HiInclusiveJetAnalyzer::analyze(const Ev
459		float genjet_pt = genjet.pt();
460
461		// threshold to reduce size of output in minbias PbPb
462	<	if(genjet_pt>20.){
462	>	if(genjet_pt>genPtMin_){
463
464		jets_.genpt[jets_.ngen] = genjet_pt;
465		jets_.geneta[jets_.ngen] = genjet.eta();
#	Line 334 | Line 472 | HiInclusiveJetAnalyzer::analyze(const Ev
472		jets_.gendrjt[jets_.ngen] = -1.0;
473		jets_.genmatchindex[jets_.ngen] = -1;
474
337	–
475		for(unsigned int ijet = 0 ; ijet < jets->size(); ++ijet){
476		const pat::Jet& jet = (*jets)[ijet];
477	<
478	<	if(jet.genJet()){
479	<	if(fabs(genjet.pt()-jet.genJet()->pt())<0.0001 &&
480	<	fabs(genjet.eta()-jet.genJet()->eta())<0.0001 &&
481	<	(fabs(genjet.phi()-jet.genJet()->phi())<0.0001 || fabs(genjet.phi()-jet.genJet()->phi() - 2.0*TMath::Pi()) < 0.0001 )){
477	>	const reco::GenJet* matchedGenJet = (*patjets)[ijet].genJet();
478	>	if(matchedGenJet){
479	>	// poor man's matching, someone fix please
480	>	if(fabs(genjet.pt()-matchedGenJet->pt())<0.0001 &&
481	>	fabs(genjet.eta()-matchedGenJet->eta())<0.0001 &&
482	>	(fabs(genjet.phi()-matchedGenJet->phi())<0.0001 || fabs(genjet.phi()-matchedGenJet->phi() - 2.0*TMath::Pi()) < 0.0001 )){
483
484		jets_.genmatchindex[jets_.ngen] = (int)ijet;
485		jets_.gendphijt[jets_.ngen] = reco::deltaPhi(jet.phi(),genjet.phi());
#	Line 358 | Line 496 | HiInclusiveJetAnalyzer::analyze(const Ev
496
497		}
498		t->Fill();
499	+	memset(&jets_,0,sizeof jets_);
500		}
501
502
#	Line 457 | Line 596 | inline bool HiInclusiveJetAnalyzer::getP
596		}
597
598
599	+	int
600	+	HiInclusiveJetAnalyzer::getPFJetMuon(const pat::Jet& pfJet, const reco::PFCandidateCollection *pfCandidateColl)
601	+	{
602	+
603	+	int pfMuonIndex = -1;
604	+	float ptMax = 0.;
605	+
606
607	+	for(unsigned icand=0;icand<pfCandidateColl->size(); icand++) {
608	+	const reco::PFCandidate pfCandidate = pfCandidateColl->at(icand);
609	+
610	+	int id = pfCandidate.particleId();
611	+	if(abs(id) != 3) continue;
612	+
613	+	if(reco::deltaR(pfJet,pfCandidate)>0.5) continue;
614	+
615	+	double pt =  pfCandidate.pt();
616	+	if(pt>ptMax){
617	+	ptMax = pt;
618	+	pfMuonIndex = (int) icand;
619	+	}
620	+	}
621	+
622	+	return pfMuonIndex;
623	+
624	+	}
625	+
626	+
627	+	double
628	+	HiInclusiveJetAnalyzer::getPtRel(const reco::PFCandidate lep, const pat::Jet& jet )
629	+	{
630	+
631	+	float lj_x = jet.p4().px();
632	+	float lj_y = jet.p4().py();
633	+	float lj_z = jet.p4().pz();
634	+
635	+	// absolute values squared
636	+	float lj2  = lj_x*lj_x+lj_y*lj_y+lj_z*lj_z;
637	+	float lep2 = lep.px()*lep.px()+lep.py()*lep.py()+lep.pz()*lep.pz();
638	+
639	+	// projection vec(mu) to lepjet axis
640	+	float lepXlj = lep.px()*lj_x+lep.py()*lj_y+lep.pz()*lj_z;
641	+
642	+	// absolute value squared and normalized
643	+	float pLrel2 = lepXlj*lepXlj/lj2;
644	+
645	+	// lep2 = pTrel2 + pLrel2
646	+	float pTrel2 = lep2-pLrel2;
647	+
648	+	return (pTrel2 > 0) ? std::sqrt(pTrel2) : 0.0;
649	+	}
650
651		DEFINE_FWK_MODULE(HiInclusiveJetAnalyzer);

Diff Legend

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