CMS 3D CMS Logo

ElectronMatchedCandidateProducer.cc
Go to the documentation of this file.
7 
8 #include <memory>
9 
11 public:
13 
14 private:
15  void produce(edm::StreamID, edm::Event &, const edm::EventSetup &) const override;
16 
17  // ----------member data ---------------------------
18 
22 };
23 
25  const edm::InputTag allelectrons("gsfElectrons");
26  electronCollectionToken_ = consumes<edm::View<reco::GsfElectron>>(
27  params.getUntrackedParameter<edm::InputTag>("ReferenceElectronCollection", allelectrons));
28  scCollectionToken_ = consumes<edm::View<reco::Candidate>>(params.getParameter<edm::InputTag>("src"));
29 
30  delRMatchingCut_ = params.getUntrackedParameter<double>("deltaR", 0.30);
31 
32  produces<edm::PtrVector<reco::Candidate>>();
33  produces<edm::RefToBaseVector<reco::Candidate>>();
34 }
35 
36 //
37 // member functions
38 //
39 
40 // ------------ method called to produce the data ------------
41 
43  // Create the output collection
44  auto outColRef = std::make_unique<edm::RefToBaseVector<reco::Candidate>>();
45  auto outColPtr = std::make_unique<edm::PtrVector<reco::Candidate>>();
46 
47  // Read electrons
49  event.getByToken(electronCollectionToken_, electrons);
50 
51  //Read candidates
53  event.getByToken(scCollectionToken_, recoCandColl);
54 
55  unsigned int counter = 0;
56 
57  // Loop over candidates
58  for (edm::View<reco::Candidate>::const_iterator scIt = recoCandColl->begin(); scIt != recoCandColl->end();
59  ++scIt, ++counter) {
60  // Now loop over electrons
61  for (edm::View<reco::GsfElectron>::const_iterator elec = electrons->begin(); elec != electrons->end(); ++elec) {
62  reco::SuperClusterRef eSC = elec->superCluster();
63 
64  double dRval = reco::deltaR((float)eSC->eta(), (float)eSC->phi(), scIt->eta(), scIt->phi());
65 
66  if (dRval < delRMatchingCut_) {
67  //outCol->push_back( *scIt );
68  outColRef->push_back(recoCandColl->refAt(counter));
69  outColPtr->push_back(recoCandColl->ptrAt(counter));
70  } // end if loop
71  } // end electron loop
72 
73  } // end candidate loop
74 
75  event.put(std::move(outColRef));
76  event.put(std::move(outColPtr));
77 }
78 
79 //define this as a plug-in
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:16
void produce(edm::StreamID, edm::Event &, const edm::EventSetup &) const override
constexpr auto deltaR(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:30
edm::EDGetTokenT< edm::View< reco::Candidate > > scCollectionToken_
ElectronMatchedCandidateProducer(const edm::ParameterSet &)
boost::indirect_iterator< typename seq_t::const_iterator > const_iterator
Definition: View.h:86
static std::atomic< unsigned int > counter
edm::EDGetTokenT< edm::View< reco::GsfElectron > > electronCollectionToken_
def move(src, dest)
Definition: eostools.py:511
Definition: event.py:1