ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/HbbAnalysis/plugins/HbbTreeMaker.cc
Revision: 1.23
Committed: Sat Sep 25 13:03:51 2010 UTC (14 years, 7 months ago) by amagnan
Content type: text/plain
Branch: MAIN
CVS Tags: v00-05-01
Changes since 1.22: +8 -7 lines
Log Message:
added gg,qg,qq separation

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