Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-06-16 03:20:14

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   ///Constructor
0014   explicit PFNuclearProducer(const edm::ParameterSet&);
0015 
0016   ///Destructor
0017   ~PFNuclearProducer() override;
0018 
0019 private:
0020   void beginRun(const edm::Run&, const edm::EventSetup&) override;
0021   void endRun(const edm::Run&, const edm::EventSetup&) override;
0022 
0023   ///Produce the PFRecTrack collection
0024   void produce(edm::Event&, const edm::EventSetup&) override;
0025 
0026   ///PFTrackTransformer
0027   PFTrackTransformer* pfTransformer_;
0028   double likelihoodCut_;
0029   std::vector<edm::EDGetTokenT<reco::NuclearInteractionCollection> > nuclearContainers_;
0030 
0031   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magneticFieldToken_;
0032 };
0033 
0034 #include "FWCore/Framework/interface/MakerMacros.h"
0035 DEFINE_FWK_MODULE(PFNuclearProducer);
0036 
0037 using namespace std;
0038 using namespace edm;
0039 PFNuclearProducer::PFNuclearProducer(const ParameterSet& iConfig)
0040     : pfTransformer_(nullptr), magneticFieldToken_(esConsumes<edm::Transition::BeginRun>()) {
0041   produces<reco::PFRecTrackCollection>();
0042   produces<reco::PFNuclearInteractionCollection>();
0043 
0044   std::vector<edm::InputTag> tags = iConfig.getParameter<vector<InputTag> >("nuclearColList");
0045 
0046   for (unsigned int i = 0; i < tags.size(); ++i)
0047     nuclearContainers_.push_back(consumes<reco::NuclearInteractionCollection>(tags[i]));
0048 
0049   likelihoodCut_ = iConfig.getParameter<double>("likelihoodCut");
0050 }
0051 
0052 PFNuclearProducer::~PFNuclearProducer() { delete pfTransformer_; }
0053 
0054 void PFNuclearProducer::produce(Event& iEvent, const EventSetup& iSetup) {
0055   typedef reco::NuclearInteraction::trackRef_iterator trackRef_iterator;
0056 
0057   //create the empty collections
0058   auto pfNuclearColl = std::make_unique<reco::PFNuclearInteractionCollection>();
0059   auto pfNuclearRecTrackColl = std::make_unique<reco::PFRecTrackCollection>();
0060 
0061   reco::PFRecTrackRefProd pfTrackRefProd = iEvent.getRefBeforePut<reco::PFRecTrackCollection>();
0062   int hid = 0;
0063 
0064   // loop on the nuclear interaction collections
0065   for (unsigned int istr = 0; istr < nuclearContainers_.size(); istr++) {
0066     Handle<reco::NuclearInteractionCollection> nuclCollH;
0067     iEvent.getByToken(nuclearContainers_[istr], nuclCollH);
0068     const reco::NuclearInteractionCollection& nuclColl = *(nuclCollH.product());
0069 
0070     // loop on all NuclearInteraction
0071     for (unsigned int icoll = 0; icoll < nuclColl.size(); icoll++) {
0072       if (nuclColl[icoll].likelihood() < likelihoodCut_)
0073         continue;
0074 
0075       reco::PFRecTrackRefVector pfRecTkcoll;
0076 
0077       // convert the secondary tracks
0078       for (trackRef_iterator it = nuclColl[icoll].secondaryTracks_begin(); it != nuclColl[icoll].secondaryTracks_end();
0079            it++) {
0080         reco::PFRecTrack pftrack(
0081             (*it)->charge(), reco::PFRecTrack::KF, it->key(), (reco::TrackRef)((*it).castTo<reco::TrackRef>()));
0082         Trajectory FakeTraj;
0083         bool valid = pfTransformer_->addPoints(pftrack, **it, FakeTraj);
0084         if (valid) {
0085           pfRecTkcoll.push_back(reco::PFRecTrackRef(pfTrackRefProd, hid++));
0086           pfNuclearRecTrackColl->push_back(pftrack);
0087         }
0088       }
0089       reco::NuclearInteractionRef niRef(nuclCollH, icoll);
0090       pfNuclearColl->push_back(reco::PFNuclearInteraction(niRef, pfRecTkcoll));
0091     }
0092   }
0093   iEvent.put(std::move(pfNuclearRecTrackColl));
0094   iEvent.put(std::move(pfNuclearColl));
0095 }
0096 
0097 // ------------ method called once each job just before starting event loop  ------------
0098 void PFNuclearProducer::beginRun(const edm::Run& run, const EventSetup& iSetup) {
0099   auto const& magneticField = &iSetup.getData(magneticFieldToken_);
0100   pfTransformer_ = new PFTrackTransformer(math::XYZVector(magneticField->inTesla(GlobalPoint(0, 0, 0))));
0101   pfTransformer_->OnlyProp();
0102 }
0103 
0104 // ------------ method called once each job just after ending the event loop  ------------
0105 void PFNuclearProducer::endRun(const edm::Run& run, const EventSetup& iSetup) {
0106   delete pfTransformer_;
0107   pfTransformer_ = nullptr;
0108 }