CMS 3D CMS Logo

IPTCorrector.cc
Go to the documentation of this file.
1 #include <vector>
2 #include <memory>
3 #include <algorithm>
4 
5 // Class header file
8 // Framework
14 //
17 // Math
18 #include "Math/GenVector/VectorUtil.h"
19 #include "Math/GenVector/PxPyPzE4D.h"
21 
23  : tok_cor_(consumes<reco::TrackCollection>(config.getParameter<edm::InputTag>("corTracksLabel"))),
24  tok_uncor_(consumes<trigger::TriggerFilterObjectWithRefs>(config.getParameter<edm::InputTag>("filterLabel"))),
25  assocCone_(config.getParameter<double>("associationCone")) {
26  // register the product
27  produces<reco::IsolatedPixelTrackCandidateCollection>();
28 }
29 
32  desc.add<edm::InputTag>("corTracksLabel", edm::InputTag("hltIter0PFlowCtfWithMaterialTracks"));
33  desc.add<edm::InputTag>("filterLabel", edm::InputTag("hltIsolPixelTrackL2Filter"));
34  desc.add<double>("associationCone", 0.2);
35  descriptions.add("iptCorrector", desc);
36 }
37 
39  auto trackCollection = std::make_unique<reco::IsolatedPixelTrackCandidateCollection>();
40 
42  theEvent.getByToken(tok_cor_, corTracks);
43 
45  theEvent.getByToken(tok_uncor_, fiCand);
46 
47  std::vector<edm::Ref<reco::IsolatedPixelTrackCandidateCollection> > isoPixTrackRefs;
48 
49  fiCand->getObjects(trigger::TriggerTrack, isoPixTrackRefs);
50 
51  int nCand = isoPixTrackRefs.size();
52 
53  //loop over input ipt
54 
55  for (int p = 0; p < nCand; p++) {
56  double iptEta = isoPixTrackRefs[p]->track()->eta();
57  double iptPhi = isoPixTrackRefs[p]->track()->phi();
58 
59  int ntrk = 0;
60  double minDR = 100;
61  reco::TrackCollection::const_iterator citSel;
62 
63  for (reco::TrackCollection::const_iterator cit = corTracks->begin(); cit != corTracks->end(); cit++) {
64  double dR = deltaR(cit->eta(), cit->phi(), iptEta, iptPhi);
65  if (dR < minDR && dR < assocCone_) {
66  minDR = dR;
67  ntrk++;
68  citSel = cit;
69  }
70  }
71 
72  if (ntrk > 0) {
73  reco::IsolatedPixelTrackCandidate newCandidate(reco::TrackRef(corTracks, citSel - corTracks->begin()),
74  isoPixTrackRefs[p]->l1tau(),
75  isoPixTrackRefs[p]->maxPtPxl(),
76  isoPixTrackRefs[p]->sumPtPxl());
77  newCandidate.setEnergyIn(isoPixTrackRefs[p]->energyIn());
78  newCandidate.setEnergyOut(isoPixTrackRefs[p]->energyOut());
79  newCandidate.setNHitIn(isoPixTrackRefs[p]->nHitIn());
80  newCandidate.setNHitOut(isoPixTrackRefs[p]->nHitOut());
81  if (isoPixTrackRefs[p]->etaPhiEcalValid()) {
82  std::pair<double, double> etaphi = (isoPixTrackRefs[p]->etaPhiEcal());
83  newCandidate.setEtaPhiEcal(etaphi.first, etaphi.second);
84  }
85  trackCollection->push_back(newCandidate);
86  }
87  }
88 
89  // put the product in the event
90  theEvent.put(std::move(trackCollection));
91 }
const edm::EDGetTokenT< reco::TrackCollection > tok_cor_
Definition: IPTCorrector.h:29
OrphanHandle< PROD > put(std::unique_ptr< PROD > product)
Put a new product.
Definition: Event.h:125
void getObjects(Vids &ids, VRphoton &refs) const
various physics-level getters:
bool getByToken(EDGetToken token, Handle< PROD > &result) const
Definition: Event.h:517
const edm::EDGetTokenT< trigger::TriggerFilterObjectWithRefs > tok_uncor_
Definition: IPTCorrector.h:30
std::vector< Track > TrackCollection
collection of Tracks
Definition: TrackFwd.h:15
Definition: config.py:1
std::pair< T, T > etaphi(T x, T y, T z)
Definition: FastMath.h:128
IPTCorrector(const edm::ParameterSet &ps)
Definition: IPTCorrector.cc:22
ParameterDescriptionBase * add(U const &iLabel, T const &value)
static void fillDescriptions(edm::ConfigurationDescriptions &descriptions)
Definition: IPTCorrector.cc:30
void add(std::string const &label, ParameterSetDescription const &psetDescription)
fixed size matrix
HLT enums.
const double assocCone_
Definition: IPTCorrector.h:31
def move(src, dest)
Definition: eostools.py:511
void produce(edm::StreamID, edm::Event &, edm::EventSetup const &) const override
Definition: IPTCorrector.cc:38