File indexing completed on 2024-06-22 02:24:09
0001 #include "DataFormats/ParticleFlowReco/interface/PFNuclearInteraction.h"
0002 #include "FWCore/Framework/interface/ESHandle.h"
0003 #include "FWCore/Framework/interface/Event.h"
0004 #include "FWCore/Framework/interface/stream/EDProducer.h"
0005 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0006 #include "MagneticField/Engine/interface/MagneticField.h"
0007 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0008 #include "RecoParticleFlow/PFTracking/interface/PFTrackTransformer.h"
0009 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0010
0011 class PFNuclearProducer : public edm::stream::EDProducer<> {
0012 public:
0013
0014 explicit PFNuclearProducer(const edm::ParameterSet&);
0015
0016
0017 ~PFNuclearProducer() override;
0018
0019 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0020
0021 private:
0022 void beginRun(const edm::Run&, const edm::EventSetup&) override;
0023 void endRun(const edm::Run&, const edm::EventSetup&) override;
0024
0025
0026 void produce(edm::Event&, const edm::EventSetup&) override;
0027
0028
0029 PFTrackTransformer* pfTransformer_;
0030 double likelihoodCut_;
0031 std::vector<edm::EDGetTokenT<reco::NuclearInteractionCollection>> nuclearContainers_;
0032
0033 const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magneticFieldToken_;
0034 };
0035
0036 #include "FWCore/Framework/interface/MakerMacros.h"
0037 DEFINE_FWK_MODULE(PFNuclearProducer);
0038
0039 void PFNuclearProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0040 edm::ParameterSetDescription desc;
0041
0042 desc.add<double>("likelihoodCut", 0.1);
0043 desc.add<std::vector<edm::InputTag>>("nuclearColList", {edm::InputTag("firstnuclearInteractionMaker")});
0044 descriptions.add("pfNuclear", desc);
0045 }
0046
0047 using namespace std;
0048 using namespace edm;
0049 PFNuclearProducer::PFNuclearProducer(const ParameterSet& iConfig)
0050 : pfTransformer_(nullptr), magneticFieldToken_(esConsumes<edm::Transition::BeginRun>()) {
0051 produces<reco::PFRecTrackCollection>();
0052 produces<reco::PFNuclearInteractionCollection>();
0053
0054 std::vector<edm::InputTag> tags = iConfig.getParameter<vector<InputTag>>("nuclearColList");
0055
0056 for (unsigned int i = 0; i < tags.size(); ++i)
0057 nuclearContainers_.push_back(consumes<reco::NuclearInteractionCollection>(tags[i]));
0058
0059 likelihoodCut_ = iConfig.getParameter<double>("likelihoodCut");
0060 }
0061
0062 PFNuclearProducer::~PFNuclearProducer() { delete pfTransformer_; }
0063
0064 void PFNuclearProducer::produce(Event& iEvent, const EventSetup& iSetup) {
0065 typedef reco::NuclearInteraction::trackRef_iterator trackRef_iterator;
0066
0067
0068 auto pfNuclearColl = std::make_unique<reco::PFNuclearInteractionCollection>();
0069 auto pfNuclearRecTrackColl = std::make_unique<reco::PFRecTrackCollection>();
0070
0071 reco::PFRecTrackRefProd pfTrackRefProd = iEvent.getRefBeforePut<reco::PFRecTrackCollection>();
0072 int hid = 0;
0073
0074
0075 for (unsigned int istr = 0; istr < nuclearContainers_.size(); istr++) {
0076 Handle<reco::NuclearInteractionCollection> nuclCollH;
0077 iEvent.getByToken(nuclearContainers_[istr], nuclCollH);
0078 const reco::NuclearInteractionCollection& nuclColl = *(nuclCollH.product());
0079
0080
0081 for (unsigned int icoll = 0; icoll < nuclColl.size(); icoll++) {
0082 if (nuclColl[icoll].likelihood() < likelihoodCut_)
0083 continue;
0084
0085 reco::PFRecTrackRefVector pfRecTkcoll;
0086
0087
0088 for (trackRef_iterator it = nuclColl[icoll].secondaryTracks_begin(); it != nuclColl[icoll].secondaryTracks_end();
0089 it++) {
0090 reco::PFRecTrack pftrack(
0091 (*it)->charge(), reco::PFRecTrack::KF, it->key(), (reco::TrackRef)((*it).castTo<reco::TrackRef>()));
0092 Trajectory FakeTraj;
0093 bool valid = pfTransformer_->addPoints(pftrack, **it, FakeTraj);
0094 if (valid) {
0095 pfRecTkcoll.push_back(reco::PFRecTrackRef(pfTrackRefProd, hid++));
0096 pfNuclearRecTrackColl->push_back(pftrack);
0097 }
0098 }
0099 reco::NuclearInteractionRef niRef(nuclCollH, icoll);
0100 pfNuclearColl->push_back(reco::PFNuclearInteraction(niRef, pfRecTkcoll));
0101 }
0102 }
0103 iEvent.put(std::move(pfNuclearRecTrackColl));
0104 iEvent.put(std::move(pfNuclearColl));
0105 }
0106
0107
0108 void PFNuclearProducer::beginRun(const edm::Run& run, const EventSetup& iSetup) {
0109 auto const& magneticField = &iSetup.getData(magneticFieldToken_);
0110 pfTransformer_ = new PFTrackTransformer(math::XYZVector(magneticField->inTesla(GlobalPoint(0, 0, 0))));
0111 pfTransformer_->OnlyProp();
0112 }
0113
0114
0115 void PFNuclearProducer::endRun(const edm::Run& run, const EventSetup& iSetup) {
0116 delete pfTransformer_;
0117 pfTransformer_ = nullptr;
0118 }