CMS 3D CMS Logo

 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
IsolatedTracksNxN.cc
Go to the documentation of this file.
1 // -*- C++ -*
2 //
3 // Package: IsolatedTracksNxN
4 // Class: IsolatedTracksNxN
5 //
13 //
14 // Original Author: Seema Sharma
15 // Created: Mon Aug 10 15:30:40 CST 2009
16 //
17 //
18 
21 
34 
36 
37  //now do what ever initialization is needed
38  doMC = iConfig.getUntrackedParameter<bool> ("DoMC", false);
39  writeAllTracks = iConfig.getUntrackedParameter<bool> ("WriteAllTracks", false );
40  myverbose_ = iConfig.getUntrackedParameter<int> ("Verbosity", 5 );
41  pvTracksPtMin_ = iConfig.getUntrackedParameter<double>("PVTracksPtMin", 1.0 );
42  debugTrks_ = iConfig.getUntrackedParameter<int> ("DebugTracks" );
43  printTrkHitPattern_ = iConfig.getUntrackedParameter<bool> ("PrintTrkHitPattern" );
44  minTrackP_ = iConfig.getUntrackedParameter<double>("minTrackP", 1.0 );
45  maxTrackEta_ = iConfig.getUntrackedParameter<double>("maxTrackEta", 5.0 );
46  debugL1Info_ = iConfig.getUntrackedParameter<bool> ("DebugL1Info", false );
47  L1TriggerAlgoInfo_ = iConfig.getUntrackedParameter<bool> ("L1TriggerAlgoInfo");
48  L1extraTauJetSource_ = iConfig.getParameter<edm::InputTag> ("L1extraTauJetSource" );
49  L1extraCenJetSource_ = iConfig.getParameter<edm::InputTag> ("L1extraCenJetSource" );
50  L1extraFwdJetSource_ = iConfig.getParameter<edm::InputTag> ("L1extraFwdJetSource" );
51  L1extraMuonSource_ = iConfig.getParameter<edm::InputTag> ("L1extraMuonSource" );
52  L1extraIsoEmSource_ = iConfig.getParameter<edm::InputTag> ("L1extraIsoEmSource" );
53  L1extraNonIsoEmSource_ = iConfig.getParameter<edm::InputTag> ("L1extraNonIsoEmSource" );
54  L1GTReadoutRcdSource_ = iConfig.getParameter<edm::InputTag> ("L1GTReadoutRcdSource" );
55  L1GTObjectMapRcdSource_= iConfig.getParameter<edm::InputTag> ("L1GTObjectMapRcdSource");
56  JetSrc_ = iConfig.getParameter<edm::InputTag> ("JetSource");
57  JetExtender_ = iConfig.getParameter<edm::InputTag> ("JetExtender");
58  HBHERecHitSource_ = iConfig.getParameter<edm::InputTag> ("HBHERecHitSource");
59  tMinE_ = iConfig.getUntrackedParameter<double>("TimeMinCutECAL", -500.);
60  tMaxE_ = iConfig.getUntrackedParameter<double>("TimeMaxCutECAL", 500.);
61  tMinH_ = iConfig.getUntrackedParameter<double>("TimeMinCutHCAL", -500.);
62  tMaxH_ = iConfig.getUntrackedParameter<double>("TimeMaxCutHCAL", 500.);
63  nbad = 0;
64 
65  // define tokens for access
66  tok_L1extTauJet_ = consumes<l1extra::L1JetParticleCollection>(L1extraTauJetSource_);
67  tok_L1extCenJet_ = consumes<l1extra::L1JetParticleCollection>(L1extraCenJetSource_);
68  tok_L1extFwdJet_ = consumes<l1extra::L1JetParticleCollection>(L1extraFwdJetSource_);
69  tok_L1extMu_ = consumes<l1extra::L1MuonParticleCollection>(L1extraMuonSource_);
70  tok_L1extIsoEm_ = consumes<l1extra::L1EmParticleCollection>(L1extraIsoEmSource_);
71  tok_L1extNoIsoEm_ = consumes<l1extra::L1EmParticleCollection>(L1extraNonIsoEmSource_);
72  tok_jets_ = consumes<reco::CaloJetCollection>(JetSrc_);
73  tok_hbhe_ = consumes<HBHERecHitCollection>(HBHERecHitSource_);
74 
75 
76  tok_genTrack_ = consumes<reco::TrackCollection>(edm::InputTag("generalTracks"));
77  tok_recVtx_ = consumes<reco::VertexCollection>(edm::InputTag("offlinePrimaryVertices"));
78  tok_bs_ = consumes<reco::BeamSpot>(edm::InputTag("offlineBeamSpot"));
79  tok_EB_ = consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit","EcalRecHitsEB"));
80  tok_EE_ = consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit","EcalRecHitsEE"));
81  tok_simTk_ = consumes<edm::SimTrackContainer>(edm::InputTag("g4SimHits"));
82  tok_simVtx_ = consumes<edm::SimVertexContainer>(edm::InputTag("g4SimHits"));
83  tok_caloEB_ = consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "EcalHitsEB"));
84  tok_caloEE_ = consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "EcalHitsEE"));
85  tok_caloHH_ = consumes<edm::PCaloHitContainer>(edm::InputTag("g4SimHits", "HcalHits"));
86 
87  if(myverbose_>=0) {
88  std::cout <<"Parameters read from config file \n"
89  <<" doMC " << doMC
90  <<"\t myverbose_ " << myverbose_
91  <<"\t minTrackP_ " << minTrackP_
92  << "\t maxTrackEta_ " << maxTrackEta_
93  << "\t tMinE_ " << tMinE_
94  << "\t tMaxE_ " << tMaxE_
95  << "\t tMinH_ " << tMinH_
96  << "\t tMaxH_ " << tMaxH_
97  << "\n debugL1Info_ " << debugL1Info_
98  << "\t L1TriggerAlgoInfo_ " << L1TriggerAlgoInfo_
99  << "\t L1extraTauJetSource_ " << L1extraTauJetSource_
100  << "\t L1extraCenJetSource_ " << L1extraCenJetSource_
101  << "\t L1extraFwdJetSource_ " << L1extraFwdJetSource_
102  << "\t L1extraMuonSource_ " << L1extraMuonSource_
103  << "\t L1extraIsoEmSource_ " << L1extraIsoEmSource_
104  << "\t L1extraNonIsoEmSource_ " << L1extraNonIsoEmSource_
105  << "\t L1GTReadoutRcdSource_ " << L1GTReadoutRcdSource_
106  << "\t L1GTObjectMapRcdSource_ " << L1GTObjectMapRcdSource_
107  << "\t JetSrc_ " << JetSrc_
108  << "\t JetExtender_ " << JetExtender_
109  << "\t HBHERecHitSource_ " << HBHERecHitSource_
110  << "\n" << std::endl;
111  }
112 
113  initL1 = false;
114 
115 }
116 
118 }
119 
121 
122  bool haveIsoTrack=false;
123 
125  iSetup.get<IdealMagneticFieldRecord>().get(bFieldH);
126  bField = bFieldH.product();
127 
129 
130  t_RunNo = iEvent.id().run();
131  t_EvtNo = iEvent.id().event();
132  t_Lumi = iEvent.luminosityBlock();
133  t_Bunch = iEvent.bunchCrossing();
134 
135  nEventProc++;
136 
138  iEvent.getByToken(tok_genTrack_, trkCollection);
139  reco::TrackCollection::const_iterator trkItr;
140  if(debugTrks_>1){
141  std::cout << "Track Collection: " << std::endl;
142  std::cout << "Number of Tracks " << trkCollection->size() << std::endl;
143  }
144  std::string theTrackQuality = "highPurity";
145  reco::TrackBase::TrackQuality trackQuality_=reco::TrackBase::qualityByName(theTrackQuality);
146 
147  //===================== save L1 Trigger information =======================
148  if( L1TriggerAlgoInfo_ ) {
149 
150  bool useL1EventSetup = true;
151  bool useL1GtTriggerMenuLite = true;
152 
153  m_l1GtUtils.getL1GtRunCache(iEvent, iSetup, useL1EventSetup, useL1GtTriggerMenuLite);
154 
155  int iErrorCode = -1;
156  int l1ConfCode = -1;
157  const bool l1Conf = m_l1GtUtils.availableL1Configuration(iErrorCode, l1ConfCode);
158  if( !l1Conf) {
159  std::cout << "\nL1 configuration code:" << l1ConfCode
160  << "\nNo valid L1 trigger configuration available."
161  << "\nSee text above for error code interpretation"
162  << "\nNo return here, in order to test each method, protected against configuration error."
163  << std::endl;
164  }
165 
166  const L1GtTriggerMenu* m_l1GtMenu = m_l1GtUtils.ptrL1TriggerMenuEventSetup(iErrorCode);
167  const AlgorithmMap& algorithmMap = m_l1GtMenu->gtAlgorithmMap();
168  const std::string& menuName = m_l1GtMenu->gtTriggerMenuName();
169 
170  if (!initL1) {
171  initL1=true;
172  std::cout << "menuName " << menuName << std::endl;
173  for (CItAlgo itAlgo = algorithmMap.begin(); itAlgo != algorithmMap.end(); itAlgo++) {
174  std::string algName = itAlgo->first;
175  int algBitNumber = ( itAlgo->second ).algoBitNumber();
176  l1AlgoMap.insert( std::pair<std::pair<unsigned int,std::string>,int>( std::pair<unsigned int,std::string>(algBitNumber, algName) , 0) ) ;
177  }
178  std::map< std::pair<unsigned int,std::string>, int>::iterator itr;
179  for(itr=l1AlgoMap.begin(); itr!=l1AlgoMap.end(); itr++) {
180  std::cout << " ********** " << (itr->first).first <<" "<<(itr->first).second <<" "<<itr->second << std::endl;
181  }
182  }
183 
184  std::vector<int> algbits;
185  for (CItAlgo itAlgo = algorithmMap.begin(); itAlgo != algorithmMap.end(); itAlgo++) {
186  std::string algName = itAlgo->first;
187  int algBitNumber = ( itAlgo->second ).algoBitNumber();
188  //bool algResultBeforeMask = m_l1GtUtils.decisionBeforeMask(iEvent, itAlgo->first, iErrorCode);
189  //bool algResultAfterMask = m_l1GtUtils.decisionAfterMask (iEvent, itAlgo->first, iErrorCode);
190  //int triggerMask = m_l1GtUtils.triggerMask (iEvent, itAlgo->first, iErrorCode);
191  bool decision = m_l1GtUtils.decision (iEvent, itAlgo->first, iErrorCode);
192  int preScale = m_l1GtUtils.prescaleFactor (iEvent, itAlgo->first, iErrorCode);
193 
194  // save the algo names which fired
195  if( decision ){
196  l1AlgoMap[std::pair<unsigned int,std::string>(algBitNumber, algName)] += 1;
197  h_L1AlgoNames->Fill( algBitNumber );
198  t_L1AlgoNames->push_back(itAlgo->first);
199  t_L1PreScale->push_back(preScale);
200  t_L1Decision[algBitNumber] = 1;
201  algbits.push_back(algBitNumber);
202  }
203  }
204 
205  if(debugL1Info_) {
206  for(unsigned int ii=0; ii<t_L1AlgoNames->size(); ii++){
207  std::cout<<ii<<" "<<(*t_L1AlgoNames)[ii]<<" "<<(*t_L1PreScale)[ii]<<" "<<algbits[ii]<<std::endl;
208  }
209  std::cout << "L1Decision: ";
210  for (int i=0; i<128; ++i) {
211  std::cout << " " << i << ":" << t_L1Decision[i];
212  if ((i+1)%32 == 0) std::cout << std::endl;
213  }
214  }
215 
216  // L1Taus
218  iEvent.getByToken(tok_L1extTauJet_,l1TauHandle);
219  l1extra::L1JetParticleCollection::const_iterator itr;
220  int iL1Obj=0;
221  for(itr = l1TauHandle->begin(),iL1Obj=0; itr != l1TauHandle->end(); ++itr,iL1Obj++) {
222  if(iL1Obj<1) {
223  t_L1TauJetPt ->push_back( itr->pt() );
224  t_L1TauJetEta ->push_back( itr->eta() );
225  t_L1TauJetPhi ->push_back( itr->phi() );
226  }
227  if(debugL1Info_) {
228  std::cout << "tauJ p/pt " << itr->momentum() << " " << itr->pt()
229  << " eta/phi " << itr->eta() << " " << itr->phi()
230  << std::endl;
231  }
232  }
233 
234  // L1 Central Jets
236  iEvent.getByToken(tok_L1extCenJet_,l1CenJetHandle);
237  for( itr = l1CenJetHandle->begin(),iL1Obj=0; itr != l1CenJetHandle->end(); ++itr,iL1Obj++ ) {
238  if(iL1Obj<1) {
239  t_L1CenJetPt ->push_back( itr->pt() );
240  t_L1CenJetEta ->push_back( itr->eta() );
241  t_L1CenJetPhi ->push_back( itr->phi() );
242  }
243  if(debugL1Info_) {
244  std::cout << "cenJ p/pt " << itr->momentum() << " " << itr->pt()
245  << " eta/phi " << itr->eta() << " " << itr->phi()
246  << std::endl;
247  }
248  }
249 
250  // L1 Forward Jets
252  iEvent.getByToken(tok_L1extFwdJet_,l1FwdJetHandle);
253  for( itr = l1FwdJetHandle->begin(),iL1Obj=0; itr != l1FwdJetHandle->end(); ++itr,iL1Obj++ ) {
254  if(iL1Obj<1) {
255  t_L1FwdJetPt ->push_back( itr->pt() );
256  t_L1FwdJetEta ->push_back( itr->eta() );
257  t_L1FwdJetPhi ->push_back( itr->phi() );
258  }
259  if(debugL1Info_) {
260  std::cout << "fwdJ p/pt " << itr->momentum() << " " << itr->pt()
261  << " eta/phi " << itr->eta() << " " << itr->phi()
262  << std::endl;
263  }
264  }
265 
266  // L1 Isolated EM onjects
267  l1extra::L1EmParticleCollection::const_iterator itrEm;
269  iEvent.getByToken(tok_L1extIsoEm_, l1IsoEmHandle);
270  for( itrEm = l1IsoEmHandle->begin(),iL1Obj=0; itrEm != l1IsoEmHandle->end(); ++itrEm,iL1Obj++ ) {
271  if(iL1Obj<1) {
272  t_L1IsoEMPt ->push_back( itrEm->pt() );
273  t_L1IsoEMEta ->push_back( itrEm->eta() );
274  t_L1IsoEMPhi ->push_back( itrEm->phi() );
275  }
276  if(debugL1Info_) {
277  std::cout << "isoEm p/pt " << itrEm->momentum() << " " << itrEm->pt()
278  << " eta/phi " << itrEm->eta() << " " << itrEm->phi()
279  << std::endl;
280  }
281  }
282 
283  // L1 Non-Isolated EM onjects
285  iEvent.getByToken(tok_L1extNoIsoEm_, l1NonIsoEmHandle);
286  for( itrEm = l1NonIsoEmHandle->begin(),iL1Obj=0; itrEm != l1NonIsoEmHandle->end(); ++itrEm,iL1Obj++ ) {
287  if(iL1Obj<1) {
288  t_L1NonIsoEMPt ->push_back( itrEm->pt() );
289  t_L1NonIsoEMEta ->push_back( itrEm->eta() );
290  t_L1NonIsoEMPhi ->push_back( itrEm->phi() );
291  }
292  if(debugL1Info_) {
293  std::cout << "nonIsoEm p/pt " << itrEm->momentum() << " " << itrEm->pt()
294  << " eta/phi " << itrEm->eta() << " " << itrEm->phi()
295  << std::endl;
296  }
297  }
298 
299  // L1 Muons
300  l1extra::L1MuonParticleCollection::const_iterator itrMu;
302  iEvent.getByToken(tok_L1extMu_, l1MuHandle);
303  for( itrMu = l1MuHandle->begin(),iL1Obj=0; itrMu != l1MuHandle->end(); ++itrMu,iL1Obj++ ) {
304  if(iL1Obj<1) {
305  t_L1MuonPt ->push_back( itrMu->pt() );
306  t_L1MuonEta ->push_back( itrMu->eta() );
307  t_L1MuonPhi ->push_back( itrMu->phi() );
308  }
309  if(debugL1Info_) {
310  std::cout << "l1muon p/pt " << itrMu->momentum() << " " << itrMu->pt()
311  << " eta/phi " << itrMu->eta() << " " << itrMu->phi()
312  << std::endl;
313  }
314  }
315  }
316 
317  //============== store the information about all the Non-Fake vertices ===============
318 
320  iEvent.getByToken(tok_recVtx_,recVtxs);
321 
322  std::vector<reco::Track> svTracks;
323  math::XYZPoint leadPV(0,0,0);
324  double sumPtMax = -1.0;
325  for(unsigned int ind=0;ind<recVtxs->size();ind++) {
326 
327  if ( !((*recVtxs)[ind].isFake()) ) {
328  double vtxTrkSumPt=0.0, vtxTrkSumPtWt=0.0;
329  int vtxTrkNWt =0;
330  double vtxTrkSumPtHP=0.0, vtxTrkSumPtHPWt=0.0;
331  int vtxTrkNHP =0, vtxTrkNHPWt =0;
332 
333  reco::Vertex::trackRef_iterator vtxTrack = (*recVtxs)[ind].tracks_begin();
334 
335  for (vtxTrack = (*recVtxs)[ind].tracks_begin(); vtxTrack!=(*recVtxs)[ind].tracks_end(); vtxTrack++) {
336 
337  if((*vtxTrack)->pt()<pvTracksPtMin_) continue;
338 
339  bool trkQuality = (*vtxTrack)->quality(trackQuality_);
340 
341  vtxTrkSumPt += (*vtxTrack)->pt();
342 
343  vtxTrkSumPt += (*vtxTrack)->pt();
344  if( trkQuality ) {
345  vtxTrkSumPtHP += (*vtxTrack)->pt();
346  vtxTrkNHP++;
347  }
348 
349  double weight = (*recVtxs)[ind].trackWeight(*vtxTrack);
350  h_PVTracksWt ->Fill(weight);
351  if( weight>0.5) {
352  vtxTrkSumPtWt += (*vtxTrack)->pt();
353  vtxTrkNWt++;
354  if( trkQuality ) {
355  vtxTrkSumPtHPWt += (*vtxTrack)->pt();
356  vtxTrkNHPWt++;
357  }
358  }
359  }
360 
361  if(vtxTrkSumPt>sumPtMax) {
362  sumPtMax = vtxTrkSumPt;
363  leadPV = math::XYZPoint( (*recVtxs)[ind].x(),(*recVtxs)[ind].y(), (*recVtxs)[ind].z() );
364  }
365 
366  t_PVx ->push_back( (*recVtxs)[ind].x() );
367  t_PVy ->push_back( (*recVtxs)[ind].y() );
368  t_PVz ->push_back( (*recVtxs)[ind].z() );
369  t_PVisValid ->push_back( (*recVtxs)[ind].isValid() );
370  t_PVNTracks ->push_back( (*recVtxs)[ind].tracksSize() );
371  t_PVndof ->push_back( (*recVtxs)[ind].ndof() );
372  t_PVNTracksWt ->push_back( vtxTrkNWt );
373  t_PVTracksSumPt ->push_back( vtxTrkSumPt );
374  t_PVTracksSumPtWt->push_back( vtxTrkSumPtWt );
375 
376  t_PVNTracksHP ->push_back( vtxTrkNHP );
377  t_PVNTracksHPWt ->push_back( vtxTrkNHPWt );
378  t_PVTracksSumPtHP ->push_back( vtxTrkSumPtHP );
379  t_PVTracksSumPtHPWt->push_back( vtxTrkSumPtHPWt );
380 
381  if(myverbose_==4) {
382  std::cout<<"PV "<<ind<<" isValid "<<(*recVtxs)[ind].isValid()<<" isFake "<<(*recVtxs)[ind].isFake()
383  <<" hasRefittedTracks() "<<ind<<" "<<(*recVtxs)[ind].hasRefittedTracks()
384  <<" refittedTrksSize "<<(*recVtxs)[ind].refittedTracks().size()
385  <<" tracksSize() "<<(*recVtxs)[ind].tracksSize()<<" sumPt "<<vtxTrkSumPt
386  <<std::endl;
387  }
388  } // if vtx is not Fake
389  } // loop over PVs
390  //===================================================================================
391 
392  // Get the beamspot
393  edm::Handle<reco::BeamSpot> beamSpotH;
394  iEvent.getByToken(tok_bs_, beamSpotH);
395  math::XYZPoint bspot;
396  bspot = ( beamSpotH.isValid() ) ? beamSpotH->position() : math::XYZPoint(0, 0, 0);
397 
398  //=====================================================================
399 
401  iEvent.getByToken(tok_jets_,jets);
402  // edm::Handle<reco::JetExtendedAssociation::Container> jetExtender;
403  // iEvent.getByLabel(JetExtender_,jetExtender);
404 
405  for(unsigned int ijet=0;ijet<(*jets).size();ijet++) {
406  t_jetPt ->push_back( (*jets)[ijet].pt() );
407  t_jetEta ->push_back( (*jets)[ijet].eta() );
408  t_jetPhi ->push_back( (*jets)[ijet].phi() );
409  t_nTrksJetVtx ->push_back( -1.0);
410  t_nTrksJetCalo->push_back( -1.0);
411  }
412 
413  //===================================================================================
414 
415  // get handles to calogeometry and calotopology
417  iSetup.get<CaloGeometryRecord>().get(pG);
418  const CaloGeometry* geo = pG.product();
419 
420  edm::ESHandle<CaloTopology> theCaloTopology;
421  iSetup.get<CaloTopologyRecord>().get(theCaloTopology);
422  const CaloTopology *caloTopology = theCaloTopology.product();
423 
425  iSetup.get<IdealGeometryRecord>().get(htopo);
426  const HcalTopology* theHBHETopology = htopo.product();
427 
428  edm::Handle<EcalRecHitCollection> barrelRecHitsHandle;
429  edm::Handle<EcalRecHitCollection> endcapRecHitsHandle;
430  iEvent.getByToken(tok_EB_,barrelRecHitsHandle);
431  iEvent.getByToken(tok_EE_,endcapRecHitsHandle);
432 
433  // Retrieve the good/bad ECAL channels from the DB
435  iSetup.get<EcalChannelStatusRcd>().get(ecalChStatus);
436  const EcalChannelStatus* theEcalChStatus = ecalChStatus.product();
437 
438  // Retrieve trigger tower map
440  iSetup.get<IdealGeometryRecord>().get(hTtmap);
441  const EcalTrigTowerConstituentsMap& ttMap = *hTtmap;
442 
444  iEvent.getByToken(tok_hbhe_, hbhe);
445  if (!hbhe.isValid()) {
446  nbad++;
447  if (nbad < 10) std::cout << "No HBHE rechit collection\n";
448  return;
449  }
450  const HBHERecHitCollection Hithbhe = *(hbhe.product());
451 
452  //get Handles to SimTracks and SimHits
454  if (doMC) iEvent.getByToken(tok_simTk_,SimTk);
455  edm::SimTrackContainer::const_iterator simTrkItr;
456 
458  if (doMC) iEvent.getByToken(tok_simVtx_,SimVtx);
459 
460  //get Handles to PCaloHitContainers of eb/ee/hbhe
462  if (doMC) iEvent.getByToken(tok_caloEB_, pcaloeb);
463 
465  if (doMC) iEvent.getByToken(tok_caloEE_, pcaloee);
466 
468  if (doMC) iEvent.getByToken(tok_caloHH_, pcalohh);
469 
470  //associates tracker rechits/simhits to a track
471  TrackerHitAssociator* associate=0;
472  if (doMC) associate = new TrackerHitAssociator(iEvent);
473 
474  //===================================================================================
475 
476  h_nTracks->Fill(trkCollection->size());
477 
478  int nTracks = 0;
479 
480  t_nTracks = trkCollection->size();
481 
482  // get the list of DetIds closest to the impact point of track on surface calorimeters
483  std::vector<spr::propagatedTrackID> trkCaloDets;
484  spr::propagateCALO(trkCollection, geo, bField, theTrackQuality, trkCaloDets, false);
485  std::vector<spr::propagatedTrackID>::const_iterator trkDetItr;
486 
487  if(myverbose_>2) {
488  for(trkDetItr = trkCaloDets.begin(); trkDetItr != trkCaloDets.end(); trkDetItr++){
489  std::cout<<trkDetItr->trkItr->p()<<" "<<trkDetItr->trkItr->eta()<<" "<<trkDetItr->okECAL<<" ";
490  if(trkDetItr->detIdECAL.subdetId() == EcalBarrel) std::cout << (EBDetId)trkDetItr->detIdECAL <<" ";
491  else std::cout << (EEDetId)trkDetItr->detIdECAL <<" ";
492  std::cout<<trkDetItr->okHCAL<<" ";
493  if(trkDetItr->okHCAL) std::cout<<(HcalDetId)trkDetItr->detIdHCAL;
494  std::cout << std::endl;
495  }
496  }
497 
498  int nvtxTracks=0;
499  for(trkDetItr = trkCaloDets.begin(),nTracks=0; trkDetItr != trkCaloDets.end(); trkDetItr++,nTracks++){
500 
501  const reco::Track* pTrack = &(*(trkDetItr->trkItr));
502 
503  // find vertex index the track is associated with
504  int pVtxTkId = -1;
505  for(unsigned int ind=0; ind<recVtxs->size(); ind++) {
506  if (!((*recVtxs)[ind].isFake())) {
507  reco::Vertex::trackRef_iterator vtxTrack = (*recVtxs)[ind].tracks_begin();
508  for (vtxTrack = (*recVtxs)[ind].tracks_begin(); vtxTrack!=(*recVtxs)[ind].tracks_end(); vtxTrack++) {
509 
510  const edm::RefToBase<reco::Track> pvtxTrack = (*vtxTrack);
511  if ( pTrack == pvtxTrack.get() ) {
512  pVtxTkId = ind;
513  break;
514  if(myverbose_>2) {
515  if( pTrack->pt()>1.0) {
516  std::cout<<"Debug the track association with vertex "<<std::endl;
517  std::cout<< pTrack << " "<< pvtxTrack.get() << std::endl;
518  std::cout<<" trkVtxIndex "<<nvtxTracks<<" vtx "<<ind<<" pt "<<pTrack->pt()
519  <<" eta "<<pTrack->eta()<<" "<<pTrack->pt()-pvtxTrack->pt()
520  <<" "<< pTrack->eta()-pvtxTrack->eta()
521  << std::endl;
522  nvtxTracks++;
523  }
524  }
525  }
526  }
527  }
528  }
529 
530  const reco::HitPattern &hitp = pTrack->hitPattern();
531 
532  int nLayersCrossed = hitp.trackerLayersWithMeasurement();
533  int nOuterHits = hitp.stripTOBLayersWithMeasurement()
535 
536  bool ifGood = pTrack->quality(trackQuality_);
537  double pt1 = pTrack->pt();
538  double p1 = pTrack->p();
539  double eta1 = pTrack->momentum().eta();
540  double phi1 = pTrack->momentum().phi();
541  double etaEcal1 = trkDetItr->etaECAL;
542  double phiEcal1 = trkDetItr->phiECAL;
543  double etaHcal1 = trkDetItr->etaHCAL;
544  double phiHcal1 = trkDetItr->phiHCAL;
545  double dxy1 = pTrack->dxy();
546  double dz1 = pTrack->dz();
547  double chisq1 = pTrack->normalizedChi2();
548  double dxybs1 = beamSpotH.isValid() ? pTrack->dxy(bspot) : pTrack->dxy();
549  double dzbs1 = beamSpotH.isValid() ? pTrack->dz(bspot) : pTrack->dz();
550  double dxypv1 = pTrack->dxy();
551  double dzpv1 = pTrack->dz();
552  if(pVtxTkId>=0) {
553  math::XYZPoint thisTkPV = math::XYZPoint( (*recVtxs)[pVtxTkId].x(),(*recVtxs)[pVtxTkId].y(), (*recVtxs)[pVtxTkId].z() );
554  dxypv1 = pTrack->dxy(thisTkPV);
555  dzpv1 = pTrack->dz (thisTkPV);
556  }
557 
558  h_recEtaPt_0->Fill(eta1, pt1);
559  h_recEtaP_0 ->Fill(eta1, p1);
560  h_recPt_0 ->Fill(pt1);
561  h_recP_0 ->Fill(p1);
562  h_recEta_0 ->Fill(eta1);
563  h_recPhi_0 ->Fill(phi1);
564 
565  if(ifGood && nLayersCrossed>7 ) {
566  h_recEtaPt_1->Fill(eta1, pt1);
567  h_recEtaP_1 ->Fill(eta1, p1);
568  h_recPt_1 ->Fill(pt1);
569  h_recP_1 ->Fill(p1);
570  h_recEta_1 ->Fill(eta1);
571  h_recPhi_1 ->Fill(phi1);
572  }
573 
574  if( ! ifGood ) continue;
575 
576  if( writeAllTracks && p1>2.0 && nLayersCrossed>7) {
577  t_trackPAll ->push_back( p1 );
578  t_trackEtaAll ->push_back( eta1 );
579  t_trackPhiAll ->push_back( phi1 );
580  t_trackPtAll ->push_back( pt1 );
581  t_trackDxyAll ->push_back( dxy1 );
582  t_trackDzAll ->push_back( dz1 );
583  t_trackDxyPVAll ->push_back( dxypv1 );
584  t_trackDzPVAll ->push_back( dzpv1 );
585  t_trackChiSqAll ->push_back( chisq1 );
586  }
587  if (doMC) {
588  edm::SimTrackContainer::const_iterator matchedSimTrkAll = spr::matchedSimTrack(iEvent, SimTk, SimVtx, pTrack, *associate, false);
589  if( writeAllTracks && matchedSimTrkAll != SimTk->end()) t_trackPdgIdAll->push_back( matchedSimTrkAll->type() );
590  }
591 
592  if( pt1>minTrackP_ && std::abs(eta1)<maxTrackEta_ && trkDetItr->okECAL) {
593 
594  double maxNearP31x31=999.0, maxNearP25x25=999.0, maxNearP21x21=999.0, maxNearP15x15=999.0;
595  maxNearP31x31 = spr::chargeIsolationEcal(nTracks, trkCaloDets, geo, caloTopology, 15,15);
596  maxNearP25x25 = spr::chargeIsolationEcal(nTracks, trkCaloDets, geo, caloTopology, 12,12);
597  maxNearP21x21 = spr::chargeIsolationEcal(nTracks, trkCaloDets, geo, caloTopology, 10,10);
598  maxNearP15x15 = spr::chargeIsolationEcal(nTracks, trkCaloDets, geo, caloTopology, 7, 7);
599 
600  int iTrkEtaBin=-1, iTrkMomBin=-1;
601  for(unsigned int ieta=0; ieta<NEtaBins; ieta++) {
602  if(std::abs(eta1)>genPartEtaBins[ieta] && std::abs(eta1)<genPartEtaBins[ieta+1] ) iTrkEtaBin = ieta;
603  }
604  for(unsigned int ipt=0; ipt<NPBins; ipt++) {
605  if( p1>genPartPBins[ipt] && p1<genPartPBins[ipt+1] ) iTrkMomBin = ipt;
606  }
607  if( iTrkMomBin>=0 && iTrkEtaBin>=0 ) {
608  h_maxNearP31x31[iTrkMomBin][iTrkEtaBin]->Fill( maxNearP31x31 );
609  h_maxNearP25x25[iTrkMomBin][iTrkEtaBin]->Fill( maxNearP25x25 );
610  h_maxNearP21x21[iTrkMomBin][iTrkEtaBin]->Fill( maxNearP21x21 );
611  h_maxNearP15x15[iTrkMomBin][iTrkEtaBin]->Fill( maxNearP15x15 );
612  }
613  if( maxNearP31x31<0.0 && nLayersCrossed>7 && nOuterHits>4) {
614  h_recEtaPt_2->Fill(eta1, pt1);
615  h_recEtaP_2 ->Fill(eta1, p1);
616  h_recPt_2 ->Fill(pt1);
617  h_recP_2 ->Fill(p1);
618  h_recEta_2 ->Fill(eta1);
619  h_recPhi_2 ->Fill(phi1);
620  }
621 
622  // if charge isolated in ECAL, store the further quantities
623  if( maxNearP31x31<1.0) {
624 
625  haveIsoTrack = true;
626 
627  // get the matching simTrack
628  double simTrackP = -1;
629  if (doMC) {
630  edm::SimTrackContainer::const_iterator matchedSimTrk = spr::matchedSimTrack(iEvent, SimTk, SimVtx, pTrack, *associate, false);
631  if( matchedSimTrk != SimTk->end() )simTrackP = matchedSimTrk->momentum().P();
632  }
633 
635  iSetup.get<EcalSeverityLevelAlgoRcd>().get(sevlv);
636 
637  // get ECal Tranverse Profile
638  std::pair<double, bool> e7x7P, e9x9P, e11x11P, e15x15P;
639  std::pair<double, bool> e7x7_10SigP, e9x9_10SigP, e11x11_10SigP, e15x15_10SigP;
640  std::pair<double, bool> e7x7_15SigP, e9x9_15SigP, e11x11_15SigP, e15x15_15SigP;
641  std::pair<double, bool> e7x7_20SigP, e9x9_20SigP, e11x11_20SigP, e15x15_20SigP;
642  std::pair<double, bool> e7x7_25SigP, e9x9_25SigP, e11x11_25SigP, e15x15_25SigP;
643  std::pair<double, bool> e7x7_30SigP, e9x9_30SigP, e11x11_30SigP, e15x15_30SigP;
644 
645  spr::caloSimInfo simInfo3x3, simInfo5x5, simInfo7x7, simInfo9x9;
646  spr::caloSimInfo simInfo11x11, simInfo13x13, simInfo15x15, simInfo21x21, simInfo25x25, simInfo31x31;
647  double trkEcalEne=0;
648 
649  const DetId isoCell = trkDetItr->detIdECAL;
650  e7x7P = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),3,3, -100.0, -100.0, tMinE_,tMaxE_);
651  e9x9P = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),4,4, -100.0, -100.0, tMinE_,tMaxE_);
652  e11x11P = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),5,5, -100.0, -100.0, tMinE_,tMaxE_);
653  e15x15P = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),7,7, -100.0, -100.0, tMinE_,tMaxE_);
654 
655  e7x7_10SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),3,3, 0.030, 0.150, tMinE_,tMaxE_);
656  e9x9_10SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),4,4, 0.030, 0.150, tMinE_,tMaxE_);
657  e11x11_10SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),5,5, 0.030, 0.150, tMinE_,tMaxE_);
658  e15x15_10SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),7,7, 0.030, 0.150, tMinE_,tMaxE_);
659 
660  e7x7_15SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology, sevlv.product(),ttMap, 3,3, 0.20,0.45, tMinE_,tMaxE_);
661  e9x9_15SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology, sevlv.product(),ttMap, 4,4, 0.20,0.45, tMinE_,tMaxE_);
662  e11x11_15SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(), ttMap, 5,5, 0.20,0.45, tMinE_,tMaxE_);
663  e15x15_15SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology, sevlv.product(),ttMap, 7,7, 0.20,0.45, tMinE_,tMaxE_, false);
664 
665  e7x7_20SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),3,3, 0.060, 0.300, tMinE_,tMaxE_);
666  e9x9_20SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),4,4, 0.060, 0.300, tMinE_,tMaxE_);
667  e11x11_20SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),5,5, 0.060, 0.300, tMinE_,tMaxE_);
668  e15x15_20SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),7,7, 0.060, 0.300, tMinE_,tMaxE_);
669 
670  e7x7_25SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),3,3, 0.075, 0.375, tMinE_,tMaxE_);
671  e9x9_25SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),4,4, 0.075, 0.375, tMinE_,tMaxE_);
672  e11x11_25SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),5,5, 0.075, 0.375, tMinE_,tMaxE_);
673  e15x15_25SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),7,7, 0.075, 0.375, tMinE_,tMaxE_);
674 
675  e7x7_30SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),3,3, 0.090, 0.450, tMinE_,tMaxE_);
676  e9x9_30SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),4,4, 0.090, 0.450, tMinE_,tMaxE_);
677  e11x11_30SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),5,5, 0.090, 0.450, tMinE_,tMaxE_);
678  e15x15_30SigP = spr::eECALmatrix(isoCell,barrelRecHitsHandle,endcapRecHitsHandle, *theEcalChStatus, geo, caloTopology,sevlv.product(),7,7, 0.090, 0.450, tMinE_,tMaxE_);
679  if(myverbose_ == 2) {
680  std::cout << "clean ecal rechit " << std::endl;
681  std::cout<<"e7x7 "<<e7x7P.first<<" e9x9 "<<e9x9P.first<<" e11x11 " << e11x11P.first << " e15x15 "<<e15x15P.first<<std::endl;
682  std::cout<<"e7x7_10Sig "<<e7x7_10SigP.first<<" e11x11_10Sig "<<e11x11_10SigP.first<<" e15x15_10Sig "<<e15x15_10SigP.first<<std::endl;
683  }
684 
685  if (doMC) {
686  // check the energy from SimHits
687  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 1,1, simInfo3x3);
688  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 2,2, simInfo5x5);
689  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 3,3, simInfo7x7);
690  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 4,4, simInfo9x9);
691  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 5,5, simInfo11x11);
692  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 6,6, simInfo13x13);
693  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 7,7, simInfo15x15, 150.0,false);
694  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 10,10, simInfo21x21);
695  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 12,12, simInfo25x25);
696  spr::eECALSimInfo(iEvent,isoCell,geo,caloTopology,pcaloeb,pcaloee,SimTk,SimVtx,pTrack, *associate, 15,15, simInfo31x31);
697 
698  trkEcalEne = spr::eCaloSimInfo(iEvent, geo, pcaloeb,pcaloee, SimTk, SimVtx, pTrack, *associate);
699  if(myverbose_ == 1) {
700  std::cout << "Track momentum " << pt1 << std::endl;
701 
702  std::cout << "ecal siminfo " << std::endl;
703  std::cout << "simInfo3x3: " << "eTotal " << simInfo3x3.eTotal << " eMatched " << simInfo3x3.eMatched << " eRest " << simInfo3x3.eRest << " eGamma "<<simInfo3x3.eGamma << " eNeutralHad " << simInfo3x3.eNeutralHad << " eChargedHad " << simInfo3x3.eChargedHad << std::endl;
704  std::cout << "simInfo5x5: " << "eTotal " << simInfo5x5.eTotal << " eMatched " << simInfo5x5.eMatched << " eRest " << simInfo5x5.eRest << " eGamma "<<simInfo5x5.eGamma << " eNeutralHad " << simInfo5x5.eNeutralHad << " eChargedHad " << simInfo5x5.eChargedHad << std::endl;
705  std::cout << "simInfo7x7: " << "eTotal " << simInfo7x7.eTotal << " eMatched " << simInfo7x7.eMatched << " eRest " << simInfo7x7.eRest << " eGamma "<<simInfo7x7.eGamma << " eNeutralHad " << simInfo7x7.eNeutralHad << " eChargedHad " << simInfo7x7.eChargedHad << std::endl;
706  std::cout << "simInfo9x9: " << "eTotal " << simInfo9x9.eTotal << " eMatched " << simInfo9x9.eMatched << " eRest " << simInfo9x9.eRest << " eGamma "<<simInfo9x9.eGamma << " eNeutralHad " << simInfo9x9.eNeutralHad << " eChargedHad " << simInfo9x9.eChargedHad << std::endl;
707  std::cout << "simInfo11x11: " << "eTotal " << simInfo11x11.eTotal << " eMatched " << simInfo11x11.eMatched << " eRest " << simInfo11x11.eRest << " eGamma "<<simInfo11x11.eGamma << " eNeutralHad " << simInfo11x11.eNeutralHad << " eChargedHad " << simInfo11x11.eChargedHad << std::endl;
708  std::cout << "simInfo15x15: " << "eTotal " << simInfo15x15.eTotal << " eMatched " << simInfo15x15.eMatched << " eRest " << simInfo15x15.eRest << " eGamma "<<simInfo15x15.eGamma << " eNeutralHad " << simInfo15x15.eNeutralHad << " eChargedHad " << simInfo15x15.eChargedHad << std::endl;
709  std::cout << "simInfo31x31: " << "eTotal " << simInfo31x31.eTotal << " eMatched " << simInfo31x31.eMatched << " eRest " << simInfo31x31.eRest << " eGamma "<<simInfo31x31.eGamma << " eNeutralHad " << simInfo31x31.eNeutralHad << " eChargedHad " << simInfo31x31.eChargedHad << std::endl;
710  std::cout << "trkEcalEne" << trkEcalEne << std::endl;
711  }
712  }
713 
714  // ======= Get HCAL information
715  double hcalScale=1.0;
716  if( std::abs(pTrack->eta())<1.4 ) {
717  hcalScale=120.0;
718  } else {
719  hcalScale=135.0;
720  }
721 
722  double maxNearHcalP3x3=-1, maxNearHcalP5x5=-1, maxNearHcalP7x7=-1;
723  maxNearHcalP3x3 = spr::chargeIsolationHcal(nTracks, trkCaloDets, theHBHETopology, 1,1);
724  maxNearHcalP5x5 = spr::chargeIsolationHcal(nTracks, trkCaloDets, theHBHETopology, 2,2);
725  maxNearHcalP7x7 = spr::chargeIsolationHcal(nTracks, trkCaloDets, theHBHETopology, 3,3);
726 
727  double h3x3=0, h5x5=0, h7x7=0;
728  double h3x3Sig=0, h5x5Sig=0, h7x7Sig=0;
729  double trkHcalEne = 0;
730  spr::caloSimInfo hsimInfo3x3, hsimInfo5x5, hsimInfo7x7;
731 
732  if(trkDetItr->okHCAL) {
733  const DetId ClosestCell(trkDetItr->detIdHCAL);
734  // bool includeHO=false, bool algoNew=true, bool debug=false
735  h3x3 = spr::eHCALmatrix(theHBHETopology, ClosestCell, hbhe,1,1, false, true, -100.0, -100.0, -100.0, -100.0, tMinH_,tMaxH_);
736  h5x5 = spr::eHCALmatrix(theHBHETopology, ClosestCell, hbhe,2,2, false, true, -100.0, -100.0, -100.0, -100.0, tMinH_,tMaxH_);
737  h7x7 = spr::eHCALmatrix(theHBHETopology, ClosestCell, hbhe,3,3, false, true, -100.0, -100.0, -100.0, -100.0, tMinH_,tMaxH_);
738  h3x3Sig = spr::eHCALmatrix(theHBHETopology, ClosestCell, hbhe,1,1, false, true, 0.7, 0.8, -100.0, -100.0, tMinH_,tMaxH_);
739  h5x5Sig = spr::eHCALmatrix(theHBHETopology, ClosestCell, hbhe,2,2, false, true, 0.7, 0.8, -100.0, -100.0, tMinH_,tMaxH_);
740  h7x7Sig = spr::eHCALmatrix(theHBHETopology, ClosestCell, hbhe,3,3, false, true, 0.7, 0.8, -100.0, -100.0, tMinH_,tMaxH_);
741 
742  if(myverbose_==2) {
743  std::cout << "HCAL 3x3 " << h3x3 << " " << h3x3Sig << " 5x5 " << h5x5 << " " << h5x5Sig << " 7x7 " << h7x7 << " " << h7x7Sig << std::endl;
744  }
745 
746  if (doMC) {
747  spr::eHCALSimInfo(iEvent, theHBHETopology, ClosestCell, geo,pcalohh, SimTk, SimVtx, pTrack, *associate, 1,1, hsimInfo3x3);
748  spr::eHCALSimInfo(iEvent, theHBHETopology, ClosestCell, geo,pcalohh, SimTk, SimVtx, pTrack, *associate, 2,2, hsimInfo5x5);
749  spr::eHCALSimInfo(iEvent, theHBHETopology, ClosestCell, geo,pcalohh, SimTk, SimVtx, pTrack, *associate, 3,3, hsimInfo7x7, 150.0, false,false);
750  trkHcalEne = spr::eCaloSimInfo(iEvent, geo,pcalohh, SimTk, SimVtx, pTrack, *associate);
751  if(myverbose_ == 1) {
752  std::cout << "Hcal siminfo " << std::endl;
753  std::cout << "hsimInfo3x3: " << "eTotal " << hsimInfo3x3.eTotal << " eMatched " << hsimInfo3x3.eMatched << " eRest " << hsimInfo3x3.eRest << " eGamma "<<hsimInfo3x3.eGamma << " eNeutralHad " << hsimInfo3x3.eNeutralHad << " eChargedHad " << hsimInfo3x3.eChargedHad << std::endl;
754  std::cout << "hsimInfo5x5: " << "eTotal " << hsimInfo5x5.eTotal << " eMatched " << hsimInfo5x5.eMatched << " eRest " << hsimInfo5x5.eRest << " eGamma "<<hsimInfo5x5.eGamma << " eNeutralHad " << hsimInfo5x5.eNeutralHad << " eChargedHad " << hsimInfo5x5.eChargedHad << std::endl;
755  std::cout << "hsimInfo7x7: " << "eTotal " << hsimInfo7x7.eTotal << " eMatched " << hsimInfo7x7.eMatched << " eRest " << hsimInfo7x7.eRest << " eGamma "<<hsimInfo7x7.eGamma << " eNeutralHad " << hsimInfo7x7.eNeutralHad << " eChargedHad " << hsimInfo7x7.eChargedHad << std::endl;
756  std::cout << "trkHcalEne " << trkHcalEne << std::endl;
757  }
758  }
759 
760  // debug the ecal and hcal matrix
761  if(myverbose_==4) {
762  std::cout<<"Run "<<iEvent.id().run()<<" Event "<<iEvent.id().event()<<std::endl;
763  std::vector<std::pair<DetId,double> > v7x7 = spr::eHCALmatrixCell(theHBHETopology, ClosestCell, hbhe,3,3, false, false);
764  double sumv=0.0;
765 
766  for(unsigned int iv=0; iv<v7x7.size(); iv++) {
767  sumv += v7x7[iv].second;
768  }
769  std::cout<<"h7x7 "<<h7x7<<" v7x7 "<<sumv << " in " << v7x7.size() <<std::endl;
770  for(unsigned int iv=0; iv<v7x7.size(); iv++) {
771  HcalDetId id = v7x7[iv].first;
772  std::cout << " Cell " << iv << " 0x" << std::hex << v7x7[iv].first() << std::dec << " " << id << " Energy " << v7x7[iv].second << std::endl;
773  }
774  }
775 
776  }
777  if (doMC) {
778  trkHcalEne = spr::eCaloSimInfo(iEvent, geo,pcalohh, SimTk, SimVtx, pTrack, *associate);
779  }
780 
781  // ====================================================================================================
782  // get diff between track outermost hit position and the propagation point at outermost surface of tracker
783  std::pair<math::XYZPoint,double> point2_TK0 = spr::propagateTrackerEnd( pTrack, bField, false);
784  math::XYZPoint diff(pTrack->outerPosition().X()-point2_TK0.first.X(),
785  pTrack->outerPosition().Y()-point2_TK0.first.Y(),
786  pTrack->outerPosition().Z()-point2_TK0.first.Z() );
787  double trackOutPosOutHitDr = diff.R();
788  double trackL = point2_TK0.second;
789  if (myverbose_==5) {
790  std::cout<<" propagted "<<point2_TK0.first<<" "<< point2_TK0.first.eta()<<" "<<point2_TK0.first.phi()<<std::endl;
791  std::cout<<" outerPosition() "<< pTrack->outerPosition() << " "<< pTrack->outerPosition().eta()<< " " << pTrack->outerPosition().phi()<< std::endl;
792  std::cout<<"diff " << diff << " diffR " <<diff.R()<<" diffR/L "<<diff.R()/point2_TK0.second <<std::endl;
793  }
794 
795  for(unsigned int ind=0;ind<recVtxs->size();ind++) {
796  if (!((*recVtxs)[ind].isFake())) {
797  reco::Vertex::trackRef_iterator vtxTrack = (*recVtxs)[ind].tracks_begin();
798  if( DeltaR(eta1,phi1, (*vtxTrack)->eta(),(*vtxTrack)->phi()) < 0.01 ) t_trackPVIdx ->push_back( ind );
799  else t_trackPVIdx ->push_back( -1 );
800  }
801  }
802 
803  // Fill the tree Branches here
804  t_trackPVIdx ->push_back( pVtxTkId );
805  t_trackP ->push_back( p1 );
806  t_trackPt ->push_back( pt1 );
807  t_trackEta ->push_back( eta1 );
808  t_trackPhi ->push_back( phi1 );
809  t_trackEcalEta ->push_back( etaEcal1 );
810  t_trackEcalPhi ->push_back( phiEcal1 );
811  t_trackHcalEta ->push_back( etaHcal1 );
812  t_trackHcalPhi ->push_back( phiHcal1 );
813  t_trackDxy ->push_back( dxy1 );
814  t_trackDz ->push_back( dz1 );
815  t_trackDxyBS ->push_back( dxybs1 );
816  t_trackDzBS ->push_back( dzbs1 );
817  t_trackDxyPV ->push_back( dxypv1 );
818  t_trackDzPV ->push_back( dzpv1 );
819  t_trackChiSq ->push_back( chisq1 );
820  t_trackNOuterHits ->push_back( nOuterHits );
821  t_NLayersCrossed ->push_back( nLayersCrossed );
822  using namespace reco;
823  t_trackHitsTOB ->push_back(hitp.stripTOBLayersWithMeasurement());
824  t_trackHitsTEC ->push_back(hitp.stripTECLayersWithMeasurement());
825 
826  t_trackHitInMissTOB ->push_back(hitp.stripTOBLayersWithoutMeasurement(HitPattern::MISSING_INNER_HITS));
827  t_trackHitInMissTEC ->push_back(hitp.stripTECLayersWithoutMeasurement(HitPattern::MISSING_INNER_HITS));
828  t_trackHitInMissTIB ->push_back(hitp.stripTIBLayersWithoutMeasurement(HitPattern::MISSING_INNER_HITS));
829  t_trackHitInMissTID ->push_back(hitp.stripTIDLayersWithoutMeasurement(HitPattern::MISSING_INNER_HITS));
830  t_trackHitInMissTIBTID ->push_back(hitp.stripTIBLayersWithoutMeasurement(HitPattern::MISSING_INNER_HITS)
831  + hitp.stripTIDLayersWithoutMeasurement(HitPattern::MISSING_INNER_HITS));
832 
833  t_trackHitOutMissTOB ->push_back(hitp.stripTOBLayersWithoutMeasurement(HitPattern::MISSING_OUTER_HITS));
834  t_trackHitOutMissTEC ->push_back(hitp.stripTECLayersWithoutMeasurement(HitPattern::MISSING_OUTER_HITS));
835  t_trackHitOutMissTIB ->push_back(hitp.stripTIBLayersWithoutMeasurement(HitPattern::MISSING_OUTER_HITS));
836  t_trackHitOutMissTID ->push_back(hitp.stripTIDLayersWithoutMeasurement(HitPattern::MISSING_OUTER_HITS));
837  t_trackHitOutMissTOBTEC ->push_back(hitp.stripTOBLayersWithoutMeasurement(HitPattern::MISSING_OUTER_HITS)
838  + hitp.stripTECLayersWithoutMeasurement(HitPattern::MISSING_OUTER_HITS));
839 
840  // The following push_back(0) where calls to methods XYZLayersWithMeasurement() on
841  // HiPatterns ExpectedMissingInner and ExpectedMissingOuter. Since those HitPatterns
842  // only contain MISSING_HITS, any call to XYZLayersWithMeasuremnt() whould return 0.
843  // As example:
844  // t_trackHitInMeasTOB ->push_back( hitpIn.stripTOBLayersWithMeasurement() );
845  // where hitpIn was expectedInnerHits.
846  t_trackHitInMeasTOB ->push_back(0);
847  t_trackHitInMeasTEC ->push_back(0);
848  t_trackHitInMeasTIB ->push_back(0);
849  t_trackHitInMeasTID ->push_back(0);
850 
851  t_trackHitOutMeasTOB ->push_back(0);
852  t_trackHitOutMeasTEC ->push_back(0);
853  t_trackHitOutMeasTIB ->push_back(0);
854  t_trackHitOutMeasTID ->push_back(0);
855 
856  t_trackOutPosOutHitDr ->push_back(trackOutPosOutHitDr );
857  t_trackL ->push_back(trackL );
858 
859  t_maxNearP31x31 ->push_back( maxNearP31x31 );
860  t_maxNearP21x21 ->push_back( maxNearP21x21 );
861 
862  t_ecalSpike11x11 ->push_back( e11x11P.second );
863  t_e7x7 ->push_back( e7x7P.first );
864  t_e9x9 ->push_back( e9x9P.first );
865  t_e11x11 ->push_back( e11x11P.first );
866  t_e15x15 ->push_back( e15x15P.first );
867 
868  t_e7x7_10Sig ->push_back( e7x7_10SigP.first );
869  t_e9x9_10Sig ->push_back( e9x9_10SigP.first );
870  t_e11x11_10Sig ->push_back( e11x11_10SigP.first );
871  t_e15x15_10Sig ->push_back( e15x15_10SigP.first );
872  t_e7x7_15Sig ->push_back( e7x7_15SigP.first );
873  t_e9x9_15Sig ->push_back( e9x9_15SigP.first );
874  t_e11x11_15Sig ->push_back( e11x11_15SigP.first );
875  t_e15x15_15Sig ->push_back( e15x15_15SigP.first );
876  t_e7x7_20Sig ->push_back( e7x7_20SigP.first );
877  t_e9x9_20Sig ->push_back( e9x9_20SigP.first );
878  t_e11x11_20Sig ->push_back( e11x11_20SigP.first );
879  t_e15x15_20Sig ->push_back( e15x15_20SigP.first );
880  t_e7x7_25Sig ->push_back( e7x7_25SigP.first );
881  t_e9x9_25Sig ->push_back( e9x9_25SigP.first );
882  t_e11x11_25Sig ->push_back( e11x11_25SigP.first );
883  t_e15x15_25Sig ->push_back( e15x15_25SigP.first );
884  t_e7x7_30Sig ->push_back( e7x7_30SigP.first );
885  t_e9x9_30Sig ->push_back( e9x9_30SigP.first );
886  t_e11x11_30Sig ->push_back( e11x11_30SigP.first );
887  t_e15x15_30Sig ->push_back( e15x15_30SigP.first );
888 
889  if (doMC) {
890  t_esim7x7 ->push_back( simInfo7x7.eTotal );
891  t_esim9x9 ->push_back( simInfo9x9.eTotal );
892  t_esim11x11 ->push_back( simInfo11x11.eTotal );
893  t_esim15x15 ->push_back( simInfo15x15.eTotal );
894 
895  t_esim7x7Matched ->push_back( simInfo7x7.eMatched );
896  t_esim9x9Matched ->push_back( simInfo9x9.eMatched );
897  t_esim11x11Matched ->push_back( simInfo11x11.eMatched );
898  t_esim15x15Matched ->push_back( simInfo15x15.eMatched );
899 
900  t_esim7x7Rest ->push_back( simInfo7x7.eRest );
901  t_esim9x9Rest ->push_back( simInfo9x9.eRest );
902  t_esim11x11Rest ->push_back( simInfo11x11.eRest );
903  t_esim15x15Rest ->push_back( simInfo15x15.eRest );
904 
905  t_esim7x7Photon ->push_back( simInfo7x7.eGamma );
906  t_esim9x9Photon ->push_back( simInfo9x9.eGamma );
907  t_esim11x11Photon ->push_back( simInfo11x11.eGamma );
908  t_esim15x15Photon ->push_back( simInfo15x15.eGamma );
909 
910  t_esim7x7NeutHad ->push_back( simInfo7x7.eNeutralHad );
911  t_esim9x9NeutHad ->push_back( simInfo9x9.eNeutralHad );
912  t_esim11x11NeutHad ->push_back( simInfo11x11.eNeutralHad );
913  t_esim15x15NeutHad ->push_back( simInfo15x15.eNeutralHad );
914 
915  t_esim7x7CharHad ->push_back( simInfo7x7.eChargedHad );
916  t_esim9x9CharHad ->push_back( simInfo9x9.eChargedHad );
917  t_esim11x11CharHad ->push_back( simInfo11x11.eChargedHad );
918  t_esim15x15CharHad ->push_back( simInfo15x15.eChargedHad );
919 
920  t_trkEcalEne ->push_back( trkEcalEne );
921  t_simTrackP ->push_back( simTrackP );
922  t_esimPdgId ->push_back( simInfo11x11.pdgMatched );
923  }
924 
925  t_maxNearHcalP3x3 ->push_back( maxNearHcalP3x3 );
926  t_maxNearHcalP5x5 ->push_back( maxNearHcalP5x5 );
927  t_maxNearHcalP7x7 ->push_back( maxNearHcalP7x7 );
928 
929  t_h3x3 ->push_back( h3x3 );
930  t_h5x5 ->push_back( h5x5 );
931  t_h7x7 ->push_back( h7x7 );
932  t_h3x3Sig ->push_back( h3x3Sig );
933  t_h5x5Sig ->push_back( h5x5Sig );
934  t_h7x7Sig ->push_back( h7x7Sig );
935 
936  t_infoHcal ->push_back( trkDetItr->okHCAL );
937  if (doMC) {
938  t_trkHcalEne ->push_back( hcalScale*trkHcalEne );
939 
940  t_hsim3x3 ->push_back( hcalScale*hsimInfo3x3.eTotal );
941  t_hsim5x5 ->push_back( hcalScale*hsimInfo5x5.eTotal );
942  t_hsim7x7 ->push_back( hcalScale*hsimInfo7x7.eTotal );
943 
944  t_hsim3x3Matched ->push_back( hcalScale*hsimInfo3x3.eMatched );
945  t_hsim5x5Matched ->push_back( hcalScale*hsimInfo5x5.eMatched );
946  t_hsim7x7Matched ->push_back( hcalScale*hsimInfo7x7.eMatched );
947 
948  t_hsim3x3Rest ->push_back( hcalScale*hsimInfo3x3.eRest );
949  t_hsim5x5Rest ->push_back( hcalScale*hsimInfo5x5.eRest );
950  t_hsim7x7Rest ->push_back( hcalScale*hsimInfo7x7.eRest );
951 
952  t_hsim3x3Photon ->push_back( hcalScale*hsimInfo3x3.eGamma );
953  t_hsim5x5Photon ->push_back( hcalScale*hsimInfo5x5.eGamma );
954  t_hsim7x7Photon ->push_back( hcalScale*hsimInfo7x7.eGamma );
955 
956  t_hsim3x3NeutHad ->push_back( hcalScale*hsimInfo3x3.eNeutralHad );
957  t_hsim5x5NeutHad ->push_back( hcalScale*hsimInfo5x5.eNeutralHad );
958  t_hsim7x7NeutHad ->push_back( hcalScale*hsimInfo7x7.eNeutralHad );
959 
960  t_hsim3x3CharHad ->push_back( hcalScale*hsimInfo3x3.eChargedHad );
961  t_hsim5x5CharHad ->push_back( hcalScale*hsimInfo5x5.eChargedHad );
962  t_hsim7x7CharHad ->push_back( hcalScale*hsimInfo7x7.eChargedHad );
963  }
964  /*
965  if(hcalScale*hsimInfo3x3["eTotal"] > 50.0) {
966 
967  std::cout << "Loosely Iso Track : eta " << eta1 << " Rec Mom " << p1 << " SimMom " << simTrackP << " h3x3 " << h3x3 << std::endl;
968 
969  std::cout <<"Closest cell Hcal (atHCAL) " << (HcalDetId)ClosestCell << std::endl;
970 
971  std::cout <<"trkHcalEne, etotal, matched, rest " <<hcalScale*trkHcalEne<<std::setw(15)<<hcalScale*hsimInfo3x3["eTotal"]
972  <<std::setw(15)<<hcalScale*hsimInfo3x3["eMatched"]<<std::setw(15)<<hcalScale*hsimInfo3x3["eRest"]
973  <<std::endl;
974  unsigned int nn = t_trkHcalEne->size();
975  std::cout <<"in Tree " << (*t_trkHcalEne)[nn-1] <<std::setw(15)<< (*t_hsim3x3)[nn-1]
976  <<std::setw(15)<< (*t_hsim3x3Matched)[nn-1] <<std::setw(15)<< (*t_hsim3x3Rest)[nn-1]
977  << std::endl;
978 
979  std::cout << "debug output \n" << std::endl;
980  std::map<std::string, double> hsimInfo3x3_debug = spr::eHCALSimInfo(iEvent, theHBHETopology, ClosestCell, geo,pcalohh, SimTk, SimVtx, pTrack, *associate, 1,1, 150.0, true);
981 
982  }
983  */
984 
985 
986  } // if loosely isolated track
987  } // check p1/eta
988  } // loop over track collection
989 
990  if(haveIsoTrack) tree->Fill();
991 }
992 
993 // ----- method called once each job just before starting event loop ----
995 
996  nEventProc=0;
997 
998  // double tempgen_TH[21] = { 1.0, 2.0, 3.0, 4.0, 5.0,
999  double tempgen_TH[16] = { 0.0, 1.0, 2.0, 3.0, 4.0,
1000  5.0, 6.0, 7.0, 9.0, 11.0,
1001  15.0, 20.0, 30.0, 50.0, 75.0, 100.0};
1002 
1003  for(int i=0; i<16; i++) genPartPBins[i] = tempgen_TH[i];
1004 
1005  double tempgen_Eta[4] = {0.0, 1.131, 1.653, 2.172};
1006 
1007  for(int i=0; i<4; i++) genPartEtaBins[i] = tempgen_Eta[i];
1008 
1009  BookHistograms();
1010 }
1011 
1012 // ----- method called once each job just after ending the event loop ----
1014 
1015  if(L1TriggerAlgoInfo_) {
1016  std::map< std::pair<unsigned int,std::string>, int>::iterator itr;
1017  for(itr=l1AlgoMap.begin(); itr!=l1AlgoMap.end(); itr++) {
1018  std::cout << " ****endjob**** " << (itr->first).first <<" "
1019  <<(itr->first).second <<" "<<itr->second
1020  << std::endl;
1021  int ibin = (itr->first).first;
1022  TString name( (itr->first).second );
1023  h_L1AlgoNames->GetXaxis()->SetBinLabel(ibin+1,name);
1024  }
1025  std::cout << "Number of Events Processed " << nEventProc << std::endl;
1026  }
1027 
1028 }
1029 
1030 //========================================================================================================
1031 
1033 
1034  t_PVx ->clear();
1035  t_PVy ->clear();
1036  t_PVz ->clear();
1037  t_PVisValid ->clear();
1038  t_PVndof ->clear();
1039  t_PVNTracks ->clear();
1040  t_PVNTracksWt ->clear();
1041  t_PVTracksSumPt ->clear();
1042  t_PVTracksSumPtWt ->clear();
1043  t_PVNTracksHP ->clear();
1044  t_PVNTracksHPWt ->clear();
1045  t_PVTracksSumPtHP ->clear();
1046  t_PVTracksSumPtHPWt ->clear();
1047 
1048  for(int i=0; i<128; i++) t_L1Decision[i]=0;
1049  t_L1AlgoNames ->clear();
1050  t_L1PreScale ->clear();
1051 
1052  t_L1CenJetPt ->clear();
1053  t_L1CenJetEta ->clear();
1054  t_L1CenJetPhi ->clear();
1055  t_L1FwdJetPt ->clear();
1056  t_L1FwdJetEta ->clear();
1057  t_L1FwdJetPhi ->clear();
1058  t_L1TauJetPt ->clear();
1059  t_L1TauJetEta ->clear();
1060  t_L1TauJetPhi ->clear();
1061  t_L1MuonPt ->clear();
1062  t_L1MuonEta ->clear();
1063  t_L1MuonPhi ->clear();
1064  t_L1IsoEMPt ->clear();
1065  t_L1IsoEMEta ->clear();
1066  t_L1IsoEMPhi ->clear();
1067  t_L1NonIsoEMPt ->clear();
1068  t_L1NonIsoEMEta ->clear();
1069  t_L1NonIsoEMPhi ->clear();
1070  t_L1METPt ->clear();
1071  t_L1METEta ->clear();
1072  t_L1METPhi ->clear();
1073 
1074  t_jetPt ->clear();
1075  t_jetEta ->clear();
1076  t_jetPhi ->clear();
1077  t_nTrksJetCalo ->clear();
1078  t_nTrksJetVtx ->clear();
1079 
1080  t_trackPAll ->clear();
1081  t_trackEtaAll ->clear();
1082  t_trackPhiAll ->clear();
1083  t_trackPdgIdAll ->clear();
1084  t_trackPtAll ->clear();
1085  t_trackDxyAll ->clear();
1086  t_trackDzAll ->clear();
1087  t_trackDxyPVAll ->clear();
1088  t_trackDzPVAll ->clear();
1089  t_trackChiSqAll ->clear();
1090 
1091  t_trackP ->clear();
1092  t_trackPt ->clear();
1093  t_trackEta ->clear();
1094  t_trackPhi ->clear();
1095  t_trackEcalEta ->clear();
1096  t_trackEcalPhi ->clear();
1097  t_trackHcalEta ->clear();
1098  t_trackHcalPhi ->clear();
1099  t_NLayersCrossed ->clear();
1100  t_trackNOuterHits ->clear();
1101  t_trackDxy ->clear();
1102  t_trackDxyBS ->clear();
1103  t_trackDz ->clear();
1104  t_trackDzBS ->clear();
1105  t_trackDxyPV ->clear();
1106  t_trackDzPV ->clear();
1107  t_trackChiSq ->clear();
1108  t_trackPVIdx ->clear();
1109  t_trackHitsTOB ->clear();
1110  t_trackHitsTEC ->clear();
1111  t_trackHitInMissTOB ->clear();
1112  t_trackHitInMissTEC ->clear();
1113  t_trackHitInMissTIB ->clear();
1114  t_trackHitInMissTID ->clear();
1115  t_trackHitInMissTIBTID ->clear();
1116  t_trackHitOutMissTOB ->clear();
1117  t_trackHitOutMissTEC ->clear();
1118  t_trackHitOutMissTIB ->clear();
1119  t_trackHitOutMissTID ->clear();
1120  t_trackHitOutMissTOBTEC ->clear();
1121  t_trackHitInMeasTOB ->clear();
1122  t_trackHitInMeasTEC ->clear();
1123  t_trackHitInMeasTIB ->clear();
1124  t_trackHitInMeasTID ->clear();
1125  t_trackHitOutMeasTOB ->clear();
1126  t_trackHitOutMeasTEC ->clear();
1127  t_trackHitOutMeasTIB ->clear();
1128  t_trackHitOutMeasTID ->clear();
1129  t_trackOutPosOutHitDr ->clear();
1130  t_trackL ->clear();
1131 
1132  t_maxNearP31x31 ->clear();
1133  t_maxNearP21x21 ->clear();
1134 
1135  t_ecalSpike11x11 ->clear();
1136  t_e7x7 ->clear();
1137  t_e9x9 ->clear();
1138  t_e11x11 ->clear();
1139  t_e15x15 ->clear();
1140 
1141  t_e7x7_10Sig ->clear();
1142  t_e9x9_10Sig ->clear();
1143  t_e11x11_10Sig ->clear();
1144  t_e15x15_10Sig ->clear();
1145  t_e7x7_15Sig ->clear();
1146  t_e9x9_15Sig ->clear();
1147  t_e11x11_15Sig ->clear();
1148  t_e15x15_15Sig ->clear();
1149  t_e7x7_20Sig ->clear();
1150  t_e9x9_20Sig ->clear();
1151  t_e11x11_20Sig ->clear();
1152  t_e15x15_20Sig ->clear();
1153  t_e7x7_25Sig ->clear();
1154  t_e9x9_25Sig ->clear();
1155  t_e11x11_25Sig ->clear();
1156  t_e15x15_25Sig ->clear();
1157  t_e7x7_30Sig ->clear();
1158  t_e9x9_30Sig ->clear();
1159  t_e11x11_30Sig ->clear();
1160  t_e15x15_30Sig ->clear();
1161 
1162  if (doMC) {
1163  t_simTrackP ->clear();
1164  t_esimPdgId ->clear();
1165  t_trkEcalEne ->clear();
1166 
1167  t_esim7x7 ->clear();
1168  t_esim9x9 ->clear();
1169  t_esim11x11 ->clear();
1170  t_esim15x15 ->clear();
1171 
1172  t_esim7x7Matched ->clear();
1173  t_esim9x9Matched ->clear();
1174  t_esim11x11Matched ->clear();
1175  t_esim15x15Matched ->clear();
1176 
1177  t_esim7x7Rest ->clear();
1178  t_esim9x9Rest ->clear();
1179  t_esim11x11Rest ->clear();
1180  t_esim15x15Rest ->clear();
1181 
1182  t_esim7x7Photon ->clear();
1183  t_esim9x9Photon ->clear();
1184  t_esim11x11Photon ->clear();
1185  t_esim15x15Photon ->clear();
1186 
1187  t_esim7x7NeutHad ->clear();
1188  t_esim9x9NeutHad ->clear();
1189  t_esim11x11NeutHad ->clear();
1190  t_esim15x15NeutHad ->clear();
1191 
1192  t_esim7x7CharHad ->clear();
1193  t_esim9x9CharHad ->clear();
1194  t_esim11x11CharHad ->clear();
1195  t_esim15x15CharHad ->clear();
1196  }
1197 
1198  t_maxNearHcalP3x3 ->clear();
1199  t_maxNearHcalP5x5 ->clear();
1200  t_maxNearHcalP7x7 ->clear();
1201 
1202  t_h3x3 ->clear();
1203  t_h5x5 ->clear();
1204  t_h7x7 ->clear();
1205  t_h3x3Sig ->clear();
1206  t_h5x5Sig ->clear();
1207  t_h7x7Sig ->clear();
1208 
1209  t_infoHcal ->clear();
1210 
1211  if (doMC) {
1212  t_trkHcalEne ->clear();
1213 
1214  t_hsim3x3 ->clear();
1215  t_hsim5x5 ->clear();
1216  t_hsim7x7 ->clear();
1217  t_hsim3x3Matched ->clear();
1218  t_hsim5x5Matched ->clear();
1219  t_hsim7x7Matched ->clear();
1220  t_hsim3x3Rest ->clear();
1221  t_hsim5x5Rest ->clear();
1222  t_hsim7x7Rest ->clear();
1223  t_hsim3x3Photon ->clear();
1224  t_hsim5x5Photon ->clear();
1225  t_hsim7x7Photon ->clear();
1226  t_hsim3x3NeutHad ->clear();
1227  t_hsim5x5NeutHad ->clear();
1228  t_hsim7x7NeutHad ->clear();
1229  t_hsim3x3CharHad ->clear();
1230  t_hsim5x5CharHad ->clear();
1231  t_hsim7x7CharHad ->clear();
1232  }
1233 }
1234 
1236 
1237  char hname[100], htit[100];
1238 
1239  TFileDirectory dir = fs->mkdir("nearMaxTrackP");
1240 
1241  for(unsigned int ieta=0; ieta<NEtaBins; ieta++) {
1242  double lowEta=-5.0, highEta= 5.0;
1243  lowEta = genPartEtaBins[ieta];
1244  highEta = genPartEtaBins[ieta+1];
1245 
1246  for(unsigned int ipt=0; ipt<NPBins; ipt++) {
1247  double lowP=0.0, highP=300.0;
1248  lowP = genPartPBins[ipt];
1249  highP = genPartPBins[ipt+1];
1250  sprintf(hname, "h_maxNearP31x31_ptBin%i_etaBin%i",ipt, ieta);
1251  sprintf(htit, "maxNearP in 31x31 (%3.2f<|#eta|<%3.2f), (%2.0f<trkP<%3.0f)", lowEta, highEta, lowP, highP );
1252  h_maxNearP31x31[ipt][ieta] = dir.make<TH1F>(hname, htit, 220, -2.0, 20.0);
1253  h_maxNearP31x31[ipt][ieta] ->Sumw2();
1254  sprintf(hname, "h_maxNearP25x25_ptBin%i_etaBin%i",ipt, ieta);
1255  sprintf(htit, "maxNearP in 25x25 (%3.2f<|#eta|<%3.2f), (%2.0f<trkP<%3.0f)", lowEta, highEta, lowP, highP );
1256  h_maxNearP25x25[ipt][ieta] = dir.make<TH1F>(hname, htit, 220, -2.0, 20.0);
1257  h_maxNearP25x25[ipt][ieta] ->Sumw2();
1258  sprintf(hname, "h_maxNearP21x21_ptBin%i_etaBin%i",ipt, ieta);
1259  sprintf(htit, "maxNearP in 21x21 (%3.2f<|#eta|<%3.2f), (%2.0f<trkP<%3.0f)", lowEta, highEta, lowP, highP );
1260  h_maxNearP21x21[ipt][ieta] = dir.make<TH1F>(hname, htit, 220, -2.0, 20.0);
1261  h_maxNearP21x21[ipt][ieta] ->Sumw2();
1262  sprintf(hname, "h_maxNearP15x15_ptBin%i_etaBin%i",ipt, ieta);
1263  sprintf(htit, "maxNearP in 15x15 (%3.2f<|#eta|<%3.2f), (%2.0f<trkP<%3.0f)", lowEta, highEta, lowP, highP );
1264  h_maxNearP15x15[ipt][ieta] = dir.make<TH1F>(hname, htit, 220, -2.0, 20.0);
1265  h_maxNearP15x15[ipt][ieta] ->Sumw2();
1266  }
1267  }
1268 
1269  h_L1AlgoNames = fs->make<TH1I>("h_L1AlgoNames", "h_L1AlgoNames:Bin Labels", 128, -0.5, 127.5);
1270 
1271  // Reconstructed Tracks
1272 
1273  h_PVTracksWt = fs->make<TH1F>("h_PVTracksWt", "h_PVTracksWt", 600, -0.1, 1.1);
1274 
1275  h_nTracks = fs->make<TH1F>("h_nTracks", "h_nTracks", 1000, -0.5, 999.5);
1276 
1277  sprintf(hname, "h_recEtaPt_0");
1278  sprintf(htit, "h_recEtaPt (all tracks Eta vs pT)");
1279  h_recEtaPt_0 = fs->make<TH2F>(hname, htit, 30, -3.0,3.0, 15, genPartPBins);
1280 
1281  sprintf(hname, "h_recEtaP_0");
1282  sprintf(htit, "h_recEtaP (all tracks Eta vs pT)");
1283  h_recEtaP_0 = fs->make<TH2F>(hname, htit, 30, -3.0,3.0, 15, genPartPBins);
1284 
1285  h_recPt_0 = fs->make<TH1F>("h_recPt_0", "Pt (all tracks)", 15, genPartPBins);
1286  h_recP_0 = fs->make<TH1F>("h_recP_0", "P (all tracks)", 15, genPartPBins);
1287  h_recEta_0 = fs->make<TH1F>("h_recEta_0", "Eta (all tracks)", 60, -3.0, 3.0);
1288  h_recPhi_0 = fs->make<TH1F>("h_recPhi_0", "Phi (all tracks)", 100, -3.2, 3.2);
1289  //-------------------------
1290  sprintf(hname, "h_recEtaPt_1");
1291  sprintf(htit, "h_recEtaPt (all good tracks Eta vs pT)");
1292  h_recEtaPt_1 = fs->make<TH2F>(hname, htit, 30, -3.0,3.0, 15, genPartPBins);
1293 
1294  sprintf(hname, "h_recEtaP_1");
1295  sprintf(htit, "h_recEtaP (all good tracks Eta vs pT)");
1296  h_recEtaP_1 = fs->make<TH2F>(hname, htit, 30, -3.0, 3.0, 15, genPartPBins);
1297 
1298  h_recPt_1 = fs->make<TH1F>("h_recPt_1", "Pt (all good tracks)", 15, genPartPBins);
1299  h_recP_1 = fs->make<TH1F>("h_recP_1", "P (all good tracks)", 15, genPartPBins);
1300  h_recEta_1 = fs->make<TH1F>("h_recEta_1", "Eta (all good tracks)", 60, -3.0, 3.0);
1301  h_recPhi_1 = fs->make<TH1F>("h_recPhi_1", "Phi (all good tracks)", 100, -3.2, 3.2);
1302  //-------------------------
1303  sprintf(hname, "h_recEtaPt_2");
1304  sprintf(htit, "h_recEtaPt (charge isolation Eta vs pT)");
1305  h_recEtaPt_2 = fs->make<TH2F>(hname, htit, 30, -3.0, 3.0, 15, genPartPBins);
1306 
1307  sprintf(hname, "h_recEtaP_2");
1308  sprintf(htit, "h_recEtaP (charge isolation Eta vs pT)");
1309  h_recEtaP_2 = fs->make<TH2F>(hname, htit, 30, -3.0, 3.0, 15, genPartPBins);
1310 
1311  h_recPt_2 = fs->make<TH1F>("h_recPt_2", "Pt (charge isolation)", 15, genPartPBins);
1312  h_recP_2 = fs->make<TH1F>("h_recP_2", "P (charge isolation)", 15, genPartPBins);
1313  h_recEta_2 = fs->make<TH1F>("h_recEta_2", "Eta (charge isolation)", 60, -3.0, 3.0);
1314  h_recPhi_2 = fs->make<TH1F>("h_recPhi_2", "Phi (charge isolation)", 100, -3.2, 3.2);
1315 
1316 
1317  tree = fs->make<TTree>("tree", "tree");
1318  tree->SetAutoSave(10000);
1319 
1320 
1321  tree->Branch("t_EvtNo" ,&t_EvtNo ,"t_EvtNo/I");
1322  tree->Branch("t_RunNo" ,&t_RunNo ,"t_RunNo/I");
1323  tree->Branch("t_Lumi" ,&t_Lumi ,"t_Lumi/I");
1324  tree->Branch("t_Bunch" ,&t_Bunch ,"t_Bunch/I");
1325 
1326 
1327  t_PVx = new std::vector<double>();
1328  t_PVy = new std::vector<double>();
1329  t_PVz = new std::vector<double>();
1330  t_PVisValid = new std::vector<int>();
1331  t_PVndof = new std::vector<int>();
1332  t_PVNTracks = new std::vector<int>();
1333  t_PVNTracksWt = new std::vector<int>();
1334  t_PVTracksSumPt = new std::vector<double>();
1335  t_PVTracksSumPtWt = new std::vector<double>();
1336  t_PVNTracksHP = new std::vector<int>();
1337  t_PVNTracksHPWt = new std::vector<int>();
1338  t_PVTracksSumPtHP = new std::vector<double>();
1339  t_PVTracksSumPtHPWt = new std::vector<double>();
1340 
1341  tree->Branch("PVx" ,"vector<double>" ,&t_PVx);
1342  tree->Branch("PVy" ,"vector<double>" ,&t_PVy);
1343  tree->Branch("PVz" ,"vector<double>" ,&t_PVz);
1344  tree->Branch("PVisValid" ,"vector<int>" ,&t_PVisValid);
1345  tree->Branch("PVndof" ,"vector<int>" ,&t_PVndof);
1346  tree->Branch("PVNTracks" ,"vector<int>" ,&t_PVNTracks);
1347  tree->Branch("PVNTracksWt" ,"vector<int>" ,&t_PVNTracksWt);
1348  tree->Branch("t_PVTracksSumPt" ,"vector<double>" ,&t_PVTracksSumPt);
1349  tree->Branch("t_PVTracksSumPtWt" ,"vector<double>" ,&t_PVTracksSumPtWt);
1350  tree->Branch("PVNTracksHP" ,"vector<int>" ,&t_PVNTracksHP);
1351  tree->Branch("PVNTracksHPWt" ,"vector<int>" ,&t_PVNTracksHPWt);
1352  tree->Branch("t_PVTracksSumPtHP" ,"vector<double>" ,&t_PVTracksSumPtHP);
1353  tree->Branch("t_PVTracksSumPtHPWt" ,"vector<double>" ,&t_PVTracksSumPtHPWt);
1354 
1355  //----- L1Trigger information
1356  for(int i=0; i<128; i++) t_L1Decision[i]=0;
1357  t_L1AlgoNames = new std::vector<std::string>();
1358  t_L1PreScale = new std::vector<int>();
1359  t_L1CenJetPt = new std::vector<double>();
1360  t_L1CenJetEta = new std::vector<double>();
1361  t_L1CenJetPhi = new std::vector<double>();
1362  t_L1FwdJetPt = new std::vector<double>();
1363  t_L1FwdJetEta = new std::vector<double>();
1364  t_L1FwdJetPhi = new std::vector<double>();
1365  t_L1TauJetPt = new std::vector<double>();
1366  t_L1TauJetEta = new std::vector<double>();
1367  t_L1TauJetPhi = new std::vector<double>();
1368  t_L1MuonPt = new std::vector<double>();
1369  t_L1MuonEta = new std::vector<double>();
1370  t_L1MuonPhi = new std::vector<double>();
1371  t_L1IsoEMPt = new std::vector<double>();
1372  t_L1IsoEMEta = new std::vector<double>();
1373  t_L1IsoEMPhi = new std::vector<double>();
1374  t_L1NonIsoEMPt = new std::vector<double>();
1375  t_L1NonIsoEMEta = new std::vector<double>();
1376  t_L1NonIsoEMPhi = new std::vector<double>();
1377  t_L1METPt = new std::vector<double>();
1378  t_L1METEta = new std::vector<double>();
1379  t_L1METPhi = new std::vector<double>();
1380 
1381  tree->Branch("t_L1Decision", t_L1Decision, "t_L1Decision[128]/I");
1382  tree->Branch("t_L1AlgoNames", "vector<string>", &t_L1AlgoNames);
1383  tree->Branch("t_L1PreScale", "vector<int>", &t_L1PreScale);
1384  tree->Branch("t_L1CenJetPt", "vector<double>", &t_L1CenJetPt);
1385  tree->Branch("t_L1CenJetEta", "vector<double>", &t_L1CenJetEta);
1386  tree->Branch("t_L1CenJetPhi", "vector<double>", &t_L1CenJetPhi);
1387  tree->Branch("t_L1FwdJetPt", "vector<double>", &t_L1FwdJetPt);
1388  tree->Branch("t_L1FwdJetEta", "vector<double>", &t_L1FwdJetEta);
1389  tree->Branch("t_L1FwdJetPhi", "vector<double>", &t_L1FwdJetPhi);
1390  tree->Branch("t_L1TauJetPt", "vector<double>", &t_L1TauJetPt);
1391  tree->Branch("t_L1TauJetEta", "vector<double>", &t_L1TauJetEta);
1392  tree->Branch("t_L1TauJetPhi", "vector<double>", &t_L1TauJetPhi);
1393  tree->Branch("t_L1MuonPt", "vector<double>", &t_L1MuonPt);
1394  tree->Branch("t_L1MuonEta", "vector<double>", &t_L1MuonEta);
1395  tree->Branch("t_L1MuonPhi", "vector<double>", &t_L1MuonPhi);
1396  tree->Branch("t_L1IsoEMPt", "vector<double>", &t_L1IsoEMPt);
1397  tree->Branch("t_L1IsoEMEta", "vector<double>", &t_L1IsoEMEta);
1398  tree->Branch("t_L1IsoEMPhi", "vector<double>", &t_L1IsoEMPhi);
1399  tree->Branch("t_L1NonIsoEMPt", "vector<double>", &t_L1NonIsoEMPt);
1400  tree->Branch("t_L1NonIsoEMEta", "vector<double>", &t_L1NonIsoEMEta);
1401  tree->Branch("t_L1NonIsoEMPhi", "vector<double>", &t_L1NonIsoEMPhi);
1402  tree->Branch("t_L1METPt", "vector<double>", &t_L1METPt);
1403  tree->Branch("t_L1METEta", "vector<double>", &t_L1METEta);
1404  tree->Branch("t_L1METPhi", "vector<double>", &t_L1METPhi);
1405 
1406  t_jetPt = new std::vector<double>();
1407  t_jetEta = new std::vector<double>();
1408  t_jetPhi = new std::vector<double>();
1409  t_nTrksJetCalo = new std::vector<double>();
1410  t_nTrksJetVtx = new std::vector<double>();
1411  tree->Branch("t_jetPt", "vector<double>",&t_jetPt);
1412  tree->Branch("t_jetEta", "vector<double>",&t_jetEta);
1413  tree->Branch("t_jetPhi", "vector<double>",&t_jetPhi);
1414  tree->Branch("t_nTrksJetCalo", "vector<double>",&t_nTrksJetCalo);
1415  tree->Branch("t_nTrksJetVtx", "vector<double>",&t_nTrksJetVtx);
1416 
1417  t_trackPAll = new std::vector<double>();
1418  t_trackEtaAll = new std::vector<double>();
1419  t_trackPhiAll = new std::vector<double>();
1420  t_trackPdgIdAll = new std::vector<double>();
1421  t_trackPtAll = new std::vector<double>();
1422  t_trackDxyAll = new std::vector<double>();
1423  t_trackDzAll = new std::vector<double>();
1424  t_trackDxyPVAll = new std::vector<double>();
1425  t_trackDzPVAll = new std::vector<double>();
1426  t_trackChiSqAll = new std::vector<double>();
1427  tree->Branch("t_trackPAll", "vector<double>", &t_trackPAll );
1428  tree->Branch("t_trackPhiAll", "vector<double>", &t_trackPhiAll );
1429  tree->Branch("t_trackEtaAll", "vector<double>", &t_trackEtaAll );
1430  tree->Branch("t_trackPtAll", "vector<double>", &t_trackPtAll );
1431  tree->Branch("t_trackDxyAll", "vector<double>", &t_trackDxyAll );
1432  tree->Branch("t_trackDzAll", "vector<double>", &t_trackDzAll );
1433  tree->Branch("t_trackDxyPVAll", "vector<double>", &t_trackDxyPVAll );
1434  tree->Branch("t_trackDzPVAll", "vector<double>", &t_trackDzPVAll );
1435  tree->Branch("t_trackChiSqAll", "vector<double>", &t_trackChiSqAll );
1436  //tree->Branch("t_trackPdgIdAll", "vector<double>", &t_trackPdgIdAll);
1437 
1438  t_trackP = new std::vector<double>();
1439  t_trackPt = new std::vector<double>();
1440  t_trackEta = new std::vector<double>();
1441  t_trackPhi = new std::vector<double>();
1442  t_trackEcalEta = new std::vector<double>();
1443  t_trackEcalPhi = new std::vector<double>();
1444  t_trackHcalEta = new std::vector<double>();
1445  t_trackHcalPhi = new std::vector<double>();
1446  t_trackNOuterHits = new std::vector<int>();
1447  t_NLayersCrossed = new std::vector<int>();
1448  t_trackDxy = new std::vector<double>();
1449  t_trackDxyBS = new std::vector<double>();
1450  t_trackDz = new std::vector<double>();
1451  t_trackDzBS = new std::vector<double>();
1452  t_trackDxyPV = new std::vector<double>();
1453  t_trackDzPV = new std::vector<double>();
1454  t_trackPVIdx = new std::vector<int>();
1455  t_trackChiSq = new std::vector<double>();
1456  t_trackHitsTOB = new std::vector<int>();
1457  t_trackHitsTEC = new std::vector<int>();
1458  t_trackHitInMissTOB = new std::vector<int>();
1459  t_trackHitInMissTEC = new std::vector<int>();
1460  t_trackHitInMissTIB = new std::vector<int>();
1461  t_trackHitInMissTID = new std::vector<int>();
1462  t_trackHitOutMissTOB = new std::vector<int>();
1463  t_trackHitOutMissTEC = new std::vector<int>();
1464  t_trackHitOutMissTIB = new std::vector<int>();
1465  t_trackHitOutMissTID = new std::vector<int>();
1466  t_trackHitInMissTIBTID = new std::vector<int>();
1467  t_trackHitOutMissTOB = new std::vector<int>();
1468  t_trackHitOutMissTEC = new std::vector<int>();
1469  t_trackHitOutMissTIB = new std::vector<int>();
1470  t_trackHitOutMissTID = new std::vector<int>();
1471  t_trackHitOutMissTOBTEC= new std::vector<int>();
1472  t_trackHitInMeasTOB = new std::vector<int>();
1473  t_trackHitInMeasTEC = new std::vector<int>();
1474  t_trackHitInMeasTIB = new std::vector<int>();
1475  t_trackHitInMeasTID = new std::vector<int>();
1476  t_trackHitOutMeasTOB = new std::vector<int>();
1477  t_trackHitOutMeasTEC = new std::vector<int>();
1478  t_trackHitOutMeasTIB = new std::vector<int>();
1479  t_trackHitOutMeasTID = new std::vector<int>();
1480  t_trackOutPosOutHitDr = new std::vector<double>();
1481  t_trackL = new std::vector<double>();
1482 
1483  tree->Branch("t_trackP", "vector<double>", &t_trackP );
1484  tree->Branch("t_trackPt", "vector<double>", &t_trackPt );
1485  tree->Branch("t_trackEta", "vector<double>", &t_trackEta );
1486  tree->Branch("t_trackPhi", "vector<double>", &t_trackPhi );
1487  tree->Branch("t_trackEcalEta", "vector<double>", &t_trackEcalEta );
1488  tree->Branch("t_trackEcalPhi", "vector<double>", &t_trackEcalPhi );
1489  tree->Branch("t_trackHcalEta", "vector<double>", &t_trackHcalEta );
1490  tree->Branch("t_trackHcalPhi", "vector<double>", &t_trackHcalPhi );
1491 
1492  tree->Branch("t_trackNOuterHits", "vector<int>", &t_trackNOuterHits );
1493  tree->Branch("t_NLayersCrossed", "vector<int>", &t_NLayersCrossed );
1494  tree->Branch("t_trackHitsTOB", "vector<int>", &t_trackHitsTOB );
1495  tree->Branch("t_trackHitsTEC", "vector<int>", &t_trackHitsTEC );
1496  tree->Branch("t_trackHitInMissTOB", "vector<int>", &t_trackHitInMissTOB );
1497  tree->Branch("t_trackHitInMissTEC", "vector<int>", &t_trackHitInMissTEC );
1498  tree->Branch("t_trackHitInMissTIB", "vector<int>", &t_trackHitInMissTIB );
1499  tree->Branch("t_trackHitInMissTID", "vector<int>", &t_trackHitInMissTID );
1500  tree->Branch("t_trackHitInMissTIBTID", "vector<int>", &t_trackHitInMissTIBTID );
1501  tree->Branch("t_trackHitOutMissTOB", "vector<int>", &t_trackHitOutMissTOB);
1502  tree->Branch("t_trackHitOutMissTEC", "vector<int>", &t_trackHitOutMissTEC);
1503  tree->Branch("t_trackHitOutMissTIB", "vector<int>", &t_trackHitOutMissTIB);
1504  tree->Branch("t_trackHitOutMissTID", "vector<int>", &t_trackHitOutMissTID);
1505  tree->Branch("t_trackHitOutMissTOBTEC","vector<int>", &t_trackHitOutMissTOBTEC);
1506  tree->Branch("t_trackHitInMeasTOB", "vector<int>", &t_trackHitInMeasTOB );
1507  tree->Branch("t_trackHitInMeasTEC", "vector<int>", &t_trackHitInMeasTEC );
1508  tree->Branch("t_trackHitInMeasTIB", "vector<int>", &t_trackHitInMeasTIB );
1509  tree->Branch("t_trackHitInMeasTID", "vector<int>", &t_trackHitInMeasTID );
1510  tree->Branch("t_trackHitOutMeasTOB", "vector<int>", &t_trackHitOutMeasTOB);
1511  tree->Branch("t_trackHitOutMeasTEC", "vector<int>", &t_trackHitOutMeasTEC);
1512  tree->Branch("t_trackHitOutMeasTIB", "vector<int>", &t_trackHitOutMeasTIB);
1513  tree->Branch("t_trackHitOutMeasTID", "vector<int>", &t_trackHitOutMeasTID);
1514  tree->Branch("t_trackOutPosOutHitDr", "vector<double>", &t_trackOutPosOutHitDr);
1515  tree->Branch("t_trackL", "vector<double>", &t_trackL);
1516 
1517  tree->Branch("t_trackDxy", "vector<double>", &t_trackDxy );
1518  tree->Branch("t_trackDxyBS", "vector<double>", &t_trackDxyBS );
1519  tree->Branch("t_trackDz", "vector<double>", &t_trackDz );
1520  tree->Branch("t_trackDzBS", "vector<double>", &t_trackDzBS );
1521  tree->Branch("t_trackDxyPV", "vector<double>", &t_trackDxyPV );
1522  tree->Branch("t_trackDzPV", "vector<double>", &t_trackDzPV );
1523  tree->Branch("t_trackChiSq", "vector<double>", &t_trackChiSq );
1524  tree->Branch("t_trackPVIdx", "vector<int>", &t_trackPVIdx );
1525 
1526  t_maxNearP31x31 = new std::vector<double>();
1527  t_maxNearP21x21 = new std::vector<double>();
1528 
1529  tree->Branch("t_maxNearP31x31", "vector<double>", &t_maxNearP31x31);
1530  tree->Branch("t_maxNearP21x21", "vector<double>", &t_maxNearP21x21);
1531 
1532  t_ecalSpike11x11 = new std::vector<int>();
1533  t_e7x7 = new std::vector<double>();
1534  t_e9x9 = new std::vector<double>();
1535  t_e11x11 = new std::vector<double>();
1536  t_e15x15 = new std::vector<double>();
1537 
1538  tree->Branch("t_ecalSpike11x11", "vector<int>", &t_ecalSpike11x11);
1539  tree->Branch("t_e7x7", "vector<double>", &t_e7x7);
1540  tree->Branch("t_e9x9", "vector<double>", &t_e9x9);
1541  tree->Branch("t_e11x11", "vector<double>", &t_e11x11);
1542  tree->Branch("t_e15x15", "vector<double>", &t_e15x15);
1543 
1544  t_e7x7_10Sig = new std::vector<double>();
1545  t_e9x9_10Sig = new std::vector<double>();
1546  t_e11x11_10Sig = new std::vector<double>();
1547  t_e15x15_10Sig = new std::vector<double>();
1548  t_e7x7_15Sig = new std::vector<double>();
1549  t_e9x9_15Sig = new std::vector<double>();
1550  t_e11x11_15Sig = new std::vector<double>();
1551  t_e15x15_15Sig = new std::vector<double>();
1552  t_e7x7_20Sig = new std::vector<double>();
1553  t_e9x9_20Sig = new std::vector<double>();
1554  t_e11x11_20Sig = new std::vector<double>();
1555  t_e15x15_20Sig = new std::vector<double>();
1556  t_e7x7_25Sig = new std::vector<double>();
1557  t_e9x9_25Sig = new std::vector<double>();
1558  t_e11x11_25Sig = new std::vector<double>();
1559  t_e15x15_25Sig = new std::vector<double>();
1560  t_e7x7_30Sig = new std::vector<double>();
1561  t_e9x9_30Sig = new std::vector<double>();
1562  t_e11x11_30Sig = new std::vector<double>();
1563  t_e15x15_30Sig = new std::vector<double>();
1564 
1565  tree->Branch("t_e7x7_10Sig" ,"vector<double>", &t_e7x7_10Sig);
1566  tree->Branch("t_e9x9_10Sig" ,"vector<double>", &t_e9x9_10Sig);
1567  tree->Branch("t_e11x11_10Sig" ,"vector<double>", &t_e11x11_10Sig);
1568  tree->Branch("t_e15x15_10Sig" ,"vector<double>", &t_e15x15_10Sig);
1569  tree->Branch("t_e7x7_15Sig" ,"vector<double>", &t_e7x7_15Sig);
1570  tree->Branch("t_e9x9_15Sig" ,"vector<double>", &t_e9x9_15Sig);
1571  tree->Branch("t_e11x11_15Sig" ,"vector<double>", &t_e11x11_15Sig);
1572  tree->Branch("t_e15x15_15Sig" ,"vector<double>", &t_e15x15_15Sig);
1573  tree->Branch("t_e7x7_20Sig" ,"vector<double>", &t_e7x7_20Sig);
1574  tree->Branch("t_e9x9_20Sig" ,"vector<double>", &t_e9x9_20Sig);
1575  tree->Branch("t_e11x11_20Sig" ,"vector<double>", &t_e11x11_20Sig);
1576  tree->Branch("t_e15x15_20Sig" ,"vector<double>", &t_e15x15_20Sig);
1577  tree->Branch("t_e7x7_25Sig" ,"vector<double>", &t_e7x7_25Sig);
1578  tree->Branch("t_e9x9_25Sig" ,"vector<double>", &t_e9x9_25Sig);
1579  tree->Branch("t_e11x11_25Sig" ,"vector<double>", &t_e11x11_25Sig);
1580  tree->Branch("t_e15x15_25Sig" ,"vector<double>", &t_e15x15_25Sig);
1581  tree->Branch("t_e7x7_30Sig" ,"vector<double>", &t_e7x7_30Sig);
1582  tree->Branch("t_e9x9_30Sig" ,"vector<double>", &t_e9x9_30Sig);
1583  tree->Branch("t_e11x11_30Sig" ,"vector<double>", &t_e11x11_30Sig);
1584  tree->Branch("t_e15x15_30Sig" ,"vector<double>", &t_e15x15_30Sig);
1585 
1586  if (doMC) {
1587  t_esim7x7 = new std::vector<double>();
1588  t_esim9x9 = new std::vector<double>();
1589  t_esim11x11 = new std::vector<double>();
1590  t_esim15x15 = new std::vector<double>();
1591 
1592  t_esim7x7Matched = new std::vector<double>();
1593  t_esim9x9Matched = new std::vector<double>();
1594  t_esim11x11Matched = new std::vector<double>();
1595  t_esim15x15Matched = new std::vector<double>();
1596 
1597  t_esim7x7Rest = new std::vector<double>();
1598  t_esim9x9Rest = new std::vector<double>();
1599  t_esim11x11Rest = new std::vector<double>();
1600  t_esim15x15Rest = new std::vector<double>();
1601 
1602  t_esim7x7Photon = new std::vector<double>();
1603  t_esim9x9Photon = new std::vector<double>();
1604  t_esim11x11Photon = new std::vector<double>();
1605  t_esim15x15Photon = new std::vector<double>();
1606 
1607  t_esim7x7NeutHad = new std::vector<double>();
1608  t_esim9x9NeutHad = new std::vector<double>();
1609  t_esim11x11NeutHad = new std::vector<double>();
1610  t_esim15x15NeutHad = new std::vector<double>();
1611 
1612  t_esim7x7CharHad = new std::vector<double>();
1613  t_esim9x9CharHad = new std::vector<double>();
1614  t_esim11x11CharHad = new std::vector<double>();
1615  t_esim15x15CharHad = new std::vector<double>();
1616 
1617  t_trkEcalEne = new std::vector<double>();
1618  t_simTrackP = new std::vector<double>();
1619  t_esimPdgId = new std::vector<double>();
1620 
1621  tree->Branch("t_esim7x7", "vector<double>", &t_esim7x7);
1622  tree->Branch("t_esim9x9", "vector<double>", &t_esim9x9);
1623  tree->Branch("t_esim11x11", "vector<double>", &t_esim11x11);
1624  tree->Branch("t_esim15x15", "vector<double>", &t_esim15x15);
1625 
1626  tree->Branch("t_esim7x7Matched", "vector<double>", &t_esim7x7Matched);
1627  tree->Branch("t_esim9x9Matched", "vector<double>", &t_esim9x9Matched);
1628  tree->Branch("t_esim11x11Matched", "vector<double>", &t_esim11x11Matched);
1629  tree->Branch("t_esim15x15Matched", "vector<double>", &t_esim15x15Matched);
1630 
1631  tree->Branch("t_esim7x7Rest", "vector<double>", &t_esim7x7Rest);
1632  tree->Branch("t_esim9x9Rest", "vector<double>", &t_esim9x9Rest);
1633  tree->Branch("t_esim11x11Rest", "vector<double>", &t_esim11x11Rest);
1634  tree->Branch("t_esim15x15Rest", "vector<double>", &t_esim15x15Rest);
1635 
1636  tree->Branch("t_esim7x7Photon", "vector<double>", &t_esim7x7Photon);
1637  tree->Branch("t_esim9x9Photon", "vector<double>", &t_esim9x9Photon);
1638  tree->Branch("t_esim11x11Photon", "vector<double>", &t_esim11x11Photon);
1639  tree->Branch("t_esim15x15Photon", "vector<double>", &t_esim15x15Photon);
1640 
1641 ;
1642  tree->Branch("t_esim7x7NeutHad", "vector<double>", &t_esim7x7NeutHad);
1643  tree->Branch("t_esim9x9NeutHad", "vector<double>", &t_esim9x9NeutHad);
1644  tree->Branch("t_esim11x11NeutHad", "vector<double>", &t_esim11x11NeutHad);
1645  tree->Branch("t_esim15x15NeutHad", "vector<double>", &t_esim15x15NeutHad);
1646 
1647  tree->Branch("t_esim7x7CharHad", "vector<double>", &t_esim7x7CharHad);
1648  tree->Branch("t_esim9x9CharHad", "vector<double>", &t_esim9x9CharHad);
1649  tree->Branch("t_esim11x11CharHad", "vector<double>", &t_esim11x11CharHad);
1650  tree->Branch("t_esim15x15CharHad", "vector<double>", &t_esim15x15CharHad);
1651 
1652  tree->Branch("t_trkEcalEne", "vector<double>", &t_trkEcalEne);
1653  tree->Branch("t_simTrackP", "vector<double>", &t_simTrackP);
1654  tree->Branch("t_esimPdgId", "vector<double>", &t_esimPdgId);
1655  }
1656 
1657  t_maxNearHcalP3x3 = new std::vector<double>();
1658  t_maxNearHcalP5x5 = new std::vector<double>();
1659  t_maxNearHcalP7x7 = new std::vector<double>();
1660  t_h3x3 = new std::vector<double>();
1661  t_h5x5 = new std::vector<double>();
1662  t_h7x7 = new std::vector<double>();
1663  t_h3x3Sig = new std::vector<double>();
1664  t_h5x5Sig = new std::vector<double>();
1665  t_h7x7Sig = new std::vector<double>();
1666  t_infoHcal = new std::vector<int>();
1667 
1668  if (doMC) {
1669  t_trkHcalEne = new std::vector<double>();
1670  t_hsim3x3 = new std::vector<double>();
1671  t_hsim5x5 = new std::vector<double>();
1672  t_hsim7x7 = new std::vector<double>();
1673  t_hsim3x3Matched = new std::vector<double>();
1674  t_hsim5x5Matched = new std::vector<double>();
1675  t_hsim7x7Matched = new std::vector<double>();
1676  t_hsim3x3Rest = new std::vector<double>();
1677  t_hsim5x5Rest = new std::vector<double>();
1678  t_hsim7x7Rest = new std::vector<double>();
1679  t_hsim3x3Photon = new std::vector<double>();
1680  t_hsim5x5Photon = new std::vector<double>();
1681  t_hsim7x7Photon = new std::vector<double>();
1682  t_hsim3x3NeutHad = new std::vector<double>();
1683  t_hsim5x5NeutHad = new std::vector<double>();
1684  t_hsim7x7NeutHad = new std::vector<double>();
1685  t_hsim3x3CharHad = new std::vector<double>();
1686  t_hsim5x5CharHad = new std::vector<double>();
1687  t_hsim7x7CharHad = new std::vector<double>();
1688  }
1689 
1690  tree->Branch("t_maxNearHcalP3x3", "vector<double>", &t_maxNearHcalP3x3);
1691  tree->Branch("t_maxNearHcalP5x5", "vector<double>", &t_maxNearHcalP5x5);
1692  tree->Branch("t_maxNearHcalP7x7", "vector<double>", &t_maxNearHcalP7x7);
1693  tree->Branch("t_h3x3", "vector<double>", &t_h3x3);
1694  tree->Branch("t_h5x5", "vector<double>", &t_h5x5);
1695  tree->Branch("t_h7x7", "vector<double>", &t_h7x7);
1696  tree->Branch("t_h3x3Sig", "vector<double>", &t_h3x3Sig);
1697  tree->Branch("t_h5x5Sig", "vector<double>", &t_h5x5Sig);
1698  tree->Branch("t_h7x7Sig", "vector<double>", &t_h7x7Sig);
1699  tree->Branch("t_infoHcal", "vector<int>", &t_infoHcal);
1700 
1701  if (doMC) {
1702  tree->Branch("t_trkHcalEne", "vector<double>", &t_trkHcalEne);
1703  tree->Branch("t_hsim3x3", "vector<double>", &t_hsim3x3);
1704  tree->Branch("t_hsim5x5", "vector<double>", &t_hsim5x5);
1705  tree->Branch("t_hsim7x7", "vector<double>", &t_hsim7x7);
1706  tree->Branch("t_hsim3x3Matched", "vector<double>", &t_hsim3x3Matched);
1707  tree->Branch("t_hsim5x5Matched", "vector<double>", &t_hsim5x5Matched);
1708  tree->Branch("t_hsim7x7Matched", "vector<double>", &t_hsim7x7Matched);
1709  tree->Branch("t_hsim3x3Rest", "vector<double>", &t_hsim3x3Rest);
1710  tree->Branch("t_hsim5x5Rest", "vector<double>", &t_hsim5x5Rest);
1711  tree->Branch("t_hsim7x7Rest", "vector<double>", &t_hsim7x7Rest);
1712  tree->Branch("t_hsim3x3Photon", "vector<double>", &t_hsim3x3Photon);
1713  tree->Branch("t_hsim5x5Photon", "vector<double>", &t_hsim5x5Photon);
1714  tree->Branch("t_hsim7x7Photon", "vector<double>", &t_hsim7x7Photon);
1715  tree->Branch("t_hsim3x3NeutHad", "vector<double>", &t_hsim3x3NeutHad);
1716  tree->Branch("t_hsim5x5NeutHad", "vector<double>", &t_hsim5x5NeutHad);
1717  tree->Branch("t_hsim7x7NeutHad", "vector<double>", &t_hsim7x7NeutHad);
1718  tree->Branch("t_hsim3x3CharHad", "vector<double>", &t_hsim3x3CharHad);
1719  tree->Branch("t_hsim5x5CharHad", "vector<double>", &t_hsim5x5CharHad);
1720  tree->Branch("t_hsim7x7CharHad", "vector<double>", &t_hsim7x7CharHad);
1721  }
1722  tree->Branch("t_nTracks", &t_nTracks, "t_nTracks/I");
1723 
1724 }
1725 
1726 
1727 double IsolatedTracksNxN::DeltaPhi(double v1, double v2) {
1728  // Computes the correctly normalized phi difference
1729  // v1, v2 = phi of object 1 and 2
1730 
1731  double pi = 3.141592654;
1732  double twopi = 6.283185307;
1733 
1734  double diff = std::abs(v2 - v1);
1735  double corr = twopi - diff;
1736  if (diff < pi){ return diff;} else { return corr;}
1737 }
1738 
1739 double IsolatedTracksNxN::DeltaR(double eta1, double phi1, double eta2, double phi2) {
1740  double deta = eta1 - eta2;
1741  double dphi = DeltaPhi(phi1, phi2);
1742  return std::sqrt(deta*deta + dphi*dphi);
1743 }
1744 
1746 {
1747  using namespace reco;
1748  std::string theTrackQuality = "highPurity";
1749  reco::TrackBase::TrackQuality trackQuality_ = reco::TrackBase::qualityByName(theTrackQuality);
1750 
1751  std::cout << " Reference Point " << pTrack->referencePoint() <<"\n"
1752  << " TrackMmentum " << pTrack->momentum()
1753  << " (pt,eta,phi)(" << pTrack->pt()<<","<<pTrack->eta()<<","<<pTrack->phi()<<")"
1754  << " p " << pTrack->p() << "\n"
1755  << " Normalized chi2 " << pTrack->normalizedChi2() <<" charge " << pTrack->charge()
1756  << " qoverp() " << pTrack->qoverp() <<"+-" << pTrack->qoverpError()
1757  << " d0 " << pTrack->d0() << "\n"
1758  << " NValidHits " << pTrack->numberOfValidHits() << " NLostHits " << pTrack->numberOfLostHits()
1759  << " TrackQuality " << pTrack->qualityName(trackQuality_) << " " << pTrack->quality(trackQuality_)
1760  << std::endl;
1761 
1762  if( printTrkHitPattern_ ) {
1763  const reco::HitPattern &p = pTrack->hitPattern();
1764 
1765  std::cout<<"default " << std::endl;
1766  for(int i = 0; i < p.numberOfHits(HitPattern::TRACK_HITS); i++){
1767  p.printHitPattern(HitPattern::TRACK_HITS, i, std::cout);
1768  }
1769 
1770  std::cout<<"trackerExpectedHitsInner() " << std::endl;
1771  for(int i = 0; i < p.numberOfHits(HitPattern::MISSING_INNER_HITS); i++){
1772  p.printHitPattern(HitPattern::MISSING_INNER_HITS, i, std::cout);
1773  }
1774 
1775  std::cout<<"trackerExpectedHitsOuter() " << std::endl;
1776  for(int i = 0; i < p.numberOfHits(HitPattern::MISSING_OUTER_HITS); i++){
1777  p.printHitPattern(HitPattern::MISSING_OUTER_HITS, i, std::cout);
1778  }
1779 
1780  std::cout << "\n \t trackerLayersWithMeasurement() "
1782  << "\n \t pixelLayersWithMeasurement() "
1784  << "\n \t stripLayersWithMeasurement() "
1786  << "\n \t pixelBarrelLayersWithMeasurement() "
1788  << "\n \t pixelEndcapLayersWithMeasurement() "
1790  << "\n \t stripTIBLayersWithMeasurement() "
1792  << "\n \t stripTIDLayersWithMeasurement() "
1794  << "\n \t stripTOBLayersWithMeasurement() "
1796  << "\n \t stripTECLayersWithMeasurement() "
1798  << std::endl;
1799  }
1800 }
1801 
1802 //define this as a plug-in
RunNumber_t run() const
Definition: EventID.h:39
std::vector< double > * t_hsim7x7CharHad
double qoverp() const
q / p
Definition: TrackBase.h:549
double DeltaPhi(double v1, double v2)
double p() const
momentum vector magnitude
Definition: TrackBase.h:591
std::vector< double > * t_hsim3x3Matched
int stripTOBLayersWithMeasurement() const
Definition: HitPattern.cc:486
std::vector< double > * t_esim7x7CharHad
edm::EDGetTokenT< l1extra::L1JetParticleCollection > tok_L1extCenJet_
const std::string & gtTriggerMenuName() const
std::vector< double > * t_e15x15
value_type const * get() const
Definition: RefToBase.h:213
std::vector< double > * t_nTrksJetCalo
T getParameter(std::string const &) const
EventNumber_t event() const
Definition: EventID.h:41
std::vector< double > * t_maxNearHcalP7x7
T getUntrackedParameter(std::string const &, T const &) const
const Point & referencePoint() const
Reference point on the track.
Definition: TrackBase.h:657
std::vector< double > * t_trackPt
int i
Definition: DBlmapReader.cc:9
std::vector< double > * t_trackPAll
std::vector< double > * t_e7x7
std::vector< double > * t_trackHcalPhi
std::vector< double > * t_trackDxyPVAll
static std::string qualityName(TrackQuality)
Definition: TrackBase.h:508
std::vector< int > * t_trackHitOutMeasTEC
std::vector< double > * t_trackPdgIdAll
edm::InputTag JetExtender_
std::vector< spr::propagatedTrackID > propagateCALO(edm::Handle< reco::TrackCollection > &trkCollection, const CaloGeometry *geo, const MagneticField *bField, std::string &theTrackQuality, bool debug=false)
void eCaloSimInfo(std::vector< DetId > vdets, const CaloGeometry *geo, edm::Handle< T > &hitsEB, edm::Handle< T > &hitsEE, edm::Handle< edm::SimTrackContainer > &SimTk, edm::Handle< edm::SimVertexContainer > &SimVtx, edm::SimTrackContainer::const_iterator trkInfo, caloSimInfo &info, double timeCut=150, bool debug=false)
edm::InputTag L1extraCenJetSource_
std::vector< double > * t_trackOutPosOutHitDr
static const bool useL1GtTriggerMenuLite(false)
std::vector< double > * t_maxNearHcalP5x5
std::vector< double > * t_L1NonIsoEMPt
std::vector< int > * t_trackHitOutMeasTIB
std::vector< double > * t_esim7x7
std::vector< int > * t_trackHitInMissTOB
std::vector< double > * t_trackDz
double d0() const
dxy parameter in perigee convention (d0 = -dxy)
Definition: TrackBase.h:573
std::map< std::string, double > eHCALSimInfo(const edm::Event &, const HcalTopology *topology, const DetId &det, const CaloGeometry *geo, edm::Handle< T > &hits, edm::Handle< edm::SimTrackContainer > &SimTk, edm::Handle< edm::SimVertexContainer > &SimVtx, const reco::Track *pTrack, TrackerHitAssociator &associate, int ieta, int iphi, double timeCut=150, bool includeHO=false, bool debug=false)
std::vector< double > * t_trackEcalEta
int stripTIBLayersWithMeasurement() const
Definition: HitPattern.cc:464
double normalizedChi2() const
chi-squared divided by n.d.o.f. (or chi-squared * 1e6 if n.d.o.f. is zero)
Definition: TrackBase.h:537
std::vector< double > * t_trackEtaAll
std::vector< double > * t_esim11x11Matched
std::vector< double > * t_L1IsoEMEta
std::vector< double > * t_PVTracksSumPt
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:446
std::vector< double > * t_e9x9_10Sig
TH1F * h_maxNearP25x25[NPBins][NEtaBins]
TrackQuality
track quality
Definition: TrackBase.h:113
std::vector< double > * t_hsim3x3CharHad
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
int stripTOBLayersWithoutMeasurement(HitCategory category) const
Definition: HitPattern.cc:554
std::vector< int > * t_PVndof
std::map< std::string, L1GtAlgorithm > AlgorithmMap
map containing the algorithms
std::vector< double > * t_hsim7x7Rest
std::vector< double > * t_L1TauJetPhi
std::vector< double > * t_PVTracksSumPtWt
const bool availableL1Configuration(int &errorCode, int &l1ConfCode) const
Definition: L1GtUtils.cc:2293
std::vector< int > * t_infoHcal
std::vector< double > * t_e11x11_15Sig
static const size_t NEtaBins
std::vector< double > * t_L1TauJetPt
std::vector< int > * t_PVNTracksWt
std::vector< double > * t_trackEcalPhi
std::vector< double > * t_esim11x11NeutHad
std::vector< double > * t_L1TauJetEta
std::vector< double > * t_e15x15_15Sig
std::vector< double > * t_L1MuonEta
std::vector< double > * t_L1MuonPt
int pixelLayersWithMeasurement() const
Definition: HitPattern.h:917
int bunchCrossing() const
Definition: EventBase.h:62
std::vector< double > * t_h3x3
edm::LuminosityBlockNumber_t luminosityBlock() const
Definition: EventBase.h:59
std::vector< double > * t_L1MuonPhi
double phi() const
azimuthal angle of momentum vector
Definition: TrackBase.h:621
unsigned short numberOfLostHits() const
number of cases where track crossed a layer without getting a hit.
Definition: TrackBase.h:802
std::vector< double > * t_trackDxyAll
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
std::vector< double > * t_maxNearHcalP3x3
edm::InputTag L1extraIsoEmSource_
IsolatedTracksNxN(const edm::ParameterSet &)
edm::EDGetTokenT< edm::SimVertexContainer > tok_simVtx_
std::vector< double > * t_hsim3x3Photon
edm::EDGetTokenT< reco::BeamSpot > tok_bs_
const Vector & momentum() const
track momentum vector
Definition: TrackBase.h:651
T eta() const
int ii
Definition: cuy.py:588
std::vector< double > * t_PVTracksSumPtHP
const math::XYZPoint & outerPosition() const
position of the outermost hit
Definition: Track.h:65
std::vector< double > * t_esim11x11Rest
int pixelEndcapLayersWithMeasurement() const
Definition: HitPattern.cc:452
edm::InputTag L1GTReadoutRcdSource_
float float float z
std::vector< int > * t_trackNOuterHits
std::vector< double > * t_h7x7Sig
std::vector< double > * t_e15x15_30Sig
std::vector< double > * t_trackPhiAll
std::vector< double > * t_e7x7_25Sig
std::vector< double > * t_trkHcalEne
edm::InputTag L1extraFwdJetSource_
std::vector< double > * t_e7x7_20Sig
std::vector< double > * t_e11x11_10Sig
const Double_t pi
AlgorithmMap::const_iterator CItAlgo
iterators through map containing the algorithms
std::vector< double > * t_trackPtAll
std::vector< double > * t_h3x3Sig
edm::Service< TFileService > fs
double eHCALmatrix(const HcalTopology *topology, const DetId &det, edm::Handle< T > &hits, int ieta, int iphi, bool includeHO=false, bool algoNew=true, double hbThr=-100, double heThr=-100, double hfThr=-100, double hoThr=-100, double tMin=-500, double tMax=500, bool debug=false)
std::vector< double > * t_esimPdgId
double eNeutralHad
Definition: CaloSimInfo.h:47
void printHitPattern(HitCategory category, int position, std::ostream &stream) const
Definition: HitPattern.cc:714
double eChargedHad
Definition: CaloSimInfo.h:48
std::vector< double > * t_e11x11_20Sig
std::vector< double > * t_jetEta
std::map< std::pair< unsigned int, std::string >, int > l1AlgoMap
std::vector< double > * t_maxNearP31x31
std::vector< int > * t_trackHitInMissTEC
double chargeIsolationEcal(unsigned int trkIndex, std::vector< spr::propagatedTrackID > &vdetIds, const CaloGeometry *geo, const CaloTopology *caloTopology, int ieta, int iphi, bool debug=false)
int iEvent
Definition: GenABIO.cc:230
std::vector< double > * t_trackEta
const bool decision(const edm::Event &iEvent, const edm::InputTag &l1GtRecordInputTag, const edm::InputTag &l1GtReadoutRecordInputTag, const std::string &nameAlgoTechTrig, int &errorCode) const
Definition: L1GtUtils.cc:1448
std::vector< double > * t_esim9x9Photon
std::vector< int > * t_trackHitInMeasTID
double eta() const
pseudorapidity of momentum vector
Definition: TrackBase.h:627
std::vector< double > * t_hsim5x5Rest
std::vector< double > * t_e11x11
std::vector< double > * t_esim15x15
edm::EDGetTokenT< edm::SimTrackContainer > tok_simTk_
std::vector< double > * t_L1CenJetPt
std::vector< double > * t_trackChiSqAll
edm::InputTag L1extraTauJetSource_
std::vector< double > * t_e15x15_25Sig
std::vector< int > * t_trackHitOutMissTIB
edm::EDGetTokenT< edm::PCaloHitContainer > tok_caloHH_
std::vector< double > * t_esim15x15Rest
std::vector< double > * t_esim11x11
std::vector< double > * t_hsim7x7
std::vector< int > * t_PVNTracksHPWt
int stripTIDLayersWithMeasurement() const
Definition: HitPattern.cc:475
int trackerLayersWithMeasurement() const
Definition: HitPattern.h:912
T sqrt(T t)
Definition: SSEVec.h:48
std::vector< int > * t_trackHitOutMissTID
std::vector< int > * t_trackHitOutMissTOBTEC
vector< PseudoJet > jets
double pt() const
track transverse momentum
Definition: TrackBase.h:597
edm::InputTag HBHERecHitSource_
std::vector< int > * t_trackHitsTEC
std::vector< double > * t_e9x9_30Sig
const MagneticField * bField
edm::EDGetTokenT< edm::PCaloHitContainer > tok_caloEB_
std::vector< double > * t_L1CenJetPhi
std::vector< double > * t_L1FwdJetEta
std::vector< double > * t_h7x7
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
std::vector< int > * t_ecalSpike11x11
std::vector< double > * t_esim15x15NeutHad
std::vector< double > * t_e9x9_20Sig
std::vector< int > * t_trackHitOutMeasTOB
std::vector< int > * t_trackHitInMeasTIB
edm::SimTrackContainer::const_iterator matchedSimTrack(const edm::Event &iEvent, edm::Handle< edm::SimTrackContainer > &SimTk, edm::Handle< edm::SimVertexContainer > &SimVtx, const reco::Track *pTrack, TrackerHitAssociator &associate, bool debug=false)
unsigned short numberOfValidHits() const
number of valid hits found
Definition: TrackBase.h:796
edm::EDGetTokenT< reco::TrackCollection > tok_genTrack_
T * make(const Args &...args) const
make new ROOT object
std::vector< double > * t_trackDzPV
std::vector< int > * t_PVNTracksHP
std::vector< double > * t_L1FwdJetPhi
edm::EDGetTokenT< l1extra::L1JetParticleCollection > tok_L1extFwdJet_
edm::EDGetTokenT< l1extra::L1EmParticleCollection > tok_L1extIsoEm_
static const bool useL1EventSetup(true)
std::vector< double > * t_nTrksJetVtx
bool first
Definition: L1TdeRCT.cc:75
bool isValid() const
Definition: HandleBase.h:76
std::vector< double > * t_esim7x7NeutHad
std::vector< double > * t_hsim3x3NeutHad
std::vector< double > * t_L1NonIsoEMPhi
std::vector< double > * t_trackDzBS
void analyze(const edm::Event &, const edm::EventSetup &)
std::vector< int > * t_trackHitOutMissTOB
std::vector< double > * t_L1CenJetEta
std::vector< double > * t_trackDzPVAll
JetCorrectorParameters corr
Definition: classes.h:5
edm::EDGetTokenT< l1extra::L1MuonParticleCollection > tok_L1extMu_
std::vector< double > * t_esim7x7Matched
std::vector< double > * t_esim9x9
std::vector< double > * t_hsim5x5
std::vector< double > * t_h5x5
double qoverpError() const
error on signed transverse curvature
Definition: TrackBase.h:733
std::vector< double > * t_trackDxyPV
std::vector< double > * t_hsim5x5NeutHad
std::vector< double > * t_trackChiSq
double dz() const
dz parameter (= dsz/cos(lambda)). This is the track z0 w.r.t (0,0,0) only if the refPoint is close to...
Definition: TrackBase.h:585
std::vector< double > * t_hsim5x5Photon
std::vector< double > * t_PVy
std::vector< int > * t_L1PreScale
void getL1GtRunCache(const edm::Run &, const edm::EventSetup &, const bool, const bool, const edm::InputTag &)
get all the run-constant quantities for L1 trigger and cache them
Definition: L1GtUtils.cc:311
std::vector< double > * t_hsim3x3Rest
int stripTIDLayersWithoutMeasurement(HitCategory category) const
Definition: HitPattern.cc:543
Definition: DetId.h:18
edm::InputTag L1extraNonIsoEmSource_
std::vector< int > * t_PVisValid
static TrackQuality qualityByName(const std::string &name)
Definition: TrackBase.cc:102
std::vector< double > * t_L1METPhi
std::vector< double > * t_esim9x9Rest
std::vector< int > * t_trackHitInMissTIB
TFileDirectory mkdir(const std::string &dir, const std::string &descr="")
create a new subdirectory
Definition: TFileService.h:69
std::vector< double > * t_e9x9_15Sig
std::vector< double > * t_e7x7_15Sig
std::vector< double > * t_trackPhi
T const * product() const
Definition: Handle.h:81
std::vector< double > * t_e7x7_10Sig
const HitPattern & hitPattern() const
Access the hit pattern, indicating in which Tracker layers the track has hits.
Definition: TrackBase.h:364
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
std::vector< double > * t_trkEcalEne
edm::EDGetTokenT< l1extra::L1EmParticleCollection > tok_L1extNoIsoEm_
std::vector< double > * t_trackDzAll
edm::EDGetTokenT< reco::VertexCollection > tok_recVtx_
TH1F * h_maxNearP31x31[NPBins][NEtaBins]
std::vector< double > * t_esim11x11CharHad
const T & get() const
Definition: EventSetup.h:55
int stripLayersWithMeasurement() const
Definition: HitPattern.h:922
int pixelBarrelLayersWithMeasurement() const
Definition: HitPattern.cc:440
T const * product() const
Definition: ESHandle.h:86
edm::EDGetTokenT< edm::PCaloHitContainer > tok_caloEE_
edm::EDGetTokenT< EcalRecHitCollection > tok_EB_
TH1F * h_maxNearP15x15[NPBins][NEtaBins]
std::vector< double > * t_e9x9_25Sig
void eECALSimInfo(const edm::Event &, const DetId &det, const CaloGeometry *geo, const CaloTopology *caloTopology, edm::Handle< T > &hitsEB, edm::Handle< T > &hitsEE, edm::Handle< edm::SimTrackContainer > &SimTk, edm::Handle< edm::SimVertexContainer > &SimVtx, const reco::Track *pTrack, TrackerHitAssociator &associate, int ieta, int iphi, caloSimInfo &info, double timeCut=150, bool debug=false)
std::vector< double > * t_esim7x7Photon
std::vector< double > * t_esim9x9NeutHad
std::vector< double > * t_L1FwdJetPt
std::vector< double > * t_simTrackP
std::vector< double > * t_hsim7x7Matched
std::vector< int > * t_trackPVIdx
bool quality(const TrackQuality) const
Track quality.
Definition: TrackBase.h:479
std::vector< int > * t_trackHitsTOB
std::vector< double > * t_esim9x9Matched
std::vector< double > * t_hsim7x7NeutHad
std::vector< double > * t_L1IsoEMPhi
std::vector< int > * t_PVNTracks
edm::EDGetTokenT< HBHERecHitCollection > tok_hbhe_
std::vector< double > * t_maxNearP21x21
edm::EventID id() const
Definition: EventBase.h:56
std::vector< double > * t_L1NonIsoEMEta
double p1[4]
Definition: TauolaWrapper.h:89
std::vector< double > * t_PVz
std::pair< math::XYZPoint, double > propagateTrackerEnd(const reco::Track *, const MagneticField *, bool debug=false)
void printTrack(const reco::Track *pTrack)
std::vector< double > * t_e15x15_20Sig
std::vector< double > * t_PVTracksSumPtHPWt
std::vector< double > * t_trackDxy
std::vector< std::string > * t_L1AlgoNames
std::vector< double > * t_hsim5x5CharHad
std::vector< double > * t_e15x15_10Sig
const AlgorithmMap & gtAlgorithmMap() const
get / set the algorithm map (by name)
int stripTECLayersWithMeasurement() const
Definition: HitPattern.cc:497
std::vector< TrackBaseRef >::const_iterator trackRef_iterator
The iteratator for the vector&lt;TrackRef&gt;
Definition: Vertex.h:37
TH1F * h_maxNearP21x21[NPBins][NEtaBins]
int stripTECLayersWithoutMeasurement(HitCategory category) const
Definition: HitPattern.cc:565
std::vector< double > * t_jetPhi
std::vector< double > * t_trackL
std::vector< double > * t_esim7x7Rest
std::vector< double > * t_trackDxyBS
tuple cout
Definition: gather_cfg.py:121
int charge() const
track electric charge
Definition: TrackBase.h:543
std::vector< double > * t_trackHcalEta
std::vector< double > * t_L1METPt
edm::InputTag L1extraMuonSource_
dbl *** dir
Definition: mlp_gen.cc:35
volatile std::atomic< bool > shutdown_flag false
std::vector< int > * t_NLayersCrossed
edm::EDGetTokenT< reco::CaloJetCollection > tok_jets_
Definition: DDAxes.h:10
double DeltaR(double eta1, double phi1, double eta2, double phi2)
edm::InputTag JetSrc_
std::vector< std::pair< DetId, double > > eHCALmatrixCell(const HcalTopology *topology, const DetId &det, edm::Handle< T > &hits, int ieta, int iphi, bool includeHO=false, double hbThr=-100, double heThr=-100, double hfThr=-100, double hoThr=-100, bool debug=false)
int weight
Definition: histoStyle.py:50
std::vector< int > * t_trackHitInMeasTOB
std::vector< double > * t_PVx
std::vector< int > * t_trackHitInMeasTEC
std::vector< double > * t_esim15x15Photon
std::vector< double > * t_e9x9
std::vector< double > * t_hsim5x5Matched
const L1GtTriggerMenu * ptrL1TriggerMenuEventSetup(int &errorCode)
return a pointer to the L1 trigger menu from event setup
Definition: L1GtUtils.cc:2228
const int prescaleFactor(const edm::Event &iEvent, const edm::InputTag &l1GtRecordInputTag, const edm::InputTag &l1GtReadoutRecordInputTag, const std::string &nameAlgoTechTrig, int &errorCode) const
return prescale factor for a given algorithm or technical trigger
Definition: L1GtUtils.cc:1485
double dxy() const
dxy parameter. (This is the transverse impact parameter w.r.t. to (0,0,0) ONLY if refPoint is close t...
Definition: TrackBase.h:567
std::vector< double > * t_esim15x15Matched
double chargeIsolationHcal(unsigned int trkIndex, std::vector< spr::propagatedTrackID > &vdetIds, const HcalTopology *topology, int ieta, int iphi, bool debug=false)
std::vector< int > * t_trackHitOutMeasTID
std::vector< double > * t_L1METEta
std::vector< double > * t_e7x7_30Sig
std::vector< int > * t_trackHitOutMissTEC
std::vector< double > * t_e11x11_25Sig
std::vector< double > * t_hsim3x3
std::vector< double > * t_trackP
static const size_t NPBins
std::vector< double > * t_hsim7x7Photon
std::vector< double > * t_jetPt
std::vector< double > * t_L1IsoEMPt
int numberOfHits(HitCategory category) const
Definition: HitPattern.h:721
std::vector< double > * t_esim11x11Photon
edm::InputTag L1GTObjectMapRcdSource_
int stripTIBLayersWithoutMeasurement(HitCategory category) const
Definition: HitPattern.cc:532
std::vector< double > * t_h5x5Sig
std::vector< int > * t_trackHitInMissTID
std::vector< double > * t_e11x11_30Sig
std::vector< int > * t_trackHitInMissTIBTID
std::vector< double > * t_esim9x9CharHad
double eECALmatrix(const DetId &detId, edm::Handle< T > &hitsEB, edm::Handle< T > &hitsEE, const CaloGeometry *geo, const CaloTopology *caloTopology, int ieta, int iphi, double ebThr=-100, double eeThr=-100, double tMin=-500, double tMax=500, bool debug=false)
edm::EDGetTokenT< EcalRecHitCollection > tok_EE_
edm::EDGetTokenT< l1extra::L1JetParticleCollection > tok_L1extTauJet_
Definition: DDAxes.h:10
std::vector< double > * t_esim15x15CharHad