ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/System8/EDModule/Analyzer/src/S8TreeMaker.cc
Revision: 1.5
Committed: Fri Dec 16 00:35:40 2011 UTC (13 years, 4 months ago) by meloam
Content type: text/plain
Branch: MAIN
CVS Tags: V00-00-04, HEAD
Changes since 1.4: +0 -2 lines
Log Message:
updating EDModule for dec 2011 run

File Contents

# Content
1 /**
2 * S8TreeMaker
3 *
4 *
5 * Created by Samvel Khalatian on Sep 29, 2010
6 * Copyright 2010, All rights reserved
7 */
8
9 #include <iostream>
10 #include <vector>
11
12 #include <boost/algorithm/string.hpp>
13 #include <boost/lexical_cast.hpp>
14 #include <boost/logic/tribool.hpp>
15 #include <boost/regex.hpp>
16
17 #include <TTree.h>
18 #include <TLorentzVector.h>
19 #include <TVector3.h>
20
21 #include "CommonTools/UtilAlgos/interface/TFileService.h"
22 #include "DataFormats/Common/interface/Handle.h"
23 #include "DataFormats/Common/interface/TriggerResults.h"
24 #include "DataFormats/PatCandidates/interface/Jet.h"
25 #include "DataFormats/PatCandidates/interface/Muon.h"
26 #include "DataFormats/PatCandidates/interface/Electron.h"
27 #include "DataFormats/VertexReco/interface/Vertex.h"
28 #include "FWCore/Common/interface/TriggerNames.h"
29 #include "FWCore/Framework/interface/MakerMacros.h"
30 #include "FWCore/Framework/interface/Event.h"
31 #include "FWCore/MessageLogger/interface/MessageLogger.h"
32 #include "FWCore/ParameterSet/interface/ParameterSet.h"
33 #include "FWCore/ServiceRegistry/interface/Service.h"
34 #include "FWCore/Utilities/interface/EDMException.h"
35 #include "FWCore/Utilities/interface/InputTag.h"
36 #include "SimDataFormats/GeneratorProducts/interface/GenEventInfoProduct.h"
37 #include "SimDataFormats/PileupSummaryInfo/interface/PileupSummaryInfo.h"
38
39 #include "Tree/System8/interface/S8TreeInfo.h"
40 #include "Tree/System8/interface/S8EventID.h"
41 #include "Tree/System8/interface/S8GenEvent.h"
42 #include "Tree/System8/interface/S8GenParticle.h"
43 #include "Tree/System8/interface/S8Jet.h"
44 #include "Tree/System8/interface/S8Lepton.h"
45 #include "Tree/System8/interface/S8PrimaryVertex.h"
46 #include "Tree/System8/interface/S8TreeInfo.h"
47 #include "Tree/System8/interface/S8Trigger.h"
48 #include "Tree/System8/interface/S8TriggerCenter.h"
49 #include "EDModule/Analyzer/interface/Tools.h"
50
51 #include "EDModule/Analyzer/interface/S8TreeMaker.h"
52
53 using std::cout;
54 using std::endl;
55 using std::string;
56 using std::vector;
57
58 using boost::lexical_cast;
59 using boost::regex;
60 using boost::smatch;
61
62 using edm::errors::ErrorCodes;
63 using edm::Handle;
64 using edm::InputTag;
65 using edm::LogWarning;
66 using edm::LogInfo;
67 using edm::ParameterSet;
68 using edm::TriggerResults;
69 using edm::TriggerNames;
70
71 using reco::Vertex;
72
73 using s8::TreeInfo;
74 using s8::TriggerCenter;
75
76 using namespace top::tools;
77
78 void fillSplit(const int &status, boost::tribool &split)
79 {
80 switch(status)
81 {
82 case 3: split = false;
83 break;
84
85 case 2: if (boost::indeterminate(split))
86 split = true;
87
88 break;
89 }
90 }
91
92 S8TreeMaker::S8TreeMaker(const edm::ParameterSet &config):
93 _didInitializeHltConfigProvider(false)
94 {
95 _primaryVertices = config.getParameter<string>("primaryVertices");
96 _jets = config.getParameter<string>("jets");
97 _muons = config.getParameter<string>("muons");
98 _electrons = config.getParameter<string>("electrons");
99 _triggers = config.getParameter<string>("triggers");
100
101 _jetSelector = config.getParameter<ParameterSet>("jetSelector");
102
103 _isPythia = config.getParameter<bool>("isPythia");
104 _saveTriggers = config.getParameter<bool>("saveTriggers");
105 }
106
107 S8TreeMaker::~S8TreeMaker()
108 {
109 }
110
111 void S8TreeMaker::beginJob()
112 {
113 cout << "Start BeginJob" << endl;
114
115 edm::Service<TFileService> fileService;
116
117 _treeInfo.reset(new TreeInfo);
118 _triggerCenter.reset(new TriggerCenter);
119
120 _tree = fileService->make<TTree>("s8", "System8 tree.");
121
122 // Prepare Branches
123 //
124 _eventID.reset(new s8::EventID());
125 _tree->Branch("eventID", _eventID.get());
126
127 _genEvent.reset(new s8::GenEvent());
128 _tree->Branch("genEvent", _genEvent.get());
129
130 _s8Electrons.reset(new s8::Leptons());
131 _tree->Branch("electrons.", _s8Electrons.get());
132
133 _s8Jets.reset(new s8::Jets());
134 _tree->Branch("jets", _s8Jets.get());
135
136 _s8Muons.reset(new s8::Leptons());
137 _tree->Branch("muons.", _s8Muons.get());
138
139 _s8PrimaryVertices.reset(new s8::PrimaryVertices());
140 _tree->Branch("primaryVertices", _s8PrimaryVertices.get());
141
142 _s8Triggers.reset(new s8::Triggers());
143 _tree->Branch("triggers", _s8Triggers.get());
144
145 _didInitializeHltConfigProvider = false;
146
147 cout << "BeginJob is over" << endl;
148 }
149
150 void S8TreeMaker::endJob()
151 {
152 if (!_tree)
153 return;
154
155 // Tree Info is disabled for the moment until Hadd is fixed
156 //
157 edm::Service<TFileService> fileService;
158 // TDirectory *dir = fileService->cd();
159 TDirectory *dir = fileService->getBareDirectory();
160 dir->WriteObject(_treeInfo.get(), "s8info");
161 dir->WriteObject(_triggerCenter.get(), "s8triggers");
162 }
163
164 void S8TreeMaker::beginRun(const edm::Run &run,
165 const edm::EventSetup &eventSetup)
166 {
167 using s8::Trigger;
168
169 // Initialize HLT Config Provider for new Run
170 //
171 bool didChange = true;
172 if (!_hltConfigProvider.init(run,
173 eventSetup,
174 InputTag(_triggers).process(),
175 didChange))
176
177 throw cms::Exception("S8TreeMaker")
178 << "Failed to initialize HLTConfig for: " << _triggers;
179
180 // Test if Trigger Menu has changed
181 //
182 if (!didChange &&
183 _didInitializeHltConfigProvider)
184
185 return;
186
187 _hlts.clear();
188
189 // Triggers are found and changed
190 //
191 typedef std::vector<std::string> Triggers;
192
193 const Triggers &triggerNames = _hltConfigProvider.triggerNames();
194
195 using s8::tools::make_hash;
196
197 TriggerCenter::TriggerMap &triggers = _triggerCenter->triggers();
198
199 for(Triggers::const_iterator trigger = triggerNames.begin();
200 triggerNames.end() != trigger;
201 ++trigger)
202 {
203 smatch matches;
204 if (!regex_match(*trigger, matches,
205 regex("^(\\w+?)(?:_[vV](\\d+))?$")))
206 {
207 cout << "Do not understand Trigger Name: " << *trigger
208 << endl;
209
210 continue;
211 }
212
213 string trigger_name = matches[1];
214 boost::to_lower(trigger_name);
215 triggers.insert(make_pair(make_hash(trigger_name), matches[1]));
216
217 // Found Trigger of the interest
218 //
219 HLT foundHLT;
220
221 foundHLT.hash = make_hash(trigger_name);
222 foundHLT.id = distance(triggerNames.begin(), trigger);
223 foundHLT.version = matches[2].matched
224 ? lexical_cast<int>(matches[2])
225 : 0;
226
227 _hlts[*trigger] = foundHLT;
228 }
229
230 if (_hlts.empty())
231 LogWarning("S8TreeMaker")
232 << "None of the searched HLT Triggers is found" << endl;
233
234 _didInitializeHltConfigProvider = true;
235 }
236
237 void S8TreeMaker::analyze(const edm::Event &event,
238 const edm::EventSetup &eventSetup)
239 {
240 if (!_tree)
241 throw cms::Exception("S8TreeMaker")
242 << "Tree does not exist";
243
244 processEventID(event);
245 processGenEvent(event);
246 processElectrons(event);
247 processJets(event);
248 processMuons(event);
249 processPrimaryVertices(event);
250
251 if (_saveTriggers)
252 processTriggers(event, eventSetup);
253
254 // Write Tree entry
255 //
256 _tree->Fill();
257
258 // Reset event
259 //
260 _eventID->reset();
261 _genEvent->reset();
262
263 for(s8::Jets::iterator jet = _s8Jets->begin();
264 _s8Jets->end() != jet;
265 ++jet)
266 {
267 delete *jet;
268 }
269 _s8Jets->clear();
270
271 for(s8::Leptons::iterator electron = _s8Electrons->begin();
272 _s8Electrons->end() != electron;
273 ++electron)
274 {
275 delete *electron;
276 }
277 _s8Electrons->clear();
278
279 for(s8::Leptons::iterator muon = _s8Muons->begin();
280 _s8Muons->end() != muon;
281 ++muon)
282 {
283 delete *muon;
284 }
285 _s8Muons->clear();
286
287 for(s8::PrimaryVertices::iterator primaryVertex = _s8PrimaryVertices->begin();
288 _s8PrimaryVertices->end() != primaryVertex;
289 ++primaryVertex)
290 {
291 delete *primaryVertex;
292 }
293 _s8PrimaryVertices->clear();
294
295 for(s8::Triggers::iterator trigger = _s8Triggers->begin();
296 _s8Triggers->end() != trigger;
297 ++trigger)
298 {
299 delete *trigger;
300 }
301 _s8Triggers->clear();
302 }
303
304 void S8TreeMaker::processEventID(const edm::Event &event)
305 {
306 _eventID->setRun(event.id().run());
307 _eventID->setLumiBlock(event.id().luminosityBlock());
308 _eventID->setEvent(event.id().event());
309 }
310
311 void S8TreeMaker::processGenEvent(const edm::Event &event)
312 {
313 using reco::GenParticleCollection;
314
315 if (event.isRealData())
316 return;
317
318 // check if Event is Pythia
319 //
320 if (_isPythia)
321 {
322 Handle<GenEventInfoProduct> generator;
323 event.getByLabel(InputTag("generator"), generator);
324
325 if (!generator.isValid())
326 {
327 LogWarning("S8TreeMaker")
328 << "failed to extract Generator";
329
330 return;
331 }
332
333 _genEvent->setPtHat(generator->qScale());
334 }
335
336 Handle<GenParticleCollection> genParticles;
337 event.getByLabel(InputTag("genParticles"), genParticles);
338
339 boost::tribool bsplit = boost::indeterminate;
340 boost::tribool csplit = boost::indeterminate;
341 for(GenParticleCollection::const_iterator particle = genParticles->begin();
342 genParticles->end() != particle;
343 ++particle)
344 {
345 const int pdgID = abs(particle->pdgId());
346
347 switch(pdgID)
348 {
349 case 5: fillSplit(particle->status(), bsplit);
350 break;
351
352 case 4: fillSplit(particle->status(), csplit);
353 break;
354
355 default: break;
356 }
357 }
358
359 if (bsplit)
360 _genEvent->setGluonSplitting(s8::GenEvent::BB, true);
361
362 if (csplit)
363 _genEvent->setGluonSplitting(s8::GenEvent::CC, true);
364
365 // add beam crossing information
366 Handle<std::vector< PileupSummaryInfo > > PupInfo;
367 if ( event.getByLabel("addPileupInfo", PupInfo) )
368 {
369 std::vector<PileupSummaryInfo>::const_iterator PVInfo;
370
371 for(PVInfo = PupInfo->begin(); PVInfo != PupInfo->end(); ++PVInfo) {
372
373 int BX = PVInfo->getBunchCrossing();
374
375 if(BX == -1) {
376 _genEvent->setBx_minus1( PVInfo->getPU_NumInteractions() );
377 }
378 if(BX == 0) {
379 _genEvent->setBx_zero( PVInfo->getPU_NumInteractions() );
380 }
381 if(BX == 1) {
382 _genEvent->setBx_plus1( PVInfo->getPU_NumInteractions() );
383 }
384
385 }
386 }
387
388
389 }
390
391 void S8TreeMaker::processElectrons(const edm::Event &event)
392 {
393 using pat::ElectronCollection;
394
395 // Extract Electrons
396 //
397 Handle<ElectronCollection> electrons;
398 event.getByLabel(InputTag(_electrons), electrons);
399
400 if (!electrons.isValid())
401 {
402 LogWarning("S8TreeMaker")
403 << "failed to extract Electrons.";
404
405 return;
406 }
407
408 // Process all Electrons
409 //
410 for(ElectronCollection::const_iterator electron = electrons->begin();
411 electrons->end() != electron;
412 ++electron)
413 {
414 s8::Lepton *s8Electron = new s8::Lepton();;
415
416 setP4(s8Electron->p4(), electron->p4());
417 setVertex(s8Electron->vertex(), electron->vertex());
418
419 s8Electron->impactParameter()->first = electron->dB();
420 s8Electron->impactParameter()->second = electron->edB();
421
422 // Extract GenParticle information
423 //
424 if (electron->genLepton())
425 {
426 s8::GenParticle *s8GenParticle = s8Electron->genParticle();
427
428 setP4(s8GenParticle->p4(), electron->genLepton()->p4());
429 setVertex(s8GenParticle->vertex(), electron->genLepton()->vertex());
430
431 s8GenParticle->setId(electron->genLepton()->pdgId());
432 s8GenParticle->setStatus(electron->genLepton()->status());
433 if (electron->genLepton()->mother())
434 s8GenParticle->setParentId(electron->genLepton()->mother()->pdgId());
435 }
436
437 _s8Electrons->push_back(s8Electron);
438 } // End loop over electrons
439 }
440
441 void S8TreeMaker::processJets(const edm::Event &event)
442 {
443 using pat::JetCollection;
444
445 // Extract Jets
446 //
447 Handle<JetCollection> jets;
448 event.getByLabel(InputTag(_jets), jets);
449
450 if (!jets.isValid())
451 {
452 LogWarning("S8TreeMaker")
453 << "failed to extract Jets.";
454
455 return;
456 }
457
458 // Process All Jets
459 //
460 //pat::strbitset jetBitset = _jetSelector.getBitTemplate();
461 for(JetCollection::const_iterator jet = jets->begin();
462 jets->end() != jet;
463 ++jet)
464 {
465 //if (!_jetSelector(*jet, jetBitset))
466 // continue;
467
468 using s8::Jet;
469 Jet *s8Jet = new Jet();
470
471 setP4(s8Jet->p4(), jet->p4());
472
473 s8Jet->setFlavour(jet->partonFlavour());
474 s8Jet->setTracks(jet->associatedTracks().size());
475
476 if (jet->genParton())
477 {
478 s8::GenParticle *s8GenParticle = s8Jet->genParticle();
479
480 setP4(s8GenParticle->p4(), jet->genParton()->p4());
481 setVertex(s8GenParticle->vertex(), jet->genParton()->vertex());
482
483 s8GenParticle->setId(jet->genParton()->pdgId());
484 s8GenParticle->setStatus(jet->genParton()->status());
485 if (jet->genParton()->mother())
486 s8GenParticle->setParentId(jet->genParton()->mother()->pdgId());
487 }
488
489 // Save b-taggers
490 //
491 s8Jet->setBTag(Jet::TCHE,
492 jet->bDiscriminator("trackCountingHighEffBJetTags"));
493
494 s8Jet->setBTag(Jet::TCHP,
495 jet->bDiscriminator("trackCountingHighPurBJetTags"));
496
497 s8Jet->setBTag(Jet::JP,
498 jet->bDiscriminator("jetProbabilityBJetTags"));
499
500 s8Jet->setBTag(Jet::SSV,
501 jet->bDiscriminator("simpleSecondaryVertexBJetTags"));
502
503 s8Jet->setBTag(Jet::SSVHE,
504 jet->bDiscriminator("simpleSecondaryVertexHighEffBJetTags"));
505
506 s8Jet->setBTag(Jet::SSVHP,
507 jet->bDiscriminator("simpleSecondaryVertexHighPurBJetTags"));
508
509 s8Jet->setBTag(Jet::CSV,
510 jet->bDiscriminator("combinedSecondaryVertexBJetTags"));
511
512 s8Jet->setBTag(Jet::CSV_MVA,
513 jet->bDiscriminator("combinedSecondaryVertexMVABJetTags"));
514
515 s8Jet->setBTag(Jet::JBP,
516 jet->bDiscriminator("jetBProbabilityBJetTags"));
517
518 _s8Jets->push_back(s8Jet);
519 }
520 }
521
522 void S8TreeMaker::processMuons(const edm::Event &event)
523 {
524 using pat::MuonCollection;
525
526 // Extract Muons
527 //
528 Handle<MuonCollection> muons;
529 event.getByLabel(InputTag(_muons), muons);
530
531 if (!muons.isValid())
532 {
533 LogWarning("S8TreeMaker")
534 << "failed to extract Muons.";
535
536 return;
537 }
538
539 // Process all Muons
540 //
541 for(MuonCollection::const_iterator muon = muons->begin();
542 muons->end() != muon;
543 ++muon)
544 {
545 // Only Muons with basic cuts are saved:
546 //
547 if (1 >= muon->numberOfMatches())
548 continue;
549
550 s8::Lepton *s8Muon = new s8::Lepton();
551
552 setP4(s8Muon->p4(), muon->p4());
553 setVertex(s8Muon->vertex(), muon->vertex());
554
555 s8Muon->impactParameter()->first = muon->dB();
556 s8Muon->impactParameter()->second = muon->edB();
557
558 // Extract GenParticle information
559 //
560 if (muon->genLepton())
561 {
562 s8::GenParticle *s8GenParticle = s8Muon->genParticle();
563
564 setP4(s8GenParticle->p4(), muon->genLepton()->p4());
565 setVertex(s8GenParticle->vertex(), muon->genLepton()->vertex());
566
567 s8GenParticle->setId(muon->genLepton()->pdgId());
568 s8GenParticle->setStatus(muon->genLepton()->status());
569 if (muon->genLepton()->mother())
570 s8GenParticle->setParentId(muon->genLepton()->mother()->pdgId());
571 }
572
573 _s8Muons->push_back(s8Muon);
574 } // End loop over muons
575 }
576
577 void S8TreeMaker::processPrimaryVertices(const edm::Event &event)
578 {
579 // Primary Vertices
580 //
581 typedef vector<Vertex> PVCollection;
582
583 Handle<PVCollection> primaryVertices;
584 event.getByLabel(InputTag(_primaryVertices), primaryVertices);
585
586 if (!primaryVertices.isValid())
587 {
588 LogWarning("S8TreeMaker")
589 << "failed to extract Primary Vertices.";
590
591 return;
592 }
593
594 if (primaryVertices->empty())
595 {
596 LogWarning("S8TreeMaker")
597 << "primary vertices collection is empty.";
598
599 return;
600 }
601 // Process Primary Vertices
602 //
603 for(PVCollection::const_iterator vertex = primaryVertices->begin();
604 primaryVertices->end() != vertex;
605 ++vertex)
606 {
607 if (!isGoodPrimaryVertex(*vertex, event.isRealData()))
608 continue;
609
610 s8::PrimaryVertex *s8Vertex = new s8::PrimaryVertex();
611
612 setVertex(s8Vertex->vertex(), vertex->position());
613 s8Vertex->setNdof(vertex->ndof());
614 s8Vertex->setRho(vertex->position().Rho());
615
616 _s8PrimaryVertices->push_back(s8Vertex);
617 }
618 }
619
620 void S8TreeMaker::processTriggers(const edm::Event &event,
621 const edm::EventSetup &eventSetup)
622 {
623 if (_hlts.empty())
624 return;
625
626 // Triggers
627 //
628 Handle<TriggerResults> triggers;
629 event.getByLabel(InputTag(_triggers), triggers);
630
631 if (!triggers.isValid())
632 {
633 LogWarning("S8TreeMaker")
634 << "failed to extract Triggers";
635
636 return;
637 }
638
639 // Process only found HLTs
640 //
641 for(HLTs::const_iterator hlt = _hlts.begin();
642 _hlts.end() != hlt;
643 ++hlt)
644 {
645 s8::Trigger *s8Trigger = new s8::Trigger();
646 s8Trigger->setHash(hlt->second.hash);
647 if (hlt->second.version)
648 s8Trigger->setVersion(hlt->second.version);
649
650 s8Trigger->setIsPass(triggers->accept(hlt->second.id));
651 s8Trigger->setPrescale(_hltConfigProvider.prescaleValue(event,
652 eventSetup, hlt->first));
653
654
655 _s8Triggers->push_back(s8Trigger);
656 }
657 }
658
659 bool S8TreeMaker::isGoodPrimaryVertex(const Vertex &vertex,
660 const bool &isData)
661 {
662 return !vertex.isFake() &&
663 4 <= vertex.ndof() &&
664 (isData ? 24 : 15) >= fabs(vertex.z()) &&
665 2 >= fabs(vertex.position().Rho());
666 }
667
668 DEFINE_FWK_MODULE(S8TreeMaker);