CMS 3D CMS Logo
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
IsoTrig.cc
Go to the documentation of this file.
1 // -*- C++ -*-//
2 // Package: IsoTrig
3 // Class: IsoTrig
4 //
12 //
13 // Original Author: Ruchi Gupta
14 // Created: Fri May 25 12:02:48 CDT 2012
15 // $Id$
16 //
17 //
18 //#define DebugLog
19 #include "IsoTrig.h"
20 
21 
22 #include "FWCore/Common/interface/TriggerNames.h"
23 #include "Calibration/IsolatedParticles/interface/ChargeIsolation.h"
24 #include "Calibration/IsolatedParticles/interface/eCone.h"
25 #include "Calibration/IsolatedParticles/interface/eECALMatrix.h"
26 #include "HLTrigger/Timer/interface/FastTimerService.h"
27 
28 IsoTrig::IsoTrig(const edm::ParameterSet& iConfig) :
29  changed(false), t_timeL2Prod(0), t_nPixCand(0), t_nPixSeed(0), t_nGoodTk(0),
30  t_TrkhCone(0), t_TrkP(0), t_TrkselTkFlag(0), t_TrkqltyFlag(0),
31  t_TrkMissFlag(0), t_TrkPVFlag(0), t_TrkNuIsolFlag(0),
32  t_PixcandP(0), t_PixcandPt(0), t_PixcandEta(0), t_PixcandPhi(0),
33  t_PixcandMaxP(0), t_PixTrkcandP(0), t_PixTrkcandPt(0), t_PixTrkcandEta(0),
34  t_PixTrkcandPhi(0), t_PixTrkcandMaxP(0), t_PixTrkcandselTk(0),
35  t_NFcandP(0), t_NFcandPt(0), t_NFcandEta(0), t_NFcandPhi(0),
36  t_NFcandEmip(0), t_NFTrkcandP(0), t_NFTrkcandPt(0), t_NFTrkcandEta(0),
37  t_NFTrkcandPhi(0), t_NFTrkcandEmip(0), t_NFTrkMinDR(0), t_NFTrkMinDP1(0),
38  t_NFTrkselTkFlag(0), t_NFTrkqltyFlag(0), t_NFTrkMissFlag(0),
39  t_NFTrkPVFlag(0), t_NFTrkPropFlag(0), t_NFTrkChgIsoFlag(0),
40  t_NFTrkNeuIsoFlag(0), t_NFTrkMipFlag(0), t_ECone(0) {
41  //now do whatever initialization is neededOA
42  trigNames = iConfig.getUntrackedParameter<std::vector<std::string> >("Triggers");
43  PixcandTag_ = iConfig.getParameter<edm::InputTag> ("PixcandTag");
44  L1candTag_ = iConfig.getParameter<edm::InputTag> ("L1candTag");
45  L2candTag_ = iConfig.getParameter<edm::InputTag> ("L2candTag");
46  doL2L3 = iConfig.getUntrackedParameter<bool>("DoL2L3",true);
47  doTiming = iConfig.getUntrackedParameter<bool>("DoTimingTree",true);
48  doMipCutTree = iConfig.getUntrackedParameter<bool>("DoMipCutTree",true);
49  doTrkResTree = iConfig.getUntrackedParameter<bool>("DoTrkResTree",true);
50  doChgIsolTree = iConfig.getUntrackedParameter<bool>("DoChgIsolTree",true);
51  doStudyIsol = iConfig.getUntrackedParameter<bool>("DoStudyIsol",true);
52  verbosity = iConfig.getUntrackedParameter<int>("Verbosity",0);
53  theTrackQuality = iConfig.getUntrackedParameter<std::string>("TrackQuality","highPurity");
54  processName = iConfig.getUntrackedParameter<std::string>("ProcessName","HLT");
55  reco::TrackBase::TrackQuality trackQuality_=reco::TrackBase::qualityByName(theTrackQuality);
56  selectionParameters.minPt = iConfig.getUntrackedParameter<double>("MinTrackPt", 10.0);
57  selectionParameters.minQuality = trackQuality_;
58  selectionParameters.maxDxyPV = iConfig.getUntrackedParameter<double>("MaxDxyPV", 0.2);
59  selectionParameters.maxDzPV = iConfig.getUntrackedParameter<double>("MaxDzPV", 5.0);
60  selectionParameters.maxChi2 = iConfig.getUntrackedParameter<double>("MaxChi2", 5.0);
61  selectionParameters.maxDpOverP = iConfig.getUntrackedParameter<double>("MaxDpOverP", 0.1);
62  selectionParameters.minOuterHit = iConfig.getUntrackedParameter<int>("MinOuterHit", 4);
63  selectionParameters.minLayerCrossed = iConfig.getUntrackedParameter<int>("MinLayerCrossed", 8);
64  selectionParameters.maxInMiss = iConfig.getUntrackedParameter<int>("MaxInMiss", 0);
65  selectionParameters.maxOutMiss = iConfig.getUntrackedParameter<int>("MaxOutMiss", 0);
66  dr_L1 = iConfig.getUntrackedParameter<double>("IsolationL1",1.0);
67  a_coneR = iConfig.getUntrackedParameter<double>("ConeRadius",34.98);
68  a_charIsoR = a_coneR + 28.9;
69  a_neutIsoR = a_charIsoR*0.726;
70  a_mipR = iConfig.getUntrackedParameter<double>("ConeRadiusMIP",14.0);
71  a_neutR1 = iConfig.getUntrackedParameter<double>("ConeRadiusNeut1",21.0);
72  a_neutR2 = iConfig.getUntrackedParameter<double>("ConeRadiusNeut2",29.0);
73  cutMip = iConfig.getUntrackedParameter<double>("MIPCut", 1.0);
74  cutCharge = iConfig.getUntrackedParameter<double>("ChargeIsolation", 2.0);
75  cutNeutral = iConfig.getUntrackedParameter<double>("NeutralIsolation", 2.0);
76  minRunNo = iConfig.getUntrackedParameter<int>("minRun");
77  maxRunNo = iConfig.getUntrackedParameter<int>("maxRun");
78  pixelTracksSources_ = iConfig.getUntrackedParameter<std::vector<edm::InputTag> >("PixelTracksSources");
79  pixelIsolationConeSizeAtEC_ = iConfig.getUntrackedParameter<std::vector<double> >("PixelIsolationConeSizeAtEC");
80  minPTrackValue_ = iConfig.getUntrackedParameter<double>("MinPTrackValue");
81  vtxCutSeed_ = iConfig.getUntrackedParameter<double>("VertexCutSeed");
82  vtxCutIsol_ = iConfig.getUntrackedParameter<double>("VertexCutIsol");
83  tauUnbiasCone_ = iConfig.getUntrackedParameter<double>("TauUnbiasCone");
84  prelimCone_ = iConfig.getUntrackedParameter<double>("PrelimCone");
85  // define tokens for access
86  tok_lumi = consumes<LumiDetails, edm::InLumi>(edm::InputTag("lumiProducer"));
87  edm::InputTag triggerEvent_ ("hltTriggerSummaryAOD","",processName);
88  tok_trigEvt = consumes<trigger::TriggerEvent>(triggerEvent_);
89  edm::InputTag theTriggerResultsLabel ("TriggerResults","",processName);
90  tok_trigRes = consumes<edm::TriggerResults>(theTriggerResultsLabel);
91  tok_genTrack_ = consumes<reco::TrackCollection>(edm::InputTag("generalTracks"));
92  tok_recVtx_ = consumes<reco::VertexCollection>(edm::InputTag("offlinePrimaryVertices"));
93  tok_bs_ = consumes<reco::BeamSpot>(edm::InputTag("offlineBeamSpot"));
94  tok_EB_ = consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit","EcalRecHitsEB"));
95  tok_EE_ = consumes<EcalRecHitCollection>(edm::InputTag("ecalRecHit","EcalRecHitsEE"));
96  tok_hbhe_ = consumes<HBHERecHitCollection>(edm::InputTag("hbhereco"));
97  tok_pixtk_ = consumes<reco::IsolatedPixelTrackCandidateCollection>(PixcandTag_);
98  tok_l1cand_ = consumes<trigger::TriggerFilterObjectWithRefs>(L1candTag_);
99  tok_l2cand_ = consumes<reco::IsolatedPixelTrackCandidateCollection>(L2candTag_);
100  if (doTiming) {
101  tok_verthb_ = consumes<reco::VertexCollection>(edm::InputTag("hltHITPixelVerticesHB"));
102  tok_verthe_ = consumes<reco::VertexCollection>(edm::InputTag("hltHITPixelVerticesHB"));
103  tok_hlt_ = consumes<trigger::TriggerFilterObjectWithRefs>(edm::InputTag("hltL1sL1SingleJet68"));
104  tok_SeedingLayerhb = consumes<SeedingLayerSetsHits>(edm::InputTag("hltPixelLayerTripletsHITHB"));
105  tok_SeedingLayerhe = consumes<SeedingLayerSetsHits>(edm::InputTag("hltPixelLayerTripletsHITHE"));
106  tok_SiPixelRecHits = consumes<SiPixelRecHitCollection>(edm::InputTag("hltSiPixelRecHits"));
107  }
108  if(doChgIsolTree) {
109  for (unsigned int k=0; k<pixelTracksSources_.size(); ++k) {
110  // edm::InputTag pix (pixelTracksSources_[k],"",processName);
111  // tok_pixtks_.push_back(consumes<reco::TrackCollection>(pix));
112  tok_pixtks_.push_back(consumes<reco::TrackCollection>(pixelTracksSources_[k]));
113  }
114  }
115 #ifdef DebugLog
116  if (verbosity>=0) {
117  std::cout <<"Parameters read from config file \n"
118  <<"\t minPt " << selectionParameters.minPt
119  <<"\t theTrackQuality " << theTrackQuality
120  <<"\t minQuality " << selectionParameters.minQuality
121  <<"\t maxDxyPV " << selectionParameters.maxDxyPV
122  <<"\t maxDzPV " << selectionParameters.maxDzPV
123  <<"\t maxChi2 " << selectionParameters.maxChi2
124  <<"\t maxDpOverP " << selectionParameters.maxDpOverP
125  <<"\t minOuterHit " << selectionParameters.minOuterHit
126  <<"\t minLayerCrossed " << selectionParameters.minLayerCrossed
127  <<"\t maxInMiss " << selectionParameters.maxInMiss
128  <<"\t maxOutMiss " << selectionParameters.maxOutMiss
129  <<"\t a_coneR " << a_coneR
130  <<"\t a_charIsoR " << a_charIsoR
131  <<"\t a_neutIsoR " << a_neutIsoR
132  <<"\t a_mipR " << a_mipR
133  <<"\t a_neutR " << a_neutR1 << ":" << a_neutR2
134  <<"\t cuts (MIP " << cutMip << " : Charge " << cutCharge
135  <<" : Neutral " << cutNeutral << ")"
136  << std::endl;
137  std::cout <<"Charge Isolation parameters:"
138  <<"\t minPTrackValue " << minPTrackValue_
139  <<"\t vtxCutSeed " << vtxCutSeed_
140  <<"\t vtxCutIsol " << vtxCutIsol_
141  <<"\t tauUnbiasCone " << tauUnbiasCone_
142  <<"\t prelimCone " << prelimCone_
143  <<"\t pixelIsolationConeSizeAtEC";
144  for (unsigned int k=0; k<pixelIsolationConeSizeAtEC_.size(); ++k)
145  std::cout << " " << pixelIsolationConeSizeAtEC_[k];
146  std::cout << std::endl;
147  }
148 #endif
149  double pl[] = {20,30,40,60,80,120};
150  for (int i=0; i<6; ++i) pLimits[i] = pl[i];
151  rEB_ = 123.8;
152  zEE_ = 317.0;
153 }
154 
155 IsoTrig::~IsoTrig() {
156  // do anything here that needs to be done at desctruction time
157  // (e.g. close files, deallocate resources etc.)
158  if (t_timeL2Prod) delete t_timeL2Prod;
159  if (t_nPixCand) delete t_nPixCand;
160  if (t_nPixSeed) delete t_nPixSeed;
161  if (t_nGoodTk) delete t_nGoodTk;
162  if (t_TrkhCone) delete t_TrkhCone;
163  if (t_TrkP) delete t_TrkP;
164  if (t_TrkselTkFlag) delete t_TrkselTkFlag;
165  if (t_TrkqltyFlag) delete t_TrkqltyFlag;
166  if (t_TrkMissFlag) delete t_TrkMissFlag;
167  if (t_TrkPVFlag) delete t_TrkPVFlag;
168  if (t_TrkNuIsolFlag) delete t_TrkNuIsolFlag;
169  if (t_PixcandP) delete t_PixcandP;
170  if (t_PixcandPt) delete t_PixcandPt;
171  if (t_PixcandEta) delete t_PixcandEta;
172  if (t_PixcandPhi) delete t_PixcandPhi;
173  if (t_PixcandMaxP) delete t_PixcandMaxP;
174  if (t_PixTrkcandP) delete t_PixTrkcandP;
175  if (t_PixTrkcandPt) delete t_PixTrkcandPt;
176  if (t_PixTrkcandEta) delete t_PixTrkcandEta;
177  if (t_PixTrkcandPhi) delete t_PixTrkcandPhi;
178  if (t_PixTrkcandMaxP) delete t_PixTrkcandMaxP;
179  if (t_PixTrkcandselTk) delete t_PixTrkcandselTk;
180  if (t_NFcandP) delete t_NFcandP;
181  if (t_NFcandPt) delete t_NFcandPt;
182  if (t_NFcandEta) delete t_NFcandEta;
183  if (t_NFcandPhi) delete t_NFcandPhi;
184  if (t_NFcandEmip) delete t_NFcandEmip;
185  if (t_NFTrkcandP) delete t_NFTrkcandP;
186  if (t_NFTrkcandPt) delete t_NFTrkcandPt;
187  if (t_NFTrkcandEta) delete t_NFTrkcandEta;
188  if (t_NFTrkcandPhi) delete t_NFTrkcandPhi;
189  if (t_NFTrkcandEmip) delete t_NFTrkcandEmip;
190  if (t_NFTrkMinDR) delete t_NFTrkMinDR;
191  if (t_NFTrkMinDP1) delete t_NFTrkMinDP1;
192  if (t_NFTrkselTkFlag) delete t_NFTrkselTkFlag;
193  if (t_NFTrkqltyFlag) delete t_NFTrkqltyFlag;
194  if (t_NFTrkMissFlag) delete t_NFTrkMissFlag;
195  if (t_NFTrkPVFlag) delete t_NFTrkPVFlag;
196  if (t_NFTrkPropFlag) delete t_NFTrkPropFlag;
197  if (t_NFTrkChgIsoFlag) delete t_NFTrkChgIsoFlag;
198  if (t_NFTrkNeuIsoFlag) delete t_NFTrkNeuIsoFlag;
199  if (t_NFTrkMipFlag) delete t_NFTrkMipFlag;
200  if (t_ECone) delete t_ECone;
201 }
202 
203 void IsoTrig::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
204  //The following says we do not know what parameters are allowed so do no validation
205  // Please change this to state exactly what you do use, even if it is no parameters
206  edm::ParameterSetDescription desc;
207  desc.setUnknown();
208  descriptions.addDefault(desc);
209 }
210 
211 void IsoTrig::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
212 #ifdef DebugLog
213  if (verbosity%10 > 0) std::cout << "Event starts====================================" << std::endl;
214 #endif
215  int RunNo = iEvent.id().run();
216 
217  iSetup.get<IdealMagneticFieldRecord>().get(bFieldH);
218  iSetup.get<CaloGeometryRecord>().get(pG);
219  const MagneticField *bField = bFieldH.product();
220  GlobalVector BField=bField->inTesla(GlobalPoint(0,0,0));
221  bfVal = BField.mag();
222 
223  trigger::TriggerEvent triggerEvent;
224  edm::Handle<trigger::TriggerEvent> triggerEventHandle;
225  iEvent.getByToken(tok_trigEvt, triggerEventHandle);
226  if (!triggerEventHandle.isValid()) {
227  edm::LogWarning("IsoTrack") << "Error! Can't get the product hltTriggerSummaryAOD";
228 
229  } else {
230  triggerEvent = *(triggerEventHandle.product());
231  }
232  const trigger::TriggerObjectCollection& TOC(triggerEvent.getObjects());
234  edm::Handle<edm::TriggerResults> triggerResults;
235  iEvent.getByToken(tok_trigRes, triggerResults);
236 
237  edm::Handle<reco::TrackCollection> trkCollection;
238  iEvent.getByToken(tok_genTrack_, trkCollection);
239 
240  iEvent.getByToken(tok_EB_, barrelRecHitsHandle);
241  iEvent.getByToken(tok_EE_, endcapRecHitsHandle);
242 
243  iEvent.getByToken(tok_hbhe_, hbhe);
244 
245  iEvent.getByToken(tok_recVtx_, recVtxs);
246  iEvent.getByToken(tok_bs_, beamSpotH);
247  if (recVtxs->size()>0 && !((*recVtxs)[0].isFake())) {
248  leadPV = math::XYZPoint( (*recVtxs)[0].x(),(*recVtxs)[0].y(), (*recVtxs)[0].z() );
249  } else if (beamSpotH.isValid()) {
250  leadPV = beamSpotH->position();
251  }
252 #ifdef DebugLog
253  if ((verbosity/100)%10>0) {
254  std::cout << "Primary Vertex " << leadPV;
255  if (beamSpotH.isValid()) std::cout << " Beam Spot "
256  << beamSpotH->position();
257  std::cout << std::endl;
258  }
259 #endif
260  pixelTrackRefsHE.clear(); pixelTrackRefsHB.clear();
261  for (unsigned int iPix=0; iPix<pixelTracksSources_.size(); iPix++) {
262  edm::Handle<reco::TrackCollection> iPixCol;
263  iEvent.getByToken(tok_pixtks_[iPix],iPixCol);
264  if(iPixCol.isValid()){
265  for (reco::TrackCollection::const_iterator pit=iPixCol->begin(); pit!=iPixCol->end(); pit++) {
266  if(iPix==0)
267  pixelTrackRefsHB.push_back(reco::TrackRef(iPixCol,pit-iPixCol->begin()));
268  pixelTrackRefsHE.push_back(reco::TrackRef(iPixCol,pit-iPixCol->begin()));
269  }
270  }
271  }
272  if (doTiming) getGoodTracks(iEvent, trkCollection);
273 
274  for (unsigned int ifilter=0; ifilter<triggerEvent.sizeFilters();
275  ++ifilter) {
276  std::string FilterNames[7] = {"hltL1sL1SingleJet68", "hltIsolPixelTrackL2FilterHE", "ecalIsolPixelTrackFilterHE", "hltIsolPixelTrackL3FilterHE",
277  "hltIsolPixelTrackL2FilterHB", "ecalIsolPixelTrackFilterHB", "hltIsolPixelTrackL3FilterHB"};
278  std::string label = triggerEvent.filterTag(ifilter).label();
279  for(int i=0; i<7; i++) {
280  if(label==FilterNames[i]) h_Filters->Fill(i);
281  }
282  }
283  edm::InputTag lumiProducer("LumiProducer", "", "RECO");
284  edm::Handle<LumiDetails> Lumid;
285  iEvent.getLuminosityBlock().getByToken(tok_lumi, Lumid);
286  float mybxlumi=-1;
287  if (Lumid.isValid())
288  mybxlumi=Lumid->lumiValue(LumiDetails::kOCC1,iEvent.bunchCrossing())*6.37;
289 #ifdef DebugLog
290  if (verbosity%10 > 0)
291  std::cout << "RunNo " << RunNo << " EvtNo " << iEvent.id().event()
292  << " Lumi " << iEvent.luminosityBlock() << " Bunch "
293  << iEvent.bunchCrossing() << " mybxlumi " << mybxlumi
294  << std::endl;
295 #endif
296  if (!triggerResults.isValid()) {
297  edm::LogWarning("IsoTrack") << "Error! Can't get the product triggerResults";
298  // boost::shared_ptr<cms::Exception> const & error = triggerResults.whyFailed();
299  // edm::LogWarning(error->category()) << error->what();
300  } else {
301  std::vector<std::string> modules;
302  h_nHLT->Fill(triggerResults->size());
303  const edm::TriggerNames & triggerNames = iEvent.triggerNames(*triggerResults);
304 
305  const std::vector<std::string> & triggerNames_ = triggerNames.triggerNames();
306 #ifdef DebugLog
307  if (verbosity%10 > 1)
308  std::cout << "number of HLTs " << triggerNames_.size() << std::endl;
309 #endif
310  int hlt(-1), preL1(-1), preHLT(-1), prescale(-1);
311  for (unsigned int i=0; i<triggerResults->size(); i++) {
312  unsigned int triggerindx = hltConfig_.triggerIndex(triggerNames_[i]);
313  const std::vector<std::string>& moduleLabels(hltConfig_.moduleLabels(triggerindx));
314 
315  for (unsigned int in=0; in<trigNames.size(); ++in) {
316  // if (triggerNames_[i].find(trigNames[in].c_str())!=std::string::npos || triggerNames_[i]==" ") {
317  if (triggerNames_[i].find(trigNames[in].c_str())!=std::string::npos) {
318 #ifdef DebugLog
319  if (verbosity%10 > 0) std::cout << "trigger that i want " << triggerNames_[i] << " accept " << triggerResults->accept(i) << std::endl;
320 #endif
321  hlt = triggerResults->accept(i);
322  h_HLT -> Fill(hlt);
323  // if (hlt>0 || triggerNames_[i]==" ") {
324  if (hlt>0) {
325  edm::Handle<reco::IsolatedPixelTrackCandidateCollection> Pixcands;
326  iEvent.getByToken(tok_pixtk_,Pixcands);
327  edm::Handle<trigger::TriggerFilterObjectWithRefs> L1cands;
328  iEvent.getByToken(tok_l1cand_, L1cands);
329 
330  const std::pair<int,int> prescales(hltConfig_.prescaleValues(iEvent,iSetup,triggerNames_[i]));
331  preL1 = prescales.first;
332  preHLT = prescales.second;
333  prescale = preL1*preHLT;
334 #ifdef DebugLog
335  if (verbosity%10 > 0)
336  std::cout << triggerNames_[i] << " accept " << hlt << " preL1 "
337  << preL1 << " preHLT " << preHLT << std::endl;
338 #endif
339  for (int iv=0; iv<3; ++iv) vec[iv].clear();
340  if (TrigList.find(RunNo) != TrigList.end() ) {
341  TrigList[RunNo] += 1;
342  } else {
343  TrigList.insert(std::pair<unsigned int, unsigned int>(RunNo,1));
344  TrigPreList.insert(std::pair<unsigned int, std::pair<int, int>>(RunNo,prescales));
345  }
346  //loop over all trigger filters in event (i.e. filters passed)
347  for (unsigned int ifilter=0; ifilter<triggerEvent.sizeFilters();
348  ++ifilter) {
349  std::vector<int> Keys;
350  std::string label = triggerEvent.filterTag(ifilter).label();
351  //loop over keys to objects passing this filter
352  for (unsigned int imodule=0; imodule<moduleLabels.size();
353  imodule++) {
354  if (label.find(moduleLabels[imodule]) != std::string::npos) {
355 #ifdef DebugLog
356  if (verbosity%10 > 0) std::cout << "FILTERNAME " << label << std::endl;
357 #endif
358  for (unsigned int ifiltrKey=0; ifiltrKey<triggerEvent.filterKeys(ifilter).size(); ++ifiltrKey) {
359  Keys.push_back(triggerEvent.filterKeys(ifilter)[ifiltrKey]);
360  const trigger::TriggerObject& TO(TOC[Keys[ifiltrKey]]);
361  math::XYZTLorentzVector v4(TO.px(), TO.py(), TO.pz(), TO.energy());
362  if (label.find("L2Filter") != std::string::npos) {
363  vec[1].push_back(v4);
364  } else if (label.find("L3Filter") != std::string::npos) {
365  vec[2].push_back(v4);
366  } else {
367  vec[0].push_back(v4);
368  h_L1ObjEnergy->Fill(TO.energy());
369  }
370 #ifdef DebugLog
371  if (verbosity%10 > 0)
372  std::cout << "key " << ifiltrKey << " : pt " << TO.pt() << " eta " << TO.eta() << " phi " << TO.phi() << " mass " << TO.mass() << " Id " << TO.id() << std::endl;
373 #endif
374  }
375  }
376  }
377  }
378  std::vector<reco::TrackCollection::const_iterator> goodTks;
379  if (doL2L3) {
380  h_nL3Objs -> Fill(vec[2].size());
381  studyTrigger(trkCollection, goodTks);
382  } else {
383  if (trkCollection.isValid()) {
384  reco::TrackCollection::const_iterator trkItr;
385  for (trkItr=trkCollection->begin();
386  trkItr!=trkCollection->end(); trkItr++)
387  goodTks.push_back(trkItr);
388  }
389  }
390  // Now study isolation etc
391  if (doStudyIsol) studyIsolation(trkCollection, goodTks);
392  if (doTrkResTree) StudyTrkEbyP(trkCollection);
393 
394  std::pair<double,double> etaphi = etaPhiTrigger();
395  edm::Handle<reco::IsolatedPixelTrackCandidateCollection> L2cands;
396  iEvent.getByToken(tok_l2cand_,L2cands);
397  if (!L2cands.isValid()) {
398 #ifdef DebugLog
399  if (verbosity%10 > 0) std::cout << " trigCand is not valid " << std::endl;
400 #endif
401  } else {
402  if(doMipCutTree) studyMipCut(trkCollection, L2cands);
403  }
404  if (pixelTracksSources_.size()>0)
405  if(doChgIsolTree && pixelTrackRefsHE.size()>0) chgIsolation(etaphi.first, etaphi.second, trkCollection, iEvent);
406  }
407  break;
408  }
409  }
410  }
411  h_PreL1 -> Fill(preL1);
412  h_PreHLT -> Fill(preHLT);
413  h_Pre -> Fill(prescale);
414  h_PreL1wt -> Fill(preL1, mybxlumi);
415  h_PreHLTwt -> Fill(preHLT, mybxlumi);
416 
417  // check if trigger names in (new) config
418  // std::cout << "changed " <<changed << std::endl;
419  if (changed) {
420  changed = false;
421 #ifdef DebugLog
422  if ((verbosity/10)%10 > 1) {
423  std::cout << "New trigger menu found !!!" << std::endl;
424  const unsigned int n(hltConfig_.size());
425  for (unsigned itrig=0; itrig<triggerNames_.size(); itrig++) {
426  unsigned int triggerindx = hltConfig_.triggerIndex(triggerNames_[itrig]);
427  std::cout << triggerNames_[itrig] << " " << triggerindx << " ";
428  if (triggerindx >= n)
429  std::cout << "does not exist in the current menu" << std::endl;
430  else
431  std::cout << "exists" << std::endl;
432  }
433  }
434 #endif
435  }
436  }
437  if (doTiming) studyTiming(iEvent);
438 }
439 
440 void IsoTrig::clearChgIsolnTreeVectors() {
441  t_PixcandP ->clear();
442  t_PixcandPt ->clear();
443  t_PixcandEta ->clear();
444  t_PixcandPhi ->clear();
445  for(unsigned int i=0; i< t_PixcandMaxP->size(); i++)
446  t_PixcandMaxP[i].clear();
447  t_PixcandMaxP ->clear();
448  t_PixTrkcandP ->clear();
449  t_PixTrkcandPt ->clear();
450  t_PixTrkcandEta ->clear();
451  t_PixTrkcandPhi ->clear();
452  t_PixTrkcandMaxP ->clear();
453  t_PixTrkcandselTk ->clear();
454 }
455 
456 void IsoTrig::clearMipCutTreeVectors() {
457  t_NFcandP ->clear();
458  t_NFcandPt ->clear();
459  t_NFcandEta ->clear();
460  t_NFcandPhi ->clear();
461  t_NFcandEmip ->clear();
462  t_NFTrkcandP ->clear();
463  t_NFTrkcandPt ->clear();
464  t_NFTrkcandEta ->clear();
465  t_NFTrkcandPhi ->clear();
466  t_NFTrkcandEmip ->clear();
467  t_NFTrkMinDR ->clear();
468  t_NFTrkMinDP1 ->clear();
469  t_NFTrkselTkFlag ->clear();
470  t_NFTrkqltyFlag ->clear();
471  t_NFTrkMissFlag ->clear();
472  t_NFTrkPVFlag ->clear();
473  t_NFTrkPropFlag ->clear();
474  t_NFTrkChgIsoFlag->clear();
475  t_NFTrkNeuIsoFlag->clear();
476  t_NFTrkMipFlag ->clear();
477  t_ECone ->clear();
478 }
479 
480 void IsoTrig::pushChgIsolnTreeVecs(math::XYZTLorentzVector &Pixcand, math::XYZTLorentzVector &Trkcand,
481  std::vector<double> &PixMaxP, double &TrkMaxP,
482  bool &selTk) {
483  t_PixcandP ->push_back(Pixcand.r());
484  t_PixcandPt ->push_back(Pixcand.pt());
485  t_PixcandEta ->push_back(Pixcand.eta());
486  t_PixcandPhi ->push_back(Pixcand.phi());
487  t_PixcandMaxP ->push_back(PixMaxP);
488  t_PixTrkcandP ->push_back(Trkcand.r());
489  t_PixTrkcandPt ->push_back(Trkcand.pt());
490  t_PixTrkcandEta ->push_back(Trkcand.eta());
491  t_PixTrkcandPhi ->push_back(Trkcand.phi());
492  t_PixTrkcandMaxP ->push_back(TrkMaxP);
493  t_PixTrkcandselTk ->push_back(selTk);
494 
495 }
496 
497 void IsoTrig::pushMipCutTreeVecs(math::XYZTLorentzVector &NFcand,
498  math::XYZTLorentzVector &Trkcand,
499  double &EmipNFcand, double &EmipTrkcand,
500  double &mindR, double &mindP1,
501  std::vector<bool> &Flags, double hCone) {
502  t_NFcandP ->push_back(NFcand.r());
503  t_NFcandPt ->push_back(NFcand.pt());
504  t_NFcandEta ->push_back(NFcand.eta());
505  t_NFcandPhi ->push_back(NFcand.phi());
506  t_NFcandEmip ->push_back(EmipNFcand);
507  t_NFTrkcandP ->push_back(Trkcand.r());
508  t_NFTrkcandPt ->push_back(Trkcand.pt());
509  t_NFTrkcandEta ->push_back(Trkcand.eta());
510  t_NFTrkcandPhi ->push_back(Trkcand.phi());
511  t_NFTrkcandEmip ->push_back(EmipTrkcand);
512  t_NFTrkMinDR ->push_back(mindR);
513  t_NFTrkMinDP1 ->push_back(mindP1);
514  t_NFTrkselTkFlag ->push_back(Flags[0]);
515  t_NFTrkqltyFlag ->push_back(Flags[1]);
516  t_NFTrkMissFlag ->push_back(Flags[2]);
517  t_NFTrkPVFlag ->push_back(Flags[3]);
518  t_NFTrkPropFlag ->push_back(Flags[4]);
519  t_NFTrkChgIsoFlag->push_back(Flags[5]);
520  t_NFTrkNeuIsoFlag->push_back(Flags[6]);
521  t_NFTrkMipFlag ->push_back(Flags[7]);
522  t_ECone ->push_back(hCone);
523 }
524 
525 void IsoTrig::beginJob() {
526  char hname[100], htit[100];
527  std::string levels[20] = {"L1", "L2", "L3",
528  "Reco", "RecoMatch", "RecoNoMatch",
529  "L2Match", "L2NoMatch", "L3Match", "L3NoMatch",
530  "HLTTrk", "HLTGoodTrk", "HLTIsoTrk", "HLTMip", "HLTSelect",
531  "nonHLTTrk", "nonHLTGoodTrk", "nonHLTIsoTrk", "nonHLTMip", "nonHLTSelect"};
532  if (doTiming) {
533  TimingTree = fs->make<TTree>("TimingTree", "TimingTree");
534  t_timeL2Prod = new std::vector<double>();
535  t_nPixCand = new std::vector<int>();
536  t_nPixSeed = new std::vector<int>();
537  t_nGoodTk = new std::vector<int>();
538 
539  TimingTree->Branch("t_timeL2Prod", "std::vector<double>", &t_timeL2Prod);
540  TimingTree->Branch("t_nPixCand", "std::vector<int>", &t_nPixCand);
541  TimingTree->Branch("t_nPixSeed", "std::vector<int>", &t_nPixSeed);
542  TimingTree->Branch("t_nGoodTk", "std::vector<int>", &t_nGoodTk);
543  }
544  if (doTrkResTree) {
545  TrkResTree = fs->make<TTree>("TrkRestree", "TrkRestree");
546  t_TrkhCone = new std::vector<double>();
547  t_TrkP = new std::vector<double>();
548  t_TrkselTkFlag = new std::vector<bool>();
549  t_TrkqltyFlag = new std::vector<bool>();
550  t_TrkMissFlag = new std::vector<bool>();
551  t_TrkPVFlag = new std::vector<bool>();
552  t_TrkNuIsolFlag= new std::vector<bool>();
553 
554  TrkResTree->Branch("t_TrkhCone", "std::vector<double>", &t_TrkhCone);
555  TrkResTree->Branch("t_TrkP", "std::vector<double>", &t_TrkP);
556  TrkResTree->Branch("t_TrkselTkFlag", "std::vector<bool>", &t_TrkselTkFlag);
557  TrkResTree->Branch("t_TrkqltyFlag", "std::vector<bool>", &t_TrkqltyFlag);
558  TrkResTree->Branch("t_TrkMissFlag", "std::vector<bool>", &t_TrkMissFlag);
559  TrkResTree->Branch("t_TrkPVFlag", "std::vector<bool>", &t_TrkPVFlag);
560  TrkResTree->Branch("t_TrkNuIsolFlag","std::vector<bool>", &t_TrkNuIsolFlag);
561  }
562  if (doChgIsolTree) {
563  ChgIsolnTree = fs->make<TTree>("ChgIsolnTree", "ChgIsolntree");
564  t_PixcandP = new std::vector<double>();
565  t_PixcandPt = new std::vector<double>();
566  t_PixcandEta = new std::vector<double>();
567  t_PixcandPhi = new std::vector<double>();
568  t_PixcandMaxP = new std::vector<std::vector<double> >();
569  t_PixTrkcandP = new std::vector<double>();
570  t_PixTrkcandPt = new std::vector<double>();
571  t_PixTrkcandEta = new std::vector<double>();
572  t_PixTrkcandPhi = new std::vector<double>();
573  t_PixTrkcandMaxP = new std::vector<double>();
574  t_PixTrkcandselTk = new std::vector<bool>();
575 
576  ChgIsolnTree->Branch("t_PixcandP", "std::vector<double>", &t_PixcandP);
577  ChgIsolnTree->Branch("t_PixcandPt", "std::vector<double>", &t_PixcandPt);
578  ChgIsolnTree->Branch("t_PixcandEta", "std::vector<double>", &t_PixcandEta);
579  ChgIsolnTree->Branch("t_PixcandPhi", "std::vector<double>", &t_PixcandPhi);
580  ChgIsolnTree->Branch("t_PixcandMaxP", "std::vector<std::vector<double> >", &t_PixcandMaxP);
581  ChgIsolnTree->Branch("t_PixTrkcandP", "std::vector<double>", &t_PixTrkcandP);
582  ChgIsolnTree->Branch("t_PixTrkcandPt", "std::vector<double>", &t_PixTrkcandPt );
583  ChgIsolnTree->Branch("t_PixTrkcandEta", "std::vector<double>", &t_PixTrkcandEta );
584  ChgIsolnTree->Branch("t_PixTrkcandPhi", "std::vector<double>", &t_PixTrkcandPhi );
585  ChgIsolnTree->Branch("t_PixTrkcandMaxP", "std::vector<double>", &t_PixTrkcandMaxP);
586  ChgIsolnTree->Branch("t_PixTrkcandselTk","std::vector<bool>", &t_PixTrkcandselTk);
587  }
588  if (doMipCutTree) {
589  MipCutTree = fs->make<TTree>("MipCutTree", "MipCuttree");
590  t_NFcandP = new std::vector<double>();
591  t_NFcandPt = new std::vector<double>();
592  t_NFcandEta = new std::vector<double>();
593  t_NFcandPhi = new std::vector<double>();
594  t_NFcandEmip = new std::vector<double>();
595  t_NFTrkcandP = new std::vector<double>();
596  t_NFTrkcandPt = new std::vector<double>();
597  t_NFTrkcandEta = new std::vector<double>();
598  t_NFTrkcandPhi = new std::vector<double>();
599  t_NFTrkcandEmip = new std::vector<double>();
600  t_NFTrkMinDR = new std::vector<double>();
601  t_NFTrkMinDP1 = new std::vector<double>();
602  t_NFTrkselTkFlag = new std::vector<bool>();
603  t_NFTrkqltyFlag = new std::vector<bool>();
604  t_NFTrkMissFlag = new std::vector<bool>();
605  t_NFTrkPVFlag = new std::vector<bool>();
606  t_NFTrkPropFlag = new std::vector<bool>();
607  t_NFTrkChgIsoFlag= new std::vector<bool>();
608  t_NFTrkNeuIsoFlag= new std::vector<bool>();
609  t_NFTrkMipFlag = new std::vector<bool>();
610  t_ECone = new std::vector<double>();
611 
612  MipCutTree->Branch("t_NFcandP", "std::vector<double>", &t_NFcandP);
613  MipCutTree->Branch("t_NFcandPt", "std::vector<double>", &t_NFcandPt);
614  MipCutTree->Branch("t_NFcandEta", "std::vector<double>", &t_NFcandEta);
615  MipCutTree->Branch("t_NFcandPhi", "std::vector<double>", &t_NFcandPhi);
616  MipCutTree->Branch("t_NFcandEmip", "std::vector<double>", &t_NFcandEmip);
617  MipCutTree->Branch("t_NFTrkcandP", "std::vector<double>", &t_NFTrkcandP);
618  MipCutTree->Branch("t_NFTrkcandPt", "std::vector<double>", &t_NFTrkcandPt );
619  MipCutTree->Branch("t_NFTrkcandEta", "std::vector<double>", &t_NFTrkcandEta );
620  MipCutTree->Branch("t_NFTrkcandPhi", "std::vector<double>", &t_NFTrkcandPhi );
621  MipCutTree->Branch("t_NFTrkcandEmip", "std::vector<double>", &t_NFTrkcandEmip);
622  MipCutTree->Branch("t_NFTrkMinDR", "std::vector<double>", &t_NFTrkMinDR);
623  MipCutTree->Branch("t_NFTrkMinDP1", "std::vector<double>", &t_NFTrkMinDP1);
624  MipCutTree->Branch("t_NFTrkselTkFlag", "std::vector<bool>", &t_NFTrkselTkFlag);
625  MipCutTree->Branch("t_NFTrkqltyFlag", "std::vector<bool>", &t_NFTrkqltyFlag);
626  MipCutTree->Branch("t_NFTrkMissFlag", "std::vector<bool>", &t_NFTrkMissFlag);
627  MipCutTree->Branch("t_NFTrkPVFlag", "std::vector<bool>", &t_NFTrkPVFlag);
628  MipCutTree->Branch("t_NFTrkPropFlag", "std::vector<bool>", &t_NFTrkPropFlag);
629  MipCutTree->Branch("t_NFTrkChgIsoFlag","std::vector<bool>", &t_NFTrkChgIsoFlag);
630  MipCutTree->Branch("t_NFTrkNeuIsoFlag","std::vector<bool>", &t_NFTrkNeuIsoFlag);
631  MipCutTree->Branch("t_NFTrkMipFlag", "std::vector<bool>", &t_NFTrkMipFlag);
632  MipCutTree->Branch("t_ECone", "std::vector<double>", &t_ECone);
633  }
634  h_Filters = fs->make<TH1I>("h_Filters", "Filter Accepts", 10, 0, 10);
635  std::string FilterNames[7] = {"hltL1sL1SingleJet68", "hltIsolPixelTrackL2FilterHE", "ecalIsolPixelTrackFilterHE", "hltIsolPixelTrackL3FilterHE",
636  "hltIsolPixelTrackL2FilterHB", "ecalIsolPixelTrackFilterHB", "hltIsolPixelTrackL3FilterHB"};
637  for(int i=0; i<7; i++) {
638  h_Filters->GetXaxis()->SetBinLabel(i+1, FilterNames[i].c_str());
639  }
640 
641  h_nHLT = fs->make<TH1I>("h_nHLT" , "Size of trigger Names", 1000, 1, 1000);
642  h_HLT = fs->make<TH1I>("h_HLT" , "HLT accept", 3, -1, 2);
643  h_PreL1 = fs->make<TH1I>("h_PreL1", "L1 Prescale", 500, 0, 500);
644  h_PreHLT = fs->make<TH1I>("h_PreHLT", "HLT Prescale", 50, 0, 50);
645  h_Pre = fs->make<TH1I>("h_Pre", "Prescale", 3000, 0, 3000);
646 
647  h_PreL1wt = fs->make<TH1D>("h_PreL1wt", "Weighted L1 Prescale", 500, 0, 500);
648  h_PreHLTwt = fs->make<TH1D>("h_PreHLTwt", "Weighted HLT Prescale", 500, 0, 100);
649  h_L1ObjEnergy = fs->make<TH1D>("h_L1ObjEnergy", "Energy of L1Object", 500, 0.0, 500.0);
650 
651  h_EnIn = fs->make<TH1D>("h_EnInEcal", "EnergyIn Ecal", 200, 0.0, 20.0);
652  h_EnOut = fs->make<TH1D>("h_EnOutEcal", "EnergyOut Ecal", 200, 0.0, 20.0);
653  h_MipEnMatch = fs->make<TH2D>("h_MipEnMatch", "MipEn: HLT level vs Reco Level (Matched)", 200, 0.0, 20.0, 200, 0.0, 20.0);
654  h_MipEnNoMatch = fs->make<TH2D>("h_MipEnNoMatch", "MipEn: HLT level vs Reco Level (No Match Found)", 200, 0.0, 20.0, 200, 0.0, 20.0);
655 
656  if (doL2L3) {
657  h_nL3Objs = fs->make<TH1I>("h_nL3Objs", "Number of L3 objects", 10, 0, 10);
658 
659  std::string pairs[9] = {"L2L3", "L2L3Match", "L2L3NoMatch", "L3Reco", "L3RecoMatch", "L3RecoNoMatch", "NewFilterReco", "NewFilterRecoMatch", "NewFilterRecoNoMatch"};
660  for (int ipair=0; ipair<9; ipair++) {
661  sprintf(hname, "h_dEta%s", pairs[ipair].c_str());
662  sprintf(htit, "#Delta#eta for %s", pairs[ipair].c_str());
663  h_dEta[ipair] = fs->make<TH1D>(hname, htit, 200, -10.0, 10.0);
664  h_dEta[ipair]->GetXaxis()->SetTitle("d#eta");
665 
666  sprintf(hname, "h_dPhi%s", pairs[ipair].c_str());
667  sprintf(htit, "#Delta#phi for %s", pairs[ipair].c_str());
668  h_dPhi[ipair] = fs->make<TH1D>(hname, htit, 140, -7.0, 7.0);
669  h_dPhi[ipair]->GetXaxis()->SetTitle("d#phi");
670 
671  sprintf(hname, "h_dPt%s", pairs[ipair].c_str());
672  sprintf(htit, "#Delta dp_{T} for %s objects", pairs[ipair].c_str());
673  h_dPt[ipair] = fs->make<TH1D>(hname, htit, 400, -200.0, 200.0);
674  h_dPt[ipair]->GetXaxis()->SetTitle("dp_{T} (GeV)");
675 
676  sprintf(hname, "h_dP%s", pairs[ipair].c_str());
677  sprintf(htit, "#Delta p for %s objects", pairs[ipair].c_str());
678  h_dP[ipair] = fs->make<TH1D>(hname, htit, 400, -200.0, 200.0);
679  h_dP[ipair]->GetXaxis()->SetTitle("dP (GeV)");
680 
681  sprintf(hname, "h_dinvPt%s", pairs[ipair].c_str());
682  sprintf(htit, "#Delta (1/p_{T}) for %s objects", pairs[ipair].c_str());
683  h_dinvPt[ipair] = fs->make<TH1D>(hname, htit, 500, -0.4, 0.1);
684  h_dinvPt[ipair]->GetXaxis()->SetTitle("d(1/p_{T})");
685  sprintf(hname, "h_mindR%s", pairs[ipair].c_str());
686  sprintf(htit, "min(#Delta R) for %s objects", pairs[ipair].c_str());
687  h_mindR[ipair] = fs->make<TH1D>(hname, htit, 500, 0.0, 1.0);
688  h_mindR[ipair]->GetXaxis()->SetTitle("dR");
689  }
690 
691  for (int lvl=0; lvl<2; lvl++) {
692  sprintf(hname, "h_dEtaL1%s", levels[lvl+1].c_str());
693  sprintf(htit, "#Delta#eta for L1 and %s objects", levels[lvl+1].c_str());
694  h_dEtaL1[lvl] = fs->make<TH1D>(hname, htit, 400, -10.0, 10.0);
695 
696  sprintf(hname, "h_dPhiL1%s", levels[lvl+1].c_str());
697  sprintf(htit, "#Delta#phi for L1 and %s objects", levels[lvl+1].c_str());
698  h_dPhiL1[lvl] = fs->make<TH1D>(hname, htit, 280, -7.0, 7.0);
699 
700  sprintf(hname, "h_dRL1%s", levels[lvl+1].c_str());
701  sprintf(htit, "#Delta R for L1 and %s objects", levels[lvl+1].c_str());
702  h_dRL1[lvl] = fs->make<TH1D>(hname, htit, 100, 0.0, 10.0);
703  }
704  }
705 
706  int levmin = (doL2L3 ? 0 : 10);
707  for (int ilevel=levmin; ilevel<20; ilevel++) {
708  sprintf(hname, "h_p%s", levels[ilevel].c_str());
709  sprintf(htit, "p for %s objects", levels[ilevel].c_str());
710  h_p[ilevel] = fs->make<TH1D>(hname, htit, 100, 0.0, 500.0);
711  h_p[ilevel]->GetXaxis()->SetTitle("p (GeV)");
712 
713  sprintf(hname, "h_pt%s", levels[ilevel].c_str());
714  sprintf(htit, "p_{T} for %s objects", levels[ilevel].c_str());
715  h_pt[ilevel] = fs->make<TH1D>(hname, htit, 100, 0.0, 500.0);
716  h_pt[ilevel]->GetXaxis()->SetTitle("p_{T} (GeV)");
717 
718  sprintf(hname, "h_eta%s", levels[ilevel].c_str());
719  sprintf(htit, "#eta for %s objects", levels[ilevel].c_str());
720  h_eta[ilevel] = fs->make<TH1D>(hname, htit, 100, -5.0, 5.0);
721  h_eta[ilevel]->GetXaxis()->SetTitle("#eta");
722 
723  sprintf(hname, "h_phi%s", levels[ilevel].c_str());
724  sprintf(htit, "#phi for %s objects", levels[ilevel].c_str());
725  h_phi[ilevel] = fs->make<TH1D>(hname, htit, 70, -3.5, 3.50);
726  h_phi[ilevel]->GetXaxis()->SetTitle("#phi");
727  }
728 
729  std::string cuts[2] = {"HLTMatched", "HLTNotMatched"};
730  std::string cuts2[2] = {"All", "Away from L1"};
731  for (int icut=0; icut<2; icut++) {
732  sprintf(hname, "h_eMip%s", cuts[icut].c_str());
733  sprintf(htit, "eMip for %s tracks", cuts[icut].c_str());
734  h_eMip[icut] =fs->make<TH1D>(hname, htit, 200, 0.0, 10.0);
735  h_eMip[icut]->GetXaxis()->SetTitle("E_{Mip} (GeV)");
736 
737  sprintf(hname, "h_eMaxNearP%s", cuts[icut].c_str());
738  sprintf(htit, "eMaxNearP for %s tracks", cuts[icut].c_str());
739  h_eMaxNearP[icut]=fs->make<TH1D>(hname, htit, 240, -2.0, 10.0);
740  h_eMaxNearP[icut]->GetXaxis()->SetTitle("E_{MaxNearP} (GeV)");
741 
742  sprintf(hname, "h_eNeutIso%s", cuts[icut].c_str());
743  sprintf(htit, "eNeutIso for %s ", cuts[icut].c_str());
744  h_eNeutIso[icut] =fs->make<TH1D>(hname, htit, 200, 0.0, 10.0);
745  h_eNeutIso[icut]->GetXaxis()->SetTitle("E_{NeutIso} (GeV)");
746 
747  for (int kcut=0; kcut<2; ++kcut) {
748  for (int lim=0; lim<5; ++lim) {
749  sprintf(hname, "h_etaCalibTracks%sCut%dLim%d", cuts[icut].c_str(), kcut, lim);
750  sprintf(htit, "#eta for %s isolated MIP tracks (%4.1f < p < %5.1f Gev/c %s)", cuts[icut].c_str(), pLimits[lim], pLimits[lim+1], cuts2[kcut].c_str());
751  h_etaCalibTracks[lim][icut][kcut]=fs->make<TH1D>(hname, htit, 60, -30.0, 30.0);
752  h_etaCalibTracks[lim][icut][kcut]->GetXaxis()->SetTitle("i#eta");
753 
754  sprintf(hname, "h_etaMipTracks%sCut%dLim%d", cuts[icut].c_str(), kcut, lim);
755  sprintf(htit, "#eta for %s charge isolated MIP tracks (%4.1f < p < %5.1f Gev/c %s)", cuts[icut].c_str(), pLimits[lim], pLimits[lim+1], cuts2[kcut].c_str());
756  h_etaMipTracks[lim][icut][kcut]=fs->make<TH1D>(hname, htit, 60, -30.0, 30.0);
757  h_etaMipTracks[lim][icut][kcut]->GetXaxis()->SetTitle("i#eta");
758  }
759  }
760  }
761 
762  std::string ecut1[3] = {"all","HLTMatched","HLTNotMatched"};
763  std::string ecut2[2] = {"without","with"};
764  int etac[48] = {-1,-2,-3,-4,-5,-6,-7,-8,-9,-10,-11,-12,-13,-14,-15,-16,-17,-18,-19,-20,-21,-22,-23,-24,
765  1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24};
766  for (int icut=0; icut<6; icut++) {
767  // int i1 = (icut>3 ? 1 : 0);
768  int i1 = (icut>2 ? 1 : 0);
769  int i2 = icut - i1*3;
770  for (int kcut=0; kcut<48; kcut++) {
771  for (int lim=0; lim<5; ++lim) {
772  sprintf(hname, "h_eta%dEnHcal%s%s%d", etac[kcut], ecut1[i2].c_str(), ecut2[i1].c_str(), lim);
773  sprintf(htit, "HCAL energy for #eta=%d for %s tracks (p=%4.1f:%5.1f Gev) %s neutral isolation", etac[kcut], ecut1[i2].c_str(), pLimits[lim], pLimits[lim+1], ecut2[i1].c_str());
774  h_eHcal[lim][icut][kcut]=fs->make<TH1D>(hname, htit, 750, 0.0, 150.0);
775  h_eHcal[lim][icut][kcut]->GetXaxis()->SetTitle("Energy (GeV)");
776  sprintf(hname, "h_eta%dEnCalo%s%s%d", etac[kcut], ecut1[i2].c_str(), ecut2[i1].c_str(), lim);
777  sprintf(htit, "Calorimter energy for #eta=%d for %s tracks (p=%4.1f:%5.1f Gev) %s neutral isolation", etac[kcut], ecut1[i2].c_str(), pLimits[lim], pLimits[lim+1], ecut2[i1].c_str());
778  h_eCalo[lim][icut][kcut]=fs->make<TH1D>(hname, htit, 750, 0.0, 150.0);
779  h_eCalo[lim][icut][kcut]->GetXaxis()->SetTitle("Energy (GeV)");
780  }
781  }
782  }
783 }
784 
785 // ------------ method called once each job just after ending the event loop ------------
786 void IsoTrig::endJob() {
787  unsigned int preL1, preHLT;
788  std::map<unsigned int, unsigned int>::iterator itr;
789  std::map<unsigned int, const std::pair<int, int>>::iterator itrPre;
790  edm::LogWarning ("IsoTrack") << trigNames.size() << "Triggers were run. RunNo vs HLT accepts for";
791  for (unsigned int i=0; i<trigNames.size(); ++i)
792  edm::LogWarning("IsoTrack") << "[" << i << "]: " << trigNames[i];
793  unsigned int n = maxRunNo - minRunNo +1;
794  g_Pre = fs->make<TH1I>("h_PrevsRN", "PreScale Vs Run Number", n, minRunNo, maxRunNo);
795  g_PreL1 = fs->make<TH1I>("h_PreL1vsRN", "L1 PreScale Vs Run Number", n, minRunNo, maxRunNo);
796  g_PreHLT = fs->make<TH1I>("h_PreHLTvsRN", "HLT PreScale Vs Run Number", n, minRunNo, maxRunNo);
797  g_Accepts = fs->make<TH1I>("h_HLTAcceptsvsRN", "HLT Accepts Vs Run Number", n, minRunNo, maxRunNo);
798 
799  for (itr=TrigList.begin(), itrPre=TrigPreList.begin(); itr!=TrigList.end(); itr++, itrPre++) {
800  preL1 = (itrPre->second).first;
801  preHLT = (itrPre->second).second;
802 #ifdef DebugLog
803  std::cout << itr->first << " " << itr->second << " " << itrPre->first << " " << preL1 << " " << preHLT << std::endl;
804 #endif
805  g_Accepts->Fill(itr->first, itr->second);
806  g_PreL1->Fill(itr->first, preL1);
807  g_PreHLT->Fill(itr->first, preHLT);
808  g_Pre->Fill(itr->first, preL1*preHLT);
809  }
810 }
811 
812 // ------------ method called when starting to processes a run ------------
813 void IsoTrig::beginRun(edm::Run const& iRun, edm::EventSetup const& iSetup) {
814  edm::LogWarning ("IsoTrack") << "Run " << iRun.run() << " hltconfig.init "
815  << hltConfig_.init(iRun,iSetup,processName,changed);
816 }
817 
818 // ------------ method called when ending the processing of a run ------------
819 void IsoTrig::endRun(edm::Run const&, edm::EventSetup const&) {}
820 
821 // ------------ method called when starting to processes a luminosity block ------------
822 void IsoTrig::beginLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&) {}
823 // ------------ method called when ending the processing of a luminosity block ------------
824 void IsoTrig::endLuminosityBlock(edm::LuminosityBlock const&, edm::EventSetup const&) {}
825 // ------------ method fills 'descriptions' with the allowed parameters for the module ------------
826 
827 void IsoTrig::StudyTrkEbyP(edm::Handle<reco::TrackCollection>& trkCollection) {
828 
829  t_TrkselTkFlag->clear();
830  t_TrkqltyFlag->clear();
831  t_TrkMissFlag->clear();
832  t_TrkPVFlag->clear();
833  t_TrkNuIsolFlag->clear();
834  t_TrkhCone->clear();
835  t_TrkP->clear();
836 
837  if (!trkCollection.isValid()) {
838 #ifdef DebugLog
839  std::cout << "trkCollection.isValid is false" << std::endl;
840 #endif
841  } else {
842  std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr;
843  const MagneticField *bField = bFieldH.product();
844  const CaloGeometry* geo = pG.product();
845  std::vector<spr::propagatedTrackDirection> trkCaloDirections1;
846  spr::propagateCALO(trkCollection, geo, bField, theTrackQuality, trkCaloDirections1, ((verbosity/100)%10>2));
847  unsigned int nTracks=0;
848  int nRH_eMipDR=0, nNearTRKs=0;
849  std::vector<bool> selFlags;
850  for (trkDetItr = trkCaloDirections1.begin(); trkDetItr != trkCaloDirections1.end(); trkDetItr++,nTracks++) {
851  double conehmaxNearP = 0, hCone=0, eMipDR=0.0;
852  const reco::Track* pTrack = &(*(trkDetItr->trkItr));
853 #ifdef DebugLog
854  if (verbosity%10>0) std::cout << "track no. " << nTracks << " p(): " << pTrack->p() << std::endl;
855 #endif
856  if (pTrack->p() > 20) {
857  math::XYZTLorentzVector v2(pTrack->px(), pTrack->py(),
858  pTrack->pz(), pTrack->p());
859  eMipDR = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
860  trkDetItr->pointHCAL, trkDetItr->pointECAL,
861  a_mipR, trkDetItr->directionECAL, nRH_eMipDR);
862  bool selectTk = spr::goodTrack(pTrack,leadPV,selectionParameters,((verbosity/100)%10>1));
863  spr::trackSelectionParameters oneCutParameters = selectionParameters;
864  oneCutParameters.maxDxyPV = 10;
865  oneCutParameters.maxDzPV = 100;
866  oneCutParameters.maxInMiss = 2;
867  oneCutParameters.maxOutMiss= 2;
868  bool qltyFlag = spr::goodTrack(pTrack,leadPV,oneCutParameters,((verbosity/100)%10>1));
869  oneCutParameters = selectionParameters;
870  oneCutParameters.maxDxyPV = 10;
871  oneCutParameters.maxDzPV = 100;
872  bool qltyMissFlag = spr::goodTrack(pTrack,leadPV,oneCutParameters,((verbosity/100)%10>1));
873  oneCutParameters = selectionParameters;
874  oneCutParameters.maxInMiss = 2;
875  oneCutParameters.maxOutMiss= 2;
876  bool qltyPVFlag = spr::goodTrack(pTrack,leadPV,oneCutParameters,((verbosity/100)%10>1));
877 #ifdef DebugLog
878  /*
879  std::cout << "sel " << selectTk << std::endl;
880  std::cout << "ntracks " << nTracks;
881  std::cout << " a_charIsoR " << a_charIsoR;
882  std::cout << " nNearTRKs " << nNearTRKs << std::endl;
883  */
884 #endif
885  conehmaxNearP = spr::chargeIsolationCone(nTracks, trkCaloDirections1, a_charIsoR, nNearTRKs, ((verbosity/100)%10>1));
886 #ifdef DebugLog
887  /*
888  std::cout << "coneh " << conehmaxNearP << std::endl;
889  std::cout << "ok " << trkDetItr->okECAL << " " << trkDetItr->okHCAL << std::endl;
890  */
891 #endif
892  double e1 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
893  trkDetItr->pointHCAL, trkDetItr->pointECAL,
894  a_neutR1, trkDetItr->directionECAL, nRH_eMipDR);
895  double e2 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
896  trkDetItr->pointHCAL, trkDetItr->pointECAL,
897  a_neutR2, trkDetItr->directionECAL, nRH_eMipDR);
898  double e_inCone = e2 - e1;
899  bool chgIsolFlag = (conehmaxNearP < cutCharge);
900  bool mipFlag = (eMipDR < cutMip);
901  bool neuIsolFlag = (e_inCone < cutNeutral);
902  bool trkpropFlag = ((trkDetItr->okECAL) && (trkDetItr->okHCAL));
903  selFlags.clear();
904  selFlags.push_back(selectTk); selFlags.push_back(qltyFlag);
905  selFlags.push_back(qltyMissFlag); selFlags.push_back(qltyPVFlag);
906 #ifdef DebugLog
907  if (verbosity%10>0) std::cout << "emip: " << eMipDR << "<" << cutMip << "(" << mipFlag << ")"
908  << " ; ok: " << trkDetItr->okECAL << "/" << trkDetItr->okHCAL
909  << " ; chgiso: " << conehmaxNearP << "<" << cutCharge << "(" << chgIsolFlag << ")" << std::endl;
910 #endif
911 
912  if(chgIsolFlag && mipFlag && trkpropFlag) {
913  double distFromHotCell=-99.0;
914  int nRecHitsCone=-99, ietaHotCell=-99, iphiHotCell=-99;
915  GlobalPoint gposHotCell(0.,0.,0.);
916  std::vector<DetId> coneRecHitDetIds;
917  hCone = spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL,
918  trkDetItr->pointECAL,
919  a_coneR, trkDetItr->directionHCAL,
920  nRecHitsCone, coneRecHitDetIds,
921  distFromHotCell, ietaHotCell, iphiHotCell,
922  gposHotCell, -1);
923  // push vectors into the Tree
924  t_TrkselTkFlag ->push_back(selFlags[0]);
925  t_TrkqltyFlag ->push_back(selFlags[1]);
926  t_TrkMissFlag ->push_back(selFlags[2]);
927  t_TrkPVFlag ->push_back(selFlags[3]);
928  t_TrkNuIsolFlag->push_back(neuIsolFlag);
929  t_TrkhCone ->push_back(hCone);
930  t_TrkP ->push_back(pTrack->p());
931  }
932  }
933  }
934 #ifdef DebugLog
935  if (verbosity%10>0) std::cout << "Filling " << t_TrkP->size() << " tracks in TrkRestree out of " << nTracks << std::endl;
936 #endif
937  }
938  TrkResTree->Fill();
939 }
940 
941 void IsoTrig::studyTiming(const edm::Event& theEvent) {
942  t_timeL2Prod->clear(); t_nPixCand->clear(); t_nPixSeed->clear();
943 
944 #ifdef DebugLog
945  edm::Handle<SeedingLayerSetsHits> hblayers, helayers;
946  theEvent.getByToken(tok_SeedingLayerhb, hblayers);
947  theEvent.getByToken(tok_SeedingLayerhe, helayers);
948  const SeedingLayerSetsHits* layershb = hblayers.product();
949  const SeedingLayerSetsHits* layershe = helayers.product();
950  std::cout << "size of Seeding TripletLayers hb/he " << layershb->size() << "/" << layershe->size() << std::endl;
951  edm::Handle<SiPixelRecHitCollection> rchts;
952  theEvent.getByToken(tok_SiPixelRecHits, rchts);
953  const SiPixelRecHitCollection* rechits = rchts.product();
954  std::cout << "size of SiPixelRechits " << rechits->size() << std::endl;;
955 #endif
956  double tHB=0.0, tHE=0.0;
957  int nCandHB=pixelTrackRefsHB.size(), nCandHE=pixelTrackRefsHE.size();
958  int nSeedHB=0, nSeedHE=0;
959 
960  if(nCandHE>0) {
961  edm::Handle<reco::VertexCollection> pVertHB, pVertHE;
962  theEvent.getByToken(tok_verthb_,pVertHB);
963  theEvent.getByToken(tok_verthe_,pVertHE);
964  edm::Handle<trigger::TriggerFilterObjectWithRefs> l1trigobj;
965  theEvent.getByToken(tok_l1cand_, l1trigobj);
966 
967  std::vector< edm::Ref<l1extra::L1JetParticleCollection> > l1tauobjref;
968  std::vector< edm::Ref<l1extra::L1JetParticleCollection> > l1jetobjref;
969  std::vector< edm::Ref<l1extra::L1JetParticleCollection> > l1forjetobjref;
970 
971  l1trigobj->getObjects(trigger::TriggerL1TauJet, l1tauobjref);
972  l1trigobj->getObjects(trigger::TriggerL1CenJet, l1jetobjref);
973  l1trigobj->getObjects(trigger::TriggerL1ForJet, l1forjetobjref);
974 
975  double ptTriggered = -10;
976  double etaTriggered = -100;
977  double phiTriggered = -100;
978  for (unsigned int p=0; p<l1tauobjref.size(); p++) {
979  if (l1tauobjref[p]->pt()>ptTriggered) {
980  ptTriggered = l1tauobjref[p]->pt();
981  phiTriggered = l1tauobjref[p]->phi();
982  etaTriggered = l1tauobjref[p]->eta();
983  }
984  }
985  for (unsigned int p=0; p<l1jetobjref.size(); p++) {
986  if (l1jetobjref[p]->pt()>ptTriggered) {
987  ptTriggered = l1jetobjref[p]->pt();
988  phiTriggered = l1jetobjref[p]->phi();
989  etaTriggered = l1jetobjref[p]->eta();
990  }
991  }
992  for (unsigned int p=0; p<l1forjetobjref.size(); p++) {
993  if (l1forjetobjref[p]->pt()>ptTriggered) {
994  ptTriggered=l1forjetobjref[p]->pt();
995  phiTriggered=l1forjetobjref[p]->phi();
996  etaTriggered=l1forjetobjref[p]->eta();
997  }
998  }
999  for(unsigned iS=0; iS<pixelTrackRefsHE.size(); iS++) {
1000  reco::VertexCollection::const_iterator vitSel;
1001  double minDZ = 100;
1002  bool vtxMatch;
1003  for (reco::VertexCollection::const_iterator vit=pVertHE->begin(); vit!=pVertHE->end(); vit++) {
1004  if (fabs(pixelTrackRefsHE[iS]->dz(vit->position()))<minDZ) {
1005  minDZ = fabs(pixelTrackRefsHE[iS]->dz(vit->position()));
1006  vitSel = vit;
1007  }
1008  }
1009  //cut on dYX:
1010  if (minDZ!=100&&fabs(pixelTrackRefsHE[iS]->dxy(vitSel->position()))<vtxCutSeed_) vtxMatch=true;
1011  if (minDZ==100) vtxMatch=true;
1012 
1013  //select tracks not matched to triggered L1 jet
1014  double R=deltaR(etaTriggered, phiTriggered, pixelTrackRefsHE[iS]->eta(), pixelTrackRefsHE[iS]->phi());
1015  if (R>tauUnbiasCone_ && vtxMatch) nSeedHE++;
1016  }
1017  for(unsigned iS=0; iS<pixelTrackRefsHB.size(); iS++) {
1018  reco::VertexCollection::const_iterator vitSel;
1019  double minDZ = 100;
1020  bool vtxMatch;
1021  for (reco::VertexCollection::const_iterator vit=pVertHB->begin(); vit!=pVertHB->end(); vit++) {
1022  if (fabs(pixelTrackRefsHB[iS]->dz(vit->position()))<minDZ) {
1023  minDZ = fabs(pixelTrackRefsHB[iS]->dz(vit->position()));
1024  vitSel = vit;
1025  }
1026  }
1027  //cut on dYX:
1028  if (minDZ!=100&&fabs(pixelTrackRefsHB[iS]->dxy(vitSel->position()))<101.0) vtxMatch=true;
1029  if (minDZ==100) vtxMatch=true;
1030 
1031  //select tracks not matched to triggered L1 jet
1032  double R=deltaR(etaTriggered, phiTriggered, pixelTrackRefsHB[iS]->eta(), pixelTrackRefsHB[iS]->phi());
1033  if (R>1.2 && vtxMatch) nSeedHB++;
1034  }
1035 
1036  edm::Service<FastTimerService> ft;
1037  tHE = ft->queryModuleTimeByLabel(theEvent.streamID(),"hltIsolPixelTrackProdHE") ;
1038  tHB = ft->queryModuleTimeByLabel(theEvent.streamID(),"hltIsolPixelTrackProdHB");
1039 #ifdef DebugLog
1040  if (verbosity%10>0) std::cout << "(HB/HE) time: " << tHB <<"/" << tHE
1041  << " nCand: " << nCandHB << "/" << nCandHE
1042  << "nSeed: " << nSeedHB << "/" << nSeedHE
1043  << std::endl;
1044 #endif
1045  }
1046  t_timeL2Prod->push_back(tHB); t_timeL2Prod->push_back(tHE);
1047  t_nPixSeed->push_back(nSeedHB); t_nPixSeed->push_back(nSeedHE);
1048  t_nPixCand->push_back(nCandHB); t_nPixCand->push_back(nCandHE);
1049 
1050  TimingTree->Fill();
1051 }
1052 void IsoTrig::studyMipCut(edm::Handle<reco::TrackCollection>& trkCollection,
1053  edm::Handle<reco::IsolatedPixelTrackCandidateCollection>& L2cands) {
1054 
1055  clearMipCutTreeVectors();
1056 #ifdef DebugLog
1057  if (verbosity%10>0) std::cout << "inside studymipcut" << std::endl;
1058 #endif
1059  if (!trkCollection.isValid()) {
1060  edm::LogWarning("IsoTrack") << "trkCollection.isValid is false";
1061  } else {
1062  std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr;
1063  const MagneticField *bField = bFieldH.product();
1064  const CaloGeometry* geo = pG.product();
1065  std::vector<spr::propagatedTrackDirection> trkCaloDirections1;
1066  spr::propagateCALO(trkCollection, geo, bField, theTrackQuality, trkCaloDirections1, ((verbosity/100)%10>2));
1067 #ifdef DebugLog
1068  if (verbosity%10>0) std::cout << "Number of L2cands:" << L2cands->size()
1069  << " to be matched to something out of "
1070  << trkCaloDirections1.size() << " reco tracks" << std::endl;
1071 #endif
1072  for (unsigned int i=0; i<L2cands->size(); i++) {
1073  edm::Ref<reco::IsolatedPixelTrackCandidateCollection> candref =
1074  edm::Ref<reco::IsolatedPixelTrackCandidateCollection>(L2cands, i);
1075  double enIn = candref->energyIn();
1076  h_EnIn->Fill(candref->energyIn());
1077  h_EnOut->Fill(candref->energyOut());
1078  math::XYZTLorentzVector v1(candref->track()->px(),candref->track()->py(),
1079  candref->track()->pz(),candref->track()->p());
1080 #ifdef DebugLog
1081  if (verbosity%10>1)
1082  std::cout << "HLT Level Candidate eta/phi/pt/enIn:"
1083  << candref->track()->eta() << "/" << candref->track()->phi()
1084  << "/" << candref->track()->pt() << "/" << candref->energyIn()
1085  << std::endl;
1086 #endif
1087  math::XYZTLorentzVector mindRvec;
1088  double mindR=999.9, mindP1=999.9, eMipDR=0.0;
1089  std::vector<bool> selFlags;
1090  unsigned int nTracks=0;
1091  double conehmaxNearP = 0, hCone=0;
1092  for (trkDetItr = trkCaloDirections1.begin(); trkDetItr != trkCaloDirections1.end(); trkDetItr++,nTracks++){
1093  const reco::Track* pTrack = &(*(trkDetItr->trkItr));
1094  math::XYZTLorentzVector v2(pTrack->px(), pTrack->py(),
1095  pTrack->pz(), pTrack->p());
1096  double dr = dR(v1,v2);
1097  double dp1 = std::abs(1./v1.r() - 1./v2.r());
1098  // std::cout <<"This recotrack(eta/phi/pt) " << pTrack->eta() << "/"
1099  // << pTrack->phi() << "/" << pTrack->pt() << " has dr/dp= "
1100  // << dr << "/" << dp << "/" << dp1 << std::endl;
1101  if (dr<mindR) {
1102  mindR = dr;
1103  mindP1= dp1;
1104  mindRvec=v2;
1105  int nRH_eMipDR=0, nNearTRKs=0;
1106  eMipDR = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
1107  trkDetItr->pointHCAL, trkDetItr->pointECAL,
1108  a_mipR, trkDetItr->directionECAL, nRH_eMipDR);
1109  bool selectTk = spr::goodTrack(pTrack,leadPV,selectionParameters,((verbosity/100)%10>1));
1110  spr::trackSelectionParameters oneCutParameters = selectionParameters;
1111  oneCutParameters.maxDxyPV = 10;
1112  oneCutParameters.maxDzPV = 100;
1113  oneCutParameters.maxInMiss = 2;
1114  oneCutParameters.maxOutMiss= 2;
1115  bool qltyFlag = spr::goodTrack(pTrack,leadPV,oneCutParameters,((verbosity/100)%10>1));
1116  oneCutParameters = selectionParameters;
1117  oneCutParameters.maxDxyPV = 10;
1118  oneCutParameters.maxDzPV = 100;
1119  bool qltyMissFlag = spr::goodTrack(pTrack,leadPV,oneCutParameters,((verbosity/100)%10>1));
1120  oneCutParameters = selectionParameters;
1121  oneCutParameters.maxInMiss = 2;
1122  oneCutParameters.maxOutMiss= 2;
1123  bool qltyPVFlag = spr::goodTrack(pTrack,leadPV,oneCutParameters,((verbosity/100)%10>1));
1124  conehmaxNearP = spr::chargeIsolationCone(nTracks, trkCaloDirections1, a_charIsoR, nNearTRKs, ((verbosity/100)%10>1));
1125  double e1 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
1126  trkDetItr->pointHCAL, trkDetItr->pointECAL,
1127  a_neutR1, trkDetItr->directionECAL, nRH_eMipDR);
1128  double e2 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
1129  trkDetItr->pointHCAL, trkDetItr->pointECAL,
1130  a_neutR2, trkDetItr->directionECAL, nRH_eMipDR);
1131  double e_inCone = e2 - e1;
1132  bool chgIsolFlag = (conehmaxNearP < cutCharge);
1133  bool mipFlag = (eMipDR < cutMip);
1134  bool neuIsolFlag = (e_inCone < cutNeutral);
1135  bool trkpropFlag = ((trkDetItr->okECAL) && (trkDetItr->okHCAL));
1136  selFlags.clear();
1137  selFlags.push_back(selectTk); selFlags.push_back(qltyFlag);
1138  selFlags.push_back(qltyMissFlag); selFlags.push_back(qltyPVFlag);
1139  selFlags.push_back(trkpropFlag); selFlags.push_back(chgIsolFlag);
1140  selFlags.push_back(neuIsolFlag); selFlags.push_back(mipFlag);
1141  double distFromHotCell=-99.0;
1142  int nRecHitsCone=-99, ietaHotCell=-99, iphiHotCell=-99;
1143  GlobalPoint gposHotCell(0.,0.,0.);
1144  std::vector<DetId> coneRecHitDetIds;
1145  hCone = spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL,
1146  trkDetItr->pointECAL,
1147  a_coneR, trkDetItr->directionHCAL,
1148  nRecHitsCone, coneRecHitDetIds,
1149  distFromHotCell, ietaHotCell, iphiHotCell,
1150  gposHotCell, -1);
1151  }
1152  }
1153  pushMipCutTreeVecs(v1, mindRvec, enIn, eMipDR, mindR, mindP1, selFlags, hCone);
1154  fillDifferences(6, v1, mindRvec, (verbosity%10 >0));
1155  if (mindR>=0.05) {
1156  fillDifferences(8, v1, mindRvec, (verbosity%10 >0));
1157  h_MipEnNoMatch->Fill(candref->energyIn(), eMipDR);
1158  } else {
1159  fillDifferences(7, v1, mindRvec, (verbosity%10 >0));
1160  h_MipEnMatch->Fill(candref->energyIn(), eMipDR);
1161  }
1162  }
1163  }
1164  MipCutTree->Fill();
1165 }
1166 
1167 void IsoTrig::studyTrigger(edm::Handle<reco::TrackCollection>& trkCollection,
1168  std::vector<reco::TrackCollection::const_iterator>& goodTks) {
1169 
1170 #ifdef DebugLog
1171  if (verbosity%10 > 0) std::cout << "Inside StudyTrigger" << std::endl;
1172 #endif
1173  for (int j=0; j<3; j++) {
1175  for (unsigned int k=0; k<vec[j].size(); k++) {
1176 #ifdef DebugLog
1177  if (verbosity%10 > 0) std::cout << "vec[" << j << "][" << k << "] pt " << vec[j][k].pt() << " eta " << vec[j][k].eta() << " phi " << vec[j][k].phi() << std::endl;
1178 #endif
1179  fillHist(j, vec[j][k]);
1180  }
1181  }
1182 
1183  double deta, dphi, dr;
1185  for (int lvl=1; lvl<3; lvl++) {
1186  for (unsigned int i=0; i<vec[lvl].size(); i++) {
1187  deta = dEta(vec[0][0],vec[lvl][i]);
1188  dphi = dPhi(vec[0][0],vec[lvl][i]);
1189  dr = dR(vec[0][0],vec[lvl][i]);
1190 #ifdef DebugLog
1191  if (verbosity%10 > 1) std::cout << "lvl " <<lvl << " i " << i << " deta " << deta << " dphi " << dphi << " dR " << dr << std::endl;
1192 #endif
1193  h_dEtaL1[lvl-1] -> Fill(deta);
1194  h_dPhiL1[lvl-1] -> Fill(dphi);
1195  h_dRL1[lvl-1] -> Fill(std::sqrt(dr));
1196  }
1197  }
1198 
1199  math::XYZTLorentzVector mindRvec;
1200  double mindR;
1201  for (unsigned int k=0; k<vec[2].size(); ++k) {
1203  mindR=999.9;
1204 #ifdef DebugLog
1205  if (verbosity%10 > 1) std::cout << "L3obj: pt " << vec[2][k].pt() << " eta " << vec[2][k].eta() << " phi " << vec[2][k].phi() << std::endl;
1206 #endif
1207  for (unsigned int j=0; j<vec[1].size(); j++) {
1208  dr = dR(vec[2][k],vec[1][j]);
1209  if (dr<mindR) {
1210  mindR=dr;
1211  mindRvec=vec[1][j];
1212  }
1213  }
1214  fillDifferences(0, vec[2][k], mindRvec, (verbosity%10 >0));
1215  if (mindR < 0.03) {
1216  fillDifferences(1, vec[2][k], mindRvec, (verbosity%10 >0));
1217  fillHist(6, mindRvec);
1218  fillHist(8, vec[2][k]);
1219  } else {
1220  fillDifferences(2, vec[2][k], mindRvec, (verbosity%10 >0));
1221  fillHist(7, mindRvec);
1222  fillHist(9, vec[2][k]);
1223  }
1224 
1226  mindR=999.9;
1227 #ifdef DebugLog
1228  if (verbosity%10 > 0)
1229  std::cout << "Now Matching L3 track with reco: L3 Track (eta, phi) "
1230  << vec[2][k].eta() << ":" << vec[2][k].phi() << " L2 Track "
1231  << mindRvec.eta() << ":" << mindRvec.phi() << " dR "
1232  << mindR << std::endl;
1233 #endif
1234  reco::TrackCollection::const_iterator goodTk = trkCollection->end();
1235  if (trkCollection.isValid()) {
1236  double mindP(9999.9);
1237  reco::TrackCollection::const_iterator trkItr;
1238  for (trkItr=trkCollection->begin();
1239  trkItr!=trkCollection->end(); trkItr++) {
1240  math::XYZTLorentzVector v4(trkItr->px(), trkItr->py(),
1241  trkItr->pz(), trkItr->p());
1242  double deltaR = dR(v4, vec[2][k]);
1243  double dp = std::abs(v4.r()/vec[2][k].r()-1.0);
1244  if (deltaR<mindR) {
1245  mindR = deltaR;
1246  mindP = dp;
1247  mindRvec = v4;
1248  goodTk = trkItr;
1249  }
1250 #ifdef DebugLog
1251  if ((verbosity/10)%10>1 && deltaR<1.0) {
1252  std::cout << "reco track: pt " << v4.pt() << " eta " << v4.eta()
1253  << " phi " << v4.phi() << " DR " << deltaR
1254  << std::endl;
1255  }
1256 #endif
1257  }
1258 #ifdef DebugLog
1259  if (verbosity%10 > 0)
1260  std::cout << "Now Matching at Reco level in step 1 DR: " << mindR
1261  << ":" << mindP << " eta:phi " << mindRvec.eta() << ":"
1262  << mindRvec.phi() << std::endl;
1263 #endif
1264  if (mindR < 0.03 && mindP > 0.1) {
1265  for (trkItr=trkCollection->begin();
1266  trkItr!=trkCollection->end(); trkItr++) {
1267  math::XYZTLorentzVector v4(trkItr->px(), trkItr->py(),
1268  trkItr->pz(), trkItr->p());
1269  double deltaR = dR(v4, vec[2][k]);
1270  double dp = std::abs(v4.r()/vec[2][k].r()-1.0);
1271  if (dp<mindP && deltaR<0.03) {
1272  mindR = deltaR;
1273  mindP = dp;
1274  mindRvec = v4;
1275  goodTk = trkItr;
1276  }
1277  }
1278 #ifdef DebugLog
1279  if (verbosity%10 > 0)
1280  std::cout << "Now Matching at Reco level in step 2 DR: " << mindR
1281  << ":" << mindP << " eta:phi " << mindRvec.eta() << ":"
1282  << mindRvec.phi() << std::endl;
1283 #endif
1284  }
1285  fillDifferences(3, vec[2][k], mindRvec, (verbosity%10 >0));
1286  fillHist(3, mindRvec);
1287  if(mindR < 0.03) {
1288  fillDifferences(4, vec[2][k], mindRvec, (verbosity%10 >0));
1289  fillHist(4, mindRvec);
1290  } else {
1291  fillDifferences(5, vec[2][k], mindRvec, (verbosity%10 >0));
1292  fillHist(5, mindRvec);
1293  }
1294  if (goodTk != trkCollection->end()) goodTks.push_back(goodTk);
1295  }
1296  }
1297 }
1298 
1299 void IsoTrig::studyIsolation(edm::Handle<reco::TrackCollection>& trkCollection,
1300  std::vector<reco::TrackCollection::const_iterator>& goodTks) {
1301 
1302  if (trkCollection.isValid()) {
1303  // get handles to calogeometry and calotopology
1304  const CaloGeometry* geo = pG.product();
1305  const MagneticField *bField = bFieldH.product();
1306  std::vector<spr::propagatedTrackDirection> trkCaloDirections;
1307  spr::propagateCALO(trkCollection, geo, bField, theTrackQuality, trkCaloDirections, ((verbosity/100)%10>2));
1308 
1309  std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr;
1310 #ifdef DebugLog
1311  if ((verbosity/1000)%10 > 1) {
1312  std::cout << "n of barrelRecHitsHandle " << barrelRecHitsHandle->size() << std::endl;
1313  for (EcalRecHitCollection::const_iterator hit = barrelRecHitsHandle->begin(); hit != barrelRecHitsHandle->end(); ++hit) {
1314  std::cout << "hit : detid(ieta,iphi) " << (EBDetId)(hit->id()) << " time " << hit->time() << " energy " << hit->energy() << std::endl;
1315  }
1316  std::cout << "n of endcapRecHitsHandle " << endcapRecHitsHandle->size() << std::endl;
1317  for (EcalRecHitCollection::const_iterator hit = endcapRecHitsHandle->begin(); hit != endcapRecHitsHandle->end(); ++hit) {
1318  std::cout << "hit : detid(ieta,iphi) " << (EEDetId)(hit->id()) << " time " << hit->time() << " energy " << hit->energy() << std::endl;
1319  }
1320  std::cout << "n of hbhe " << hbhe->size() << std::endl;
1321  for (HBHERecHitCollection::const_iterator hit = hbhe->begin(); hit != hbhe->end(); ++hit) {
1322  std::cout << "hit : detid(ieta,iphi) " << hit->id() << " time " << hit->time() << " energy " << hit->energy() << std::endl;
1323  }
1324  }
1325 #endif
1326  unsigned int nTracks=0, ngoodTk=0, nselTk=0;
1327  int ieta=999;
1328  for (trkDetItr = trkCaloDirections.begin(); trkDetItr != trkCaloDirections.end(); trkDetItr++,nTracks++){
1329  bool l3Track = (std::find(goodTks.begin(),goodTks.end(),trkDetItr->trkItr) != goodTks.end());
1330  const reco::Track* pTrack = &(*(trkDetItr->trkItr));
1331  math::XYZTLorentzVector v4(pTrack->px(), pTrack->py(),
1332  pTrack->pz(), pTrack->p());
1333  bool selectTk = spr::goodTrack(pTrack,leadPV,selectionParameters,((verbosity/100)%10>1));
1334  double eMipDR=9999., e_inCone=0, conehmaxNearP=0, mindR=999.9, hCone=0;
1335  if (trkDetItr->okHCAL) {
1336  HcalDetId detId = (HcalDetId)(trkDetItr->detIdHCAL);
1337  ieta = detId.ieta();
1338  }
1339  for (unsigned k=0; k<vec[0].size(); ++k) {
1340  double deltaR = dR(v4, vec[0][k]);
1341  if (deltaR<mindR) mindR = deltaR;
1342  }
1343 #ifdef DebugLog
1344  if ((verbosity/100)%10 > 1) std::cout << "Track ECAL " << trkDetItr->okECAL << " HCAL " << trkDetItr->okHCAL << " Flag " << selectTk << std::endl;
1345 #endif
1346  if (selectTk && trkDetItr->okECAL && trkDetItr->okHCAL) {
1347  ngoodTk++;
1348  int nRH_eMipDR=0, nNearTRKs=0;
1349  double e1 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
1350  trkDetItr->pointHCAL, trkDetItr->pointECAL,
1351  a_neutR1, trkDetItr->directionECAL, nRH_eMipDR);
1352  double e2 = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
1353  trkDetItr->pointHCAL, trkDetItr->pointECAL,
1354  a_neutR2, trkDetItr->directionECAL, nRH_eMipDR);
1355  eMipDR = spr::eCone_ecal(geo, barrelRecHitsHandle, endcapRecHitsHandle,
1356  trkDetItr->pointHCAL, trkDetItr->pointECAL,
1357  a_mipR, trkDetItr->directionECAL, nRH_eMipDR);
1358  conehmaxNearP = spr::chargeIsolationCone(nTracks, trkCaloDirections, a_charIsoR, nNearTRKs, ((verbosity/100)%10>1));
1359  e_inCone = e2 - e1;
1360  double distFromHotCell=-99.0;
1361  int nRecHitsCone=-99, ietaHotCell=-99, iphiHotCell=-99;
1362  GlobalPoint gposHotCell(0.,0.,0.);
1363  std::vector<DetId> coneRecHitDetIds;
1364  hCone = spr::eCone_hcal(geo, hbhe, trkDetItr->pointHCAL,
1365  trkDetItr->pointECAL,
1366  a_coneR, trkDetItr->directionHCAL,
1367  nRecHitsCone, coneRecHitDetIds,
1368  distFromHotCell, ietaHotCell, iphiHotCell,
1369  gposHotCell, -1);
1370  if (eMipDR<1.0) nselTk++;
1371  }
1372  if (l3Track) {
1373  fillHist(10,v4);
1374  if (selectTk) {
1375  fillHist(11,v4);
1376  fillCuts(0, eMipDR, conehmaxNearP, e_inCone, v4, ieta, (mindR>dr_L1));
1377  if (conehmaxNearP < cutCharge) {
1378  fillHist(12,v4);
1379  if (eMipDR < cutMip) {
1380  fillHist(13,v4);
1381  fillEnergy(1, ieta, hCone, eMipDR, v4);
1382  fillEnergy(0, ieta, hCone, eMipDR, v4);
1383  if (e_inCone < cutNeutral) {
1384  fillHist(14, v4);
1385  fillEnergy(4, ieta, hCone, eMipDR, v4);
1386  fillEnergy(3, ieta, hCone, eMipDR, v4);
1387  }
1388  }
1389  }
1390  }
1391  } else if (doL2L3) {
1392  fillHist(15,v4);
1393  if (selectTk) {
1394  fillHist(16,v4);
1395  fillCuts(1, eMipDR, conehmaxNearP, e_inCone, v4, ieta, (mindR>dr_L1));
1396  if (conehmaxNearP < cutCharge) {
1397  fillHist(17,v4);
1398  if (eMipDR < cutMip) {
1399  fillHist(18,v4);
1400  fillEnergy(2, ieta, hCone, eMipDR, v4);
1401  fillEnergy(0, ieta, hCone, eMipDR, v4);
1402  if (e_inCone < cutNeutral) {
1403  fillHist(19, v4);
1404  fillEnergy(5, ieta, hCone, eMipDR, v4);
1405  fillEnergy(3, ieta, hCone, eMipDR, v4);
1406  }
1407  }
1408  }
1409  }
1410  }
1411  }
1412  // std::cout << "Number of tracks selected offline " << nselTk << std::endl;
1413  }
1414 }
1415 
1416 void IsoTrig::chgIsolation(double& etaTriggered, double& phiTriggered,
1417  edm::Handle<reco::TrackCollection>& trkCollection,
1418  const edm::Event& theEvent) {
1419  clearChgIsolnTreeVectors();
1420 #ifdef DebugLog
1421  if (verbosity%10>0) std::cout << "Inside chgIsolation() with eta/phi Triggered: " << etaTriggered << "/" << phiTriggered << std::endl;
1422 #endif
1423  std::vector<double> maxP;
1424 
1425  std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr;
1426  const MagneticField *bField = bFieldH.product();
1427  const CaloGeometry* geo = pG.product();
1428  std::vector<spr::propagatedTrackDirection> trkCaloDirections1;
1429  spr::propagateCALO(trkCollection, geo, bField, theTrackQuality, trkCaloDirections1, ((verbosity/100)%10>2));
1430 #ifdef DebugLog
1431  if (verbosity%10>0) std::cout << "Propagated TrkCollection" << std::endl;
1432 #endif
1433  for (unsigned int k=0; k<pixelIsolationConeSizeAtEC_.size(); ++k)
1434  maxP.push_back(0);
1435  unsigned i = pixelTrackRefsHE.size();
1436  std::vector<std::pair<unsigned int, std::pair<double, double>>> VecSeedsatEC;
1437  //loop to select isolated tracks
1438  for (unsigned iS=0; iS<pixelTrackRefsHE.size(); iS++) {
1439  if (pixelTrackRefsHE[iS]->p() > minPTrackValue_) {
1440  bool vtxMatch = false;
1441  //associate to vertex (in Z)
1442  unsigned int ivSel = recVtxs->size();
1443  double minDZ = 100;
1444  for (unsigned int iv = 0; iv < recVtxs->size(); ++iv) {
1445  if (fabs(pixelTrackRefsHE[iS]->dz((*recVtxs)[iv].position()))<minDZ) {
1446  minDZ = fabs(pixelTrackRefsHE[iS]->dz((*recVtxs)[iv].position()));
1447  ivSel = iv;
1448  }
1449  }
1450  //cut on dYX:
1451  if (ivSel == recVtxs->size()) {
1452  vtxMatch = true;
1453  } else if (fabs(pixelTrackRefsHE[iS]->dxy((*recVtxs)[ivSel].position()))<vtxCutSeed_){
1454  vtxMatch = true;
1455  }
1456  //select tracks not matched to triggered L1 jet
1457  double R = deltaR(etaTriggered, phiTriggered, pixelTrackRefsHE[iS]->eta(),
1458  pixelTrackRefsHE[iS]->phi());
1459  if (R > tauUnbiasCone_ && vtxMatch) {
1460  //propagate seed track to ECAL surface:
1461  std::pair<double,double> seedCooAtEC;
1462  // in case vertex is found:
1463  if (minDZ!=100) seedCooAtEC=GetEtaPhiAtEcal(pixelTrackRefsHE[iS]->eta(), pixelTrackRefsHE[iS]->phi(), pixelTrackRefsHE[iS]->pt(), pixelTrackRefsHE[iS]->charge(), (*recVtxs)[ivSel].z());
1464  //in case vertex is not found:
1465  else seedCooAtEC=GetEtaPhiAtEcal(pixelTrackRefsHE[iS]->eta(), pixelTrackRefsHE[iS]->phi(), pixelTrackRefsHE[iS]->pt(), pixelTrackRefsHE[iS]->charge(), 0);
1466  VecSeedsatEC.push_back(std::make_pair(iS, seedCooAtEC));
1467  }
1468  }
1469  }
1470  for (unsigned int l=0; l<VecSeedsatEC.size(); l++) {
1471  unsigned int iSeed = VecSeedsatEC[l].first;
1472  math::XYZTLorentzVector v1(pixelTrackRefsHE[iSeed]->px(),pixelTrackRefsHE[iSeed]->py(),
1473  pixelTrackRefsHE[iSeed]->pz(),pixelTrackRefsHE[iSeed]->p());
1474 
1475  for (unsigned int j=0; j<VecSeedsatEC.size(); j++) {
1476  unsigned int iSurr = VecSeedsatEC[j].first;
1477  if (iSeed != iSurr) {
1478  //define preliminary cone around seed track impact point from which tracks will be extrapolated:
1479  // edm::Ref<reco::IsolatedPixelTrackCandidateCollection> cand2ref =
1480  // edm::Ref<reco::IsolatedPixelTrackCandidateCollection>(L2cands, iSurr);
1481  if (deltaR(pixelTrackRefsHE[iSeed]->eta(), pixelTrackRefsHE[iSeed]->phi(), pixelTrackRefsHE[iSurr]->eta(),
1482  pixelTrackRefsHE[iSurr]->phi()) < prelimCone_) {
1483  unsigned int ivSel = recVtxs->size();
1484  double minDZ2=100;
1485  for (unsigned int iv = 0; iv < recVtxs->size(); ++iv) {
1486  if (fabs(pixelTrackRefsHE[iSurr]->dz((*recVtxs)[iv].position()))<minDZ2) {
1487  minDZ2 = fabs(pixelTrackRefsHE[iSurr]->dz((*recVtxs)[iv].position()));
1488  ivSel = iv;
1489  }
1490  }
1491  //cut ot dXY:
1492  if (minDZ2==100 || fabs(pixelTrackRefsHE[iSurr]->dxy((*recVtxs)[ivSel].position()))<vtxCutIsol_) {
1493  //calculate distance at ECAL surface and update isolation:
1494  double dist = getDistInCM(VecSeedsatEC[i].second.first, VecSeedsatEC[i].second.second, VecSeedsatEC[j].second.first, VecSeedsatEC[j].second.second);
1495  for (unsigned int k=0; k<pixelIsolationConeSizeAtEC_.size(); ++k) {
1496  if (dist<pixelIsolationConeSizeAtEC_[k]) {
1497  if (pixelTrackRefsHE[iSurr]->p() > maxP[k])
1498  maxP[k] = pixelTrackRefsHE[iSurr]->p();
1499  }
1500  }
1501  }
1502  }
1503  }
1504  }
1505 
1506  double conehmaxNearP = -1; bool selectTk=false;
1507  double mindR=999.9; int nTracks=0;
1508  math::XYZTLorentzVector mindRvec;
1509  for (trkDetItr = trkCaloDirections1.begin(); trkDetItr != trkCaloDirections1.end(); trkDetItr++, nTracks++){
1510  int nNearTRKs=0;
1511  const reco::Track* pTrack = &(*(trkDetItr->trkItr));
1512  math::XYZTLorentzVector v2(pTrack->px(), pTrack->py(),
1513  pTrack->pz(), pTrack->p());
1514  double dr = dR(v1,v2);
1515  if (dr<mindR) {
1516  selectTk = spr::goodTrack(pTrack,leadPV,selectionParameters,((verbosity/100)%10>1));
1517  conehmaxNearP = spr::chargeIsolationCone(nTracks, trkCaloDirections1, a_charIsoR, nNearTRKs, ((verbosity/100)%10>1));
1518  mindR = dr;
1519  mindRvec = v2;
1520  }
1521  }
1522  pushChgIsolnTreeVecs(v1, mindRvec, maxP, conehmaxNearP, selectTk);
1523  }
1524  ChgIsolnTree->Fill();
1525 }
1526 
1527 void IsoTrig::getGoodTracks(const edm::Event& iEvent,
1528  edm::Handle<reco::TrackCollection>& trkCollection){
1529 
1530  t_nGoodTk->clear();
1531  std::vector<int> nGood(4,0);
1532  if (trkCollection.isValid()) {
1533  // get handles to calogeometry and calotopology
1534  const CaloGeometry* geo = pG.product();
1535  const MagneticField *bField = bFieldH.product();
1536  std::vector<spr::propagatedTrackDirection> trkCaloDirections;
1537  spr::propagateCALO(trkCollection, geo, bField, theTrackQuality, trkCaloDirections, ((verbosity/100)%10>2));
1538 
1539  // get the trigger jet
1540  edm::Handle<trigger::TriggerFilterObjectWithRefs> l1trigobj;
1541  iEvent.getByToken(tok_l1cand_, l1trigobj);
1542 
1543  std::vector< edm::Ref<l1extra::L1JetParticleCollection> > l1tauobjref;
1544  l1trigobj->getObjects(trigger::TriggerL1TauJet, l1tauobjref);
1545  std::vector< edm::Ref<l1extra::L1JetParticleCollection> > l1jetobjref;
1546  l1trigobj->getObjects(trigger::TriggerL1CenJet, l1jetobjref);
1547  std::vector< edm::Ref<l1extra::L1JetParticleCollection> > l1forjetobjref;
1548  l1trigobj->getObjects(trigger::TriggerL1ForJet, l1forjetobjref);
1549 
1550  double ptTriggered(-10), etaTriggered(-100), phiTriggered(-100);
1551  for (unsigned int p=0; p<l1tauobjref.size(); p++) {
1552  if (l1tauobjref[p]->pt()>ptTriggered) {
1553  ptTriggered = l1tauobjref[p]->pt();
1554  phiTriggered = l1tauobjref[p]->phi();
1555  etaTriggered = l1tauobjref[p]->eta();
1556  }
1557  }
1558  for (unsigned int p=0; p<l1jetobjref.size(); p++) {
1559  if (l1jetobjref[p]->pt()>ptTriggered) {
1560  ptTriggered = l1jetobjref[p]->pt();
1561  phiTriggered = l1jetobjref[p]->phi();
1562  etaTriggered = l1jetobjref[p]->eta();
1563  }
1564  }
1565  for (unsigned int p=0; p<l1forjetobjref.size(); p++) {
1566  if (l1forjetobjref[p]->pt()>ptTriggered) {
1567  ptTriggered=l1forjetobjref[p]->pt();
1568  phiTriggered=l1forjetobjref[p]->phi();
1569  etaTriggered=l1forjetobjref[p]->eta();
1570  }
1571  }
1572  double pTriggered = ptTriggered*cosh(etaTriggered);
1573  math::XYZTLorentzVector pTrigger(ptTriggered*cos(phiTriggered),
1574  ptTriggered*sin(phiTriggered),
1575  pTriggered*tanh(etaTriggered), pTriggered);
1576 
1577  std::vector<spr::propagatedTrackDirection>::const_iterator trkDetItr;
1578  unsigned int nTracks(0);
1579  for (trkDetItr = trkCaloDirections.begin(); trkDetItr != trkCaloDirections.end(); trkDetItr++,nTracks++){
1580  const reco::Track* pTrack = &(*(trkDetItr->trkItr));
1581  math::XYZTLorentzVector v4(pTrack->px(), pTrack->py(),
1582  pTrack->pz(), pTrack->p());
1583  bool selectTk = spr::goodTrack(pTrack,leadPV,selectionParameters,((verbosity/100)%10>1));
1584  double mindR = dR(v4, pTrigger);
1585 #ifdef DebugLog
1586  if ((verbosity/100)%10 > 1) std::cout << "Track ECAL " << trkDetItr->okECAL << " HCAL " << trkDetItr->okHCAL << " Flag " << selectTk << std::endl;
1587 #endif
1588  if (selectTk && trkDetItr->okECAL && trkDetItr->okHCAL && mindR > 1.0) {
1589  int nRH_eMipDR(0), nNearTRKs(0);
1590  double eMipDR = spr::eCone_ecal(geo, barrelRecHitsHandle,
1591  endcapRecHitsHandle,
1592  trkDetItr->pointHCAL,
1593  trkDetItr->pointECAL, a_mipR,
1594  trkDetItr->directionECAL, nRH_eMipDR);
1595  double conehmaxNearP = spr::chargeIsolationCone(nTracks,
1596  trkCaloDirections,
1597  a_charIsoR, nNearTRKs,
1598  ((verbosity/100)%10>1));
1599  if (conehmaxNearP < 2.0 && eMipDR<1.0) {
1600  if (pTrack->p() >= 20 && pTrack->p() < 30) {
1601  ++nGood[0];
1602  } else if (pTrack->p() >= 30 && pTrack->p() < 40) {
1603  ++nGood[1];
1604  } else if (pTrack->p() >= 40 && pTrack->p() < 60) {
1605  ++nGood[2];
1606  } else if (pTrack->p() >= 60 && pTrack->p() < 100) {
1607  ++nGood[3];
1608  }
1609  }
1610  }
1611  }
1612  }
1613 
1614  for (unsigned int ii=0; ii<nGood.size(); ++ii)
1615  t_nGoodTk->push_back(nGood[ii]);
1616 }
1617 
1618 void IsoTrig::fillHist(int indx, math::XYZTLorentzVector& vec) {
1619  h_p[indx]->Fill(vec.r());
1620  h_pt[indx]->Fill(vec.pt());
1621  h_eta[indx]->Fill(vec.eta());
1622  h_phi[indx]->Fill(vec.phi());
1623 }
1624 
1625 void IsoTrig::fillDifferences(int indx, math::XYZTLorentzVector& vec1,
1626  math::XYZTLorentzVector& vec2, bool debug) {
1627  double dr = dR(vec1,vec2);
1628  double deta = dEta(vec1, vec2);
1629  double dphi = dPhi(vec1, vec2);
1630  double dpt = dPt(vec1, vec2);
1631  double dp = dP(vec1, vec2);
1632  double dinvpt = dinvPt(vec1, vec2);
1633  h_dEta[indx] ->Fill(deta);
1634  h_dPhi[indx] ->Fill(dphi);
1635  h_dPt[indx] ->Fill(dpt);
1636  h_dP[indx] ->Fill(dp);
1637  h_dinvPt[indx]->Fill(dinvpt);
1638  h_mindR[indx] ->Fill(dr);
1639 #ifdef DebugLog
1640  if (debug) std::cout << "mindR for index " << indx << " is " << dr << " deta " << deta << " dphi " << dphi << " dpt " << dpt << " dinvpt " << dinvpt <<std::endl;
1641 #endif
1642 }
1643 
1644 void IsoTrig::fillCuts(int indx, double eMipDR, double conehmaxNearP, double e_inCone, math::XYZTLorentzVector& vec, int ieta, bool cut) {
1645  h_eMip[indx] ->Fill(eMipDR);
1646  h_eMaxNearP[indx]->Fill(conehmaxNearP);
1647  h_eNeutIso[indx] ->Fill(e_inCone);
1648  if ((conehmaxNearP < cutCharge) && (eMipDR < cutMip)) {
1649  for (int lim=0; lim<5; ++lim) {
1650  if ((vec.r()>pLimits[lim]) && (vec.r()<=pLimits[lim+1])) {
1651  h_etaMipTracks[lim][indx][0]->Fill((double)(ieta));
1652  if (cut) h_etaMipTracks[lim][indx][1]->Fill((double)(ieta));
1653  if (e_inCone < cutNeutral) {
1654  h_etaCalibTracks[lim][indx][0]->Fill((double)(ieta));
1655  if (cut) h_etaCalibTracks[lim][indx][1]->Fill((double)(ieta));
1656  }
1657  }
1658  }
1659  }
1660 }
1661 
1662 void IsoTrig::fillEnergy(int indx, int ieta, double hCone, double eMipDR, math::XYZTLorentzVector& vec) {
1663  int kk(-1);
1664  if (ieta > 0 && ieta < 25) kk = 23 + ieta;
1665  else if (ieta > -25 && ieta < 0) kk = -(ieta + 1);
1666  if (kk >= 0 && eMipDR > 0.01 && hCone > 1.0) {
1667  for (int lim=0; lim<5; ++lim) {
1668  if ((vec.r()>pLimits[lim]) && (vec.r()<=pLimits[lim+1])) {
1669  h_eHcal[lim][indx][kk] ->Fill(hCone);
1670  h_eCalo[lim][indx][kk] ->Fill(hCone+eMipDR);
1671  }
1672  }
1673  }
1674 }
1675 
1676 double IsoTrig::dEta(math::XYZTLorentzVector& vec1, math::XYZTLorentzVector& vec2) {
1677  return (vec1.eta()-vec2.eta());
1678 }
1679 
1680 double IsoTrig::dPhi(math::XYZTLorentzVector& vec1, math::XYZTLorentzVector& vec2) {
1681 
1682  double phi1 = vec1.phi();
1683  if (phi1 < 0) phi1 += 2.0*M_PI;
1684  double phi2 = vec2.phi();
1685  if (phi2 < 0) phi2 += 2.0*M_PI;
1686  double dphi = phi1-phi2;
1687  if (dphi > M_PI) dphi -= 2.*M_PI;
1688  else if (dphi < -M_PI) dphi += 2.*M_PI;
1689  return dphi;
1690 }
1691 
1692 double IsoTrig::dR(math::XYZTLorentzVector& vec1, math::XYZTLorentzVector& vec2) {
1693  double deta = dEta(vec1,vec2);
1694  double dphi = dPhi(vec1,vec2);
1695  return std::sqrt(deta*deta + dphi*dphi);
1696 }
1697 
1698 double IsoTrig::dPt(math::XYZTLorentzVector& vec1, math::XYZTLorentzVector& vec2) {
1699  return (vec1.pt()-vec2.pt());
1700 }
1701 
1702 double IsoTrig::dP(math::XYZTLorentzVector& vec1, math::XYZTLorentzVector& vec2) {
1703  return (std::abs(vec1.r()-vec2.r()));
1704 }
1705 
1706 double IsoTrig::dinvPt(math::XYZTLorentzVector& vec1, math::XYZTLorentzVector& vec2) {
1707  return ((1/vec1.pt())-(1/vec2.pt()));
1708 }
1709 
1710 std::pair<double,double> IsoTrig::etaPhiTrigger() {
1711  double eta(0), phi(0), ptmax(0);
1712  for (unsigned int k=0; k<vec[0].size(); ++k) {
1713  if (k == 0) {
1714  eta = vec[0][k].eta();
1715  phi = vec[0][k].phi();
1716  ptmax = vec[0][k].pt();
1717  } else if (vec[0][k].pt() > ptmax) {
1718  eta = vec[0][k].eta();
1719  phi = vec[0][k].phi();
1720  ptmax = vec[0][k].pt();
1721  }
1722  }
1723  return std::pair<double,double>(eta,phi);
1724 }
1725 
1726 std::pair<double,double> IsoTrig::GetEtaPhiAtEcal(double etaIP, double phiIP,
1727  double pT, int charge,
1728  double vtxZ) {
1729 
1730  double deltaPhi=0;
1731  double etaEC = 100;
1732  double phiEC = 100;
1733 
1734  double Rcurv = 9999999;
1735  if (bfVal!=0) Rcurv=pT*33.3*100/(bfVal*10); //r(m)=pT(GeV)*33.3/B(kG)
1736 
1737  double ecDist = zEE_;
1738  double ecRad = rEB_; //radius of ECAL barrel (cm)
1739  double theta = 2*atan(exp(-etaIP));
1740  double zNew = 0;
1741  if (theta>0.5*M_PI) theta=M_PI-theta;
1742  if (fabs(etaIP)<1.479) {
1743  if ((0.5*ecRad/Rcurv)>1) {
1744  etaEC = 10000;
1745  deltaPhi = 0;
1746  } else {
1747  deltaPhi =-charge*asin(0.5*ecRad/Rcurv);
1748  double alpha1 = 2*asin(0.5*ecRad/Rcurv);
1749  double z = ecRad/tan(theta);
1750  if (etaIP>0) zNew = z*(Rcurv*alpha1)/ecRad+vtxZ; //new z-coordinate of track
1751  else zNew =-z*(Rcurv*alpha1)/ecRad+vtxZ; //new z-coordinate of track
1752  double zAbs=fabs(zNew);
1753  if (zAbs<ecDist) {
1754  etaEC = -log(tan(0.5*atan(ecRad/zAbs)));
1755  deltaPhi = -charge*asin(0.5*ecRad/Rcurv);
1756  }
1757  if (zAbs>ecDist) {
1758  zAbs = (fabs(etaIP)/etaIP)*ecDist;
1759  double Zflight = fabs(zAbs-vtxZ);
1760  double alpha = (Zflight*ecRad)/(z*Rcurv);
1761  double Rec = 2*Rcurv*sin(alpha/2);
1762  deltaPhi =-charge*alpha/2;
1763  etaEC =-log(tan(0.5*atan(Rec/ecDist)));
1764  }
1765  }
1766  } else {
1767  zNew = (fabs(etaIP)/etaIP)*ecDist;
1768  double Zflight = fabs(zNew-vtxZ);
1769  double Rvirt = fabs(Zflight*tan(theta));
1770  double Rec = 2*Rcurv*sin(Rvirt/(2*Rcurv));
1771  deltaPhi =-(charge)*(Rvirt/(2*Rcurv));
1772  etaEC =-log(tan(0.5*atan(Rec/ecDist)));
1773  }
1774 
1775  if (zNew<0) etaEC=-etaEC;
1776  phiEC = phiIP+deltaPhi;
1777 
1778  if (phiEC<-M_PI) phiEC += 2*M_PI;
1779  if (phiEC>M_PI) phiEC -= 2*M_PI;
1780 
1781  std::pair<double,double> retVal(etaEC,phiEC);
1782  return retVal;
1783 }
1784 
1785 double IsoTrig::getDistInCM(double eta1,double phi1, double eta2,double phi2) {
1786  double Rec;
1787  double theta1=2*atan(exp(-eta1));
1788  double theta2=2*atan(exp(-eta2));
1789  if (fabs(eta1)<1.479) Rec=rEB_; //radius of ECAL barrel
1790  else if (fabs(eta1)>1.479&&fabs(eta1)<7.0) Rec=tan(theta1)*zEE_; //distance from IP to ECAL endcap
1791  else return 1000;
1792 
1793  //|vect| times tg of acos(scalar product)
1794  double angle=acos((sin(theta1)*sin(theta2)*(sin(phi1)*sin(phi2)+cos(phi1)*cos(phi2))+cos(theta1)*cos(theta2)));
1795  if (angle<0.5*M_PI) return fabs((Rec/sin(theta1))*tan(angle));
1796  else return 1000;
1797 }
1798 
1799 //define this as a plug-in
1800 DEFINE_FWK_MODULE(IsoTrig);
edm::EventID::run
RunNumber_t run() const
Definition: EventID.h:39
HLTConfigProvider::size
unsigned int size() const
number of trigger paths in trigger table
Definition: HLTConfigProvider.h:80
IsoTrig::prelimCone_
double prelimCone_
Definition: IsoTrig.h:133
reco::TrackBase::p
double p() const
momentum vector magnitude
Definition: TrackBase.h:663
edm::ParameterSet::getParameter
T getParameter(std::string const &) const
edm::EventID::event
EventNumber_t event() const
Definition: EventID.h:41
IsoTrig::t_NFTrkPVFlag
std::vector< bool > * t_NFTrkPVFlag
Definition: IsoTrig.h:213
IsoTrig::t_PixTrkcandPhi
std::vector< double > * t_PixTrkcandPhi
Definition: IsoTrig.h:194
IsoTrig::t_PixTrkcandP
std::vector< double > * t_PixTrkcandP
Definition: IsoTrig.h:191
edm::ParameterSet::getUntrackedParameter
T getUntrackedParameter(std::string const &, T const &) const
IsoTrig::tauUnbiasCone_
double tauUnbiasCone_
Definition: IsoTrig.h:133
i
int i
Definition: DBlmapReader.cc:9
IsoTrig::tok_l2cand_
edm::EDGetTokenT< reco::IsolatedPixelTrackCandidateCollection > tok_l2cand_
Definition: IsoTrig.h:155
alpha
float alpha
Definition: AMPTWrapper.h:95
edm::Event::triggerNames
virtual edm::TriggerNames const & triggerNames(edm::TriggerResults const &triggerResults) const
Definition: Event.cc:208
IsoTrig::a_coneR
double a_coneR
Definition: IsoTrig.h:136
IsoTrig::L1candTag_
edm::InputTag L1candTag_
Definition: IsoTrig.h:125
create_public_lumi_plots.exp
tuple exp
Definition: create_public_lumi_plots.py:1087
spr::propagateCALO
std::vector< spr::propagatedTrackID > propagateCALO(edm::Handle< reco::TrackCollection > &trkCollection, const CaloGeometry *geo, const MagneticField *bField, std::string &theTrackQuality, bool debug=false)
Definition: CaloPropagateTrack.cc:109
IsoTrig::t_NFTrkMinDR
std::vector< double > * t_NFTrkMinDR
Definition: IsoTrig.h:208
trigger::TriggerObject::id
int id() const
getters
Definition: TriggerObject.h:55
FastTimerService_cff.lumiProducer
tuple lumiProducer
Definition: FastTimerService_cff.py:11
edm::RunBase::run
RunNumber_t run() const
Definition: RunBase.h:42
trigger::TriggerEvent
The single EDProduct to be saved for each event (AOD case)
Definition: TriggerEvent.h:25
trigger::TriggerEvent::sizeFilters
trigger::size_type sizeFilters() const
Definition: TriggerEvent.h:135
IsoTrig::chgIsolation
void chgIsolation(double &etaTriggered, double &phiTriggered, edm::Handle< reco::TrackCollection > &trkCollection, const edm::Event &theEvent)
Definition: IsoTrig.cc:1416
IsoTrig::endLuminosityBlock
virtual void endLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &)
Definition: IsoTrig.cc:824
IsoTrig::t_NFTrkMissFlag
std::vector< bool > * t_NFTrkMissFlag
Definition: IsoTrig.h:212
IsoTrig::pixelTrackRefsHE
std::vector< reco::TrackRef > pixelTrackRefsHE
Definition: IsoTrig.h:158
IsoTrig::dr_L1
double dr_L1
Definition: IsoTrig.h:136
IsoTrig::tok_lumi
edm::EDGetTokenT< LumiDetails > tok_lumi
Definition: IsoTrig.h:140
IsoTrig::h_phi
TH1D * h_phi[20]
Definition: IsoTrig.h:225
IsoTrig::t_PixTrkcandselTk
std::vector< bool > * t_PixTrkcandselTk
Definition: IsoTrig.h:196
IsoTrig::h_EnIn
TH1D * h_EnIn
Definition: IsoTrig.h:220
IsoTrig::doStudyIsol
bool doStudyIsol
Definition: IsoTrig.h:128
edm::Event::getByToken
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:449
IsoTrig::h_dEta
TH1D * h_dEta[9]
Definition: IsoTrig.h:227
IsoTrig::tok_SiPixelRecHits
edm::EDGetTokenT< SiPixelRecHitCollection > tok_SiPixelRecHits
Definition: IsoTrig.h:152
trigger::TriggerObject::phi
float phi() const
Definition: TriggerObject.h:58
edm::LuminosityBlock::getByToken
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: LuminosityBlock.h:242
reco::TrackBase::TrackQuality
TrackQuality
track quality
Definition: TrackBase.h:133
IsoTrig::studyIsolation
void studyIsolation(edm::Handle< reco::TrackCollection > &, std::vector< reco::TrackCollection::const_iterator > &)
Definition: IsoTrig.cc:1299
MagneticField::inTesla
virtual GlobalVector inTesla(const GlobalPoint &gp) const =0
Field value ad specified global point, in Tesla.
DEFINE_FWK_MODULE
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
IsoTrig::t_NFTrkcandEta
std::vector< double > * t_NFTrkcandEta
Definition: IsoTrig.h:205
IsoTrig::t_NFcandP
std::vector< double > * t_NFcandP
Definition: IsoTrig.h:198
IsoTrig::barrelRecHitsHandle
edm::Handle< EcalRecHitCollection > barrelRecHitsHandle
Definition: IsoTrig.h:162
IsoTrig::h_dEtaL1
TH1D * h_dEtaL1[2]
Definition: IsoTrig.h:226
IsoTrig::doTiming
bool doTiming
Definition: IsoTrig.h:127
IsoTrig::theTrackQuality
std::string theTrackQuality
Definition: IsoTrig.h:135
funct::sin
Sin< T >::type sin(const T &t)
Definition: Sin.h:22
IsoTrig::tok_EB_
edm::EDGetTokenT< EcalRecHitCollection > tok_EB_
Definition: IsoTrig.h:147
trigger::TriggerEvent::filterKeys
const Keys & filterKeys(trigger::size_type index) const
Definition: TriggerEvent.h:111
IsoTrig::g_PreHLT
TH1I * g_PreHLT
Definition: IsoTrig.h:232
GlobalPoint
Global3DPoint GlobalPoint
Definition: GlobalPoint.h:10
IsoTrig::leadPV
math::XYZPoint leadPV
Definition: IsoTrig.h:166
create_public_lumi_plots.log
log
Definition: create_public_lumi_plots.py:1102
edm::SortedCollection< EcalRecHit >::const_iterator
std::vector< EcalRecHit >::const_iterator const_iterator
Definition: SortedCollection.h:82
theta
Geom::Theta< T > theta() const
Definition: Basic3DVectorLD.h:170
IsoTrig::pixelTrackRefsHB
std::vector< reco::TrackRef > pixelTrackRefsHB
Definition: IsoTrig.h:158
edm::EventBase::bunchCrossing
int bunchCrossing() const
Definition: EventBase.h:66
IsoTrig::changed
bool changed
Definition: IsoTrig.h:170
edm::EventBase::luminosityBlock
edm::LuminosityBlockNumber_t luminosityBlock() const
Definition: EventBase.h:63
IsoTrig::vec
std::vector< math::XYZTLorentzVector > vec[3]
Definition: IsoTrig.h:233
IsoTrig::t_nPixSeed
std::vector< int > * t_nPixSeed
Definition: IsoTrig.h:175
IsoTrig::t_TrkhCone
std::vector< double > * t_TrkhCone
Definition: IsoTrig.h:178
IsoTrig::rEB_
double rEB_
Definition: IsoTrig.h:130
IsoTrig::tok_pixtk_
edm::EDGetTokenT< reco::IsolatedPixelTrackCandidateCollection > tok_pixtk_
Definition: IsoTrig.h:153
IsoTrig::t_NFTrkqltyFlag
std::vector< bool > * t_NFTrkqltyFlag
Definition: IsoTrig.h:211
TFileService::make
T * make(const Args &...args) const
make new ROOT object
Definition: TFileService.h:64
IsoTrig::doTrkResTree
bool doTrkResTree
Definition: IsoTrig.h:128
IsoTrig::GetEtaPhiAtEcal
std::pair< double, double > GetEtaPhiAtEcal(double etaIP, double phiIP, double pT, int charge, double vtxZ)
Definition: IsoTrig.cc:1726
IsoTrig::t_PixTrkcandMaxP
std::vector< double > * t_PixTrkcandMaxP
Definition: IsoTrig.h:195
trigger::TriggerObject::energy
float energy() const
Definition: TriggerObject.h:65
spr::chargeIsolationCone
double chargeIsolationCone(unsigned int trkIndex, std::vector< spr::propagatedTrackDirection > &trkDirs, double dR, int &nNearTRKs, bool debug=false)
Definition: ChargeIsolation.cc:235
trigger::TriggerObject::eta
float eta() const
Definition: TriggerObject.h:57
reco::TrackBase::px
double px() const
x coordinate of momentum vector
Definition: TrackBase.h:675
eta
T eta() const
Definition: Basic3DVectorLD.h:177
IsoTrig::t_NFTrkselTkFlag
std::vector< bool > * t_NFTrkselTkFlag
Definition: IsoTrig.h:210
IsoTrig::h_HLT
TH1I * h_HLT
Definition: IsoTrig.h:222
cuy.ii
int ii
Definition: cuy.py:588
spr::find
void find(edm::Handle< EcalRecHitCollection > &hits, DetId thisDet, std::vector< EcalRecHitCollection::const_iterator > &hit, bool debug=false)
Definition: FindCaloHit.cc:7
IsoTrig::dP
double dP(math::XYZTLorentzVector &, math::XYZTLorentzVector &)
Definition: IsoTrig.cc:1702
IsoTrig::t_PixcandP
std::vector< double > * t_PixcandP
Definition: IsoTrig.h:186
IsoTrig::trigNames
std::vector< std::string > trigNames
Definition: IsoTrig.h:124
IsoTrig::etaPhiTrigger
std::pair< double, double > etaPhiTrigger()
Definition: IsoTrig.cc:1710
edm::TriggerNames::triggerNames
Strings const & triggerNames() const
Definition: TriggerNames.cc:24
HLT_25ns14e33_v1_cff.InputTag
tuple InputTag
Definition: HLT_25ns14e33_v1_cff.py:45333
IsoTrig::h_dP
TH1D * h_dP[9]
Definition: IsoTrig.h:227
detailsBasic3DVector::z
float float float z
Definition: extBasic3DVector.h:15
IsoTrig::beginJob
virtual void beginJob()
Definition: IsoTrig.cc:525
AlCaHLTBitMon_QueryRunRegistry.string
string string
Definition: AlCaHLTBitMon_QueryRunRegistry.py:255
IsoTrig::hltConfig_
HLTConfigProvider hltConfig_
Definition: IsoTrig.h:123
IsoTrig::vtxCutIsol_
double vtxCutIsol_
Definition: IsoTrig.h:132
IsoTrig::t_NFTrkMipFlag
std::vector< bool > * t_NFTrkMipFlag
Definition: IsoTrig.h:217
IsoTrig::a_neutIsoR
double a_neutIsoR
Definition: IsoTrig.h:136
IsoTrig::pushMipCutTreeVecs
void pushMipCutTreeVecs(math::XYZTLorentzVector &NFcand, math::XYZTLorentzVector &Trkcand, double &EmipNFcand, double &EmipTrkcand, double &mindR, double &mindP1, std::vector< bool > &Flags, double hCone)
Definition: IsoTrig.cc:497
HLTConfigProvider::triggerIndex
unsigned int triggerIndex(const std::string &triggerName) const
slot position of trigger path in trigger table (0 to size-1)
Definition: HLTConfigProvider.h:100
trigger::TriggerObject
Single trigger physics object (e.g., an isolated muon)
Definition: TriggerObject.h:22
edm::second
U second(std::pair< T, U > const &p)
Definition: ParameterSet.cc:248
IsoTrig::t_NFTrkMinDP1
std::vector< double > * t_NFTrkMinDP1
Definition: IsoTrig.h:209
IsoTrig::h_EnOut
TH1D * h_EnOut
Definition: IsoTrig.h:220
IsoTrig::h_dPhiL1
TH1D * h_dPhiL1[2]
Definition: IsoTrig.h:226
IsoTrig::h_L1ObjEnergy
TH1D * h_L1ObjEnergy
Definition: IsoTrig.h:224
math::XYZTLorentzVector
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
IsoTrig::PixcandTag_
edm::InputTag PixcandTag_
Definition: IsoTrig.h:125
IsoTrig::studyMipCut
void studyMipCut(edm::Handle< reco::TrackCollection > &trkCollection, edm::Handle< reco::IsolatedPixelTrackCandidateCollection > &L2cands)
Definition: IsoTrig.cc:1052
IsoTrig::dinvPt
double dinvPt(math::XYZTLorentzVector &, math::XYZTLorentzVector &)
Definition: IsoTrig.cc:1706
IsoTrig::maxRunNo
int maxRunNo
Definition: IsoTrig.h:138
IsoTrig::bfVal
double bfVal
Definition: IsoTrig.h:130
IsoTrig::verbosity
int verbosity
Definition: IsoTrig.h:129
IsoTrig::h_dPhi
TH1D * h_dPhi[9]
Definition: IsoTrig.h:227
IsoTrig::pG
edm::ESHandle< CaloGeometry > pG
Definition: IsoTrig.h:160
iEvent
int iEvent
Definition: GenABIO.cc:230
PV3DBase::mag
T mag() const
Definition: PV3DBase.h:67
IsoTrig::studyTrigger
void studyTrigger(edm::Handle< reco::TrackCollection > &, std::vector< reco::TrackCollection::const_iterator > &)
Definition: IsoTrig.cc:1167
IsoTrig::fillHist
void fillHist(int, math::XYZTLorentzVector &)
Definition: IsoTrig.cc:1618
edm::ConfigurationDescriptions::addDefault
void addDefault(ParameterSetDescription const &psetDescription)
Definition: ConfigurationDescriptions.cc:99
IsoTrig::t_NFTrkcandP
std::vector< double > * t_NFTrkcandP
Definition: IsoTrig.h:203
IsoTrig::hbhe
edm::Handle< HBHERecHitCollection > hbhe
Definition: IsoTrig.h:161
IsoTrig::t_TrkNuIsolFlag
std::vector< bool > * t_TrkNuIsolFlag
Definition: IsoTrig.h:184
GetRecoTauVFromDQM_MC_cff.kk
tuple kk
Definition: GetRecoTauVFromDQM_MC_cff.py:83
IsoTrig::h_nHLT
TH1I * h_nHLT
Definition: IsoTrig.h:222
IsoTrig::a_mipR
double a_mipR
Definition: IsoTrig.h:136
IsoTrig::tok_bs_
edm::EDGetTokenT< reco::BeamSpot > tok_bs_
Definition: IsoTrig.h:146
spr::goodTrack
bool goodTrack(const reco::Track *pTrack, math::XYZPoint leadPV, trackSelectionParameters parameters, bool debug=false)
Definition: TrackSelection.cc:8
IsoTrig::t_TrkPVFlag
std::vector< bool > * t_TrkPVFlag
Definition: IsoTrig.h:183
trigger::TriggerEvent::getObjects
const TriggerObjectCollection & getObjects() const
Definition: TriggerEvent.h:98
HcalObjRepresent::Fill
void Fill(HcalDetId &id, double val, std::vector< TH2F > &depth)
Definition: HcalObjRepresent.h:611
IsoTrig::t_ECone
std::vector< double > * t_ECone
Definition: IsoTrig.h:218
vec1
std::vector< double > vec1
Definition: HCALResponse.h:15
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:48
hitfit::clear
void clear(CLHEP::HepGenMatrix &m)
Helper function: Reset all elements of a matrix to 0.
Definition: matutil.cc:167
IsoTrig::tok_genTrack_
edm::EDGetTokenT< reco::TrackCollection > tok_genTrack_
Definition: IsoTrig.h:144
IsoTrig::cutNeutral
double cutNeutral
Definition: IsoTrig.h:137
IsoTrig::ChgIsolnTree
TTree * ChgIsolnTree
Definition: IsoTrig.h:172
IsoTrig::h_Pre
TH1I * h_Pre
Definition: IsoTrig.h:223
IsoTrig::h_PreHLT
TH1I * h_PreHLT
Definition: IsoTrig.h:222
IsoTrig::h_dinvPt
TH1D * h_dinvPt[9]
Definition: IsoTrig.h:228
IsoTrig::fillEnergy
void fillEnergy(int, int, double, double, math::XYZTLorentzVector &)
Definition: IsoTrig.cc:1662
FastTimerService::queryModuleTimeByLabel
double queryModuleTimeByLabel(edm::StreamID, const std::string &) const
Definition: FastTimerService.cc:1328
IsoTrig::a_charIsoR
double a_charIsoR
Definition: IsoTrig.h:136
IsoTrig::recVtxs
edm::Handle< reco::VertexCollection > recVtxs
Definition: IsoTrig.h:165
IsoTrig::zEE_
double zEE_
Definition: IsoTrig.h:130
funct::cos
Cos< T >::type cos(const T &t)
Definition: Cos.h:22
IsoTrig::t_NFTrkChgIsoFlag
std::vector< bool > * t_NFTrkChgIsoFlag
Definition: IsoTrig.h:215
HcalDetId::ieta
int ieta() const
get the cell ieta
Definition: HcalDetId.h:36
IsoTrig::tok_pixtks_
std::vector< edm::EDGetTokenT< reco::TrackCollection > > tok_pixtks_
Definition: IsoTrig.h:156
IsoTrig::t_TrkP
std::vector< double > * t_TrkP
Definition: IsoTrig.h:179
IsoTrig::minPTrackValue_
double minPTrackValue_
Definition: IsoTrig.h:132
IsoTrig::h_eta
TH1D * h_eta[20]
Definition: IsoTrig.h:225
funct::tan
Tan< T >::type tan(const T &t)
Definition: Tan.h:22
IsoTrig::cutCharge
double cutCharge
Definition: IsoTrig.h:137
fastmath::etaphi
std::pair< T, T > etaphi(T x, T y, T z)
Definition: FastMath.h:128
IsoTrig::h_Filters
TH1I * h_Filters
Definition: IsoTrig.h:223
funct::abs
Abs< T >::type abs(const T &t)
Definition: Abs.h:22
j
int j
Definition: DBlmapReader.cc:9
IsoTrig::t_NFcandEta
std::vector< double > * t_NFcandEta
Definition: IsoTrig.h:200
IsoTrig::t_timeL2Prod
std::vector< double > * t_timeL2Prod
Definition: IsoTrig.h:173
ecalTB2006H4_GenSimDigiReco_cfg.bField
tuple bField
Definition: ecalTB2006H4_GenSimDigiReco_cfg.py:235
spr::eCone_ecal
double eCone_ecal(const CaloGeometry *geo, edm::Handle< T > &barrelhits, edm::Handle< T > &endcaphits, const GlobalPoint &hpoint1, const GlobalPoint &point1, double dR, const GlobalVector &trackMom, int &nRecHits, double ebThr=-100, double eeThr=-100, double tMin=-500, double tMax=500, bool debug=false)
IsoTrig::analyze
virtual void analyze(const edm::Event &, const edm::EventSetup &)
Definition: IsoTrig.cc:211
IsoTrig::h_eHcal
TH1D * h_eHcal[5][6][48]
Definition: IsoTrig.h:231
edm::Event::getLuminosityBlock
LuminosityBlock const & getLuminosityBlock() const
Definition: Event.h:77
IsoTrig::tok_EE_
edm::EDGetTokenT< EcalRecHitCollection > tok_EE_
Definition: IsoTrig.h:148
IsoTrig::MipCutTree
TTree * MipCutTree
Definition: IsoTrig.h:172
IsoTrig::doL2L3
bool doL2L3
Definition: IsoTrig.h:127
triggerResults
static std::string const triggerResults
Definition: EdmProvDump.cc:40
IsoTrig::IsoTrig
IsoTrig(const edm::ParameterSet &)
Definition: IsoTrig.cc:28
IsoTrig::endJob
virtual void endJob()
Definition: IsoTrig.cc:786
IsoTrig::t_PixcandPt
std::vector< double > * t_PixcandPt
Definition: IsoTrig.h:187
edm::HandleBase::isValid
bool isValid() const
Definition: HandleBase.h:75
IsoTrig::a_neutR2
double a_neutR2
Definition: IsoTrig.h:137
IsoTrig::h_PreL1wt
TH1D * h_PreL1wt
Definition: IsoTrig.h:224
IsoTrig::TrigPreList
std::map< unsigned int, const std::pair< int, int > > TrigPreList
Definition: IsoTrig.h:169
IsoTrig::t_NFcandEmip
std::vector< double > * t_NFcandEmip
Definition: IsoTrig.h:202
IsoTrig::h_eMaxNearP
TH1D * h_eMaxNearP[2]
Definition: IsoTrig.h:229
IsoTrig::dPhi
double dPhi(math::XYZTLorentzVector &, math::XYZTLorentzVector &)
Definition: IsoTrig.cc:1680
IsoTrig::studyTiming
void studyTiming(const edm::Event &theEvent)
Definition: IsoTrig.cc:941
M_PI
#define M_PI
Definition: readTestVector.cc:25
IsoTrig::clearChgIsolnTreeVectors
void clearChgIsolnTreeVectors()
Definition: IsoTrig.cc:440
IsoTrig::t_NFTrkPropFlag
std::vector< bool > * t_NFTrkPropFlag
Definition: IsoTrig.h:214
reco::TrackBase::pz
double pz() const
z coordinate of momentum vector
Definition: TrackBase.h:687
IsoTrig::t_NFTrkcandPhi
std::vector< double > * t_NFTrkcandPhi
Definition: IsoTrig.h:206
HLTConfigProvider::moduleLabels
const std::vector< std::string > & moduleLabels(unsigned int trigger) const
label(s) of module(s) on a trigger path
Definition: HLTConfigProvider.h:105
hit::id
unsigned int id
Definition: SiStripHitEffFromCalibTree.cc:89
IsoTrig::h_MipEnMatch
TH2D * h_MipEnMatch
Definition: IsoTrig.h:221
deltaR
double deltaR(double eta1, double eta2, double phi1, double phi2)
Definition: TreeUtility.cc:17
IsoTrig::h_p
TH1D * h_p[20]
Definition: IsoTrig.h:225
IsoTrig::h_eNeutIso
TH1D * h_eNeutIso[2]
Definition: IsoTrig.h:229
IsoTrig::h_PreL1
TH1I * h_PreL1
Definition: IsoTrig.h:222
IsoTrig::tok_trigRes
edm::EDGetTokenT< edm::TriggerResults > tok_trigRes
Definition: IsoTrig.h:142
trigger::TriggerObjectCollection
std::vector< TriggerObject > TriggerObjectCollection
collection of trigger physics objects (e.g., all isolated muons)
Definition: TriggerObject.h:81
trigger::TriggerEvent::filterTag
const edm::InputTag filterTag(trigger::size_type index) const
Definition: TriggerEvent.h:103
IsoTrig::minRunNo
int minRunNo
Definition: IsoTrig.h:138
IsoTrig::L2candTag_
edm::InputTag L2candTag_
Definition: IsoTrig.h:125
IsoTrig::t_TrkMissFlag
std::vector< bool > * t_TrkMissFlag
Definition: IsoTrig.h:182
IsoTrig::getDistInCM
double getDistInCM(double eta1, double phi1, double eta2, double phi2)
Definition: IsoTrig.cc:1785
IsoTrig::tok_recVtx_
edm::EDGetTokenT< reco::VertexCollection > tok_recVtx_
Definition: IsoTrig.h:145
IsoTrig::t_nPixCand
std::vector< int > * t_nPixCand
Definition: IsoTrig.h:174
reco::TrackBase::qualityByName
static TrackQuality qualityByName(const std::string &name)
Definition: TrackBase.cc:108
SiPixelRawToDigiRegional_cfi.deltaPhi
tuple deltaPhi
Definition: SiPixelRawToDigiRegional_cfi.py:9
debug
#define debug
Definition: HDRShower.cc:19
IsoTrig::StudyTrkEbyP
void StudyTrkEbyP(edm::Handle< reco::TrackCollection > &trkCollection)
Definition: IsoTrig.cc:827
IsoTrig::t_NFTrkcandPt
std::vector< double > * t_NFTrkcandPt
Definition: IsoTrig.h:204
IsoTrig::h_eCalo
TH1D * h_eCalo[5][6][48]
Definition: IsoTrig.h:231
edm::Handle::product
T const * product() const
Definition: Handle.h:81
IsoTrig::g_Pre
TH1I * g_Pre
Definition: IsoTrig.h:232
spr::eCone_hcal
double eCone_hcal(const CaloGeometry *geo, edm::Handle< T > &hits, const GlobalPoint &hpoint1, const GlobalPoint &point1, double dR, const GlobalVector &trackMom, int &nRecHits, double hbThr=-100, double heThr=-100, double hfThr=-100, double hoThr=-100, double tMin=-500, double tMax=500, int detOnly=-1, bool debug=false)
IsoTrig::endRun
virtual void endRun(edm::Run const &, edm::EventSetup const &)
Definition: IsoTrig.cc:819
math::XYZPoint
XYZPointD XYZPoint
point in space with cartesian internal representation
Definition: Point3D.h:12
trigger::Keys
std::vector< size_type > Keys
Definition: TriggerTypeDefs.h:20
IsoTrig::beamSpotH
edm::Handle< reco::BeamSpot > beamSpotH
Definition: IsoTrig.h:164
IsoTrig::t_PixTrkcandPt
std::vector< double > * t_PixTrkcandPt
Definition: IsoTrig.h:192
IsoTrig::dEta
double dEta(math::XYZTLorentzVector &, math::XYZTLorentzVector &)
Definition: IsoTrig.cc:1676
HLTConfigProvider::init
bool init(const edm::Run &iRun, const edm::EventSetup &iSetup, const std::string &processName, bool &changed)
d&#39;tor
Definition: HLTConfigProvider.cc:42
edm::EventSetup::get
const T & get() const
Definition: EventSetup.h:55
edm::ESHandle::product
T const * product() const
Definition: ESHandle.h:86
IsoTrig::tok_verthb_
edm::EDGetTokenT< reco::VertexCollection > tok_verthb_
Definition: IsoTrig.h:150
HLTConfigProvider::prescaleValues
std::pair< int, int > prescaleValues(const edm::Event &iEvent, const edm::EventSetup &iSetup, const std::string &trigger) const
Combined L1T (pair.first) and HLT (pair.second) prescales per HLT path.
Definition: HLTConfigProvider.cc:282
IsoTrig::doChgIsolTree
bool doChgIsolTree
Definition: IsoTrig.h:128
IsoTrig::t_NFcandPt
std::vector< double > * t_NFcandPt
Definition: IsoTrig.h:199
IsoTrig::fillDescriptions
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
Definition: IsoTrig.cc:203
IsoTrig::selectionParameters
spr::trackSelectionParameters selectionParameters
Definition: IsoTrig.h:134
IsoTrig::h_pt
TH1D * h_pt[20]
Definition: IsoTrig.h:225
IsoTrig::TrigList
std::map< unsigned int, unsigned int > TrigList
Definition: IsoTrig.h:168
IsoTrig::t_TrkqltyFlag
std::vector< bool > * t_TrkqltyFlag
Definition: IsoTrig.h:181
IsoTrig::h_etaMipTracks
TH1D * h_etaMipTracks[5][2][2]
Definition: IsoTrig.h:230
edm::InputTag::label
std::string const & label() const
Definition: InputTag.h:42
AlCaHLTBitMon_ParallelJobs.p
tuple p
Definition: AlCaHLTBitMon_ParallelJobs.py:152
IsoTrig::t_PixcandMaxP
std::vector< std::vector< double > > * t_PixcandMaxP
Definition: IsoTrig.h:190
IsoTrig::t_NFTrkcandEmip
std::vector< double > * t_NFTrkcandEmip
Definition: IsoTrig.h:207
IsoTrig::tok_l1cand_
edm::EDGetTokenT< trigger::TriggerFilterObjectWithRefs > tok_l1cand_
Definition: IsoTrig.h:154
IsoTrig::pixelTracksSources_
std::vector< edm::InputTag > pixelTracksSources_
Definition: IsoTrig.h:126
edmNew::DetSetVector::size
size_type size() const
Definition: DetSetVectorNew.h:363
IsoTrig::tok_SeedingLayerhb
edm::EDGetTokenT< SeedingLayerSetsHits > tok_SeedingLayerhb
Definition: IsoTrig.h:151
edm::EventBase::id
edm::EventID id() const
Definition: EventBase.h:60
HI_PhotonSkim_cff.rechits
tuple rechits
Definition: HI_PhotonSkim_cff.py:76
IsoTrig::dR
double dR(math::XYZTLorentzVector &, math::XYZTLorentzVector &)
Definition: IsoTrig.cc:1692
IsoTrig::h_mindR
TH1D * h_mindR[9]
Definition: IsoTrig.h:228
IsoTrig::a_neutR1
double a_neutR1
Definition: IsoTrig.h:137
Reconstruction_Data_cff.cuts
tuple cuts
Definition: Reconstruction_Data_cff.py:22
IsoTrig::clearMipCutTreeVectors
void clearMipCutTreeVectors()
Definition: IsoTrig.cc:456
IsoTrig::t_PixcandEta
std::vector< double > * t_PixcandEta
Definition: IsoTrig.h:188
IsoTrig::tok_hbhe_
edm::EDGetTokenT< HBHERecHitCollection > tok_hbhe_
Definition: IsoTrig.h:149
position
static int position[264][3]
Definition: ReadPGInfo.cc:509
spr::trackSelectionParameters::minQuality
reco::TrackBase::TrackQuality minQuality
Definition: TrackSelection.h:22
IsoTrig::vtxCutSeed_
double vtxCutSeed_
Definition: IsoTrig.h:132
edm::Event::streamID
StreamID streamID() const
Definition: Event.h:72
IsoTrig::tok_SeedingLayerhe
edm::EDGetTokenT< SeedingLayerSetsHits > tok_SeedingLayerhe
Definition: IsoTrig.h:151
IsoTrig::fillDifferences
void fillDifferences(int, math::XYZTLorentzVector &, math::XYZTLorentzVector &, bool)
Definition: IsoTrig.cc:1625
IsoTrig::h_dPt
TH1D * h_dPt[9]
Definition: IsoTrig.h:227
IsoTrig::h_MipEnNoMatch
TH2D * h_MipEnNoMatch
Definition: IsoTrig.h:221
IsoTrig::dPt
double dPt(math::XYZTLorentzVector &, math::XYZTLorentzVector &)
Definition: IsoTrig.cc:1698
IsoTrig::t_NFTrkNeuIsoFlag
std::vector< bool > * t_NFTrkNeuIsoFlag
Definition: IsoTrig.h:216
IsoTrig::endcapRecHitsHandle
edm::Handle< EcalRecHitCollection > endcapRecHitsHandle
Definition: IsoTrig.h:163
IsoTrig::TrkResTree
TTree * TrkResTree
Definition: IsoTrig.h:172
gather_cfg.cout
tuple cout
Definition: gather_cfg.py:121
IsoTrig::t_PixcandPhi
std::vector< double > * t_PixcandPhi
Definition: IsoTrig.h:189
IsoTrig::t_PixTrkcandEta
std::vector< double > * t_PixTrkcandEta
Definition: IsoTrig.h:193
IsoTrig::doMipCutTree
bool doMipCutTree
Definition: IsoTrig.h:127
IsoTrig::fs
edm::Service< TFileService > fs
Definition: IsoTrig.h:171
IsoTrig::t_nGoodTk
std::vector< int > * t_nGoodTk
Definition: IsoTrig.h:176
IsoTrig::TimingTree
TTree * TimingTree
Definition: IsoTrig.h:172
IsoTrig::pLimits
double pLimits[6]
Definition: IsoTrig.h:139
edm::false
volatile std::atomic< bool > shutdown_flag false
Definition: UnixSignalHandlers.cc:22
IsoTrig::g_PreL1
TH1I * g_PreL1
Definition: IsoTrig.h:232
x
Definition: DDAxes.h:10
IsoTrig::bFieldH
edm::ESHandle< MagneticField > bFieldH
Definition: IsoTrig.h:159
SQLiteEnsembleGenerator_cfg.BField
float BField
Definition: SQLiteEnsembleGenerator_cfg.py:49
trigger::TriggerObject::mass
float mass() const
Definition: TriggerObject.h:59
IsoTrig::tok_verthe_
edm::EDGetTokenT< reco::VertexCollection > tok_verthe_
Definition: IsoTrig.h:150
relval_steps.k
list k
Definition: relval_steps.py:1469
filesave_online.pairs
list pairs
sort tag files by run number
Definition: filesave_online.py:88
IsoTrig::h_nL3Objs
TH1I * h_nL3Objs
Definition: IsoTrig.h:223
IsoTrig::tok_hlt_
edm::EDGetTokenT< trigger::TriggerFilterObjectWithRefs > tok_hlt_
Definition: IsoTrig.h:143
IsoTrig::~IsoTrig
~IsoTrig()
Definition: IsoTrig.cc:155
IsoTrig::getGoodTracks
void getGoodTracks(const edm::Event &, edm::Handle< reco::TrackCollection > &)
Definition: IsoTrig.cc:1527
vec2
std::vector< vec1 > vec2
Definition: HCALResponse.h:16
IsoTrig::h_dRL1
TH1D * h_dRL1[2]
Definition: IsoTrig.h:226
IsoTrig::h_etaCalibTracks
TH1D * h_etaCalibTracks[5][2][2]
Definition: IsoTrig.h:230
SeedingLayerSetsHits::size
unsigned short size() const
Get the number of SeedingLayerSets.
Definition: SeedingLayerSetsHits.h:182
IsoTrig::t_TrkselTkFlag
std::vector< bool > * t_TrkselTkFlag
Definition: IsoTrig.h:180
IsoTrig::cutMip
double cutMip
Definition: IsoTrig.h:137
IsoTrig::processName
std::string processName
Definition: IsoTrig.h:135
findQualityFiles.size
tuple size
Write out results.
Definition: findQualityFiles.py:442
IsoTrig::fillCuts
void fillCuts(int, double, double, double, math::XYZTLorentzVector &, int, bool)
Definition: IsoTrig.cc:1644
reco::TrackBase::py
double py() const
y coordinate of momentum vector
Definition: TrackBase.h:681
edm::Run
Definition: Run.h:41
IsoTrig::beginLuminosityBlock
virtual void beginLuminosityBlock(edm::LuminosityBlock const &, edm::EventSetup const &)
Definition: IsoTrig.cc:822
GOODCOLL_filter_cfg.cut
tuple cut
Definition: GOODCOLL_filter_cfg.py:65
IsoTrig::h_PreHLTwt
TH1D * h_PreHLTwt
Definition: IsoTrig.h:224
IsoTrig::tok_trigEvt
edm::EDGetTokenT< trigger::TriggerEvent > tok_trigEvt
Definition: IsoTrig.h:141
IsoTrig::pixelIsolationConeSizeAtEC_
std::vector< double > pixelIsolationConeSizeAtEC_
Definition: IsoTrig.h:131
IsoTrig::beginRun
virtual void beginRun(edm::Run const &, edm::EventSetup const &)
Definition: IsoTrig.cc:813
IsoTrig::pushChgIsolnTreeVecs
void pushChgIsolnTreeVecs(math::XYZTLorentzVector &Pixcand, math::XYZTLorentzVector &Trkcand, std::vector< double > &PixMaxP, double &TrkMaxP, bool &selTk)
Definition: IsoTrig.cc:480
IsoTrig::g_Accepts
TH1I * g_Accepts
Definition: IsoTrig.h:232
angle
T angle(T x1, T y1, T z1, T x2, T y2, T z2)
Definition: angle.h:11
IsoTrig::h_eMip
TH1D * h_eMip[2]
Definition: IsoTrig.h:228
phi
Definition: DDAxes.h:10
IsoTrig::t_NFcandPhi
std::vector< double > * t_NFcandPhi
Definition: IsoTrig.h:201
Generated for CMSSW Reference Manual by doxygen 1.8.5