File indexing completed on 2023-03-17 11:25:51
0001
0002
0003
0004
0005
0006
0007
0008
0009 #include "DataFormats/Common/interface/Association.h"
0010 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0011 #include "FWCore/Framework/interface/stream/EDProducer.h"
0012 #include "SimTracker/TrackHistory/interface/TrackHistory.h"
0013
0014 namespace edm {
0015 class ParameterSet;
0016 }
0017
0018 using namespace edm;
0019 using namespace std;
0020 using namespace reco;
0021
0022 class GenTrackMatcher : public edm::stream::EDProducer<> {
0023 public:
0024
0025 GenTrackMatcher(const edm::ParameterSet &);
0026
0027 private:
0028 void produce(edm::Event &evt, const edm::EventSetup &es) override;
0029 TrackHistory tracer_;
0030 edm::EDGetTokenT<View<Track>> tracks_;
0031 edm::EDGetTokenT<GenParticleCollection> genParticles_;
0032 edm::EDGetTokenT<vector<int>> genParticleInts_;
0033 typedef edm::Association<reco::GenParticleCollection> GenParticleMatch;
0034 };
0035
0036 #include "DataFormats/Common/interface/Handle.h"
0037 #include "FWCore/Framework/interface/ESHandle.h"
0038 #include "FWCore/Framework/interface/Event.h"
0039 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0040 #include "FWCore/Utilities/interface/EDMException.h"
0041
0042 GenTrackMatcher::GenTrackMatcher(const ParameterSet &p)
0043 : tracer_(p, consumesCollector()),
0044 tracks_(consumes<View<Track>>(p.getUntrackedParameter<edm::InputTag>("trackProducer"))),
0045 genParticles_(consumes<GenParticleCollection>(p.getUntrackedParameter<edm::InputTag>("genParticles"))),
0046 genParticleInts_(consumes<vector<int>>(p.getUntrackedParameter<edm::InputTag>("genParticles"))) {
0047 produces<GenParticleMatch>();
0048 }
0049
0050 void GenTrackMatcher::produce(Event &evt, const EventSetup &es) {
0051 Handle<View<Track>> tracks;
0052 evt.getByToken(tracks_, tracks);
0053 Handle<vector<int>> barCodes;
0054 evt.getByToken(genParticles_, barCodes);
0055 Handle<GenParticleCollection> genParticles;
0056 evt.getByToken(genParticles_, genParticles);
0057 unique_ptr<GenParticleMatch> match(new GenParticleMatch(GenParticleRefProd(genParticles)));
0058 GenParticleMatch::Filler filler(*match);
0059 size_t n = tracks->size();
0060 vector<int> indices(n, -1);
0061 tracer_.newEvent(evt, es);
0062 for (size_t i = 0; i < n; ++i) {
0063 RefToBase<Track> track(tracks, i);
0064 if (tracer_.evaluate(track)) {
0065 const HepMC::GenParticle *particle = tracer_.genParticle();
0066 if (particle) {
0067 int barCode = particle->barcode();
0068 vector<int>::const_iterator b = barCodes->begin(), e = barCodes->end(), f = find(b, e, barCode);
0069 if (f == e) {
0070 edm::EDConsumerBase::Labels labels;
0071 labelsForToken(genParticles_, labels);
0072 throw edm::Exception(errors::InvalidReference)
0073 << "found matching particle with barcode" << *f << " which has not been found in " << labels.module;
0074 }
0075 indices[i] = *f;
0076 }
0077 }
0078 }
0079 filler.insert(tracks, indices.begin(), indices.end());
0080 filler.fill();
0081 evt.put(std::move(match));
0082 }
0083
0084 #include "FWCore/Framework/interface/MakerMacros.h"
0085
0086 DEFINE_FWK_MODULE(GenTrackMatcher);