ViewVC Help
View File | Revision Log | Show Annotations | Root Listing
root/cvsroot/UserCode/HbbAnalysis/plugins/HbbTreeMaker.cc
Revision: 1.28
Committed: Thu Mar 17 12:47:39 2011 UTC (14 years, 2 months ago) by amagnan
Content type: text/plain
Branch: MAIN
CVS Tags: v00-05-02
Changes since 1.27: +26 -20 lines
Log Message:
JEC

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