ViewVC Help

View File | Revision Log | Show Annotations | Root Listing

root/cvsroot/UserCode/Morgan/src/JetAnalyzer.cc

Comparing UserCode/Morgan/src/JetAnalyzer.cc (file contents):
Revision 1.3 by lethuill, Mon Dec 15 19:08:52 2008 UTC vs.
Revision 1.13 by lethuill, Fri Nov 6 08:03:06 2009 UTC

#	Line 4 | Line 4 | using namespace std;
4		using namespace reco;
5		using namespace edm;
6
7	–	JetAnalyzer::JetAnalyzer(const edm::ParameterSet& producersNames):verbosity_(0),useMC_(false)
8	–	{
9	–	dataType_ = producersNames.getUntrackedParameter<string>("dataType","unknown");
10	–	jetProducer_ = producersNames.getParameter<edm::InputTag>("jetProducer");
11	–	}
12	–
13	–	JetAnalyzer::JetAnalyzer(const edm::ParameterSet& producersNames, int verbosity):verbosity_(verbosity),useMC_(false)
14	–	{
15	–	dataType_ = producersNames.getUntrackedParameter<string>("dataType","unknown");
16	–	jetProducer_ = producersNames.getParameter<edm::InputTag>("jetProducer");
17	–	}
18	–
7		JetAnalyzer::JetAnalyzer(const edm::ParameterSet& producersNames, const edm::ParameterSet& myConfig, int verbosity):verbosity_(verbosity)
8		{
9		dataType_ = producersNames.getUntrackedParameter<string>("dataType","unknown");
10		jetProducer_ = producersNames.getParameter<edm::InputTag>("jetProducer");
11		useMC_ = myConfig.getUntrackedParameter<bool>("doJetMC");
12	+	allowMissingCollection_ = producersNames.getUntrackedParameter<bool>("allowMissingCollection", false);
13		}
14
15		JetAnalyzer::~JetAnalyzer()
#	Line 29 | Line 18 | JetAnalyzer::~JetAnalyzer()
18
19		bool Rsortrule (std::pair <double,double> p1, std::pair <double,double> p2 )
20		{
21	<	return p1.second<p2.second;
21	>	return p1.second<p2.second;
22		}
23
24	<	void JetAnalyzer::Process(const edm::Event& iEvent, TClonesArray* rootJets)
24	>	bool JetAnalyzer::process(const edm::Event& iEvent, TClonesArray* rootJets)
25		{
26	<
26	>
27		unsigned int nJets=0;
28	<
28	>
29		// reco::CaloJet or reco::PFJet ?
30		std::string jetType = "BASIC";
31		if( jetProducer_.label()=="kt4CaloJets"
#	Line 44 | Line 33 | void JetAnalyzer::Process(const edm::Eve
33		|| jetProducer_.label()=="iterativeCone5CaloJets"
34		|| jetProducer_.label()=="sisCone5CaloJets"
35		|| jetProducer_.label()=="sisCone7CaloJets"
36	<	) jetType="CALO";
37	<
36	>	|| jetProducer_.label()=="antikt5CaloJets"
37	>	) jetType="CALO";
38	>
39		if( jetProducer_.label()=="kt4PFJets"
40		|| jetProducer_.label()=="kt6PFJets"
41		|| jetProducer_.label()=="iterativeCone5PFJets"
42		|| jetProducer_.label()=="sisCone5PFJets"
43		|| jetProducer_.label()=="sisCone7PFJets"
44	<	) jetType="PF";
45	<
44	>	|| jetProducer_.label()=="antikt5PFJets"
45	>	) jetType="PF";
46	>
47		edm::Handle < std::vector <reco::CaloJet> > recoCaloJets;
48	<	if( (dataType_=="RECO" || dataType_=="AOD") && jetType=="CALO" )
48	>	if( dataType_=="RECO" && jetType=="CALO" )
49		{
50	<	iEvent.getByLabel(jetProducer_, recoCaloJets);
51	<	nJets = recoCaloJets->size();
50	>	try
51	>	{
52	>	iEvent.getByLabel(jetProducer_, recoCaloJets);
53	>	nJets = recoCaloJets->size();
54	>	}
55	>	catch (cms::Exception& exception)
56	>	{
57	>	if ( !allowMissingCollection_ )
58	>	{
59	>	cout << "  ##### ERROR IN  JetAnalyzer::process => reco::CaloJet collection is missing #####"<<endl;
60	>	throw exception;
61	>	}
62	>	if(verbosity_>1) cout <<  "   ===> No reco::CaloJet collection, skip jet info" << endl;
63	>	return false;
64	>	}
65		}
66	<
66	>
67		edm::Handle < std::vector <reco::PFJet> > recoPFJets;
68	<	if( (dataType_=="RECO" || dataType_=="AOD") && jetType=="PF" )
68	>	if( dataType_=="RECO" && jetType=="PF" )
69		{
70	<	iEvent.getByLabel(jetProducer_, recoPFJets);
71	<	nJets = recoPFJets->size();
70	>	try
71	>	{
72	>	iEvent.getByLabel(jetProducer_, recoPFJets);
73	>	nJets = recoPFJets->size();
74	>	}
75	>	catch (cms::Exception& exception)
76	>	{
77	>	if ( !allowMissingCollection_ )
78	>	{
79	>	cout << "  ##### ERROR IN  JetAnalyzer::process => reco::PFJet collection is missing #####"<<endl;
80	>	throw exception;
81	>	}
82	>	if(verbosity_>1) cout <<  "   ===> No reco::PFJet collection, skip jet info" << endl;
83	>	return false;
84	>	}
85		}
86	<
86	>
87		edm::Handle < std::vector <pat::Jet> > patJets;
88	<	if( dataType_=="PAT" || dataType_=="PATAOD" )
88	>	if( dataType_=="PAT" )
89		{
90	<	iEvent.getByLabel(jetProducer_, patJets);
91	<	nJets = patJets->size();
90	>	try
91	>	{
92	>	iEvent.getByLabel(jetProducer_, patJets);
93	>	nJets = patJets->size();
94	>	}
95	>	catch (cms::Exception& exception)
96	>	{
97	>	if ( !allowMissingCollection_ )
98	>	{
99	>	cout << "  ##### ERROR IN  JetAnalyzer::process => pat::Jet collection is missing #####"<<endl;
100	>	throw exception;
101	>	}
102	>	if(verbosity_>1) cout <<  "   ===> No pat::Jet collection, skip jet info" << endl;
103	>	return false;
104	>	}
105		}
106	<
106	>
107		if(verbosity_>1) std::cout << "   Number of jets = " << nJets << "   Label: " << jetProducer_.label() << "   Instance: " << jetProducer_.instance() << std::endl;
108	<
108	>
109		for (unsigned int j=0; j<nJets; j++)
110		{
111	<	const reco::Jet* jet = 0;
112	<	if( (dataType_=="RECO" || dataType_=="AOD") && jetType=="CALO" ) jet = (const reco::Jet*) ( & ((*recoCaloJets)[j]) );
113	<	if( (dataType_=="RECO" || dataType_=="AOD") && jetType=="PF" ) jet = (const reco::Jet*) ( & ((*recoPFJets)[j]) );
114	<	if( dataType_=="PAT" || dataType_=="PATAOD" )
111	>	const reco::Jet* jet = 0;
112	>	if( dataType_=="RECO" && jetType=="CALO" ) jet = (const reco::Jet*) ( & ((*recoCaloJets)[j]) );
113	>	if( dataType_=="RECO" && jetType=="PF" ) jet = (const reco::Jet*) ( & ((*recoPFJets)[j]) );
114	>	if( dataType_=="PAT" )
115		{
116		jet = (const reco::Jet*) ( & ((*patJets)[j]) );
117		if( (*patJets)[j].isCaloJet() ) jetType="CALO";
118		if( (*patJets)[j].isPFJet() ) jetType="PF";
119		}
120	<
120	>
121		TRootJet localJet(
122	<	jet->px()
123	<	,jet->py()
124	<	,jet->pz()
125	<	,jet->energy()
126	<	,jet->vx()
127	<	,jet->vy()
128	<	,jet->vz()
129	<	,abs(jet->pdgId())
130	<	,jet->charge()
131	<	);
132	<
122	>	jet->px()
123	>	,jet->py()
124	>	,jet->pz()
125	>	,jet->energy()
126	>	,jet->vx()
127	>	,jet->vy()
128	>	,jet->vz()
129	>	,jet->pdgId()
130	>	,jet->charge()
131	>	);
132	>
133		localJet.setNConstituents(jet->nConstituents());
134	+	//std::vector< edm::Ptr <Candidate> > constituents = jet->getJetConstituents();
135		localJet.setN90(jet->nCarrying(0.9));
136		localJet.setN60(jet->nCarrying(0.6));
137		localJet.setJetArea(jet->jetArea());
#	Line 115 | Line 146 | void JetAnalyzer::Process(const edm::Eve
146		localJet.setDR04EnergyFraction( jet->etInAnnulus(0.,0.4) / totalEnergy );
147		localJet.setDR05EnergyFraction( jet->etInAnnulus(0.,0.5) / totalEnergy );
148		}
149	<
150	<
151	<	if( dataType_=="RECO" || dataType_=="AOD" )
149	>
150	>
151	>	if( dataType_=="RECO" )
152		{
153		// Some specific methods to reco::Jet
154		// Do association to genParticle ?
155	<
155	>
156		if( jetType=="CALO" )
157		{
158		// Variables from reco::CaloJet
#	Line 131 | Line 162 | void JetAnalyzer::Process(const edm::Eve
162		localJet.setHcalEnergyFraction(caloJet->energyFractionHadronic());
163		//std::vector<CaloTowerPtr> getCaloConstituents () const;
164		}
165	<
165	>
166		if( jetType=="PF" )
167		{
168		// Variables from reco::PFJet
#	Line 143 | Line 174 | void JetAnalyzer::Process(const edm::Eve
174		localJet.setChargedMultiplicity((int)pfJet->chargedMultiplicity()) ;
175		//std::vector <const reco::PFCandidate*> getPFConstituents () const;
176		}
177	<
177	>
178		}
179	<
180	<
181	<	if( dataType_=="PATAOD" || dataType_=="PAT" )
179	>
180	>
181	>	if( dataType_=="PAT" )
182		{
183		// Some specific methods to pat::Jet
184		const pat::Jet *patJet = dynamic_cast<const pat::Jet*>(&*jet);
185	<
185	>
186		// Variables from pat::Jet (Basic)
187		localJet.setBtag_trackCountingHighEff(patJet->bDiscriminator("trackCountingHighEffBJetTags"));
188		localJet.setBtag_trackCountingHighPur(patJet->bDiscriminator("trackCountingHighPurBJetTags"));
189		localJet.setBtag_jetProbability(patJet->bDiscriminator("jetProbabilityBJetTags"));
190	<	localJet.setBCorrection(patJet->correctionFactor(pat::Jet::bCorrection));
191	<	localJet.setCCorrection(patJet->correctionFactor(pat::Jet::cCorrection));
192	<	localJet.setUDSCorrection(patJet->correctionFactor(pat::Jet::udsCorrection));
193	<	localJet.setGCorrection(patJet->correctionFactor(pat::Jet::gCorrection));
194	<
190	>
191	>	// see DataFormats/PatCandidates/interface/JetCorrFactors.h
192	>	localJet.setBCorrection(patJet->corrFactor("part","b"));
193	>	localJet.setCCorrection(patJet->corrFactor("part","c"));
194	>	localJet.setUDSCorrection(patJet->corrFactor("part","uds"));
195	>	localJet.setGCorrection(patJet->corrFactor("part","glu"));
196	>
197		// Use  associated tracks to calculate charged broadness of the jet
198		// FIXME - Check generalTracks collection is present
199		reco::TrackRefVector tracks =  patJet->associatedTracks() ;
#	Line 174 | Line 207 | void JetAnalyzer::Process(const edm::Eve
207		// TODO - Use std::map....
208		vector < pair < Float_t , Float_t > > tracksVPair ;
209		pair < Float_t , Float_t > tracksPair;
210	<
210	>
211		for(unsigned int iTracks = 0; iTracks< patJet->associatedTracks().size() ; iTracks++)
212		{
213		deltaR = localJet.DeltaR(TLorentzVector(tracks[iTracks]->px(),tracks[iTracks]->py(),tracks[iTracks]->pz(),0));
#	Line 210 | Line 243 | void JetAnalyzer::Process(const edm::Eve
243		localJet.setChargedBroadnessDR04(pTrackerCone04/pTrackerTot);
244		localJet.setChargedBroadnessDR05(pTrackerCone05/pTrackerTot);
245		}
246	<
247	<
246	>
247	>
248		if ( patJet->isCaloJet() ) {
249		// Variables from pat::Jet (CaloSpecific)
250		localJet.setJetType(1);
#	Line 220 | Line 253 | void JetAnalyzer::Process(const edm::Eve
253		localJet.setChargedMultiplicity(patJet->associatedTracks().size()) ;
254		//std::vector<CaloTowerPtr> getCaloConstituents () const;
255		}
256	<
256	>
257		if(patJet->isPFJet()) {
258		// Variables from pat::Jet (PFSpecific)
259		localJet.setJetType(2);
#	Line 230 | Line 263 | void JetAnalyzer::Process(const edm::Eve
263		localJet.setChargedMultiplicity((int)patJet->chargedMultiplicity()) ;
264		//std::vector <const reco::PFCandidate*> getPFConstituents () const;
265		}
266	<
266	>
267		// Matched genParticle
268		if (useMC_)
269		{
270	<	if ( patJet->genParton() != NULL )
270	>	// MC truth associator index
271	>	if ((patJet->genParticleRef()).isNonnull()) {
272	>	localJet.setGenParticleIndex((patJet->genParticleRef()).index());
273	>	} else {
274	>	localJet.setGenParticleIndex(-1);
275	>	}
276	>
277	>	// check if jet comes from a top
278	>	bool IsTopJet =  false;
279	>	if(patJet->genParton())
280		{
281	<	localJet.setMomentumMCParton(TLorentzVector(patJet->genParton()->px(),patJet->genParton()->py(),patJet->genParton()->pz(),patJet->genParton()->energy()));
282	<	localJet.setVertexMCParton(TVector3(patJet->genParton()->vx(),patJet->genParton()->vy(),patJet->genParton()->vz() ) );
283	<	localJet.setPdgIdMCParton(patJet->genParton()->pdgId());
281	>	const reco::Candidate* gen = patJet->genParton();
282	>	while(gen->mother())
283	>	{
284	>	if(abs((gen->mother())->pdgId()) == 6)
285	>	{
286	>	IsTopJet =  true;
287	>	break;
288	>	}
289	>	else
290	>	{
291	>	gen = (gen->mother() );
292	>	}
293	>	}
294		}
295	<
295	>	localJet.setIsTopJet(IsTopJet);
296	>
297	>	// Matched generated jet
298		if ( patJet->genJet() != NULL )
299		{
300	<	localJet.setMomentumMCJet(TLorentzVector(patJet->genJet()->px(),patJet->genJet()->py(),patJet->genJet()->pz(),patJet->genJet()->energy()));
247	<	localJet.setVertexMCJet(TVector3(patJet->genJet()->vx(),patJet->genJet()->vy(),patJet->genJet()->vz() ));
248	<	localJet.setPdgIdMCJet(patJet->genJet()->pdgId());
300	>	localJet.setGenJetEnergy( patJet->genJet()->energy() );
301		}
302	+
303		}
304		}
305	<
305	>
306		new( (*rootJets)[j] ) TRootJet(localJet);
307		if(verbosity_>2) cout << "   ["<< setw(3) << j << "] " << localJet << endl;
308		}
309	+
310	+	return true;
311		}

Diff Legend

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