CMS 3D CMS Logo

All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Pages
MCTrackMatcher.cc
Go to the documentation of this file.
1 
12 
13 namespace edm { class ParameterSet; }
14 
16  public:
19 
20  private:
21  void produce( edm::Event& evt, const edm::EventSetup& es );
25 };
26 
35 using namespace edm;
36 using namespace std;
37 using namespace reco;
38 
40  associator_(p.getParameter<string>("associator")),
41  tracks_(p.getParameter<InputTag>("tracks")),
42  genParticles_( p.getParameter<InputTag>("genParticles")),
43  trackingParticles_( p.getParameter<InputTag>("trackingParticles")) {
44  produces<GenParticleMatch>();
45 }
46 
47 void MCTrackMatcher::produce(Event& evt, const EventSetup& es) {
49  es.get<TrackAssociatorRecord>().get(associator_,assoc);
50  const TrackAssociatorBase * associator = assoc.product();
52  evt.getByLabel(tracks_, tracks);
54  evt.getByLabel(trackingParticles_,trackingParticles);
55  Handle<vector<int> > barCodes;
56  evt.getByLabel(genParticles_,barCodes );
58  evt.getByLabel(genParticles_, genParticles );
59  RecoToSimCollection associations = associator->associateRecoToSim ( tracks, trackingParticles, & evt, &es );
60  auto_ptr<GenParticleMatch> match(new GenParticleMatch(GenParticleRefProd(genParticles)));
61  GenParticleMatch::Filler filler(*match);
62  size_t n = tracks->size();
63  vector<int> indices(n,-1);
64  for (size_t i = 0; i < n; ++ i ) {
65  RefToBase<Track> track(tracks, i);
66  RecoToSimCollection::const_iterator f = associations.find(track);
67  if ( f != associations.end() ) {
68  TrackingParticleRef tp = f->val.front().first;
69  TrackingParticle::genp_iterator j, b = tp->genParticle_begin(), e = tp->genParticle_end();
70  for( j = b; j != e; ++ j ) {
71  const reco::GenParticle * p = j->get();
72  if (p->status() == 1) {
73  indices[i] = j->key();
74  break;
75  }
76  }
77  }
78  }
79  filler.insert(tracks, indices.begin(), indices.end());
80  filler.fill();
81  evt.put(match);
82 }
83 
85 
87 
edm::Association< reco::GenParticleCollection > GenParticleMatch
int i
Definition: DBlmapReader.cc:9
const_iterator end() const
last iterator over the map (read only)
#define DEFINE_FWK_MODULE(type)
Definition: MakerMacros.h:17
const_iterator find(const key_type &k) const
find element with specified reference key
void insert(const H &h, I begin, I end)
Definition: ValueMap.h:53
edm::RefProd< GenParticleCollection > GenParticleRefProd
persistent reference to a GenParticle collection
void produce(edm::Event &evt, const edm::EventSetup &es)
std::string associator_
virtual int status() const GCC11_FINAL
status word
OrphanHandle< PROD > put(std::auto_ptr< PROD > product)
Put a new product.
Definition: Event.h:94
edm::InputTag genParticles_
int j
Definition: DBlmapReader.cc:9
double f[11][100]
bool getByLabel(InputTag const &tag, Handle< PROD > &result) const
Definition: Event.h:361
virtual reco::RecoToSimCollection associateRecoToSim(edm::Handle< edm::View< reco::Track > > &tCH, edm::Handle< TrackingParticleCollection > &tPCH, const edm::Event *event, const edm::EventSetup *setup) const
compare reco to sim the handle of reco::Track and TrackingParticle collections
tuple tracks
Definition: testEve_cfg.py:39
edm::InputTag trackingParticles_
const T & get() const
Definition: EventSetup.h:55
T const * product() const
Definition: ESHandle.h:62
double b
Definition: hdecay.h:120
MCTrackMatcher(const edm::ParameterSet &)
constructor
edm::InputTag tracks_
std::pair< typename Association::data_type::first_type, double > match(Reference key, Association association, bool bestMatchByMaxValue)
Generic matching function.
Definition: Utils.h:6
const reco::PFCandidateRefVector & tracks_