ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/HbbAnalysis/plugins/HbbTreeMaker.cc
Revision: 1.2
Committed: Fri Oct 2 11:05:53 2009 UTC (15 years, 7 months ago) by amagnan
Content type: text/plain
Branch: MAIN
Changes since 1.1: +14 -1 lines
Log Message:
add histos classes to fill from tree

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/VertexReco/interface/Vertex.h"
10 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
11 #include "DataFormats/TrackReco/interface/TrackFwd.h"
12
13 #include "DataFormats/PatCandidates/interface/Lepton.h"
14 #include "DataFormats/PatCandidates/interface/Muon.h"
15 #include "DataFormats/PatCandidates/interface/Electron.h"
16 #include "DataFormats/PatCandidates/interface/Tau.h"
17 #include "DataFormats/PatCandidates/interface/Jet.h"
18 #include "DataFormats/PatCandidates/interface/MET.h"
19
20 #include "FWCore/ServiceRegistry/interface/Service.h"
21 #include "PhysicsTools/UtilAlgos/interface/TFileService.h"
22
23
24 #include "UserCode/HbbAnalysis/interface/Electron.hh"
25 #include "UserCode/HbbAnalysis/interface/Muon.hh"
26 #include "UserCode/HbbAnalysis/interface/Tau.hh"
27 #include "UserCode/HbbAnalysis/interface/Jet.hh"
28 #include "UserCode/HbbAnalysis/interface/Met.hh"
29 #include "UserCode/HbbAnalysis/interface/Trigger.hh"
30
31 #include "UserCode/HbbAnalysis/plugins/HbbTreeMaker.hh"
32
33 using namespace HbbAnalysis;
34
35 HbbTreeMaker::HbbTreeMaker(const edm::ParameterSet & pset):
36 debug_(pset.getParameter<int>("DEBUG")),
37 flavour_(pset.getParameter<unsigned int>("JetFlavour")),
38 doGen_(pset.getParameter<bool>("DOGEN")),
39 genParticleSrc_(pset.getParameter<edm::InputTag>("GenParticles")),
40 electronSrc_(pset.getParameter<edm::InputTag>("Electrons")),
41 muonSrc_(pset.getParameter<edm::InputTag>("Muons")),
42 caloTauSrc_(pset.getParameter<edm::InputTag>("CaloTaus")),
43 pfTauSrc_(pset.getParameter<edm::InputTag>("PFTaus")),
44 caloJetSrc_(pset.getParameter<edm::InputTag>("CaloJets")),
45 jptJetSrc_(pset.getParameter<edm::InputTag>("JPTJets")),
46 pfJetSrc_(pset.getParameter<edm::InputTag>("PFJets")),
47 caloMetSrc_(pset.getParameter<edm::InputTag>("CaloMET")),
48 tcMetSrc_(pset.getParameter<edm::InputTag>("TCMET")),
49 pfMetSrc_(pset.getParameter<edm::InputTag>("PFMET")),
50 //pairSrc_(pset.getParameter<edm::InputTag>("Pair")),
51 //mmPairSrc_(pset.getParameter<edm::InputTag>("MuMuPair")),
52 //etPairSrc_(pset.getParameter<edm::InputTag>("ETauPair")),
53 //mtPairSrc_(pset.getParameter<edm::InputTag>("MuTauPair")),
54 vertexSrc_(pset.getParameter<edm::InputTag>("Vertex")),
55 triggerSrc_(pset.getParameter<edm::InputTag>("Trigger")),
56 hltPaths_(pset.getParameter<std::vector<std::string> >("HLTPaths")),
57 event_(0),
58 tree_(0)
59 //processVec_(pset.getParameter<std::vector<std::string> >("ProcessVec"))
60 {//constructor
61
62
63 }//constructor
64
65 HbbTreeMaker::~HbbTreeMaker(){//destructor
66 }//destructor
67
68
69
70 void HbbTreeMaker::beginJob(const edm::EventSetup&){//beginJob
71
72
73 edm::Service<TFileService> lFileService;
74
75 TFileDirectory lDir = lFileService->mkdir("HbbAnalysis");
76
77 tree_ = lDir.make<TTree>("Tree","Tree");
78 tree_->Branch("HbbEvent","HbbAnalysis::HbbEvent",&event_);
79 event_ = new HbbEvent();
80
81 if (debug_) std::cout << "Initialising JetFlavour : " << std::endl;
82
83 lDir = lFileService->mkdir("JetFlavours");
84
85 jetFlav_.Initialise(lDir, debug_, flavour_);
86
87 }//beginJob
88
89 void HbbTreeMaker::endJob(){//endJob
90 jetFlav_.printSummary();
91
92 //tree_->Write();
93 //delete tree_;
94 //delete event_;
95
96 }//endJob
97
98 void HbbTreeMaker::analyze(const edm::Event& aEvt, const edm::EventSetup& aEvtSetup){//analyze
99
100 //dumpContent(aEvt,"","MuTauPAT");
101 if (debug_) std::cout << "Processing event #" << aEvt.id().event() << std::endl;
102
103
104 static bool firstEvent = true;
105
106 event_->Clear();
107 event_->event(aEvt.id().event());
108
109 edm::Handle<reco::GenParticleCollection> lGenParticles;
110 try {
111 aEvt.getByLabel(genParticleSrc_,lGenParticles);
112 if (debug_) std::cout << "** ngenParticles = " << lGenParticles->size() << std::endl;
113 } catch(cms::Exception& e) {
114 std::cout << "AMM: Collection genParticles not available! Exception : " << e.what() << ". " << std::endl;
115 }
116
117 if (doGen_) HbbParticles(lGenParticles,event_->particles());
118
119 unsigned int lNPartons = jetFlav_.fillPartons(lGenParticles);
120
121 if (debug_) std::cout << "--- Number of partons = " << lNPartons << "." << std::endl;
122
123
124 edm::Handle<std::vector<reco::Vertex> > lRecoVertices;
125 try {
126 aEvt.getByLabel(vertexSrc_, lRecoVertices);
127 if (!lRecoVertices.isValid()){
128 edm::LogInfo("ERROR")<< "Error! Can't get vertex by label. ";
129 }
130 if (debug_) std::cout << "** vertexcollection = " << lRecoVertices->size() << " elements." << std::endl;
131 } catch(cms::Exception& e) {
132 std::cout << "AMM: Collection " << vertexSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
133 }
134
135 edm::Handle<std::vector<pat::Electron> > lElectronCollection;
136
137 try {
138 aEvt.getByLabel(electronSrc_,lElectronCollection);
139 if (!lElectronCollection.isValid()){
140 edm::LogInfo("ERROR")<< "Error! Can't get electron by label. ";
141 }
142 if (debug_) std::cout << "** electroncollection = " << lElectronCollection->size() << " elements." << std::endl;
143 } catch(cms::Exception& e) {
144 std::cout << "AMM: Collection " << electronSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
145 }
146
147 HbbElectrons(lElectronCollection,event_->electrons());
148
149
150 edm::Handle<std::vector<pat::Muon> > lMuonCollection;
151
152 try {
153 aEvt.getByLabel(muonSrc_,lMuonCollection);
154 if (!lMuonCollection.isValid()){
155 edm::LogInfo("ERROR")<< "Error! Can't get muon by label. ";
156 }
157 if (debug_) std::cout << "** muoncollection = " << lMuonCollection->size() << " elements." << std::endl;
158 } catch(cms::Exception& e) {
159 std::cout << "AMM: Collection " << muonSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
160 }
161
162 HbbMuons(lMuonCollection,event_->muons());
163
164
165 edm::Handle<std::vector<pat::Tau> > lTauCollection;
166
167 try {
168 aEvt.getByLabel(caloTauSrc_,lTauCollection);
169 if (!lTauCollection.isValid()){
170 edm::LogInfo("ERROR")<< "Error! Can't get caloTau by label. ";
171 }
172 if (debug_) std::cout << "** caloTaucollection = " << lTauCollection->size() << " elements." << std::endl;
173 } catch(cms::Exception& e) {
174 std::cout << "AMM: Collection " << caloTauSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
175 }
176
177 HbbTaus(lTauCollection,lRecoVertices,event_->caloTaus());
178
179 edm::Handle<std::vector<pat::Tau> > lPFTauCollection;
180
181 try {
182 aEvt.getByLabel(pfTauSrc_,lPFTauCollection);
183 if (!lPFTauCollection.isValid()){
184 edm::LogInfo("ERROR")<< "Error! Can't get pfTau by label. ";
185 }
186 if (debug_) std::cout << "** pfTaucollection = " << lPFTauCollection->size() << " elements." << std::endl;
187 } catch(cms::Exception& e) {
188 std::cout << "AMM: Collection " << pfTauSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
189 }
190
191 HbbTaus(lPFTauCollection,lRecoVertices,event_->pfTaus());
192
193 edm::Handle<std::vector<pat::Jet> > lCaloJetCollection;
194
195 try {
196 aEvt.getByLabel(caloJetSrc_,lCaloJetCollection);
197 if (!lCaloJetCollection.isValid()){
198 edm::LogInfo("ERROR")<< "Error! Can't get caloJets by label. ";
199 }
200 if (debug_) std::cout << "** caloJetcollection = " << lCaloJetCollection->size() << " elements." << std::endl;
201 } catch(cms::Exception& e) {
202 std::cout << "AMM: Collection " << caloJetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
203 }
204
205 HbbJets(lCaloJetCollection,jetFlav_,lGenParticles,event_->caloJets());
206
207 edm::Handle<std::vector<pat::Jet> > lJptJetCollection;
208
209 try {
210 aEvt.getByLabel(jptJetSrc_,lJptJetCollection);
211 if (!lJptJetCollection.isValid()){
212 edm::LogInfo("ERROR")<< "Error! Can't get jptJets by label. ";
213 }
214 if (debug_) std::cout << "** jptJetcollection = " << lJptJetCollection->size() << " elements." << std::endl;
215 } catch(cms::Exception& e) {
216 std::cout << "AMM: Collection " << jptJetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
217 }
218
219 HbbJets(lJptJetCollection,jetFlav_,lGenParticles,event_->jptJets());
220
221 edm::Handle<std::vector<pat::Jet> > lPfJetCollection;
222
223 try {
224 aEvt.getByLabel(pfJetSrc_,lPfJetCollection);
225 if (!lPfJetCollection.isValid()){
226 edm::LogInfo("ERROR")<< "Error! Can't get pfJets by label. ";
227 }
228 if (debug_) std::cout << "** pfJetcollection = " << lPfJetCollection->size() << " elements." << std::endl;
229 } catch(cms::Exception& e) {
230 std::cout << "AMM: Collection " << pfJetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
231 }
232
233 HbbJets(lPfJetCollection,jetFlav_,lGenParticles,event_->pfJets());
234
235 edm::Handle<std::vector<pat::MET> > lCaloMetCol;
236
237 try {
238 aEvt.getByLabel(caloMetSrc_,lCaloMetCol);
239 if (!lCaloMetCol.isValid()){
240 edm::LogInfo("ERROR")<< "Error! Can't get caloMet by label. ";
241 }
242 if (debug_) std::cout << "** caloMet = " << lCaloMetCol->size() << " elements." << std::endl;
243 } catch(cms::Exception& e) {
244 std::cout << "AMM: Collection " << caloMetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
245 }
246
247 HbbMet(lCaloMetCol,event_->caloMet());
248
249 edm::Handle<std::vector<pat::MET> > lTcMetCol;
250
251 try {
252 aEvt.getByLabel(tcMetSrc_,lTcMetCol);
253 if (!lTcMetCol.isValid()){
254 edm::LogInfo("ERROR")<< "Error! Can't get tcMet by label. ";
255 }
256 if (debug_) std::cout << "** tcMet = " << lTcMetCol->size() << " elements." << std::endl;
257 } catch(cms::Exception& e) {
258 std::cout << "AMM: Collection " << tcMetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
259 }
260
261 HbbMet(lTcMetCol,event_->tcMet());
262
263 edm::Handle<std::vector<pat::MET> > lPfMetCol;
264
265 try {
266 aEvt.getByLabel(pfMetSrc_,lPfMetCol);
267 if (!lPfMetCol.isValid()){
268 edm::LogInfo("ERROR")<< "Error! Can't get pfMet by label. ";
269 }
270 if (debug_) std::cout << "** pfMet = " << lPfMetCol->size() << " elements." << std::endl;
271 } catch(cms::Exception& e) {
272 std::cout << "AMM: Collection " << pfMetSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
273 }
274
275 HbbMet(lPfMetCol,event_->pfMet());
276
277
278
279 //triggers
280 edm::Handle<edm::TriggerResults> lTrigCol;
281 try {
282 aEvt.getByLabel(triggerSrc_,lTrigCol);
283 if (!lTrigCol.isValid()){
284 edm::LogInfo("ERROR")<< "Error! Can't get triggers by label. ";
285 }
286 if (debug_) std::cout << "** triggers = " << lTrigCol->size() << std::endl;
287 } catch(cms::Exception& e) {
288 std::cout << "AMM: Collection " << triggerSrc_ << " not available! Exception : " << e.what() << ". " << std::endl;
289 }
290
291 HbbTrigger(lTrigCol,event_->triggers());
292
293 //event_->order();
294 tree_->Fill();
295 firstEvent = false;
296
297 }//analyze
298
299
300 void HbbTreeMaker::HbbElectrons(const edm::Handle<std::vector<pat::Electron> > & aCol,
301 std::vector<HbbAnalysis::Electron> & aVec)
302 {//HbbElectrons
303
304 if (aCol.isValid()){
305 if (aCol->size() > 0) {
306
307 unsigned int iEle = 0;
308 for (std::vector<pat::Electron>::const_iterator iter = aCol->begin();
309 iter != aCol->end();
310 iter++)
311 {
312 //if (debug_ > 1) std::cout << "**** Ele #" << iEle << ", pT,eta,phi = " << (*iter).pt() << " " << (*iter).eta() << " " << (*iter).phi() << std::endl;
313
314 HbbAnalysis::GenVars lGen;
315 if ((*iter).genLepton()){
316 lGen.valid = true;
317 lGen.E = (*iter).genLepton()->energy();
318 lGen.pT = (*iter).genLepton()->pt();
319 lGen.eta = (*iter).genLepton()->eta();
320 lGen.phi = (*iter).genLepton()->phi();
321 lGen.charge = (*iter).genLepton()->charge();
322 lGen.pdgId = (*iter).genLepton()->pdgId();
323 lGen.status = (*iter).genLepton()->status();
324 lGen.mass = (*iter).genLepton()->mass();
325 lGen.vx = (*iter).genLepton()->vx();
326 lGen.vy = (*iter).genLepton()->vy();
327 lGen.vz = (*iter).genLepton()->vz();
328 }
329 else {
330 lGen.valid = false;
331 lGen.E = 0;
332 lGen.pT = 0;
333 lGen.eta = 0;
334 lGen.phi = 0;
335 lGen.charge = 0;
336 lGen.pdgId = 0;
337 lGen.status = 0;
338 lGen.mass = 0;
339 lGen.vx = 0;
340 lGen.vy = 0;
341 lGen.vz = 0;
342 }
343
344 HbbAnalysis::BaseVars lReco;
345 lReco.E = (*iter).energy();
346 lReco.pT = (*iter).pt();
347 lReco.eta = (*iter).eta();
348 lReco.phi = (*iter).phi();
349 lReco.charge = (*iter).charge();
350 lReco.vx = (*iter).vx();
351 lReco.vy = (*iter).vy();
352 lReco.vz = (*iter).vz();
353
354 HbbAnalysis::SCVars lSC;
355 lSC.sigmaEtaEta = (*iter).scSigmaEtaEta();
356 lSC.sigmaIEtaIEta = (*iter).scSigmaIEtaIEta();
357 lSC.e1x5 = (*iter).scE1x5();
358 lSC.e2x5Max = (*iter).scE2x5Max();
359 lSC.e5x5 = (*iter).scE5x5();
360 lSC.eOverP = (*iter).eSuperClusterOverP();
361
362 HbbAnalysis::EleIsoVars lIso;
363 lIso.calo = (*iter).caloIso();
364 lIso.track = (*iter).trackIso();
365 lIso.ecal = (*iter).ecalIso();
366 lIso.hcal = (*iter).hcalIso();
367
368 HbbAnalysis::EleIDVars lID;
369 lID.electronIDs = (*iter).electronIDs();
370 lID.hOverE = (*iter).hadronicOverEm();
371 lID.deltaPhiIn = (*iter).deltaPhiSuperClusterTrackAtVtx();
372 lID.deltaEtaIn = (*iter).deltaEtaSuperClusterTrackAtVtx();
373
374 HbbAnalysis::Electron lObj(lGen,lReco,lSC,lIso,lID);
375 aVec.push_back(lObj);
376 iEle++;
377 }
378 }
379 }
380
381 }//HbbElectrons
382
383
384 void HbbTreeMaker::HbbMuons(const edm::Handle<std::vector<pat::Muon> > & aCol,
385 std::vector<HbbAnalysis::Muon> & aVec)
386 {//HbbMuons
387
388 if (aCol.isValid()){
389 if (aCol->size() > 0) {
390
391 unsigned int iEle = 0;
392 for (std::vector<pat::Muon>::const_iterator iter = aCol->begin();
393 iter != aCol->end();
394 iter++)
395 {
396 const reco::Muon* recoMuon = dynamic_cast<const reco::Muon*>((*iter).originalObject());
397
398 HbbAnalysis::GenVars lGen;
399 if ((*iter).genLepton()){
400 lGen.valid = true;
401 lGen.E = (*iter).genLepton()->energy();
402 lGen.pT = (*iter).genLepton()->pt();
403 lGen.eta = (*iter).genLepton()->eta();
404 lGen.phi = (*iter).genLepton()->phi();
405 lGen.charge = (*iter).genLepton()->charge();
406 lGen.pdgId = (*iter).genLepton()->pdgId();
407 lGen.status = (*iter).genLepton()->status();
408 lGen.mass = (*iter).genLepton()->mass();
409 lGen.vx = (*iter).genLepton()->vx();
410 lGen.vy = (*iter).genLepton()->vy();
411 lGen.vz = (*iter).genLepton()->vz();
412 }
413 else {
414 lGen.valid = false;
415 lGen.E = 0;
416 lGen.pT = 0;
417 lGen.eta = 0;
418 lGen.phi = 0;
419 lGen.charge = 0;
420 lGen.pdgId = 0;
421 lGen.status = 0;
422 lGen.mass = 0;
423 lGen.vx = 0;
424 lGen.vy = 0;
425 lGen.vz = 0;
426 }
427
428 HbbAnalysis::BaseVars lReco;
429 lReco.E = (*iter).energy();
430 lReco.pT = (*iter).pt();
431 lReco.eta = (*iter).eta();
432 lReco.phi = (*iter).phi();
433 lReco.charge = (*iter).charge();
434 lReco.vx = (*iter).vx();
435 lReco.vy = (*iter).vy();
436 lReco.vz = (*iter).vz();
437
438 HbbAnalysis::MuIsoVars lIsoR03;
439 lIsoR03.sumPt = recoMuon->isolationR03().sumPt;
440 lIsoR03.emEt = recoMuon->isolationR03().emEt;
441 lIsoR03.hadEt = recoMuon->isolationR03().hadEt;
442 lIsoR03.nTracks = recoMuon->isolationR03().nTracks;
443 lIsoR03.nJets = recoMuon->isolationR03().nJets;
444
445 HbbAnalysis::MuIsoVars lIsoR05;
446 lIsoR05.sumPt = recoMuon->isolationR05().sumPt;
447 lIsoR05.emEt = recoMuon->isolationR05().emEt;
448 lIsoR05.hadEt = recoMuon->isolationR05().hadEt;
449 lIsoR05.nTracks = recoMuon->isolationR05().nTracks;
450 lIsoR05.nJets = recoMuon->isolationR05().nJets;
451
452 HbbAnalysis::MuIDVars lID;
453 if (recoMuon->isGlobalMuon()) lID.type = 1;
454 else if (recoMuon->isTrackerMuon()) lID.type = 2;
455 else if (recoMuon->isStandAloneMuon()) lID.type = 3;
456 else if (recoMuon->isCaloMuon()) lID.type = 4;
457 else lID.type = 0;
458
459 if (recoMuon->isGood(reco::Muon::AllGlobalMuons)) lID.ids.push_back(1);
460 if (recoMuon->isGood(reco::Muon::AllStandAloneMuons)) lID.ids.push_back(2);
461 if (recoMuon->isGood(reco::Muon::AllTrackerMuons)) lID.ids.push_back(3);
462 if (recoMuon->isGood(reco::Muon::TrackerMuonArbitrated)) lID.ids.push_back(4);
463 if (recoMuon->isGood(reco::Muon::AllArbitrated)) lID.ids.push_back(5);
464 if (recoMuon->isGood(reco::Muon::GlobalMuonPromptTight)) lID.ids.push_back(6);
465 if (recoMuon->isGood(reco::Muon::TMLastStationLoose)) lID.ids.push_back(7);
466 if (recoMuon->isGood(reco::Muon::TMLastStationTight)) lID.ids.push_back(8);
467 if (recoMuon->isGood(reco::Muon::TM2DCompatibilityLoose)) lID.ids.push_back(9);
468 if (recoMuon->isGood(reco::Muon::TM2DCompatibilityTight)) lID.ids.push_back(10);
469 if (recoMuon->isGood(reco::Muon::TMOneStationLoose)) lID.ids.push_back(11);
470 if (recoMuon->isGood(reco::Muon::TMOneStationTight)) lID.ids.push_back(12);
471 if (recoMuon->isGood(reco::Muon::TMLastStationOptimizedLowPtLoose)) lID.ids.push_back(13);
472 if (recoMuon->isGood(reco::Muon::TMLastStationOptimizedLowPtTight)) lID.ids.push_back(14);
473
474 lID.caloCompat = recoMuon->caloCompatibility();
475 lID.segCompat = recoMuon->segmentCompatibility();
476 lID.nChambers = recoMuon->numberOfChambers();
477 lID.nMatchesLoose = recoMuon->numberOfMatches(reco::Muon::NoArbitration);
478 lID.nMatchesMedium = recoMuon->numberOfMatches(reco::Muon::SegmentArbitration);
479 lID.nMatchesTight = recoMuon->numberOfMatches(reco::Muon::SegmentAndTrackArbitration);
480
481 HbbAnalysis::Muon lObj(lGen,lReco,lIsoR03,lIsoR05,lID);
482 aVec.push_back(lObj);
483 iEle++;
484 }
485 }
486 }
487
488 }//HbbMuons
489
490 void HbbTreeMaker::HbbTaus(const edm::Handle<std::vector<pat::Tau> > & aCol,
491 const edm::Handle<std::vector<reco::Vertex> > & aRecoVertices,
492 std::vector<HbbAnalysis::Tau> & aVec)
493 {//HbbTaus
494
495 if (aCol.isValid()){
496 if (aCol->size() > 0) {
497
498 unsigned int iEle = 0;
499 for (std::vector<pat::Tau>::const_iterator iter = aCol->begin();
500 iter != aCol->end();
501 iter++)
502 {
503 //if (debug_ > 1) std::cout << "**** Ele #" << iEle << ", pT,eta,phi = " << (*iter).pt() << " " << (*iter).eta() << " " << (*iter).phi() << std::endl;
504
505 HbbAnalysis::GenVars lGen;
506 if ((*iter).genLepton()){
507 lGen.valid = true;
508 lGen.E = (*iter).genLepton()->energy();
509 lGen.pT = (*iter).genLepton()->pt();
510 lGen.eta = (*iter).genLepton()->eta();
511 lGen.phi = (*iter).genLepton()->phi();
512 lGen.charge = (*iter).genLepton()->charge();
513 lGen.pdgId = (*iter).genLepton()->pdgId();
514 lGen.status = (*iter).genLepton()->status();
515 lGen.mass = (*iter).genLepton()->mass();
516 lGen.vx = (*iter).genLepton()->vx();
517 lGen.vy = (*iter).genLepton()->vy();
518 lGen.vz = (*iter).genLepton()->vz();
519 }
520 else {
521 lGen.valid = false;
522 lGen.E = 0;
523 lGen.pT = 0;
524 lGen.eta = 0;
525 lGen.phi = 0;
526 lGen.charge = 0;
527 lGen.pdgId = 0;
528 lGen.status = 0;
529 lGen.mass = 0;
530 lGen.vx = 0;
531 lGen.vy = 0;
532 lGen.vz = 0;
533 }
534
535
536 HbbAnalysis::GenVars lGenJet;
537 if ((*iter).genJet()){
538 lGenJet.valid = true;
539 lGenJet.pT = (*iter).genJet()->pt();
540 lGenJet.eta = (*iter).genJet()->eta();
541 lGenJet.phi = (*iter).genJet()->phi();
542 lGenJet.charge = (*iter).genJet()->charge();
543 lGenJet.pdgId = (*iter).genJet()->pdgId();
544 lGenJet.status = (*iter).genJet()->status();
545 lGenJet.mass = (*iter).genJet()->mass();
546 lGenJet.vx = (*iter).genJet()->vx();
547 lGenJet.vy = (*iter).genJet()->vy();
548 lGenJet.vz = (*iter).genJet()->vz();
549 }
550 else {
551 lGenJet.valid = false;
552 lGenJet.pT = 0;
553 lGenJet.eta = 0;
554 lGenJet.phi = 0;
555 lGenJet.charge = 0;
556 lGenJet.pdgId = 0;
557 lGenJet.status = 0;
558 lGenJet.mass = 0;
559 lGenJet.vx = 0;
560 lGenJet.vy = 0;
561 lGenJet.vz = 0;
562 }
563
564 HbbAnalysis::BaseVars lReco;
565 lReco.E = (*iter).energy();
566 lReco.pT = (*iter).pt();
567 lReco.eta = (*iter).eta();
568 lReco.phi = (*iter).phi();
569 lReco.charge = (*iter).charge();
570 lReco.vx = (*iter).vx();
571 lReco.vy = (*iter).vy();
572 lReco.vz = (*iter).vz();
573
574 bool isCalo = (*iter).isCaloTau();
575 bool isPF = (*iter).isPFTau();
576
577 assert (isCalo == !isPF);
578
579 HbbAnalysis::TauLeadTrkVars lLead;
580 if (isCalo) {
581 lLead.signedSipt = (*iter).leadTracksignedSipt();
582
583 reco::TrackRef lCaloTrk = (*iter).leadTrack();
584
585 if ( lCaloTrk.isAvailable() && lCaloTrk.isNonnull() ) {
586 lLead.pT = lCaloTrk->pt();
587 lLead.eta = lCaloTrk->eta();
588 lLead.phi = lCaloTrk->phi();
589
590 lLead.matchDist = reco::deltaR(lCaloTrk->momentum(), (*iter).p4());
591
592 if ( aRecoVertices->size() >= 1 ) {
593 const reco::Vertex& thePrimaryEventVertex = (*aRecoVertices->begin());
594 lLead.IPxy = lCaloTrk->dxy(thePrimaryEventVertex.position());
595 lLead.IPz = lCaloTrk->dz(thePrimaryEventVertex.position());
596 }
597 else {
598 lLead.IPxy = 0;
599 lLead.IPz = 0;
600 }
601 }
602 else {
603 lLead.pT = 0;
604 lLead.eta = 0;
605 lLead.phi = 0;
606 lLead.matchDist = 0;
607 lLead.IPxy = 0;
608 lLead.IPz = 0;
609 lLead.signedSipt = 0;
610 }
611 }//calo
612 else {
613 lLead.signedSipt = (*iter).leadPFChargedHadrCandsignedSipt();
614
615 reco::PFCandidateRef lHadrCand = (*iter).leadPFChargedHadrCand();
616
617 if ( lHadrCand.isAvailable() && lHadrCand.isNonnull() ) {
618 lLead.pT = lHadrCand->pt();
619 lLead.eta = lHadrCand->eta();
620 lLead.phi = lHadrCand->phi();
621
622 lLead.matchDist = reco::deltaR(lHadrCand->p4(), (*iter).p4());
623
624 reco::TrackRef lTrk = lHadrCand->trackRef();
625
626 if ( aRecoVertices->size() >= 1) {
627 const reco::Vertex& thePrimaryEventVertex = (*aRecoVertices->begin());
628 lLead.IPxy = lTrk->dxy(thePrimaryEventVertex.position());
629 lLead.IPz = lTrk->dz(thePrimaryEventVertex.position());
630 }
631 else {
632 lLead.IPxy = 0;
633 lLead.IPz = 0;
634 }
635 }
636 else {
637 lLead.pT = 0;
638 lLead.eta = 0;
639 lLead.phi = 0;
640 lLead.matchDist = 0;
641 lLead.IPxy = 0;
642 lLead.IPz = 0;
643 lLead.signedSipt = 0;
644 }
645
646 }//pf
647
648
649
650 const std::vector<std::pair<std::string,float> > & lIDs = (*iter).tauIDs();
651 bool lPrint = false;
652 //a few warning if additional discriminants are found:
653 if ((isPF && lIDs.size() != 7) ||
654 (isCalo && lIDs.size() != 3))
655 lPrint = true;
656
657 if (lPrint) {
658 std::cout << "!!!!!!! Discriminants changed, please update histograms !!!!!!!" << std::endl;
659 std::cout << "--- isCaloTau = " << isCalo << ", isPFTau = " << isPF << std::endl;
660 std::cout << "------ ID names = " << std::endl;
661 }
662
663
664 HbbAnalysis::PFTauIDVars lpfId;
665 HbbAnalysis::CaloTauIDVars lcaloId;
666
667 for (unsigned int id(0); id<lIDs.size(); id++){
668
669 std::string lName = lIDs.at(id).first;
670 float lDiscri = lIDs.at(id).second;
671
672 if (lPrint) std::cout << "--------- " << lName << " = " << lDiscri << std::endl;
673 //pf
674
675 if (isPF) {
676 if (lName.find("leadingTrackFinding") != lName.npos) lpfId.byLeadingTrackFinding = lDiscri;
677 if (lName.find("leadingTrackPtCut") != lName.npos) lpfId.byLeadingTrackPtCut = lDiscri;
678 if (lName.find("trackIsolation") != lName.npos) lpfId.byTrackIsolation = lDiscri;
679 if (lName.find("ecalIsolation") != lName.npos) lpfId.byECALIsolation = lDiscri;
680 if (lName.find("byIsolation") != lName.npos) lpfId.byIsolation = lDiscri;
681 if (lName.find("againstElectron") != lName.npos) lpfId.againstElectron = lDiscri;
682 if (lName.find("againstMuon") != lName.npos) lpfId.againstMuon = lDiscri;
683 }
684 if (isCalo){
685 if (lName.find("byIsolation") != lName.npos) lcaloId.byIsolation = lDiscri;
686 if (lName.find("leadingTrackFinding") != lName.npos) lcaloId.byLeadingTrackFinding = lDiscri;
687 if (lName.find("leadingTrackPtCut") != lName.npos) lcaloId.byLeadingTrackPtCut = lDiscri;
688 }
689
690
691 }
692
693
694 if (isCalo){
695
696 HbbAnalysis::CaloTauIsoVars lcaloIso;
697 lcaloIso.nIsoTracks = (*iter).isolationTracks().size();
698 lcaloIso.nSigTracks = (*iter).signalTracks().size();
699 lcaloIso.leadTrackHCAL3x3hitsEtSum = (*iter).leadTrackHCAL3x3hitsEtSum();
700 lcaloIso.leadTrackHCAL3x3hottesthitDEta = (*iter).leadTrackHCAL3x3hottesthitDEta();
701 lcaloIso.signalTracksInvariantMass = (*iter).signalTracksInvariantMass();
702 lcaloIso.tracksInvariantMass = (*iter).TracksInvariantMass();
703 lcaloIso.isolationTracksPtSum = (*iter).isolationTracksPtSum();
704 lcaloIso.isolationECALhitsEtSum = (*iter).isolationECALhitsEtSum();
705 lcaloIso.maximumHCALhitEt = (*iter).maximumHCALhitEt();
706
707 HbbAnalysis::Tau lObj(lGen,lGenJet,lReco,lLead,lcaloIso,lcaloId);
708 aVec.push_back(lObj);
709 }
710 else {
711
712 HbbAnalysis::PFTauIsoVars lpfIso;
713 lpfIso.nSigCands = (*iter).signalPFCands().size();
714 lpfIso.nIsoCands = (*iter).isolationPFCands().size();
715 lpfIso.maximumHCALPFClusterEt = (*iter).maximumHCALPFClusterEt();
716 //lpfIso.emFraction = (*iter).emFraction();
717 lpfIso.hcalTotOverPLead = (*iter).hcalTotOverPLead();
718 lpfIso.hcalMaxOverPLead = (*iter).hcalMaxOverPLead();
719 lpfIso.hcal3x3OverPLead = (*iter).hcal3x3OverPLead();
720 lpfIso.ecalStripSumEOverPLead = (*iter).ecalStripSumEOverPLead();
721 lpfIso.bremsRecoveryEOverPLead = (*iter).bremsRecoveryEOverPLead();
722
723 HbbAnalysis::PFTauCandVars lHadr;
724 lHadr.nSigCands = (*iter).signalPFChargedHadrCands().size();
725 lHadr.nIsoCands = (*iter).isolationPFChargedHadrCands().size();
726 lHadr.isolationPtSum = (*iter).isolationPFChargedHadrCandsPtSum();
727
728 HbbAnalysis::PFTauCandVars lNeutr;
729 lNeutr.nSigCands = (*iter).signalPFNeutrHadrCands().size();
730 lNeutr.nIsoCands = (*iter).isolationPFNeutrHadrCands().size();
731 lNeutr.isolationPtSum = 0;
732
733 HbbAnalysis::PFTauCandVars lGamma;
734 lGamma.nSigCands = (*iter).signalPFGammaCands().size();
735 lGamma.nIsoCands = (*iter).isolationPFGammaCands().size();
736 lGamma.isolationPtSum = (*iter).isolationPFGammaCandsEtSum();
737
738 HbbAnalysis::PFTauEleIDVars lEleId;
739 if ( (*iter).electronPreIDTrack().isAvailable() && (*iter).electronPreIDTrack().isNonnull() ) {
740 lEleId.pT = (*iter).electronPreIDTrack()->pt();
741 lEleId.eta = (*iter).electronPreIDTrack()->eta();
742 lEleId.phi = (*iter).electronPreIDTrack()->phi();
743 }
744 else {
745 lEleId.pT = 0;
746 lEleId.eta = 0;
747 lEleId.phi = 0;
748 }
749
750 lEleId.output = (*iter).electronPreIDOutput();
751 lEleId.decision = (*iter).electronPreIDDecision();
752
753 HbbAnalysis::PFTauMuIDVars lMuId;
754 lMuId.caloCompat = (*iter).caloComp();
755 lMuId.segCompat = (*iter).segComp();
756 lMuId.decision = (*iter).muonDecision();
757
758 HbbAnalysis::Tau lObj(lGen,lGenJet,lReco,lLead,lpfIso,lHadr,lNeutr,lGamma,lpfId,lEleId,lMuId);
759 aVec.push_back(lObj);
760 }
761 iEle++;
762 }
763 }
764 }
765
766 }//HbbTaus
767
768 void HbbTreeMaker::HbbJets(const edm::Handle<std::vector<pat::Jet> > & aCol,
769 const HbbAnalysis::JetFlavour & aJetFlav,
770 const edm::Handle<reco::GenParticleCollection> & aGenParticles,
771 std::vector<HbbAnalysis::Jet> & aVec)
772 {//HbbJets
773
774 if (aCol.isValid()){//valid
775 if (aCol->size() > 0) {//non empty
776
777 unsigned int iEle = 0;
778 for (std::vector<pat::Jet>::const_iterator iter = aCol->begin();
779 iter != aCol->end();
780 iter++)
781 {//loop on element
782 //if (debug_ > 1) std::cout << "**** Ele #" << iEle << ", pT,eta,phi = " << (*iter).pt() << " " << (*iter).eta() << " " << (*iter).phi() << std::endl;
783
784 HbbAnalysis::GenVars lGen;
785 if ((*iter).genParton()){
786 lGen.valid = true;
787 lGen.E = (*iter).genParton()->energy();
788 lGen.pT = (*iter).genParton()->pt();
789 lGen.eta = (*iter).genParton()->eta();
790 lGen.phi = (*iter).genParton()->phi();
791 lGen.charge = (*iter).genParton()->charge();
792 lGen.pdgId = (*iter).genParton()->pdgId();
793 lGen.status = (*iter).genParton()->status();
794 lGen.mass = (*iter).genParton()->mass();
795 lGen.vx = (*iter).genParton()->vx();
796 lGen.vy = (*iter).genParton()->vy();
797 lGen.vz = (*iter).genParton()->vz();
798 }
799 else {
800 lGen.valid = false;
801 lGen.E = 0;
802 lGen.pT = 0;
803 lGen.eta = 0;
804 lGen.phi = 0;
805 lGen.charge = 0;
806 lGen.pdgId = 0;
807 lGen.status = 0;
808 lGen.mass = 0;
809 lGen.vx = 0;
810 lGen.vy = 0;
811 lGen.vz = 0;
812 }
813
814
815 HbbAnalysis::GenVars lGenJet;
816 if ((*iter).genJet()){
817 lGenJet.valid = true;
818 lGenJet.pT = (*iter).genJet()->pt();
819 lGenJet.eta = (*iter).genJet()->eta();
820 lGenJet.phi = (*iter).genJet()->phi();
821 lGenJet.charge = (*iter).genJet()->charge();
822 lGenJet.pdgId = (*iter).genJet()->pdgId();
823 lGenJet.status = (*iter).genJet()->status();
824 lGenJet.mass = (*iter).genJet()->mass();
825 lGenJet.vx = (*iter).genJet()->vx();
826 lGenJet.vy = (*iter).genJet()->vy();
827 lGenJet.vz = (*iter).genJet()->vz();
828 }
829 else {
830 lGenJet.valid = false;
831 lGenJet.pT = 0;
832 lGenJet.eta = 0;
833 lGenJet.phi = 0;
834 lGenJet.charge = 0;
835 lGenJet.pdgId = 0;
836 lGenJet.status = 0;
837 lGenJet.mass = 0;
838 lGenJet.vx = 0;
839 lGenJet.vy = 0;
840 lGenJet.vz = 0;
841 }
842
843 HbbAnalysis::BaseVars lReco;
844 lReco.E = (*iter).energy();
845 lReco.pT = (*iter).pt();
846 lReco.eta = (*iter).eta();
847 lReco.phi = (*iter).phi();
848 lReco.charge = (*iter).jetCharge();
849 lReco.vx = (*iter).vx();
850 lReco.vy = (*iter).vy();
851 lReco.vz = (*iter).vz();
852
853 unsigned int lFlav = 0;
854 if (aJetFlav.nPartons()) {
855 lFlav = aJetFlav.partonMatchingGenJet((*iter),aGenParticles,0.4).second;
856 }
857
858 HbbAnalysis::JetVars lCommon;
859 lCommon.flavour = lFlav;
860 lCommon.partonFlavour = (*iter).partonFlavour();
861 lCommon.nAssociatedTracks = ((*iter).associatedTracks()).size();
862 lCommon.rawpT = (*iter).pt();
863 if ((*iter).hasJetCorrFactors())
864 lCommon.rawpT = (*iter).pt()/((*iter).jetCorrFactors().scaleDefault());
865 //lCommon.rawEta = (*iter).;
866 //lCommon.rawPhi = (*iter).;
867 // p_hasJetCorrFactors->Fill(aJet.hasJetCorrFactors());
868
869 HbbAnalysis::JetBtagVars lBtag;
870 lBtag.cSV = (*iter).bDiscriminator("combinedSecondaryVertexBJetTags");
871 lBtag.cSVMVA = (*iter).bDiscriminator("combinedSecondaryVertexMVABJetTags");
872 lBtag.iPMVA = (*iter).bDiscriminator("impactParameterMVABJetTags");
873 lBtag.bProba = (*iter).bDiscriminator("jetBProbabilityBJetTags");
874 lBtag.probability = (*iter).bDiscriminator("jetProbabilityBJetTags");
875 lBtag.sSV = (*iter).bDiscriminator("simpleSecondaryVertexBJetTags");
876 lBtag.softElectron = (*iter).bDiscriminator("softElectronBJetTags");
877 lBtag.softMuon = (*iter).bDiscriminator("softMuonBJetTags");
878 lBtag.softMuonNoIP = (*iter).bDiscriminator("softMuonNoIPBJetTags");
879 lBtag.tCHE = (*iter).bDiscriminator("trackCountingHighEffBJetTags");
880 lBtag.tCHP = (*iter).bDiscriminator("trackCountingHighPurBJetTags");
881
882 bool isCalo = (*iter).isCaloJet();
883 bool isPF = (*iter).isPFJet();
884
885 assert (isCalo == !isPF);
886 if (isCalo) {
887 HbbAnalysis::CaloJetVars lCalo;
888 lCalo.maxEInEmTowers = (*iter).maxEInEmTowers();
889 lCalo.maxEInHadTowers = (*iter).maxEInHadTowers();
890 lCalo.energyFractionHadronic = (*iter).energyFractionHadronic();
891 lCalo.emEnergyFraction = (*iter).emEnergyFraction();
892 lCalo.hadEnergyInHB = (*iter).hadEnergyInHB();
893 lCalo.hadEnergyInHO = (*iter).hadEnergyInHO();
894 lCalo.hadEnergyInHE = (*iter).hadEnergyInHE();
895 lCalo.hadEnergyInHF = (*iter).hadEnergyInHF();
896 lCalo.emEnergyInEB = (*iter).emEnergyInEB();
897 lCalo.emEnergyInEE = (*iter).emEnergyInEE();
898 lCalo.emEnergyInHF = (*iter).emEnergyInHF();
899 lCalo.towersArea = (*iter).towersArea();
900 lCalo.n90 = (*iter).n90();
901 lCalo.n60 = (*iter).n60();
902
903 HbbAnalysis::Jet lObj(lGen,lGenJet,lReco,lCommon,lCalo,lBtag);
904 aVec.push_back(lObj);
905
906 }
907
908 if (isPF){
909 HbbAnalysis::PFJetVars lPF;
910 lPF.chargedHadronEnergy = (*iter).chargedHadronEnergy();
911 lPF.chargedHadronEnergyFraction = (*iter).chargedHadronEnergyFraction();
912 lPF.neutralHadronEnergy = (*iter).neutralHadronEnergy();
913 lPF.neutralHadronEnergyFraction = (*iter).neutralHadronEnergyFraction();
914 lPF.chargedEmEnergy = (*iter).chargedEmEnergy();
915 lPF.chargedEmEnergyFraction = (*iter).chargedEmEnergyFraction();
916 lPF.chargedMuEnergy = (*iter).chargedMuEnergy();
917 lPF.chargedMuEnergyFraction = (*iter).chargedMuEnergyFraction();
918 lPF.neutralEmEnergy = (*iter).neutralEmEnergy();
919 lPF.neutralEmEnergyFraction = (*iter).neutralEmEnergyFraction();
920 lPF.chargedMultiplicity = (*iter).chargedMultiplicity();
921 lPF.neutralMultiplicity = (*iter).neutralMultiplicity();
922 lPF.muonMultiplicity = (*iter).muonMultiplicity();
923
924 HbbAnalysis::Jet lObj(lGen,lGenJet,lReco,lCommon,lPF,lBtag);
925 aVec.push_back(lObj);
926
927 }
928
929 iEle++;
930 }//loop on element
931 }//non empty
932 }//valid
933
934 }//HbbJets
935
936 void HbbTreeMaker::HbbMet(const edm::Handle<std::vector<pat::MET> > & aCol,
937 HbbAnalysis::Met & aMet)
938 {//HbbMet
939 HbbAnalysis::MetVars lGen;
940 HbbAnalysis::MetVars lReco;
941
942 assert(aCol->size() == 1);
943 const pat::MET & lMet = *(aCol->begin());
944
945 lReco.mET = lMet.pt();
946 lReco.mEx = lMet.px();
947 lReco.mEy = lMet.py();
948 lReco.sumET = lMet.sumEt();
949 lReco.phi = lMet.phi();
950 lReco.mEtSig = lMet.mEtSig();
951
952 const reco::GenMET *lGenMET = lMet.genMET();
953
954 if (lGenMET){
955 lGen.mET = lGenMET->pt();
956 lGen.mEx = lGenMET->px();
957 lGen.mEy = lGenMET->py();
958 lGen.sumET = lGenMET->sumEt();
959 lGen.phi = lGenMET->phi();
960 lGen.mEtSig = lGenMET->mEtSig();
961 }
962 else {
963 lGen.mET = 0;
964 lGen.mEx = 0;
965 lGen.mEy = 0;
966 lGen.sumET = 0;
967 lGen.phi = 0;
968 lGen.mEtSig = 0;
969 }
970
971 aMet.genVars(lGen);
972 aMet.recoVars(lReco);
973
974 }//HbbMet
975
976 void HbbTreeMaker::HbbTrigger(const edm::Handle<edm::TriggerResults> & aCol,
977 std::vector<HbbAnalysis::Trigger> & aVec)
978 {//HbbTrigger
979
980 if (aCol.isValid()){
981 edm::TriggerNames lNames;
982 lNames.init(*aCol);
983
984 //for (unsigned int j=0; j<hltPaths_.size(); j++) {
985 //bool valid=false;
986 for (unsigned int i=0; i<aCol->size(); i++){
987 std::string trueName = lNames.triggerNames().at(i);
988
989 HbbAnalysis::TriggerVars lTrigVars;
990 lTrigVars.name = trueName;
991 lTrigVars.index = i;
992 lTrigVars.accept = aCol->accept(i);
993
994 HbbAnalysis::Trigger lObj(lTrigVars);
995
996 aVec.push_back(lObj);
997 }
998
999 }
1000
1001 }
1002
1003 void HbbTreeMaker::HbbParticles(const edm::Handle<reco::GenParticleCollection> & aCol,
1004 std::vector<HbbAnalysis::GenParticle> & aVec)
1005 {//genparticles
1006
1007 //add genParticle inside vector
1008 for (unsigned int mccount = 0;mccount<aCol->size();mccount++)
1009 {//loop on particles
1010 const reco::Candidate & p = (*aCol)[mccount];
1011
1012 HbbAnalysis::MCVars lMC;
1013 lMC.index = mccount;
1014 lMC.E = p.energy();
1015 lMC.pT = p.pt();
1016 lMC.eta = p.eta();
1017 lMC.phi = p.phi();
1018 lMC.pdgId = p.pdgId();
1019 lMC.status = p.status();
1020 HbbAnalysis::GenParticle lObj(lMC);
1021 aVec.push_back(lObj);
1022
1023 }//loop on particles
1024
1025 //now need to get the list of parents....
1026
1027 for (unsigned int mccount = 0;mccount<aCol->size();mccount++)
1028 {//loop on particles
1029 const reco::Candidate & p = (*aCol)[mccount];
1030
1031 unsigned int nD = p.numberOfDaughters();
1032 //std::cout << mccount << " has " << nD << " daughter(s) : " ;
1033 for(unsigned int dau=0;dau<nD;dau++){//loop on daughters
1034 const reco::Candidate * pDau = p.daughter(dau);
1035 //std::cout << pDau << " ("<<std::flush;
1036 for (unsigned int mccount2 = 0;mccount2<aCol->size();mccount2++)
1037 {//loop on particles
1038 const reco::Candidate & p2 = (*aCol)[mccount2];
1039 //std::cout<<"DBG: mccount2 = "<<mccount2<<" gen size = "<<data_->gen().size()<<std::endl;
1040
1041 HbbAnalysis::GenParticle & part2 = aVec.at(mccount2);
1042 if (pDau == &p2) {
1043 part2.setParent(mccount);
1044 //std::cout << &p2 << ", index = " << mccount2 << "), "<<std::flush;
1045 break;
1046 }
1047 //else std::cout << std::endl << "****no match " << mccount2 << " " << &p2 << std::endl;
1048 //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){
1049 //part2.parent(mccount);
1050 //}
1051 }//loop on particles
1052 }//loop on daughters
1053 //std::cout << std::endl;
1054
1055 }//loop on particles
1056
1057 if (debug_){
1058 for (unsigned int imc(0); imc < aVec.size(); imc++){
1059 HbbAnalysis::GenParticle & p = aVec.at(imc);
1060 p.print();
1061 }
1062 }
1063
1064
1065 }//genparticles