Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2024-04-06 12:24:50 

0001 // -*- C++ -*-
0002 //
0003 // Package:    ElectronProducers
0004 // Class:      LowPtGsfElectronSeedProducer
0005 //
0006 /**\class LowPtGsfElectronSeedProducer RecoEgamma/EgammaElectronProducers/plugins/LowPtGsfElectronSeedProducer.cc
0007  Description: EDProducer of ElectronSeed objects
0008  Implementation:
0009      <Notes on implementation>
0010 */
0011 // Original Author:  Robert Bainbridge
0012 
0013 #include "DataFormats/GeometryVector/interface/GlobalPoint.h"
0014 #include "DataFormats/GeometryVector/interface/GlobalVector.h"
0015 #include "DataFormats/Common/interface/ValueMap.h"
0016 #include "DataFormats/Math/interface/deltaPhi.h"
0017 #include "DataFormats/Math/interface/deltaR.h"
0018 #include "DataFormats/Math/interface/LorentzVector.h"
0019 #include "DataFormats/ParticleFlowReco/interface/PFTrajectoryPoint.h"
0020 #include "DataFormats/TrackReco/interface/TrackBase.h"
0021 #include "CommonTools/BaseParticlePropagator/interface/BaseParticlePropagator.h"
0022 #include "CommonTools/BaseParticlePropagator/interface/RawParticle.h"
0023 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0024 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0025 #include "FWCore/Utilities/interface/InputTag.h"
0026 #include "MagneticField/Records/interface/IdealMagneticFieldRecord.h"
0027 #include "RecoTracker/TransientTrackingRecHit/interface/TkClonerImpl.h"
0028 #include "RecoTracker/TransientTrackingRecHit/interface/TkTransientTrackingRecHitBuilder.h"
0029 #include "TrackingTools/PatternTools/interface/Trajectory.h"
0030 #include "TrackingTools/Records/interface/TransientRecHitRecord.h"
0031 #include "TrackingTools/TrajectoryParametrization/interface/GlobalTrajectoryParameters.h"
0032 #include "TrackingTools/TrajectoryState/interface/TrajectoryStateOnSurface.h"
0033 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0034 #include "DataFormats/Common/interface/Handle.h"
0035 #include "DataFormats/Common/interface/RefProd.h"
0036 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
0037 #include "DataFormats/EgammaReco/interface/ElectronSeed.h"
0038 #include "DataFormats/EgammaReco/interface/ElectronSeedFwd.h"
0039 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0040 #include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
0041 #include "DataFormats/ParticleFlowReco/interface/PFRecTrack.h"
0042 #include "DataFormats/ParticleFlowReco/interface/PFRecTrackFwd.h"
0043 #include "DataFormats/ParticleFlowReco/interface/PreId.h"
0044 #include "DataFormats/ParticleFlowReco/interface/PreIdFwd.h"
0045 #include "DataFormats/TrackReco/interface/Track.h"
0046 #include "DataFormats/TrackReco/interface/TrackFwd.h"
0047 #include "FWCore/Framework/interface/ESHandle.h"
0048 #include "FWCore/Framework/interface/stream/EDProducer.h"
0049 #include "FWCore/ParameterSet/interface/ConfigurationDescriptions.h"
0050 #include "MagneticField/Engine/interface/MagneticField.h"
0051 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterLazyTools.h"
0052 #include "TrackingTools/PatternTools/interface/TrajectorySmoother.h"
0053 #include "TrackingTools/TrackFitters/interface/TrajectoryFitter.h"
0054 
0055 #include "LowPtGsfElectronSeedHeavyObjectCache.h"
0056 
0057 class LowPtGsfElectronSeedProducer final
0058     : public edm::stream::EDProducer<edm::GlobalCache<lowptgsfeleseed::HeavyObjectCache> > {
0059 public:
0060   using TrackIndxMap = std::unordered_map<reco::TrackRef::key_type, size_t>;
0061   explicit LowPtGsfElectronSeedProducer(const edm::ParameterSet&, const lowptgsfeleseed::HeavyObjectCache*);
0062 
0063   static std::unique_ptr<lowptgsfeleseed::HeavyObjectCache> initializeGlobalCache(const edm::ParameterSet& conf) {
0064     return std::make_unique<lowptgsfeleseed::HeavyObjectCache>(lowptgsfeleseed::HeavyObjectCache(conf));
0065   }
0066 
0067   static void globalEndJob(lowptgsfeleseed::HeavyObjectCache const*) {}
0068 
0069   void produce(edm::Event&, const edm::EventSetup&) override;
0070 
0071   static void fillDescriptions(edm::ConfigurationDescriptions&);
0072 
0073 private:  // member functions
0074   template <typename T>
0075   void loop(const edm::Handle<std::vector<T> >& handle,
0076             edm::Handle<reco::PFClusterCollection>& hcalClusters,
0077             reco::ElectronSeedCollection& seeds,
0078             reco::PreIdCollection& ecalPreIds,
0079             reco::PreIdCollection& hcalPreIds,
0080             TrackIndxMap& trksToPreIdIndx,
0081             edm::Event&,
0082             const edm::EventSetup&);
0083 
0084   // Overloaded methods to retrieve reco::TrackRef
0085 
0086   reco::TrackRef getBaseRef(edm::Handle<std::vector<reco::Track> > handle, int idx) const;
0087   reco::TrackRef getBaseRef(edm::Handle<std::vector<reco::PFRecTrack> > handle, int idx) const;
0088 
0089   // Overloaded methods to populate PreIds (using PF or KF tracks)
0090 
0091   void propagateTrackToCalo(const reco::PFRecTrackRef& pfTrackRef,
0092                             const edm::Handle<reco::PFClusterCollection>& ecalClusters,
0093                             const edm::Handle<reco::PFClusterCollection>& hcalClusters,
0094                             std::vector<int>& matchedEcalClusters,
0095                             std::vector<int>& matchedHcalClusters,
0096                             reco::PreId& ecalPreId,
0097                             reco::PreId& hcalPreId);
0098 
0099   void propagateTrackToCalo(const reco::PFRecTrackRef& pfTrackRef,
0100                             const edm::Handle<reco::PFClusterCollection>& clusters,
0101                             std::vector<int>& matchedClusters,
0102                             reco::PreId& preId,
0103                             bool ecal);
0104 
0105   void propagateTrackToCalo(const reco::TrackRef& pfTrack,
0106                             const edm::Handle<reco::PFClusterCollection>& ecalClusters,
0107                             const edm::Handle<reco::PFClusterCollection>& hcalClusters,
0108                             std::vector<int>& matchedEcalClusters,
0109                             std::vector<int>& matchedHcalClusters,
0110                             reco::PreId& ecalPreId,
0111                             reco::PreId& hcalPreId);
0112   template <typename CollType>
0113   void fillPreIdRefValueMap(edm::Handle<CollType> tracksHandle,
0114                             const TrackIndxMap& trksToPreIdIndx,
0115                             const edm::OrphanHandle<reco::PreIdCollection>& preIdHandle,
0116                             edm::ValueMap<reco::PreIdRef>::Filler& filler);
0117 
0118   // Overloaded methods to evaluate BDTs (using PF or KF tracks)
0119 
0120   bool decision(const reco::PFRecTrackRef& pfTrackRef,
0121                 reco::PreId& ecal,
0122                 reco::PreId& hcal,
0123                 double rho,
0124                 const reco::BeamSpot& spot,
0125                 noZS::EcalClusterLazyTools& ecalTools);
0126 
0127   bool decision(const reco::TrackRef& kfTrackRef,
0128                 reco::PreId& ecal,
0129                 reco::PreId& hcal,
0130                 double rho,
0131                 const reco::BeamSpot& spot,
0132                 noZS::EcalClusterLazyTools& ecalTools);
0133 
0134   // Perform lightweight GSF tracking
0135   bool lightGsfTracking(reco::PreId&, const reco::TrackRef&, const reco::ElectronSeed&);
0136 
0137 private:  // member data
0138   edm::ESHandle<MagneticField> field_;
0139   std::unique_ptr<TrajectoryFitter> fitterPtr_;
0140   std::unique_ptr<TrajectorySmoother> smootherPtr_;
0141   edm::EDGetTokenT<reco::TrackCollection> kfTracks_;
0142   edm::EDGetTokenT<reco::PFRecTrackCollection> pfTracks_;
0143   const edm::EDGetTokenT<reco::PFClusterCollection> ecalClusters_;
0144   edm::EDGetTokenT<reco::PFClusterCollection> hcalClusters_;
0145   const edm::EDGetTokenT<EcalRecHitCollection> ebRecHits_;
0146   const edm::EDGetTokenT<EcalRecHitCollection> eeRecHits_;
0147   const edm::EDGetTokenT<double> rho_;
0148   const edm::EDGetTokenT<reco::BeamSpot> beamSpot_;
0149 
0150   const edm::ESGetToken<TrajectoryFitter, TrajectoryFitter::Record> trajectoryFitterToken_;
0151   const edm::ESGetToken<TrajectorySmoother, TrajectoryFitter::Record> trajectorySmootherToken_;
0152   const edm::ESGetToken<TransientTrackingRecHitBuilder, TransientRecHitRecord> builderToken_;
0153   const edm::ESGetToken<MagneticField, IdealMagneticFieldRecord> magToken_;
0154   const noZS::EcalClusterLazyTools::ESGetTokens ecalClusterToolsESGetTokens_;
0155 
0156   const bool passThrough_;
0157   const bool usePfTracks_;
0158   const double minPtThreshold_;
0159   const double maxPtThreshold_;
0160 
0161   // pow( sinh(1.65), 2. )
0162   static constexpr double boundary_ = 2.50746495928 * 2.50746495928;
0163   // pow( ele_mass, 2. )
0164   static constexpr double mass_ = 0.000511 * 0.000511;
0165 };
0166 
0167 //////////////////////////////////////////////////////////////////////////////////////////
0168 //
0169 LowPtGsfElectronSeedProducer::LowPtGsfElectronSeedProducer(const edm::ParameterSet& conf,
0170                                                            const lowptgsfeleseed::HeavyObjectCache*)
0171     : field_(),
0172       fitterPtr_(),
0173       smootherPtr_(),
0174       kfTracks_(),
0175       pfTracks_(),
0176       ecalClusters_{consumes(conf.getParameter<edm::InputTag>("ecalClusters"))},
0177       hcalClusters_(),
0178       ebRecHits_{consumes(conf.getParameter<edm::InputTag>("EBRecHits"))},
0179       eeRecHits_{consumes(conf.getParameter<edm::InputTag>("EERecHits"))},
0180       rho_(consumes(conf.getParameter<edm::InputTag>("rho"))),
0181       beamSpot_(consumes(conf.getParameter<edm::InputTag>("BeamSpot"))),
0182       trajectoryFitterToken_{esConsumes(conf.getParameter<edm::ESInputTag>("Fitter"))},
0183       trajectorySmootherToken_{esConsumes(conf.getParameter<edm::ESInputTag>("Smoother"))},
0184       builderToken_{esConsumes(conf.getParameter<edm::ESInputTag>("TTRHBuilder"))},
0185       magToken_{esConsumes()},
0186       ecalClusterToolsESGetTokens_{consumesCollector()},
0187       passThrough_(conf.getParameter<bool>("PassThrough")),
0188       usePfTracks_(conf.getParameter<bool>("UsePfTracks")),
0189       minPtThreshold_(conf.getParameter<double>("MinPtThreshold")),
0190       maxPtThreshold_(conf.getParameter<double>("MaxPtThreshold")) {
0191   if (usePfTracks_) {
0192     pfTracks_ = consumes(conf.getParameter<edm::InputTag>("pfTracks"));
0193     hcalClusters_ = consumes(conf.getParameter<edm::InputTag>("hcalClusters"));
0194   }
0195   kfTracks_ = consumes(conf.getParameter<edm::InputTag>("tracks"));
0196 
0197   produces<reco::ElectronSeedCollection>();
0198   produces<reco::PreIdCollection>();
0199   produces<reco::PreIdCollection>("HCAL");
0200   produces<edm::ValueMap<reco::PreIdRef> >();  // indexed by edm::Ref<ElectronSeed>.index()
0201 }
0202 
0203 //////////////////////////////////////////////////////////////////////////////////////////
0204 //
0205 void LowPtGsfElectronSeedProducer::produce(edm::Event& event, const edm::EventSetup& setup) {
0206   field_ = setup.getHandle(magToken_);
0207   // Products
0208   auto seeds = std::make_unique<reco::ElectronSeedCollection>();
0209   auto ecalPreIds = std::make_unique<reco::PreIdCollection>();
0210   auto hcalPreIds = std::make_unique<reco::PreIdCollection>();
0211 
0212   const edm::RefProd<reco::PreIdCollection> preIdsRefProd = event.getRefBeforePut<reco::PreIdCollection>();
0213 
0214   // HCAL clusters (only used with PF tracks)
0215   edm::Handle<reco::PFClusterCollection> hcalClusters;
0216 
0217   //we always need kftracks as we link the preids to them
0218   edm::Handle<reco::TrackCollection> kfTracks;
0219   event.getByToken(kfTracks_, kfTracks);
0220 
0221   TrackIndxMap trksToPreIdIndx;
0222   if (usePfTracks_) {
0223     edm::Handle<reco::PFRecTrackCollection> pfTracks;
0224     event.getByToken(pfTracks_, pfTracks);
0225     event.getByToken(hcalClusters_, hcalClusters);
0226 
0227     //check consistency between kfTracks and pfTracks collection
0228     for (auto& trk : *pfTracks) {
0229       if (trk.trackRef().isNonnull()) {
0230         if (trk.trackRef().id() != kfTracks.id()) {
0231           throw cms::Exception("ConfigError")
0232               << "kfTracks is not the collection that pfTracks was built from, please fix this";
0233         }
0234         break;  //we only need one valid track ref to check this so end the loop
0235       }
0236     }
0237 
0238     loop(pfTracks,  // PF tracks
0239          hcalClusters,
0240          *seeds,
0241          *ecalPreIds,
0242          *hcalPreIds,
0243          trksToPreIdIndx,
0244          event,
0245          setup);
0246 
0247   } else {
0248     loop(kfTracks,  // KF tracks
0249          hcalClusters,
0250          *seeds,
0251          *ecalPreIds,
0252          *hcalPreIds,
0253          trksToPreIdIndx,
0254          event,
0255          setup);
0256   }
0257 
0258   auto ecalPreIdsHandle = event.put(std::move(ecalPreIds));
0259   event.put(std::move(hcalPreIds), "HCAL");
0260   event.put(std::move(seeds));
0261 
0262   auto preIdVMOut = std::make_unique<edm::ValueMap<reco::PreIdRef> >();
0263   edm::ValueMap<reco::PreIdRef>::Filler mapFiller(*preIdVMOut);
0264   fillPreIdRefValueMap(kfTracks, trksToPreIdIndx, ecalPreIdsHandle, mapFiller);
0265   mapFiller.fill();
0266   event.put(std::move(preIdVMOut));
0267 }
0268 
0269 //////////////////////////////////////////////////////////////////////////////////////////
0270 // Return reco::Track from edm::Ref<T>
0271 
0272 reco::TrackRef LowPtGsfElectronSeedProducer::getBaseRef(edm::Handle<std::vector<reco::Track> > handle, int idx) const {
0273   return reco::TrackRef(handle, idx);
0274 }
0275 
0276 reco::TrackRef LowPtGsfElectronSeedProducer::getBaseRef(edm::Handle<std::vector<reco::PFRecTrack> > handle,
0277                                                         int idx) const {
0278   return reco::PFRecTrackRef(handle, idx)->trackRef();
0279 }
0280 
0281 //////////////////////////////////////////////////////////////////////////////////////////
0282 // Template function, instantiated for both reco::Tracks and reco::PFRecTracks
0283 template <typename T>
0284 void LowPtGsfElectronSeedProducer::loop(const edm::Handle<std::vector<T> >& handle,  // PF or KF tracks
0285                                         edm::Handle<reco::PFClusterCollection>& hcalClusters,
0286                                         reco::ElectronSeedCollection& seeds,
0287                                         reco::PreIdCollection& ecalPreIds,
0288                                         reco::PreIdCollection& hcalPreIds,
0289                                         TrackIndxMap& trksToPreIdIndx,
0290                                         edm::Event& event,
0291                                         const edm::EventSetup& setup) {
0292   // Pileup
0293   auto const& rho = event.get(rho_);
0294 
0295   // Beam spot
0296   auto const& spot = event.get(beamSpot_);
0297 
0298   // Track fitter
0299   fitterPtr_ = setup.getData(trajectoryFitterToken_).clone();
0300 
0301   // Track smoother
0302   smootherPtr_.reset(setup.getData(trajectorySmootherToken_).clone());
0303 
0304   // RecHit cloner
0305   TkClonerImpl hitCloner = static_cast<TkTransientTrackingRecHitBuilder const&>(setup.getData(builderToken_)).cloner();
0306   fitterPtr_->setHitCloner(&hitCloner);
0307   smootherPtr_->setHitCloner(&hitCloner);
0308 
0309   // ECAL clusters
0310   auto ecalClusters = event.getHandle(ecalClusters_);
0311 
0312   // Utility to access to shower shape vars
0313   noZS::EcalClusterLazyTools ecalTools(event, ecalClusterToolsESGetTokens_.get(setup), ebRecHits_, eeRecHits_);
0314 
0315   // Ensure each cluster is only matched once to a track
0316   std::vector<int> matchedEcalClusters;
0317   std::vector<int> matchedHcalClusters;
0318 
0319   // Reserve
0320   seeds.reserve(handle->size());
0321   ecalPreIds.reserve(handle->size());
0322   hcalPreIds.reserve(handle->size());
0323 
0324   // Iterate through (PF or KF) tracks
0325   for (unsigned int itrk = 0; itrk < handle.product()->size(); itrk++) {
0326     edm::Ref<std::vector<T> > templatedRef(handle, itrk);  // TrackRef or PFRecTrackRef
0327     reco::TrackRef trackRef = getBaseRef(handle, itrk);
0328 
0329     if (!(trackRef->quality(reco::TrackBase::qualityByName("highPurity")))) {
0330       continue;
0331     }
0332     if (!passThrough_ && (trackRef->pt() < minPtThreshold_)) {
0333       continue;
0334     }
0335 
0336     // Create ElectronSeed
0337     reco::ElectronSeed seed(*(trackRef->seedRef()));
0338     if (seed.nHits() == 0) {  //if DeepCore is used in jetCore iteration the seed are hitless, in case skip
0339       continue;
0340     }
0341     seed.setCtfTrack(trackRef);
0342 
0343     // Create PreIds
0344     unsigned int nModels = globalCache()->modelNames().size();
0345     reco::PreId ecalPreId(nModels);
0346     reco::PreId hcalPreId(nModels);
0347 
0348     // Add track ref to PreId
0349     ecalPreId.setTrack(trackRef);
0350     hcalPreId.setTrack(trackRef);
0351 
0352     // Add Track-Calo matching variables to PreIds
0353     propagateTrackToCalo(
0354         templatedRef, ecalClusters, hcalClusters, matchedEcalClusters, matchedHcalClusters, ecalPreId, hcalPreId);
0355 
0356     // Add variables related to GSF tracks to PreId
0357     lightGsfTracking(ecalPreId, trackRef, seed);
0358 
0359     // Decision based on BDT
0360     bool result = decision(templatedRef, ecalPreId, hcalPreId, rho, spot, ecalTools);
0361 
0362     // If fails BDT, do not store seed
0363     if (!result) {
0364       continue;
0365     }
0366 
0367     // Store PreId
0368     ecalPreIds.push_back(ecalPreId);
0369     hcalPreIds.push_back(hcalPreId);
0370     trksToPreIdIndx[trackRef.index()] = ecalPreIds.size() - 1;
0371 
0372     // Store ElectronSeed and corresponding edm::Ref<PreId>.index()
0373     seeds.push_back(seed);
0374   }
0375 }
0376 
0377 //////////////////////////////////////////////////////////////////////////////////////////
0378 // Template instantiation for reco::Tracks
0379 template void LowPtGsfElectronSeedProducer::loop<reco::Track>(const edm::Handle<std::vector<reco::Track> >&,
0380                                                               edm::Handle<reco::PFClusterCollection>& hcalClusters,
0381                                                               reco::ElectronSeedCollection& seeds,
0382                                                               reco::PreIdCollection& ecalPreIds,
0383                                                               reco::PreIdCollection& hcalPreIds,
0384                                                               TrackIndxMap& trksToPreIdIndx,
0385                                                               edm::Event&,
0386                                                               const edm::EventSetup&);
0387 
0388 //////////////////////////////////////////////////////////////////////////////////////////
0389 // Template instantiation for reco::PFRecTracks
0390 template void LowPtGsfElectronSeedProducer::loop<reco::PFRecTrack>(const edm::Handle<std::vector<reco::PFRecTrack> >&,
0391                                                                    edm::Handle<reco::PFClusterCollection>& hcalClusters,
0392                                                                    reco::ElectronSeedCollection& seeds,
0393                                                                    reco::PreIdCollection& ecalPreIds,
0394                                                                    reco::PreIdCollection& hcalPreIds,
0395                                                                    TrackIndxMap& trksToPreIdIndx,
0396                                                                    edm::Event&,
0397                                                                    const edm::EventSetup&);
0398 
0399 //////////////////////////////////////////////////////////////////////////////////////////
0400 // Loops through both ECAL and HCAL clusters
0401 void LowPtGsfElectronSeedProducer::propagateTrackToCalo(const reco::PFRecTrackRef& pfTrackRef,
0402                                                         const edm::Handle<reco::PFClusterCollection>& ecalClusters,
0403                                                         const edm::Handle<reco::PFClusterCollection>& hcalClusters,
0404                                                         std::vector<int>& matchedEcalClusters,
0405                                                         std::vector<int>& matchedHcalClusters,
0406                                                         reco::PreId& ecalPreId,
0407                                                         reco::PreId& hcalPreId) {
0408   propagateTrackToCalo(pfTrackRef, ecalClusters, matchedEcalClusters, ecalPreId, true);
0409   propagateTrackToCalo(pfTrackRef, hcalClusters, matchedHcalClusters, hcalPreId, false);
0410 }
0411 
0412 //////////////////////////////////////////////////////////////////////////////////////////
0413 // Loops through ECAL or HCAL clusters (called twice)
0414 void LowPtGsfElectronSeedProducer::propagateTrackToCalo(const reco::PFRecTrackRef& pfTrackRef,
0415                                                         const edm::Handle<reco::PFClusterCollection>& clusters,
0416                                                         std::vector<int>& matched,
0417                                                         reco::PreId& preId,
0418                                                         bool ecal) {
0419   // Store info for PreId
0420   struct Info {
0421     reco::PFClusterRef cluRef = reco::PFClusterRef();
0422     float dr2min = 1.e6;
0423     float deta = 1.e6;
0424     float dphi = 1.e6;
0425     math::XYZPoint showerPos = math::XYZPoint(0., 0., 0.);
0426   } info;
0427 
0428   // Find closest "seed cluster" to KF track extrapolated to ECAL (or HCAL)
0429   reco::PFTrajectoryPoint point;
0430   if (ecal) {
0431     point = pfTrackRef->extrapolatedPoint(reco::PFTrajectoryPoint::LayerType::ECALShowerMax);
0432   } else {
0433     point = pfTrackRef->extrapolatedPoint(reco::PFTrajectoryPoint::LayerType::HCALEntrance);
0434   }
0435 
0436   if (point.isValid()) {
0437     Info info;
0438     for (unsigned int iclu = 0; iclu < clusters.product()->size(); iclu++) {
0439       if (std::find(matched.begin(), matched.end(), iclu) == matched.end()) {
0440         reco::PFClusterRef cluRef(clusters, iclu);
0441 
0442         // Determine dR squared
0443         float dr2 = reco::deltaR2(cluRef->positionREP(), point.positionREP());
0444 
0445         if (dr2 < info.dr2min) {
0446           info.dr2min = dr2;
0447           info.cluRef = cluRef;
0448           info.deta = cluRef->positionREP().eta() - point.positionREP().eta();
0449           info.dphi =
0450               reco::deltaPhi(cluRef->positionREP().phi(), point.positionREP().phi()) * pfTrackRef->trackRef()->charge();
0451           info.showerPos = point.position();
0452         }
0453       }
0454     }
0455 
0456     // Set PreId content if match found
0457     if (info.dr2min < 1.e5) {
0458       float ep = info.cluRef->correctedEnergy() / std::sqrt(pfTrackRef->trackRef()->innerMomentum().mag2());
0459       preId.setECALMatchingProperties(info.cluRef,
0460                                       point.position(),  // ECAL or HCAL surface
0461                                       info.showerPos,    //
0462                                       info.deta,
0463                                       info.dphi,
0464                                       0.f,                                       // chieta
0465                                       0.f,                                       // chiphi
0466                                       pfTrackRef->trackRef()->normalizedChi2(),  // chi2
0467                                       ep);
0468     }
0469 
0470   }  // clusters
0471 }
0472 
0473 //////////////////////////////////////////////////////////////////////////////////////////
0474 // Original implementation in GoodSeedProducer, loops over ECAL clusters only
0475 void LowPtGsfElectronSeedProducer::propagateTrackToCalo(
0476     const reco::TrackRef& kfTrackRef,
0477     const edm::Handle<reco::PFClusterCollection>& ecalClusters,
0478     const edm::Handle<reco::PFClusterCollection>& hcalClusters,  // not used
0479     std::vector<int>& matchedEcalClusters,
0480     std::vector<int>& matchedHcalClusters,  // not used
0481     reco::PreId& ecalPreId,
0482     reco::PreId& hcalPreId /* not used */) {
0483   // Store info for PreId
0484   struct Info {
0485     reco::PFClusterRef cluRef = reco::PFClusterRef();
0486     float dr2min = 1.e6;
0487     float deta = 1.e6;
0488     float dphi = 1.e6;
0489     math::XYZPoint showerPos = math::XYZPoint(0., 0., 0.);
0490   } info;
0491 
0492   // Propagate 'electron' to ECAL surface
0493   float energy = sqrt(mass_ + kfTrackRef->outerMomentum().Mag2());
0494   XYZTLorentzVector mom = XYZTLorentzVector(
0495       kfTrackRef->outerMomentum().x(), kfTrackRef->outerMomentum().y(), kfTrackRef->outerMomentum().z(), energy);
0496   XYZTLorentzVector pos = XYZTLorentzVector(
0497       kfTrackRef->outerPosition().x(), kfTrackRef->outerPosition().y(), kfTrackRef->outerPosition().z(), 0.);
0498   math::XYZVector field(field_->inTesla(GlobalPoint(0, 0, 0)));
0499   BaseParticlePropagator particle(RawParticle(mom, pos, kfTrackRef->charge()), 0, 0, field.z());
0500   particle.propagateToEcalEntrance(false);
0501   if (particle.getSuccess() == 0) {
0502     return;
0503   }
0504 
0505   // ECAL entry point for track
0506   GlobalPoint ecal_pos(
0507       particle.particle().vertex().x(), particle.particle().vertex().y(), particle.particle().vertex().z());
0508 
0509   // Preshower limit
0510   bool below_ps = pow(ecal_pos.z(), 2.) > boundary_ * ecal_pos.perp2();
0511 
0512   // Iterate through ECAL clusters
0513   for (unsigned int iclu = 0; iclu < ecalClusters.product()->size(); iclu++) {
0514     reco::PFClusterRef cluRef(ecalClusters, iclu);
0515 
0516     // Correct ecal_pos for shower depth
0517     double shower_depth = reco::PFCluster::getDepthCorrection(cluRef->correctedEnergy(), below_ps, false);
0518     GlobalPoint showerPos = ecal_pos + GlobalVector(particle.particle().momentum().x(),
0519                                                     particle.particle().momentum().y(),
0520                                                     particle.particle().momentum().z())
0521                                                .unit() *
0522                                            shower_depth;
0523 
0524     // Determine dR squared
0525     float dr2 = reco::deltaR2(cluRef->positionREP(), showerPos);
0526 
0527     // Find nearest ECAL cluster
0528     if (dr2 < info.dr2min) {
0529       info.dr2min = dr2;
0530       info.cluRef = cluRef;
0531       info.deta = std::abs(cluRef->positionREP().eta() - showerPos.eta());
0532       info.dphi = std::abs(reco::deltaPhi(cluRef->positionREP().phi(), showerPos.phi())) * kfTrackRef->charge();
0533       info.showerPos = showerPos;
0534     }
0535   }
0536 
0537   // Populate PreId object
0538   math::XYZPoint point(ecal_pos.x(), ecal_pos.y(), ecal_pos.z());
0539 
0540   // Set PreId content
0541   ecalPreId.setECALMatchingProperties(
0542       info.cluRef,
0543       point,
0544       info.showerPos,
0545       info.deta,
0546       info.dphi,
0547       0.f,                                                                              // chieta
0548       0.f,                                                                              // chiphi
0549       kfTrackRef->normalizedChi2(),                                                     // chi2
0550       info.cluRef->correctedEnergy() / std::sqrt(kfTrackRef->innerMomentum().mag2()));  // E/p
0551 }
0552 
0553 //////////////////////////////////////////////////////////////////////////////////////////
0554 // Original implementation for "lightweight" GSF tracking
0555 bool LowPtGsfElectronSeedProducer::lightGsfTracking(reco::PreId& preId,
0556                                                     const reco::TrackRef& trackRef,
0557                                                     const reco::ElectronSeed& seed) {
0558   Trajectory::ConstRecHitContainer hits;
0559   for (unsigned int ihit = 0; ihit < trackRef->recHitsSize(); ++ihit) {
0560     hits.push_back(trackRef->recHit(ihit)->cloneSH());
0561   }
0562 
0563   GlobalVector gv(trackRef->innerMomentum().x(), trackRef->innerMomentum().y(), trackRef->innerMomentum().z());
0564   GlobalPoint gp(trackRef->innerPosition().x(), trackRef->innerPosition().y(), trackRef->innerPosition().z());
0565 
0566   GlobalTrajectoryParameters gtps(gp, gv, trackRef->charge(), &*field_);
0567 
0568   TrajectoryStateOnSurface tsos(gtps, trackRef->innerStateCovariance(), *hits[0]->surface());
0569 
0570   // Track fitted and smoothed under electron hypothesis
0571   Trajectory traj1 = fitterPtr_->fitOne(seed, hits, tsos);
0572   if (!traj1.isValid()) {
0573     return false;
0574   }
0575   Trajectory traj2 = smootherPtr_->trajectory(traj1);
0576   if (!traj2.isValid()) {
0577     return false;
0578   }
0579 
0580   // Set PreId content
0581   float chi2Ratio = trackRef->chi2() > 0. ? traj2.chiSquared() / trackRef->chi2() : -1.;
0582   float gsfReducedChi2 = chi2Ratio > -1. ? chi2Ratio * trackRef->normalizedChi2() : -1.;
0583   float ptOut = traj2.firstMeasurement().updatedState().globalMomentum().perp();
0584   float ptIn = traj2.lastMeasurement().updatedState().globalMomentum().perp();
0585   float gsfDpt = (ptIn > 0) ? fabs(ptOut - ptIn) / ptIn : 0.;
0586   preId.setTrackProperties(gsfReducedChi2, chi2Ratio, gsfDpt);
0587 
0588   return true;
0589 }
0590 
0591 //////////////////////////////////////////////////////////////////////////////////////////
0592 // Decision based on OR of outputs from list of models
0593 bool LowPtGsfElectronSeedProducer::decision(const reco::PFRecTrackRef& pfTrackRef,
0594                                             reco::PreId& ecalPreId,
0595                                             reco::PreId& hcalPreId,
0596                                             double rho,
0597                                             const reco::BeamSpot& spot,
0598                                             noZS::EcalClusterLazyTools& ecalTools) {
0599   bool result = false;
0600   for (auto& name : globalCache()->modelNames()) {
0601     result |= globalCache()->eval(name, ecalPreId, hcalPreId, rho, spot, ecalTools);
0602   }
0603   return passThrough_ || (pfTrackRef->trackRef()->pt() > maxPtThreshold_) || result;
0604 }
0605 
0606 //////////////////////////////////////////////////////////////////////////////////////////
0607 //
0608 bool LowPtGsfElectronSeedProducer::decision(const reco::TrackRef& kfTrackRef,
0609                                             reco::PreId& ecalPreId,
0610                                             reco::PreId& hcalPreId,
0611                                             double rho,
0612                                             const reco::BeamSpot& spot,
0613                                             noZS::EcalClusterLazyTools& ecalTools) {
0614   // No implementation currently
0615   return passThrough_;
0616 }
0617 
0618 template <typename CollType>
0619 void LowPtGsfElectronSeedProducer::fillPreIdRefValueMap(edm::Handle<CollType> tracksHandle,
0620                                                         const TrackIndxMap& trksToPreIdIndx,
0621                                                         const edm::OrphanHandle<reco::PreIdCollection>& preIdHandle,
0622                                                         edm::ValueMap<reco::PreIdRef>::Filler& filler) {
0623   std::vector<reco::PreIdRef> values;
0624 
0625   unsigned ntracks = tracksHandle->size();
0626   for (unsigned itrack = 0; itrack < ntracks; ++itrack) {
0627     edm::Ref<CollType> trackRef(tracksHandle, itrack);
0628     auto preIdRefIt = trksToPreIdIndx.find(trackRef.index());
0629     if (preIdRefIt == trksToPreIdIndx.end()) {
0630       values.push_back(reco::PreIdRef());
0631     } else {
0632       edm::Ref<reco::PreIdCollection> preIdRef(preIdHandle, preIdRefIt->second);
0633       values.push_back(preIdRef);
0634     }
0635   }
0636   filler.insert(tracksHandle, values.begin(), values.end());
0637 }
0638 
0639 //////////////////////////////////////////////////////////////////////////////////////////
0640 //
0641 void LowPtGsfElectronSeedProducer::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0642   edm::ParameterSetDescription desc;
0643   desc.add<edm::InputTag>("tracks", edm::InputTag("generalTracks"));
0644   desc.add<edm::InputTag>("pfTracks", edm::InputTag("lowPtGsfElePfTracks"));
0645   desc.add<edm::InputTag>("ecalClusters", edm::InputTag("particleFlowClusterECAL"));
0646   desc.add<edm::InputTag>("hcalClusters", edm::InputTag("particleFlowClusterHCAL"));
0647   desc.add<edm::InputTag>("EBRecHits", edm::InputTag("ecalRecHit", "EcalRecHitsEB"));
0648   desc.add<edm::InputTag>("EERecHits", edm::InputTag("ecalRecHit", "EcalRecHitsEE"));
0649   desc.add<edm::InputTag>("rho", edm::InputTag("fixedGridRhoFastjetAllTmp"));
0650   desc.add<edm::InputTag>("BeamSpot", edm::InputTag("offlineBeamSpot"));
0651   desc.add<edm::ESInputTag>("Fitter", edm::ESInputTag{"", "GsfTrajectoryFitter_forPreId"});
0652   desc.add<edm::ESInputTag>("Smoother", edm::ESInputTag{"", "GsfTrajectorySmoother_forPreId"});
0653   desc.add<edm::ESInputTag>("TTRHBuilder", edm::ESInputTag{"", "WithAngleAndTemplate"});
0654   desc.add<std::vector<std::string> >("ModelNames", {});
0655   desc.add<std::vector<std::string> >("ModelWeights", {});
0656   desc.add<std::vector<double> >("ModelThresholds", {});
0657   desc.add<bool>("PassThrough", false);
0658   desc.add<bool>("UsePfTracks", true);
0659   desc.add<double>("MinPtThreshold", 1.0);
0660   desc.add<double>("MaxPtThreshold", 15.);
0661   descriptions.add("lowPtGsfElectronSeeds", desc);
0662 }
0663 
0664 //////////////////////////////////////////////////////////////////////////////////////////
0665 //
0666 #include "FWCore/Framework/interface/MakerMacros.h"
0667 DEFINE_FWK_MODULE(LowPtGsfElectronSeedProducer);
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!