File indexing completed on 2025-01-08 03:36:27
0001
0002
0003
0004
0005
0006
0007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0008
0009
0010 #include "FWCore/Framework/interface/Event.h"
0011 #include "DataFormats/Common/interface/View.h"
0012 #include "DataFormats/JetReco/interface/Jet.h"
0013 #include "DataFormats/TrackReco/interface/Track.h"
0014 #include "DataFormats/JetReco/interface/JetTracksAssociation.h"
0015
0016 #include "JetTracksAssociatorAtVertex.h"
0017
0018 JetTracksAssociatorAtVertex::JetTracksAssociatorAtVertex(const edm::ParameterSet& fConfig)
0019 : mJets{consumes<edm::View<reco::Jet> >(fConfig.getParameter<edm::InputTag>("jets"))},
0020 mTracks{consumes<reco::TrackCollection>(fConfig.getParameter<edm::InputTag>("tracks"))},
0021 mAssociator(fConfig.getParameter<double>("coneSize")),
0022 mAssociatorAssigned(fConfig.getParameter<double>("coneSize")),
0023 useAssigned(fConfig.getParameter<bool>("useAssigned")),
0024 pvSrc() {
0025 if (useAssigned) {
0026 pvSrc = consumes<reco::VertexCollection>(fConfig.getParameter<edm::InputTag>("pvSrc"));
0027 }
0028 produces<reco::JetTracksAssociation::Container>();
0029 }
0030
0031 void JetTracksAssociatorAtVertex::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0032 edm::ParameterSetDescription desc;
0033 desc.add<edm::InputTag>("jets", edm::InputTag(""));
0034 desc.add<edm::InputTag>("tracks", edm::InputTag("generalTracks"));
0035 desc.add<double>("coneSize", 0.4);
0036 desc.add<bool>("useAssigned", false);
0037 desc.add<edm::InputTag>("pvSrc", edm::InputTag("offlinePrimaryVertices"));
0038 descriptions.addWithDefaultLabel(desc);
0039 }
0040
0041 void JetTracksAssociatorAtVertex::produce(edm::Event& fEvent, const edm::EventSetup& fSetup) {
0042 edm::Handle<edm::View<reco::Jet> > jets_h;
0043 fEvent.getByToken(mJets, jets_h);
0044 edm::Handle<reco::TrackCollection> tracks_h;
0045 fEvent.getByToken(mTracks, tracks_h);
0046
0047 auto jetTracks = std::make_unique<reco::JetTracksAssociation::Container>(reco::JetRefBaseProd(jets_h));
0048
0049
0050 std::vector<edm::RefToBase<reco::Jet> > allJets;
0051 allJets.reserve(jets_h->size());
0052 for (unsigned i = 0; i < jets_h->size(); ++i)
0053 allJets.push_back(jets_h->refAt(i));
0054 std::vector<reco::TrackRef> allTracks;
0055 allTracks.reserve(tracks_h->size());
0056
0057 for (unsigned i = 0; i < tracks_h->size(); ++i) {
0058 allTracks.push_back(reco::TrackRef(tracks_h, i));
0059 }
0060 if (!useAssigned) {
0061 mAssociator.produce(&*jetTracks, allJets, allTracks);
0062 } else {
0063 edm::Handle<reco::VertexCollection> pvHandle;
0064 fEvent.getByToken(pvSrc, pvHandle);
0065 const reco::VertexCollection& vertices = *pvHandle.product();
0066 mAssociatorAssigned.produce(&*jetTracks, allJets, allTracks, vertices);
0067 }
0068
0069
0070 fEvent.put(std::move(jetTracks));
0071 }