ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/HbbAnalysis/plugins/HbbTreeMaker.cc
Revision: 1.18
Committed: Fri Jun 4 17:11:48 2010 UTC (14 years, 11 months ago) by amagnan
Content type: text/plain
Branch: MAIN
Changes since 1.17: +2 -0 lines
Log Message:
bug if no file read

File Contents

# Content
1 #include "DataFormats/Common/interface/Handle.h"
2 #include "DataFormats/Common/interface/TriggerResults.h"
3 #include "DataFormats/Common/interface/HLTenums.h"
4 #include "DataFormats/Common/interface/ValueMap.h"
5 #include "DataFormats/Candidate/interface/Candidate.h"
6 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
7 #include "DataFormats/MuonReco/interface/Muon.h"
8 #include "DataFormats/MuonReco/interface/MuonFwd.h"
9 #include "DataFormats/MuonReco/interface/MuonSelectors.h"
10 #include "DataFormats/VertexReco/interface/Vertex.h"
11 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
12 #include "DataFormats/TrackReco/interface/TrackFwd.h"
13
14 #include "DataFormats/PatCandidates/interface/Lepton.h"
15 #include "DataFormats/PatCandidates/interface/Muon.h"
16 #include "DataFormats/PatCandidates/interface/Electron.h"
17 #include "DataFormats/PatCandidates/interface/Tau.h"
18 #include "DataFormats/PatCandidates/interface/Jet.h"
19 #include "DataFormats/PatCandidates/interface/MET.h"
20
21 #include "FWCore/ServiceRegistry/interface/Service.h"
22 #include "CommonTools/UtilAlgos/interface/TFileService.h"
23
24
25 #include "UserCode/HbbAnalysis/interface/Electron.hh"
26 #include "UserCode/HbbAnalysis/interface/Muon.hh"
27 #include "UserCode/HbbAnalysis/interface/Tau.hh"
28 #include "UserCode/HbbAnalysis/interface/Jet.hh"
29 #include "UserCode/HbbAnalysis/interface/Met.hh"
30 #include "UserCode/HbbAnalysis/interface/Trigger.hh"
31 #include "UserCode/HbbAnalysis/interface/Vertex.hh"
32 #include "UserCode/HbbAnalysis/interface/L1Object.hh"
33 #include "UserCode/HbbAnalysis/interface/HLTObject.hh"
34
35 #include "UserCode/HbbAnalysis/plugins/HbbTreeMaker.hh"
36
37 using namespace HbbAnalysis;
38
39 HbbTreeMaker::HbbTreeMaker(const edm::ParameterSet & pset):
40 debug_(pset.getParameter<int>("DEBUG")),
41 processData_(false),
42 flavour_(pset.getParameter<unsigned int>("JetFlavour")),
43 doGen_(pset.getParameter<bool>("DOGEN")),
44 genParticleSrc_(pset.getParameter<edm::InputTag>("GenParticles")),
45 electronSrc_(pset.getParameter<edm::InputTag>("Electrons")),
46 muonSrc_(pset.getParameter<edm::InputTag>("Muons")),
47 caloTauSrc_(pset.getParameter<edm::InputTag>("CaloTaus")),
48 pfTauSrc_(pset.getParameter<edm::InputTag>("PFTaus")),
49 caloJetSrc_(pset.getParameter<edm::InputTag>("CaloJets")),
50 jptJetSrc_(pset.getParameter<edm::InputTag>("JPTJets")),
51 pfJetSrc_(pset.getParameter<edm::InputTag>("PFJets")),
52 caloMetSrc_(pset.getParameter<edm::InputTag>("CaloMET")),
53 tcMetSrc_(pset.getParameter<edm::InputTag>("TCMET")),
54 pfMetSrc_(pset.getParameter<edm::InputTag>("PFMET")),
55 //pairSrc_(pset.getParameter<edm::InputTag>("Pair")),
56 //mmPairSrc_(pset.getParameter<edm::InputTag>("MuMuPair")),
57 //etPairSrc_(pset.getParameter<edm::InputTag>("ETauPair")),
58 //mtPairSrc_(pset.getParameter<edm::InputTag>("MuTauPair")),
59 vertexSrc_(pset.getParameter<edm::InputTag>("Vertex")),
60 triggerSrc_(pset.getParameter<edm::InputTag>("Trigger")),
61 hltPaths_(pset.getParameter<std::vector<std::string> >("HLTPaths")),
62 l1CenJetSrc_(pset.getParameter<edm::InputTag>("L1CenJet")),
63 l1TauJetSrc_(pset.getParameter<edm::InputTag>("L1TauJet")),
64 l1FwdJetSrc_(pset.getParameter<edm::InputTag>("L1FwdJet")),
65 hltSummarySrc_(pset.getParameter<edm::InputTag>("HLTSummary")),
66 event_(0),
67 tree_(0),
68 hltTags_(pset.getParameter<std::vector<edm::InputTag> >("HLTTags"))
69 //processVec_(pset.getParameter<std::vector<std::string> >("ProcessVec"))
70 {//constructor
71
72
73 }//constructor
74
75 HbbTreeMaker::~HbbTreeMaker(){//destructor
76 }//destructor
77
78
79
80 void HbbTreeMaker::beginJob(){//beginJob
81
82
83 edm::Service<TFileService> lFileService;
84
85 TFileDirectory lDir = lFileService->mkdir("HbbAnalysis");
86
87 tree_ = lDir.make<TTree>("Tree","Tree");
88 tree_->Branch("HbbEvent","HbbAnalysis::HbbEvent",&event_);
89 event_ = new HbbEvent();
90
91 if (debug_) std::cout << "Initialising JetFlavour : " << std::endl;
92
93 lDir = lFileService->mkdir("JetFlavours");
94 jetFlav_.Initialise(lDir, debug_, flavour_);
95
96
97 }//beginJob
98
99 void HbbTreeMaker::endJob(){//endJob
100 if (!processData_) jetFlav_.printSummary();
101
102 //tree_->Write();
103 //delete tree_;
104 //delete event_;
105
106 }//endJob
107
108 void HbbTreeMaker::analyze(const edm::Event& aEvt, const edm::EventSetup& aEvtSetup){//analyze
109
110 //dumpContent(aEvt,"","MuTauPAT");
111 if (debug_) std::cout << "Processing event #" << aEvt.id().event() << std::endl;
112
113
114 static bool firstEvent = true;
115
116 event_->Clear();
117 event_->event(aEvt.id().event());
118 processData_ = aEvt.isRealData();
119
120 event_->run(aEvt.run());
121 event_->lumiBlock(aEvt.luminosityBlock());
122 event_->isRealData(aEvt.isRealData());
123
124 edm::Handle<reco::GenParticleCollection> lGenParticles;
125 try {
126 aEvt.getByLabel(genParticleSrc_,lGenParticles);
127 if (debug_) std::cout << "** ngenParticles = " << lGenParticles->size() << std::endl;
128 } catch(cms::Exception& e) {
129 if (!processData_) std::cout << "AMM: Collection genParticles not available! Exception : " << e.what() << ". " << std::endl;
130 }
131
132 if (doGen_ && !processData_) HbbParticles(lGenParticles,event_->particles());
133
134 unsigned int lNPartons = 0;
135 if (!processData_) lNPartons = jetFlav_.fillPartons(lGenParticles);
136
137 if (debug_ && !processData_) std::cout << "--- Number of partons = " << lNPartons << "." << std::endl;
138
139
140 edm::Handle<std::vector<reco::Vertex> > lRecoVertices;
141 try {
142 aEvt.getByLabel(vertexSrc_, lRecoVertices);
143 if (!lRecoVertices.isValid()){
144 edm::LogInfo("ERROR")<< "Error! Can't get vertex by label. ";
145 }
146 if (debug_) std::cout << "** vertexcollection = " << lRecoVertices->size() << " elements." << std::endl;
147 } catch(cms::Exception& e) {
148 std::cout << "AMM: Collection " << vertexSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
149 }
150
151 HbbVertices(lRecoVertices,event_->vertices());
152
153 edm::Handle<std::vector<pat::Electron> > lElectronCollection;
154
155 try {
156 aEvt.getByLabel(electronSrc_,lElectronCollection);
157 if (!lElectronCollection.isValid()){
158 edm::LogInfo("ERROR")<< "Error! Can't get electron by label. ";
159 }
160 if (debug_) std::cout << "** electroncollection = " << lElectronCollection->size() << " elements." << std::endl;
161 } catch(cms::Exception& e) {
162 std::cout << "AMM: Collection " << electronSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
163 }
164
165 HbbElectrons(lElectronCollection,event_->electrons());
166
167
168 edm::Handle<std::vector<pat::Muon> > lMuonCollection;
169
170 try {
171 aEvt.getByLabel(muonSrc_,lMuonCollection);
172 if (!lMuonCollection.isValid()){
173 edm::LogInfo("ERROR")<< "Error! Can't get muon by label. ";
174 }
175 if (debug_) std::cout << "** muoncollection = " << lMuonCollection->size() << " elements." << std::endl;
176 } catch(cms::Exception& e) {
177 std::cout << "AMM: Collection " << muonSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
178 }
179
180 HbbMuons(lMuonCollection,lRecoVertices,event_->muons());
181
182
183 if (!( caloTauSrc_.label()=="" && caloTauSrc_.instance()=="" )) {
184 edm::Handle<std::vector<pat::Tau> > lTauCollection;
185
186 try {
187 aEvt.getByLabel(caloTauSrc_,lTauCollection);
188 if (!lTauCollection.isValid()){
189 edm::LogInfo("ERROR")<< "Error! Can't get caloTau by label. ";
190 }
191 if (debug_) std::cout << "** caloTaucollection = " << lTauCollection->size() << " elements." << std::endl;
192 } catch(cms::Exception& e) {
193 std::cout << "AMM: Collection " << caloTauSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
194 }
195
196 HbbTaus(lTauCollection,lRecoVertices,event_->caloTaus());
197 }
198
199 edm::Handle<std::vector<pat::Tau> > lPFTauCollection;
200
201 try {
202 aEvt.getByLabel(pfTauSrc_,lPFTauCollection);
203 if (!lPFTauCollection.isValid()){
204 edm::LogInfo("ERROR")<< "Error! Can't get pfTau by label. ";
205 }
206 if (debug_) std::cout << "** pfTaucollection = " << lPFTauCollection->size() << " elements." << std::endl;
207 } catch(cms::Exception& e) {
208 std::cout << "AMM: Collection " << pfTauSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
209 }
210
211 HbbTaus(lPFTauCollection,lRecoVertices,event_->pfTaus());
212
213 edm::Handle<std::vector<pat::Jet> > lCaloJetCollection;
214
215 try {
216 aEvt.getByLabel(caloJetSrc_,lCaloJetCollection);
217 if (!lCaloJetCollection.isValid()){
218 edm::LogInfo("ERROR")<< "Error! Can't get caloJets by label. ";
219 }
220 if (debug_) std::cout << "** caloJetcollection = " << lCaloJetCollection->size() << " elements." << std::endl;
221 } catch(cms::Exception& e) {
222 std::cout << "AMM: Collection " << caloJetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
223 }
224
225 //std::cout << "Processing calo jets:" << std::endl;
226 HbbJets(lCaloJetCollection,jetFlav_,lGenParticles,event_->caloJets());
227
228 edm::Handle<std::vector<pat::Jet> > lJptJetCollection;
229
230 try {
231 aEvt.getByLabel(jptJetSrc_,lJptJetCollection);
232 if (!lJptJetCollection.isValid()){
233 edm::LogInfo("ERROR")<< "Error! Can't get jptJets by label. ";
234 }
235 if (debug_) std::cout << "** jptJetcollection = " << lJptJetCollection->size() << " elements." << std::endl;
236 } catch(cms::Exception& e) {
237 std::cout << "AMM: Collection " << jptJetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
238 }
239
240 //std::cout << "Processing JPT jets:" << std::endl;
241 HbbJets(lJptJetCollection,jetFlav_,lGenParticles,event_->jptJets());
242
243 edm::Handle<std::vector<pat::Jet> > lPfJetCollection;
244
245 try {
246 aEvt.getByLabel(pfJetSrc_,lPfJetCollection);
247 if (!lPfJetCollection.isValid()){
248 edm::LogInfo("ERROR")<< "Error! Can't get pfJets by label. ";
249 }
250 if (debug_) std::cout << "** pfJetcollection = " << lPfJetCollection->size() << " elements." << std::endl;
251 } catch(cms::Exception& e) {
252 std::cout << "AMM: Collection " << pfJetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
253 }
254
255 //std::cout << "Processing PF jets:" << std::endl;
256 HbbJets(lPfJetCollection,jetFlav_,lGenParticles,event_->pfJets());
257
258 edm::Handle<std::vector<pat::MET> > lCaloMetCol;
259
260 try {
261 aEvt.getByLabel(caloMetSrc_,lCaloMetCol);
262 if (!lCaloMetCol.isValid()){
263 edm::LogInfo("ERROR")<< "Error! Can't get caloMet by label. ";
264 }
265 if (debug_) std::cout << "** caloMet = " << lCaloMetCol->size() << " elements." << std::endl;
266 } catch(cms::Exception& e) {
267 std::cout << "AMM: Collection " << caloMetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
268 }
269
270 HbbMet(lCaloMetCol,event_->caloMet());
271
272 edm::Handle<std::vector<pat::MET> > lTcMetCol;
273
274 try {
275 aEvt.getByLabel(tcMetSrc_,lTcMetCol);
276 if (!lTcMetCol.isValid()){
277 edm::LogInfo("ERROR")<< "Error! Can't get tcMet by label. ";
278 }
279 if (debug_) std::cout << "** tcMet = " << lTcMetCol->size() << " elements." << std::endl;
280 } catch(cms::Exception& e) {
281 std::cout << "AMM: Collection " << tcMetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
282 }
283
284 HbbMet(lTcMetCol,event_->tcMet());
285
286 edm::Handle<std::vector<pat::MET> > lPfMetCol;
287
288 try {
289 aEvt.getByLabel(pfMetSrc_,lPfMetCol);
290 if (!lPfMetCol.isValid()){
291 edm::LogInfo("ERROR")<< "Error! Can't get pfMet by label. ";
292 }
293 if (debug_) std::cout << "** pfMet = " << lPfMetCol->size() << " elements." << std::endl;
294 } catch(cms::Exception& e) {
295 std::cout << "AMM: Collection " << pfMetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
296 }
297
298 HbbMet(lPfMetCol,event_->pfMet());
299
300
301
302 //triggers
303 edm::Handle<edm::TriggerResults> lTrigCol;
304 try {
305 aEvt.getByLabel(triggerSrc_,lTrigCol);
306 if (!lTrigCol.isValid()){
307 edm::LogInfo("ERROR")<< "Error! Can't get triggers by label. ";
308 }
309 if (debug_) std::cout << "** triggers = " << lTrigCol->size() << std::endl;
310 } catch(cms::Exception& e) {
311 std::cout << "AMM: Collection " << triggerSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
312 }
313 if (lTrigCol.isValid()) {
314 const edm::TriggerNames & lNames = aEvt.triggerNames(*lTrigCol);
315 HbbTrigger(lTrigCol,lNames,event_->triggers());
316 }
317
318 //L1 jet objects
319
320 edm::Handle<l1extra::L1JetParticleCollection> lCenJet;
321 edm::Handle<l1extra::L1JetParticleCollection> lFwdJet;
322 edm::Handle<l1extra::L1JetParticleCollection> lTauJet;
323 try {
324 aEvt.getByLabel(l1CenJetSrc_,lCenJet);
325 aEvt.getByLabel(l1TauJetSrc_,lTauJet);
326 aEvt.getByLabel(l1FwdJetSrc_,lFwdJet);
327 } catch(cms::Exception& e) {
328 std::cout << "AMM: L1 Collections " << l1CenJetSrc_ << ", " << l1TauJetSrc_ << ", " << l1FwdJetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
329 }
330 if (lCenJet.isValid()) {
331 HbbL1Objects(lCenJet,event_->l1Jets(),1);
332 }
333 if (lTauJet.isValid()) {
334 HbbL1Objects(lTauJet,event_->l1Jets(),2);
335 }
336 if (lFwdJet.isValid()) {
337 HbbL1Objects(lFwdJet,event_->l1Jets(),3);
338 }
339
340
341 edm::Handle<trigger::TriggerEvent> lHltSummary;
342 try {
343 aEvt.getByLabel(hltSummarySrc_,lHltSummary);
344 } catch(cms::Exception& e) {
345 std::cout << "AMM: collection " << hltSummarySrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
346 }
347 if (lHltSummary.isValid()) {
348 HbbHLTObjects(lHltSummary,event_->hltJets());
349 }
350
351
352
353 //event_->order();
354 tree_->Fill();
355 firstEvent = false;
356
357 }//analyze
358
359
360 void HbbTreeMaker::HbbElectrons(const edm::Handle<std::vector<pat::Electron> > & aCol,
361 std::vector<HbbAnalysis::Electron> & aVec)
362 {//HbbElectrons
363
364 if (aCol.isValid()){
365 if (aCol->size() > 0) {
366
367 unsigned int iEle = 0;
368 for (std::vector<pat::Electron>::const_iterator iter = aCol->begin();
369 iter != aCol->end();
370 iter++)
371 {
372 //if (debug_ > 1) std::cout << "**** Ele #" << iEle << ", pT,eta,phi = " << (*iter).pt() << " " << (*iter).eta() << " " << (*iter).phi() << std::endl;
373
374 HbbAnalysis::GenVars lGen;
375 if (!processData_ && (*iter).genLepton()){
376 lGen.valid = true;
377 lGen.E = (*iter).genLepton()->energy();
378 lGen.pT = (*iter).genLepton()->pt();
379 lGen.eta = (*iter).genLepton()->eta();
380 lGen.y = (*iter).genLepton()->y();
381 lGen.phi = (*iter).genLepton()->phi();
382 lGen.charge = (*iter).genLepton()->charge();
383 lGen.pdgId = (*iter).genLepton()->pdgId();
384 lGen.status = (*iter).genLepton()->status();
385 lGen.mass = (*iter).genLepton()->mass();
386 lGen.vx = (*iter).genLepton()->vx();
387 lGen.vy = (*iter).genLepton()->vy();
388 lGen.vz = (*iter).genLepton()->vz();
389 }
390 else {
391 lGen.valid = false;
392 lGen.E = 0;
393 lGen.pT = 0;
394 lGen.eta = 0;
395 lGen.y = 0;
396 lGen.phi = 0;
397 lGen.charge = 0;
398 lGen.pdgId = 0;
399 lGen.status = 0;
400 lGen.mass = 0;
401 lGen.vx = 0;
402 lGen.vy = 0;
403 lGen.vz = 0;
404 }
405
406 HbbAnalysis::BaseVars lReco;
407 lReco.E = (*iter).energy();
408 lReco.pT = (*iter).pt();
409 lReco.eta = (*iter).eta();
410 lReco.y = (*iter).y();
411 lReco.phi = (*iter).phi();
412 lReco.charge = (*iter).charge();
413 lReco.vx = (*iter).vx();
414 lReco.vy = (*iter).vy();
415 lReco.vz = (*iter).vz();
416
417 HbbAnalysis::SCVars lSC;
418 lSC.sigmaEtaEta = (*iter).scSigmaEtaEta();
419 lSC.sigmaIEtaIEta = (*iter).scSigmaIEtaIEta();
420 lSC.e1x5 = (*iter).scE1x5();
421 lSC.e2x5Max = (*iter).scE2x5Max();
422 lSC.e5x5 = (*iter).scE5x5();
423 lSC.eOverP = (*iter).eSuperClusterOverP();
424
425 HbbAnalysis::EleIsoVars lIso;
426 lIso.calo = (*iter).caloIso();
427 lIso.track = (*iter).trackIso();
428 lIso.ecal = (*iter).ecalIso();
429 lIso.hcal = (*iter).hcalIso();
430
431 HbbAnalysis::EleIDVars lID;
432 lID.electronIDs = (*iter).electronIDs();
433 lID.hOverE = (*iter).hadronicOverEm();
434 lID.deltaPhiIn = (*iter).deltaPhiSuperClusterTrackAtVtx();
435 lID.deltaEtaIn = (*iter).deltaEtaSuperClusterTrackAtVtx();
436 lID.ecalDrivenSeed = (*iter).ecalDrivenSeed();
437 lID.trackerDrivenSeed = (*iter).trackerDrivenSeed();
438
439 HbbAnalysis::Electron lObj(lGen,lReco,lSC,lIso,lID);
440 aVec.push_back(lObj);
441 iEle++;
442 }
443 }
444 }
445
446 }//HbbElectrons
447
448
449 void HbbTreeMaker::HbbMuons(const edm::Handle<std::vector<pat::Muon> > & aCol,
450 const edm::Handle<std::vector<reco::Vertex> > & aRecoVertices,
451 std::vector<HbbAnalysis::Muon> & aVec)
452 {//HbbMuons
453
454 if (aCol.isValid()){
455 if (aCol->size() > 0) {
456
457 unsigned int iEle = 0;
458 for (std::vector<pat::Muon>::const_iterator iter = aCol->begin();
459 iter != aCol->end();
460 iter++)
461 {
462 const reco::Muon* recoMuon = dynamic_cast<const reco::Muon*>((*iter).originalObject());
463
464 HbbAnalysis::GenVars lGen;
465 if (!processData_ && (*iter).genLepton()){
466 lGen.valid = true;
467 lGen.E = (*iter).genLepton()->energy();
468 lGen.pT = (*iter).genLepton()->pt();
469 lGen.eta = (*iter).genLepton()->eta();
470 lGen.y = (*iter).genLepton()->y();
471 lGen.phi = (*iter).genLepton()->phi();
472 lGen.charge = (*iter).genLepton()->charge();
473 lGen.pdgId = (*iter).genLepton()->pdgId();
474 lGen.status = (*iter).genLepton()->status();
475 lGen.mass = (*iter).genLepton()->mass();
476 lGen.vx = (*iter).genLepton()->vx();
477 lGen.vy = (*iter).genLepton()->vy();
478 lGen.vz = (*iter).genLepton()->vz();
479 }
480 else {
481 lGen.valid = false;
482 lGen.E = 0;
483 lGen.pT = 0;
484 lGen.eta = 0;
485 lGen.y = 0;
486 lGen.phi = 0;
487 lGen.charge = 0;
488 lGen.pdgId = 0;
489 lGen.status = 0;
490 lGen.mass = 0;
491 lGen.vx = 0;
492 lGen.vy = 0;
493 lGen.vz = 0;
494 }
495
496 HbbAnalysis::BaseVars lReco;
497 lReco.E = (*iter).energy();
498 lReco.pT = (*iter).pt();
499 lReco.eta = (*iter).eta();
500 lReco.y = (*iter).y();
501 lReco.phi = (*iter).phi();
502 lReco.charge = (*iter).charge();
503 lReco.vx = (*iter).vx();
504 lReco.vy = (*iter).vy();
505 lReco.vz = (*iter).vz();
506
507 HbbAnalysis::MuTrkVars lTrk;
508
509 reco::TrackRef lTrackerTrk = (*iter).innerTrack();
510 if ( lTrackerTrk.isAvailable() && lTrackerTrk.isNonnull() ) {
511 if ( aRecoVertices->size() >= 1 ) {
512 const reco::Vertex& thePrimaryEventVertex = (*aRecoVertices->begin());
513 lTrk.IPd0 = -lTrackerTrk->dxy(thePrimaryEventVertex.position());
514 lTrk.IPdz = lTrackerTrk->dz(thePrimaryEventVertex.position());
515 }
516 else {
517 lTrk.IPd0 = 0;
518 lTrk.IPdz = 0;
519 }
520
521 lTrk.nHits = lTrackerTrk->numberOfValidHits();
522 }
523 else {
524 lTrk.IPd0 = 0;
525 lTrk.IPdz = 0;
526 lTrk.nHits = 0;
527 }
528
529 HbbAnalysis::MuIsoVars lIsoR03;
530 lIsoR03.sumPt = recoMuon->isolationR03().sumPt;
531 lIsoR03.emEt = recoMuon->isolationR03().emEt;
532 lIsoR03.hadEt = recoMuon->isolationR03().hadEt;
533 lIsoR03.nTracks = recoMuon->isolationR03().nTracks;
534 lIsoR03.nJets = recoMuon->isolationR03().nJets;
535
536 HbbAnalysis::MuIsoVars lIsoR05;
537 lIsoR05.sumPt = recoMuon->isolationR05().sumPt;
538 lIsoR05.emEt = recoMuon->isolationR05().emEt;
539 lIsoR05.hadEt = recoMuon->isolationR05().hadEt;
540 lIsoR05.nTracks = recoMuon->isolationR05().nTracks;
541 lIsoR05.nJets = recoMuon->isolationR05().nJets;
542
543 HbbAnalysis::MuIDVars lID;
544 if (recoMuon->isGlobalMuon()) lID.type = 1;
545 else if (recoMuon->isTrackerMuon()) lID.type = 2;
546 else if (recoMuon->isStandAloneMuon()) lID.type = 3;
547 else if (recoMuon->isCaloMuon()) lID.type = 4;
548 else lID.type = 0;
549
550 if (muon::isGoodMuon(*recoMuon,muon::AllGlobalMuons)) lID.ids.push_back(1);
551 if (muon::isGoodMuon(*recoMuon,muon::AllStandAloneMuons)) lID.ids.push_back(2);
552 if (muon::isGoodMuon(*recoMuon,muon::AllTrackerMuons)) lID.ids.push_back(3);
553 if (muon::isGoodMuon(*recoMuon,muon::TrackerMuonArbitrated)) lID.ids.push_back(4);
554 if (muon::isGoodMuon(*recoMuon,muon::AllArbitrated)) lID.ids.push_back(5);
555 if (muon::isGoodMuon(*recoMuon,muon::GlobalMuonPromptTight)) lID.ids.push_back(6);
556 if (muon::isGoodMuon(*recoMuon,muon::TMLastStationLoose)) lID.ids.push_back(7);
557 if (muon::isGoodMuon(*recoMuon,muon::TMLastStationTight)) lID.ids.push_back(8);
558 if (muon::isGoodMuon(*recoMuon,muon::TM2DCompatibilityLoose)) lID.ids.push_back(9);
559 if (muon::isGoodMuon(*recoMuon,muon::TM2DCompatibilityTight)) lID.ids.push_back(10);
560 if (muon::isGoodMuon(*recoMuon,muon::TMOneStationLoose)) lID.ids.push_back(11);
561 if (muon::isGoodMuon(*recoMuon,muon::TMOneStationTight)) lID.ids.push_back(12);
562 if (muon::isGoodMuon(*recoMuon,muon::TMLastStationOptimizedLowPtLoose)) lID.ids.push_back(13);
563 if (muon::isGoodMuon(*recoMuon,muon::TMLastStationOptimizedLowPtTight)) lID.ids.push_back(14);
564 if (muon::isGoodMuon(*recoMuon,muon::GMTkChiCompatibility)) lID.ids.push_back(15);
565 if (muon::isGoodMuon(*recoMuon,muon::GMStaChiCompatibility)) lID.ids.push_back(16);
566 if (muon::isGoodMuon(*recoMuon,muon::GMTkKinkTight)) lID.ids.push_back(17);
567 if (muon::isGoodMuon(*recoMuon,muon::TMLastStationAngLoose)) lID.ids.push_back(18);
568 if (muon::isGoodMuon(*recoMuon,muon::TMLastStationAngTight)) lID.ids.push_back(19);
569 if (muon::isGoodMuon(*recoMuon,muon::TMOneStationAngLoose)) lID.ids.push_back(20);
570 if (muon::isGoodMuon(*recoMuon,muon::TMOneStationAngTight)) lID.ids.push_back(21);
571 if (muon::isGoodMuon(*recoMuon,muon::TMLastStationOptimizedBarrelLowPtLoose)) lID.ids.push_back(22);
572 if (muon::isGoodMuon(*recoMuon,muon::TMLastStationOptimizedBarrelLowPtTight)) lID.ids.push_back(23);
573
574 lID.caloCompat = recoMuon->caloCompatibility();
575 lID.segCompat = muon::segmentCompatibility(*recoMuon);
576 lID.nChambers = recoMuon->numberOfChambers();
577 lID.nMatchesLoose = recoMuon->numberOfMatches(reco::Muon::NoArbitration);
578 lID.nMatchesMedium = recoMuon->numberOfMatches(reco::Muon::SegmentArbitration);
579 lID.nMatchesTight = recoMuon->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
580
581 HbbAnalysis::Muon lObj(lGen,lReco,lTrk,lIsoR03,lIsoR05,lID);
582 aVec.push_back(lObj);
583 iEle++;
584 }
585 }
586 }
587
588 }//HbbMuons
589
590 void HbbTreeMaker::HbbTaus(const edm::Handle<std::vector<pat::Tau> > & aCol,
591 const edm::Handle<std::vector<reco::Vertex> > & aRecoVertices,
592 std::vector<HbbAnalysis::Tau> & aVec)
593 {//HbbTaus
594
595 if (aCol.isValid()){
596 if (aCol->size() > 0) {
597
598 unsigned int iEle = 0;
599 for (std::vector<pat::Tau>::const_iterator iter = aCol->begin();
600 iter != aCol->end();
601 iter++)
602 {
603 //if (debug_ > 1) std::cout << "**** Ele #" << iEle << ", pT,eta,phi = " << (*iter).pt() << " " << (*iter).eta() << " " << (*iter).phi() << std::endl;
604
605 HbbAnalysis::GenVars lGen;
606 if (!processData_ && (*iter).genLepton()){
607 lGen.valid = true;
608 lGen.E = (*iter).genLepton()->energy();
609 lGen.pT = (*iter).genLepton()->pt();
610 lGen.eta = (*iter).genLepton()->eta();
611 lGen.y = (*iter).genLepton()->y();
612 lGen.phi = (*iter).genLepton()->phi();
613 lGen.charge = (*iter).genLepton()->charge();
614 lGen.pdgId = (*iter).genLepton()->pdgId();
615 lGen.status = (*iter).genLepton()->status();
616 lGen.mass = (*iter).genLepton()->mass();
617 lGen.vx = (*iter).genLepton()->vx();
618 lGen.vy = (*iter).genLepton()->vy();
619 lGen.vz = (*iter).genLepton()->vz();
620 }
621 else {
622 lGen.valid = false;
623 lGen.E = 0;
624 lGen.pT = 0;
625 lGen.eta = 0;
626 lGen.y = 0;
627 lGen.phi = 0;
628 lGen.charge = 0;
629 lGen.pdgId = 0;
630 lGen.status = 0;
631 lGen.mass = 0;
632 lGen.vx = 0;
633 lGen.vy = 0;
634 lGen.vz = 0;
635 }
636
637
638 HbbAnalysis::GenVars lGenJet;
639 if (!processData_ && (*iter).genJet()){
640 lGenJet.valid = true;
641 lGenJet.E = (*iter).genJet()->energy();
642 lGenJet.pT = (*iter).genJet()->pt();
643 lGenJet.eta = (*iter).genJet()->eta();
644 lGenJet.y = (*iter).genJet()->y();
645 lGenJet.phi = (*iter).genJet()->phi();
646 lGenJet.charge = (*iter).genJet()->charge();
647 lGenJet.pdgId = (*iter).genJet()->pdgId();
648 lGenJet.status = (*iter).genJet()->status();
649 lGenJet.mass = (*iter).genJet()->mass();
650 lGenJet.vx = (*iter).genJet()->vx();
651 lGenJet.vy = (*iter).genJet()->vy();
652 lGenJet.vz = (*iter).genJet()->vz();
653 }
654 else {
655 lGenJet.valid = false;
656 lGenJet.E = 0;
657 lGenJet.pT = 0;
658 lGenJet.eta = 0;
659 lGenJet.y = 0;
660 lGenJet.phi = 0;
661 lGenJet.charge = 0;
662 lGenJet.pdgId = 0;
663 lGenJet.status = 0;
664 lGenJet.mass = 0;
665 lGenJet.vx = 0;
666 lGenJet.vy = 0;
667 lGenJet.vz = 0;
668 }
669
670 HbbAnalysis::BaseVars lReco;
671 lReco.E = (*iter).energy();
672 lReco.pT = (*iter).pt();
673 lReco.eta = (*iter).eta();
674 lReco.y = (*iter).y();
675 lReco.phi = (*iter).phi();
676 lReco.charge = (*iter).charge();
677 lReco.vx = (*iter).vx();
678 lReco.vy = (*iter).vy();
679 lReco.vz = (*iter).vz();
680
681 bool isCalo = (*iter).isCaloTau();
682 bool isPF = (*iter).isPFTau();
683
684 assert (isCalo == !isPF);
685
686 HbbAnalysis::TrkVars lLead;
687 if (isCalo) {
688 lLead.signedSipt = (*iter).leadTracksignedSipt();
689
690 reco::TrackRef lCaloTrk = (*iter).leadTrack();
691
692 if ( lCaloTrk.isAvailable() && lCaloTrk.isNonnull() ) {
693 lLead.pT = lCaloTrk->pt();
694 lLead.eta = lCaloTrk->eta();
695 lLead.phi = lCaloTrk->phi();
696
697 lLead.matchDist = reco::deltaR(lCaloTrk->momentum(), (*iter).p4());
698
699 if ( aRecoVertices->size() >= 1 ) {
700 const reco::Vertex& thePrimaryEventVertex = (*aRecoVertices->begin());
701 lLead.IPxy = lCaloTrk->dxy(thePrimaryEventVertex.position());
702 lLead.IPz = lCaloTrk->dz(thePrimaryEventVertex.position());
703 }
704 else {
705 lLead.IPxy = 0;
706 lLead.IPz = 0;
707 }
708 }
709 else {
710 lLead.pT = 0;
711 lLead.eta = 0;
712 lLead.phi = 0;
713 lLead.matchDist = 0;
714 lLead.IPxy = 0;
715 lLead.IPz = 0;
716 lLead.signedSipt = 0;
717 }
718 }//calo
719 else {
720 lLead.signedSipt = (*iter).leadPFChargedHadrCandsignedSipt();
721
722 reco::PFCandidateRef lHadrCand = (*iter).leadPFChargedHadrCand();
723
724 if ( lHadrCand.isAvailable() && lHadrCand.isNonnull() ) {
725 lLead.pT = lHadrCand->pt();
726 lLead.eta = lHadrCand->eta();
727 lLead.phi = lHadrCand->phi();
728
729 lLead.matchDist = reco::deltaR(lHadrCand->p4(), (*iter).p4());
730
731 reco::TrackRef lTrk = lHadrCand->trackRef();
732
733 if ( aRecoVertices->size() >= 1) {
734 const reco::Vertex& thePrimaryEventVertex = (*aRecoVertices->begin());
735 lLead.IPxy = lTrk->dxy(thePrimaryEventVertex.position());
736 lLead.IPz = lTrk->dz(thePrimaryEventVertex.position());
737 }
738 else {
739 lLead.IPxy = 0;
740 lLead.IPz = 0;
741 }
742 }
743 else {
744 lLead.pT = 0;
745 lLead.eta = 0;
746 lLead.phi = 0;
747 lLead.matchDist = 0;
748 lLead.IPxy = 0;
749 lLead.IPz = 0;
750 lLead.signedSipt = 0;
751 }
752
753 }//pf
754
755
756
757 const std::vector<std::pair<std::string,float> > & lIDs = (*iter).tauIDs();
758 bool lPrint = false;
759 //a few warning if additional discriminants are found:
760 if ((isPF && lIDs.size() != 16) ||
761 (isCalo && lIDs.size() != 3))
762 lPrint = true;
763
764 if (lPrint) {
765 std::cout << "!!!!!!! Discriminants changed, please update tree !!!!!!!" << std::endl;
766 std::cout << "--- isCaloTau = " << isCalo << ", isPFTau = " << isPF << std::endl;
767 std::cout << "------ ID names = " << std::endl;
768 }
769
770
771 HbbAnalysis::CaloTauIDVars lcaloId;
772 lcaloId.byIsolation = 0;
773 lcaloId.byLeadingTrackFinding = 0;
774 lcaloId.byLeadingTrackPtCut = 0;
775
776 HbbAnalysis::PFTauIDVars lpfId;
777 lpfId.byLeadingTrackFinding = 0;
778 lpfId.byLeadingTrackPtCut = 0;
779 lpfId.byTrackIsolation = 0;
780 lpfId.byECALIsolation = 0;
781 lpfId.byIsolation = 0;
782 lpfId.againstElectron = 0;
783 lpfId.againstMuon = 0;
784 lpfId.byIsolationUsingLeadingPion = 0;
785 lpfId.byTaNC = 0;
786 lpfId.byTaNCfrHalfPercent = 0;
787 lpfId.byTaNCfrOnePercent = 0;
788 lpfId.byTaNCfrQuarterPercent = 0;
789 lpfId.byTaNCfrTenthPercent = 0;
790 lpfId.ecalIsolationUsingLeadingPion = 0;
791 lpfId.leadingPionPtCut = 0;
792 lpfId.trackIsolationUsingLeadingPion = 0;
793
794
795
796 for (unsigned int id(0); id<lIDs.size(); id++){
797
798 std::string lName = lIDs.at(id).first;
799 float lDiscri = lIDs.at(id).second;
800
801 if (lPrint) std::cout << "--------- " << lName << " = " << lDiscri << std::endl;
802 //pf
803
804 if (isPF) {
805 if (lName.find("leadingTrackFinding") != lName.npos) lpfId.byLeadingTrackFinding = lDiscri;
806 if (lName.find("leadingTrackPtCut") != lName.npos) lpfId.byLeadingTrackPtCut = lDiscri;
807 if (lName.find("trackIsolationUsingLeadingPion") != lName.npos) lpfId.trackIsolationUsingLeadingPion = lDiscri;
808 if (lName.find("trackIsolation") != lName.npos &&
809 lName.find("trackIsolationUsingLeadingPion") == lName.npos) lpfId.byTrackIsolation = lDiscri;
810 if (lName.find("ecalIsolationUsingLeadingPion") != lName.npos) lpfId.ecalIsolationUsingLeadingPion = lDiscri;
811 if (lName.find("ecalIsolation") != lName.npos &&
812 lName.find("ecalIsolationUsingLeadingPion") == lName.npos) lpfId.byECALIsolation = lDiscri;
813 if (lName.find("byIsolationUsingLeadingPion") != lName.npos) lpfId.byIsolationUsingLeadingPion = lDiscri;
814 if (lName.find("byIsolation") != lName.npos &&
815 lName.find("byIsolationUsingLeadingPion") == lName.npos) lpfId.byIsolation = lDiscri;
816 if (lName.find("againstElectron") != lName.npos) lpfId.againstElectron = lDiscri;
817 if (lName.find("againstMuon") != lName.npos) lpfId.againstMuon = lDiscri;
818 if (lName.find("byTaNCfrHalfPercent") != lName.npos) lpfId.byTaNCfrHalfPercent = lDiscri;
819 if (lName.find("byTaNCfrOnePercent") != lName.npos) lpfId.byTaNCfrOnePercent = lDiscri;
820 if (lName.find("byTaNCfrQuarterPercent") != lName.npos) lpfId.byTaNCfrQuarterPercent = lDiscri;
821 if (lName.find("byTaNCfrTenthPercent") != lName.npos) lpfId.byTaNCfrTenthPercent = lDiscri;
822 if (lName.find("byTaNC") != lName.npos &&
823 lName.find("byTaNCfr") == lName.npos) lpfId.byTaNC = lDiscri;
824 if (lName.find("leadingPionPtCut") != lName.npos) lpfId.leadingPionPtCut = lDiscri;
825 }
826 if (isCalo){
827 if (lName.find("byIsolation") != lName.npos) lcaloId.byIsolation = lDiscri;
828 if (lName.find("leadingTrackFinding") != lName.npos) lcaloId.byLeadingTrackFinding = lDiscri;
829 if (lName.find("leadingTrackPtCut") != lName.npos) lcaloId.byLeadingTrackPtCut = lDiscri;
830 }
831
832
833 }
834
835
836 if (isCalo){
837
838 HbbAnalysis::CaloTauIsoVars lcaloIso;
839 lcaloIso.nIsoTracks = (*iter).isolationTracks().size();
840 lcaloIso.nSigTracks = (*iter).signalTracks().size();
841 lcaloIso.leadTrackHCAL3x3hitsEtSum = (*iter).leadTrackHCAL3x3hitsEtSum();
842 lcaloIso.leadTrackHCAL3x3hottesthitDEta = (*iter).leadTrackHCAL3x3hottesthitDEta();
843 lcaloIso.signalTracksInvariantMass = (*iter).signalTracksInvariantMass();
844 lcaloIso.tracksInvariantMass = (*iter).TracksInvariantMass();
845 lcaloIso.isolationTracksPtSum = (*iter).isolationTracksPtSum();
846 lcaloIso.isolationECALhitsEtSum = (*iter).isolationECALhitsEtSum();
847 lcaloIso.maximumHCALhitEt = (*iter).maximumHCALhitEt();
848
849 HbbAnalysis::Tau lObj(lGen,lGenJet,lReco,lLead,lcaloIso,lcaloId);
850 aVec.push_back(lObj);
851 }
852 else {
853
854 HbbAnalysis::PFTauIsoVars lpfIso;
855 lpfIso.nSigCands = (*iter).signalPFCands().size();
856 lpfIso.nIsoCands = (*iter).isolationPFCands().size();
857 lpfIso.maximumHCALPFClusterEt = (*iter).maximumHCALPFClusterEt();
858 //lpfIso.emFraction = (*iter).emFraction();
859 lpfIso.hcalTotOverPLead = (*iter).hcalTotOverPLead();
860 lpfIso.hcalMaxOverPLead = (*iter).hcalMaxOverPLead();
861 lpfIso.hcal3x3OverPLead = (*iter).hcal3x3OverPLead();
862 lpfIso.ecalStripSumEOverPLead = (*iter).ecalStripSumEOverPLead();
863 lpfIso.bremsRecoveryEOverPLead = (*iter).bremsRecoveryEOverPLead();
864
865 HbbAnalysis::PFTauCandVars lHadr;
866 lHadr.nSigCands = (*iter).signalPFChargedHadrCands().size();
867 lHadr.nIsoCands = (*iter).isolationPFChargedHadrCands().size();
868 lHadr.isolationPtSum = (*iter).isolationPFChargedHadrCandsPtSum();
869
870 HbbAnalysis::PFTauCandVars lNeutr;
871 lNeutr.nSigCands = (*iter).signalPFNeutrHadrCands().size();
872 lNeutr.nIsoCands = (*iter).isolationPFNeutrHadrCands().size();
873 lNeutr.isolationPtSum = (*iter).neutralHadronIso();
874
875 HbbAnalysis::PFTauCandVars lGamma;
876 lGamma.nSigCands = (*iter).signalPFGammaCands().size();
877 lGamma.nIsoCands = (*iter).isolationPFGammaCands().size();
878 lGamma.isolationPtSum = (*iter).isolationPFGammaCandsEtSum();
879
880 HbbAnalysis::PFTauEleIDVars lEleId;
881 if ( (*iter).electronPreIDTrack().isAvailable() && (*iter).electronPreIDTrack().isNonnull() ) {
882 lEleId.pT = (*iter).electronPreIDTrack()->pt();
883 lEleId.eta = (*iter).electronPreIDTrack()->eta();
884 lEleId.phi = (*iter).electronPreIDTrack()->phi();
885 }
886 else {
887 lEleId.pT = 0;
888 lEleId.eta = 0;
889 lEleId.phi = 0;
890 }
891
892 lEleId.output = (*iter).electronPreIDOutput();
893 lEleId.decision = (*iter).electronPreIDDecision();
894
895 HbbAnalysis::PFTauMuIDVars lMuId;
896 lMuId.caloCompat = (*iter).caloComp();
897 lMuId.segCompat = (*iter).segComp();
898 lMuId.decision = (*iter).muonDecision();
899
900 HbbAnalysis::Tau lObj(lGen,lGenJet,lReco,lLead,lpfIso,lHadr,lNeutr,lGamma,lpfId,lEleId,lMuId);
901 aVec.push_back(lObj);
902 }
903 iEle++;
904 }
905 }
906 }
907
908 }//HbbTaus
909
910 void HbbTreeMaker::HbbJets(const edm::Handle<std::vector<pat::Jet> > & aCol,
911 const HbbAnalysis::JetFlavour & aJetFlav,
912 const edm::Handle<reco::GenParticleCollection> & aGenParticles,
913 std::vector<HbbAnalysis::Jet> & aVec)
914 {//HbbJets
915
916 if (aCol.isValid()){//valid
917 if (aCol->size() > 0) {//non empty
918
919 unsigned int iEle = 0;
920 for (std::vector<pat::Jet>::const_iterator iter = aCol->begin();
921 iter != aCol->end();
922 iter++)
923 {//loop on element
924 //if (debug_ > 1) std::cout << "**** Ele #" << iEle << ", pT,eta,phi = " << (*iter).pt() << " " << (*iter).eta() << " " << (*iter).phi() << std::endl;
925
926 HbbAnalysis::GenVars lGen;
927 if (!processData_ && (*iter).genParton()){
928 lGen.valid = true;
929 lGen.E = (*iter).genParton()->energy();
930 lGen.pT = (*iter).genParton()->pt();
931 lGen.eta = (*iter).genParton()->eta();
932 lGen.y = (*iter).genParton()->y();
933 lGen.phi = (*iter).genParton()->phi();
934 lGen.charge = (*iter).genParton()->charge();
935 lGen.pdgId = (*iter).genParton()->pdgId();
936 lGen.status = (*iter).genParton()->status();
937 lGen.mass = (*iter).genParton()->mass();
938 lGen.vx = (*iter).genParton()->vx();
939 lGen.vy = (*iter).genParton()->vy();
940 lGen.vz = (*iter).genParton()->vz();
941 }
942 else {
943 lGen.valid = false;
944 lGen.E = 0;
945 lGen.pT = 0;
946 lGen.eta = 0;
947 lGen.y = 0;
948 lGen.phi = 0;
949 lGen.charge = 0;
950 lGen.pdgId = 0;
951 lGen.status = 0;
952 lGen.mass = 0;
953 lGen.vx = 0;
954 lGen.vy = 0;
955 lGen.vz = 0;
956 }
957
958
959 HbbAnalysis::GenVars lGenJet;
960 if (!processData_ && (*iter).genJet()){
961 lGenJet.valid = true;
962 lGenJet.E = (*iter).genJet()->energy();
963 lGenJet.pT = (*iter).genJet()->pt();
964 lGenJet.eta = (*iter).genJet()->eta();
965 lGenJet.y = (*iter).genJet()->y();
966 lGenJet.phi = (*iter).genJet()->phi();
967 lGenJet.charge = (*iter).genJet()->charge();
968 lGenJet.pdgId = (*iter).genJet()->pdgId();
969 lGenJet.status = (*iter).genJet()->status();
970 lGenJet.mass = (*iter).genJet()->mass();
971 lGenJet.vx = (*iter).genJet()->vx();
972 lGenJet.vy = (*iter).genJet()->vy();
973 lGenJet.vz = (*iter).genJet()->vz();
974 }
975 else {
976 lGenJet.valid = false;
977 lGenJet.E = 0;
978 lGenJet.pT = 0;
979 lGenJet.eta = 0;
980 lGenJet.y = 0;
981 lGenJet.phi = 0;
982 lGenJet.charge = 0;
983 lGenJet.pdgId = 0;
984 lGenJet.status = 0;
985 lGenJet.mass = 0;
986 lGenJet.vx = 0;
987 lGenJet.vy = 0;
988 lGenJet.vz = 0;
989 }
990
991 HbbAnalysis::BaseVars lReco;
992 lReco.E = (*iter).energy();
993 lReco.pT = (*iter).pt();
994 lReco.eta = (*iter).eta();
995 //lReco.etaDet = (*iter).detectorEta((*iter).vz(),(*iter).eta());
996 lReco.y = (*iter).y();
997 lReco.phi = (*iter).phi();
998 lReco.charge = (*iter).jetCharge();
999 lReco.vx = (*iter).vx();
1000 lReco.vy = (*iter).vy();
1001 lReco.vz = (*iter).vz();
1002 //lReco.myEtaDet = HbbAnalysis::EtaDetector(lReco);
1003
1004 unsigned int lFlav = 0;
1005 if (!processData_ && aJetFlav.nPartons()) {
1006 lFlav = aJetFlav.partonMatchingGenJet((*iter),aGenParticles,0.4).second;
1007 }
1008
1009 HbbAnalysis::JetVars lCommon;
1010 lCommon.flavour = lFlav;
1011 lCommon.partonFlavour = (*iter).partonFlavour();
1012 lCommon.nAssociatedTracks = ((*iter).associatedTracks()).size();
1013 lCommon.rawpT = (*iter).pt();
1014 if ((*iter).hasCorrFactors() && (*iter).corrFactor("raw")>0)
1015 lCommon.rawpT = (*iter).pt()/((*iter).corrFactor("raw"));
1016
1017 lCommon.etaMean = (*iter).etaPhiStatistics().etaMean;
1018 lCommon.phiMean = (*iter).etaPhiStatistics().phiMean;
1019 lCommon.etaEtaMoment = (*iter).etaPhiStatistics().etaEtaMoment;
1020 lCommon.phiPhiMoment = (*iter).etaPhiStatistics().phiPhiMoment;
1021 lCommon.etaPhiMoment = (*iter).etaPhiStatistics().etaPhiMoment;
1022
1023 lCommon.l1Match = false;//TODO
1024 lCommon.hltMatch = false;//TODO
1025
1026 //lCommon.rawEta = (*iter).;
1027 //lCommon.rawPhi = (*iter).;
1028 // p_hasJetCorrFactors->Fill(aJet.hasJetCorrFactors());
1029
1030 HbbAnalysis::JetBtagVars lBtag;
1031 lBtag.cSV = (*iter).bDiscriminator("combinedSecondaryVertexBJetTags");
1032 lBtag.cSVMVA = (*iter).bDiscriminator("combinedSecondaryVertexMVABJetTags");
1033 lBtag.iPMVA = (*iter).bDiscriminator("impactParameterMVABJetTags");
1034 lBtag.bProba = (*iter).bDiscriminator("jetBProbabilityBJetTags");
1035 lBtag.probability = (*iter).bDiscriminator("jetProbabilityBJetTags");
1036 lBtag.sSVHP = (*iter).bDiscriminator("simpleSecondaryVertexHighPurBJetTags");
1037 lBtag.sSVHE = (*iter).bDiscriminator("simpleSecondaryVertexHighEffBJetTags");
1038 lBtag.softElectronByPt = ((*iter).bDiscriminator("softElectronByPtBJetTags") < -1000) ? -1000 : (*iter).bDiscriminator("softElectronByPtBJetTags");
1039 lBtag.softElectronByIP3d = ((*iter).bDiscriminator("softElectronByIP3dBJetTags") < -1000) ? -1000 : (*iter).bDiscriminator("softElectronByIP3dBJetTags");
1040 lBtag.softMuon = ((*iter).bDiscriminator("softMuonBJetTags") < -1000) ? -1000 : (*iter).bDiscriminator("softMuonBJetTags");
1041 lBtag.softMuonByPt = ((*iter).bDiscriminator("softMuonByPtBJetTags") < -1000) ? -1000 : (*iter).bDiscriminator("softMuonByPtBJetTags");
1042 lBtag.softMuonByIP3d = ((*iter).bDiscriminator("softMuonByIP3dBJetTags") < -1000) ? -1000 : (*iter).bDiscriminator("softMuonByIP3dBJetTags");
1043
1044 //std::cout << " -- New values of b-discri for jet : " << iEle << std::endl
1045 // << " ---- softElecs: " << lBtag.softElectronByPt << " " << lBtag.softElectronByIP3d << std::endl
1046 // << " ---- softMus: " << lBtag.softMuon << " " << lBtag.softMuonByPt << " " << lBtag.softMuonByIP3d << std::endl;
1047
1048
1049 lBtag.tCHE = (*iter).bDiscriminator("trackCountingHighEffBJetTags");
1050 lBtag.tCHP = (*iter).bDiscriminator("trackCountingHighPurBJetTags");
1051
1052
1053 bool isCalo = (*iter).isCaloJet();
1054 bool isPF = (*iter).isPFJet();
1055 bool isJPT = (*iter).isJPTJet();
1056
1057 //std::cout << " -- isCalo = " << isCalo << ", isPF = " << isPF << std::endl;
1058
1059 //assert (isCalo == !isPF);
1060 HbbAnalysis::JetIDVars lId;
1061 lId.fHPD = 0;
1062 lId.fRBX = 0;
1063 lId.n90Hits = 0;
1064 lId.hitsInN90 = 0;
1065 lId.restrictedEMF = 0;
1066 lId.nHCALTowers = 0;
1067 lId.nECALTowers = 0;
1068 if (isCalo || isJPT) {
1069 lId.fHPD = (*iter).jetID().fHPD;
1070 lId.fRBX = (*iter).jetID().fRBX;
1071 lId.n90Hits = (*iter).jetID().n90Hits;
1072 lId.hitsInN90 = (*iter).jetID().hitsInN90;
1073 lId.restrictedEMF = (*iter).jetID().restrictedEMF;
1074 lId.nHCALTowers = (*iter).jetID().nHCALTowers;
1075 lId.nECALTowers = (*iter).jetID().nECALTowers;
1076 }
1077 if (isCalo){
1078 HbbAnalysis::CaloJetVars lCalo;
1079 lCalo.maxEInEmTowers = (*iter).maxEInEmTowers();
1080 lCalo.maxEInHadTowers = (*iter).maxEInHadTowers();
1081 lCalo.energyFractionHadronic = (*iter).energyFractionHadronic();
1082 lCalo.emEnergyFraction = (*iter).emEnergyFraction();
1083 lCalo.hadEnergyInHB = (*iter).hadEnergyInHB();
1084 lCalo.hadEnergyInHO = (*iter).hadEnergyInHO();
1085 lCalo.hadEnergyInHE = (*iter).hadEnergyInHE();
1086 lCalo.hadEnergyInHF = (*iter).hadEnergyInHF();
1087 lCalo.emEnergyInEB = (*iter).emEnergyInEB();
1088 lCalo.emEnergyInEE = (*iter).emEnergyInEE();
1089 lCalo.emEnergyInHF = (*iter).emEnergyInHF();
1090 lCalo.towersArea = (*iter).towersArea();
1091 lCalo.n90 = (*iter).n90();
1092 lCalo.n60 = (*iter).n60();
1093
1094 HbbAnalysis::Jet lObj(lGen,lGenJet,lReco,lCommon,lCalo,lBtag,lId);
1095 aVec.push_back(lObj);
1096 }
1097 HbbAnalysis::JPTPFJetVars lJPTPF;
1098 lJPTPF.chargedHadronEnergy = 0;
1099 lJPTPF.chargedHadronEnergyFraction = 0;
1100 lJPTPF.neutralHadronEnergy = 0;
1101 lJPTPF.neutralHadronEnergyFraction = 0;
1102 lJPTPF.chargedEmEnergy = 0;
1103 lJPTPF.chargedEmEnergyFraction = 0;
1104 lJPTPF.neutralEmEnergy = 0;
1105 lJPTPF.neutralEmEnergyFraction = 0;
1106 lJPTPF.chargedMultiplicity = 0;
1107 lJPTPF.muonMultiplicity = 0;
1108
1109 if (isPF || isJPT){
1110 lJPTPF.chargedEmEnergy = (*iter).chargedEmEnergy();
1111 lJPTPF.neutralEmEnergy = (*iter).neutralEmEnergy();
1112 lJPTPF.muonMultiplicity = (*iter).muonMultiplicity();
1113 lJPTPF.chargedHadronEnergy = (*iter).chargedHadronEnergy();
1114 lJPTPF.neutralHadronEnergy = (*iter).neutralHadronEnergy();
1115 lJPTPF.chargedMultiplicity = (*iter).chargedMultiplicity();
1116
1117 lJPTPF.chargedHadronEnergyFraction = (*iter).chargedHadronEnergyFraction();
1118 lJPTPF.neutralHadronEnergyFraction = (*iter).neutralHadronEnergyFraction();
1119 lJPTPF.chargedEmEnergyFraction = (*iter).chargedEmEnergyFraction();
1120 lJPTPF.neutralEmEnergyFraction = (*iter).neutralEmEnergyFraction();
1121 }
1122 if (isJPT) {
1123 HbbAnalysis::JPTJetVars lJPT;
1124 // const reco::TrackRefVector& pionsInVertexInCalo () const{return jptSpecific().pionsInVertexInCalo; }
1125 // const reco::TrackRefVector& pionsInVertexOutCalo() const{return jptSpecific().pionsInVertexOutCalo;}
1126 // const reco::TrackRefVector& pionsOutVertexInCalo() const{return jptSpecific().pionsOutVertexInCalo;}
1127 // const reco::TrackRefVector& muonsInVertexInCalo () const{return jptSpecific().muonsInVertexInCalo; }
1128 // const reco::TrackRefVector& muonsInVertexOutCalo() const{return jptSpecific().muonsInVertexOutCalo;}
1129 // const reco::TrackRefVector& muonsOutVertexInCalo() const{return jptSpecific().muonsOutVertexInCalo;}
1130 // const reco::TrackRefVector& elecsInVertexInCalo () const{return jptSpecific().elecsInVertexInCalo; }
1131 // const reco::TrackRefVector& elecsInVertexOutCalo() const{return jptSpecific().elecsInVertexOutCalo;}
1132 // const reco::TrackRefVector& elecsOutVertexInCalo() const{return jptSpecific().elecsOutVertexInCalo;}
1133 lJPT.zspCorrection = (*iter).zspCorrection();
1134 lJPT.elecMultiplicity = (*iter).elecMultiplicity();
1135
1136
1137 HbbAnalysis::CaloJetVars lCalo;
1138 const reco::JPTJet * lJetRef = dynamic_cast<const reco::JPTJet *>((*iter).originalObject());
1139 edm::RefToBase<reco::Jet> jptjetRef = lJetRef->getCaloJetRef();
1140
1141 // Dynamic_cast to CaloJet object
1142 reco::CaloJet const * rawcalojet =
1143 dynamic_cast<reco::CaloJet const *>( &* jptjetRef);
1144
1145 lCalo.maxEInEmTowers = rawcalojet->maxEInEmTowers();
1146 lCalo.maxEInHadTowers = rawcalojet->maxEInHadTowers();
1147 lCalo.energyFractionHadronic = rawcalojet->energyFractionHadronic();
1148 lCalo.emEnergyFraction = rawcalojet->emEnergyFraction();
1149 lCalo.hadEnergyInHB = rawcalojet->hadEnergyInHB();
1150 lCalo.hadEnergyInHO = rawcalojet->hadEnergyInHO();
1151 lCalo.hadEnergyInHE = rawcalojet->hadEnergyInHE();
1152 lCalo.hadEnergyInHF = rawcalojet->hadEnergyInHF();
1153 lCalo.emEnergyInEB = rawcalojet->emEnergyInEB();
1154 lCalo.emEnergyInEE = rawcalojet->emEnergyInEE();
1155 lCalo.emEnergyInHF = rawcalojet->emEnergyInHF();
1156 lCalo.towersArea = rawcalojet->towersArea();
1157 lCalo.n90 = rawcalojet->n90();
1158 lCalo.n60 = rawcalojet->n60();
1159
1160 lCommon.etaMean = rawcalojet->etaPhiStatistics().etaMean;
1161 lCommon.phiMean = rawcalojet->etaPhiStatistics().phiMean;
1162 lCommon.etaEtaMoment = rawcalojet->etaPhiStatistics().etaEtaMoment;
1163 lCommon.phiPhiMoment = rawcalojet->etaPhiStatistics().phiPhiMoment;
1164 lCommon.etaPhiMoment = rawcalojet->etaPhiStatistics().etaPhiMoment;
1165
1166 HbbAnalysis::Jet lObj(lGen,lGenJet,lReco,lCommon,lCalo,lJPT,lJPTPF,lBtag,lId);
1167 aVec.push_back(lObj);
1168 }
1169 if (isPF) {
1170 HbbAnalysis::PFJetVars lPF;
1171 lPF.chargedMuEnergy = (*iter).chargedMuEnergy();
1172 lPF.chargedMuEnergyFraction = (*iter).chargedMuEnergyFraction();
1173 lPF.neutralMultiplicity = (*iter).neutralMultiplicity();
1174
1175 HbbAnalysis::Jet lObj(lGen,lGenJet,lReco,lCommon,lPF,lJPTPF,lBtag);
1176 aVec.push_back(lObj);
1177 }
1178
1179 iEle++;
1180 }//loop on element
1181 }//non empty
1182 }//valid
1183
1184 }//HbbJets
1185
1186 void HbbTreeMaker::HbbMet(const edm::Handle<std::vector<pat::MET> > & aCol,
1187 HbbAnalysis::Met & aMet)
1188 {//HbbMet
1189 HbbAnalysis::MetVars lGen;
1190 HbbAnalysis::MetVars lReco;
1191
1192 assert(aCol->size() == 1);
1193 const pat::MET & lMet = *(aCol->begin());
1194
1195 lReco.mET = lMet.pt();
1196 lReco.mEx = lMet.px();
1197 lReco.mEy = lMet.py();
1198 lReco.sumET = lMet.sumEt();
1199 lReco.phi = lMet.phi();
1200 lReco.mEtSig = lMet.mEtSig();
1201
1202 const reco::GenMET *lGenMET = lMet.genMET();
1203
1204 if (!processData_ && lGenMET){
1205 lGen.mET = lGenMET->pt();
1206 lGen.mEx = lGenMET->px();
1207 lGen.mEy = lGenMET->py();
1208 lGen.sumET = lGenMET->sumEt();
1209 lGen.phi = lGenMET->phi();
1210 lGen.mEtSig = lGenMET->mEtSig();
1211 }
1212 else {
1213 lGen.mET = 0;
1214 lGen.mEx = 0;
1215 lGen.mEy = 0;
1216 lGen.sumET = 0;
1217 lGen.phi = 0;
1218 lGen.mEtSig = 0;
1219 }
1220
1221 aMet.genVars(lGen);
1222 aMet.recoVars(lReco);
1223
1224 }//HbbMet
1225
1226 void HbbTreeMaker::HbbTrigger(const edm::Handle<edm::TriggerResults> & aCol,
1227 const edm::TriggerNames & aNames,
1228 std::vector<HbbAnalysis::Trigger> & aVec)
1229 {//HbbTrigger
1230
1231 if (aCol.isValid()){
1232 //edm::TriggerNames lNames;
1233 //lNames.init(*aCol);
1234
1235 //for (unsigned int j=0; j<hltPaths_.size(); j++) {
1236 //bool valid=false;
1237 for (unsigned int i=0; i<aCol->size(); i++){
1238 std::string trueName = aNames.triggerName(i);//.at(i);
1239
1240 HbbAnalysis::TriggerVars lTrigVars;
1241 lTrigVars.name = trueName;
1242 lTrigVars.index = i;
1243 lTrigVars.accept = aCol->accept(i);
1244
1245 HbbAnalysis::Trigger lObj(lTrigVars);
1246
1247 aVec.push_back(lObj);
1248 }
1249
1250 }
1251
1252 }
1253
1254 void HbbTreeMaker::HbbParticles(const edm::Handle<reco::GenParticleCollection> & aCol,
1255 std::vector<HbbAnalysis::GenParticle> & aVec)
1256 {//genparticles
1257
1258 //add genParticle inside vector
1259 for (unsigned int mccount = 0;mccount<aCol->size();mccount++)
1260 {//loop on particles
1261 const reco::Candidate & p = (*aCol)[mccount];
1262
1263 HbbAnalysis::MCVars lMC;
1264 lMC.index = mccount;
1265 lMC.E = p.energy();
1266 lMC.pT = p.pt();
1267 lMC.eta = p.eta();
1268 lMC.y = p.y();
1269 lMC.phi = p.phi();
1270 lMC.pdgId = p.pdgId();
1271 lMC.status = p.status();
1272 HbbAnalysis::GenParticle lObj(lMC);
1273 aVec.push_back(lObj);
1274
1275 }//loop on particles
1276
1277 //now need to get the list of parents....
1278
1279 for (unsigned int mccount = 0;mccount<aCol->size();mccount++)
1280 {//loop on particles
1281 const reco::Candidate & p = (*aCol)[mccount];
1282
1283 unsigned int nD = p.numberOfDaughters();
1284 //std::cout << mccount << " has " << nD << " daughter(s) : " ;
1285 for(unsigned int dau=0;dau<nD;dau++){//loop on daughters
1286 const reco::Candidate * pDau = p.daughter(dau);
1287 //std::cout << pDau << " ("<<std::flush;
1288 for (unsigned int mccount2 = 0;mccount2<aCol->size();mccount2++)
1289 {//loop on particles
1290 const reco::Candidate & p2 = (*aCol)[mccount2];
1291 //std::cout<<"DBG: mccount2 = "<<mccount2<<" gen size = "<<data_->gen().size()<<std::endl;
1292
1293 HbbAnalysis::GenParticle & part2 = aVec.at(mccount2);
1294 if (pDau == &p2) {
1295 part2.setParent(mccount);
1296 //std::cout << &p2 << ", index = " << mccount2 << "), "<<std::flush;
1297 break;
1298 }
1299 //else std::cout << std::endl << "****no match " << mccount2 << " " << &p2 << std::endl;
1300 //if (p2.pdgId() == pDau->pdgId() && fabs(p2.energy()-pDau->energy())<0.000001 && fabs(p2.eta()-pDau->eta())<0.000001 && fabs(p2.phi()-pDau->phi())<0.000001){
1301 //part2.parent(mccount);
1302 //}
1303 }//loop on particles
1304 }//loop on daughters
1305 //std::cout << std::endl;
1306
1307 }//loop on particles
1308
1309 if (debug_){
1310 for (unsigned int imc(0); imc < aVec.size(); imc++){
1311 HbbAnalysis::GenParticle & p = aVec.at(imc);
1312 p.print();
1313 }
1314 }
1315
1316
1317 }//genparticles
1318
1319
1320
1321 void HbbTreeMaker::HbbVertices(const edm::Handle<std::vector<reco::Vertex> > & aCol,
1322 std::vector<HbbAnalysis::Vertex> & aVec)
1323 {
1324
1325 for (std::vector<reco::Vertex>::const_iterator iter = aCol->begin();
1326 iter != aCol->end();
1327 iter++)
1328 {
1329
1330 if (!((*iter).isValid()) || (*iter).isFake()) continue;
1331
1332 HbbAnalysis::VertexVars lVtx;
1333
1334 if ((*iter).tracksSize() > 0) {
1335 for (std::vector<reco::TrackBaseRef>::const_iterator lTrk = (*iter).tracks_begin();
1336 lTrk != (*iter).tracks_end();
1337 lTrk++) {
1338 lVtx.trackWeights.push_back((*iter).trackWeight(*lTrk));
1339 }
1340 }
1341
1342 if (lVtx.trackWeights.size() != (*iter).tracksSize()) {
1343 std::cout<< " -- Problem with tracks, size is not as expected ! "
1344 << std::endl
1345 << " --- Size of recoVertex tracks = " << (*iter).tracksSize()
1346 << std::endl
1347 << " --- Size of trackWeights = " << lVtx.trackWeights.size()
1348 << std::endl;
1349 }
1350
1351 lVtx.chi2 = (*iter).chi2();
1352 lVtx.ndof = (*iter).ndof();
1353 lVtx.x = (*iter).x();
1354 lVtx.y = (*iter).y();
1355 lVtx.z = (*iter).z();
1356 lVtx.xError = (*iter).xError();
1357 lVtx.yError = (*iter).yError();
1358 lVtx.zError = (*iter).zError();
1359 lVtx.cov01 = (*iter).covariance(0,1);
1360 lVtx.cov02 = (*iter).covariance(0,2);
1361 lVtx.cov12 = (*iter).covariance(1,2);
1362
1363 HbbAnalysis::Vertex lObj(lVtx);
1364 aVec.push_back(lObj);
1365 }
1366
1367 }
1368
1369 void HbbTreeMaker::HbbL1Objects(const edm::Handle<l1extra::L1JetParticleCollection> & aCol,
1370 std::vector<HbbAnalysis::L1Object> & aVec,
1371 const unsigned int aType)
1372 {
1373
1374 for (l1extra::L1JetParticleCollection::const_iterator iter = aCol->begin();
1375 iter != aCol->end();
1376 iter++)
1377 {
1378
1379 HbbAnalysis::L1Vars lVars;
1380 lVars.pT = (*iter).pt();
1381 lVars.ET = (*iter).et();
1382 lVars.eta = (*iter).eta();
1383 lVars.phi = (*iter).phi();
1384 lVars.bx = (*iter).bx();
1385 lVars.type = aType;
1386
1387 HbbAnalysis::L1Object lObj(lVars);
1388 aVec.push_back(lObj);
1389
1390 }
1391
1392 }
1393 void HbbTreeMaker::HbbHLTObjects(const edm::Handle<trigger::TriggerEvent> & aCol,
1394 std::vector<HbbAnalysis::HLTObject> & aVec)
1395 {
1396
1397 const trigger::TriggerObjectCollection & toc(aCol->getObjects());
1398
1399 //cout
1400 //for (unsigned i=0; i<aCol->sizeFilters(); i++){
1401 //const edm::InputTag & lTag = aCol->filterTag(i);
1402 //std::cout << "ele " << i << ": " << lTag << std::endl;
1403 //}
1404
1405 for (unsigned i(0); i<hltTags_.size(); i++){
1406
1407 const int index = aCol->filterIndex(hltTags_[i]);
1408 if (index >= aCol->sizeFilters()) {
1409 //std::cout << " -- Filter " << hltTags_[i] << " not found !" << std::endl;
1410 continue;
1411 }
1412
1413 const trigger::Keys & k = aCol->filterKeys(index);
1414
1415 for (trigger::Keys::const_iterator ki = k.begin(); ki !=k.end(); ++ki ) {
1416 HbbAnalysis::HLTVars lVars;
1417 lVars.pT = toc[*ki].pt();
1418 lVars.eta = toc[*ki].eta();
1419 lVars.phi = toc[*ki].phi();
1420 lVars.id = toc[*ki].id();
1421 lVars.mass = toc[*ki].mass();
1422 HbbAnalysis::HLTObject lObj(lVars);
1423 aVec.push_back(lObj);
1424 }
1425 }
1426
1427 }
1428
1429
1430 #include "FWCore/Framework/interface/MakerMacros.h"
1431 DEFINE_FWK_MODULE(HbbTreeMaker);
1432
1433