CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
GsfElectronLinker.cc
Go to the documentation of this file.
2 
7 
8 
11  = iConfig.getParameter<edm::InputTag>("PFCandidate");
13  = iConfig.getParameter<edm::InputTag>("GsfElectrons");
15  = iConfig.getParameter<std::string>("OutputPF");
16 
17  produces<reco::PFCandidateCollection>(nameOutputPF_);
18  produces<edm::ValueMap<reco::PFCandidateRef> > (nameOutputPF_);
19 }
20 
22 
24 
26  electronCandidateMap_.clear();
27  std::auto_ptr<reco::PFCandidateCollection>
28  pfCandidates_p(new reco::PFCandidateCollection);
29 
30  std::auto_ptr<edm::ValueMap<reco::PFCandidateRef> >
32  edm::ValueMap<reco::PFCandidateRef>::Filler pfMapFiller(*pfMap_p);
33 
35  bool status=fetchCandidateCollection(pfCandidates,
37  iEvent );
38 
40  status=fetchGsfElectronCollection(gsfElectrons,
42  iEvent );
43 
44  unsigned ncand=(status)?pfCandidates->size():0;
45 
46  for( unsigned i=0; i<ncand; ++i ) {
47  edm::Ptr<reco::PFCandidate> candPtr(pfCandidates,i);
48  reco::PFCandidate cand = (*pfCandidates)[i];
49  cand.setSourceCandidatePtr(candPtr);
50 
51  // if not an electron or not GsfTrackRef
52  if( (cand.particleId()!=reco::PFCandidate::e) ) {
53  pfCandidates_p->push_back(cand);
54  continue; // Watch out ! Continue
55  }
56 
57  // if it is an electron. Find the GsfElectron with the same GsfTrack
58  const reco::GsfTrackRef & gsfTrackRef(cand.gsfTrackRef());
59  GsfElectronEqual myEqual(gsfTrackRef);
60  std::vector<reco::GsfElectron>::const_iterator itcheck=find_if(gsfElectrons->begin(),gsfElectrons->end(),myEqual);
61  if(itcheck==gsfElectrons->end()) {
62  std::ostringstream err;
63  err << " Problem in GsfElectronLinker: no GsfElectron " << std::endl;
64  edm::LogError("GsfElectronLinker") << err.str();
65  continue; // Watch out ! Continue
66  }
67  reco::GsfElectronRef electronRef(gsfElectrons,itcheck-gsfElectrons->begin());
68  cand.setGsfElectronRef(electronRef);
69  electronCandidateMap_[electronRef]=i;
70  pfCandidates_p->push_back(cand);
71  }
72 
73  const edm::OrphanHandle<reco::PFCandidateCollection> pfCandidateRefProd =
74  iEvent.put(pfCandidates_p,nameOutputPF_);
75 
76  fillValueMap(gsfElectrons,pfCandidateRefProd,pfMapFiller);
77  iEvent.put(pfMap_p,nameOutputPF_);
78 }
79 
80 
82  const edm::InputTag& tag,
83  const edm::Event& iEvent) const {
84  bool found = iEvent.getByLabel(tag, c);
85 
86  if(!found )
87  {
88  std::ostringstream err;
89  err<<" cannot get PFCandidates: "
90  <<tag<<std::endl;
91  edm::LogError("GsfElectronLinker")<<err.str();
92  }
93  return found;
94 
95 }
96 
98  const edm::InputTag& tag,
99  const edm::Event& iEvent) const {
100  bool found = iEvent.getByLabel(tag, c);
101 
102  if(!found )
103  {
104  std::ostringstream err;
105  err<<" cannot get GsfElectrons: "
106  <<tag<<std::endl;
107  edm::LogError("GsfElectronLinker")<<err.str();
108  }
109  return found;
110 
111 }
112 
116  unsigned nElectrons=electrons->size();
117  std::vector<reco::PFCandidateRef> values;
118  for(unsigned ielec=0;ielec<nElectrons;++ielec) {
119  reco::GsfElectronRef gsfElecRef(electrons,ielec);
120  std::map<reco::GsfElectronRef,unsigned>::const_iterator itcheck=electronCandidateMap_.find(gsfElecRef);
121  if(itcheck==electronCandidateMap_.end()) {
122  values.push_back(reco::PFCandidateRef());
123  } else {
124  values.push_back(reco::PFCandidateRef(pfHandle,itcheck->second));
125  }
126  }
127  filler.insert(electrons,values.begin(),values.end());
128 }
edm::InputTag inputTagPFCandidates_
Input PFCandidates.
T getParameter(std::string const &) const
int i
Definition: DBlmapReader.cc:9
void setGsfElectronRef(const reco::GsfElectronRef &ref)
set GsfElectronRef
Definition: PFCandidate.cc:466
void insert(const H &h, I begin, I end)
Definition: ValueMap.h:53
void setSourceCandidatePtr(const PFCandidatePtr &ptr)
Definition: PFCandidate.h:111
std::vector< PFCandidatePtr > pfCandidates(const PFJet &jet, int particleId, bool sort=true)
void fillValueMap(edm::Handle< reco::GsfElectronCollection > &c, const edm::OrphanHandle< reco::PFCandidateCollection > &pfHandle, edm::ValueMap< reco::PFCandidateRef >::Filler &filler) const
int iEvent
Definition: GenABIO.cc:243
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:84
std::map< reco::GsfElectronRef, unsigned > electronCandidateMap_
map GsfElectron PFCandidate (index)
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:355
bool fetchGsfElectronCollection(edm::Handle< reco::GsfElectronCollection > &c, const edm::InputTag &tag, const edm::Event &iEvent) const
std::vector< reco::PFCandidate > PFCandidateCollection
collection of PFCandidates
virtual void produce(edm::Event &, const edm::EventSetup &)
Particle reconstructed by the particle flow algorithm.
Definition: PFCandidate.h:33
virtual void beginRun(edm::Run &run, const edm::EventSetup &es)
GsfElectronLinker(const edm::ParameterSet &)
reco::GsfTrackRef gsfTrackRef() const
Definition: PFCandidate.cc:368
bool fetchCandidateCollection(edm::Handle< reco::PFCandidateCollection > &c, const edm::InputTag &tag, const edm::Event &iEvent) const
tuple status
Definition: ntuplemaker.py:245
virtual ParticleType particleId() const
Definition: PFCandidate.h:314
edm::InputTag inputTagGsfElectrons_
Input GsfElectrons.
std::string nameOutputPF_
name of output collection of PFCandidate
Definition: Run.h:32