plugins/linkers/TrackAndHCALLinker.cc

0001 #include "RecoParticleFlow/PFProducer/interface/BlockElementLinkerBase.h"
0002 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
0003 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementCluster.h"
0004 #include "DataFormats/ParticleFlowReco/interface/PFBlockElementTrack.h"
0005 #include "DataFormats/ParticleFlowReco/interface/PFTrajectoryPoint.h"
0006 #include "RecoParticleFlow/PFClusterTools/interface/LinkByRecHit.h"
0007 #include "DataFormats/Math/interface/deltaPhi.h"
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009
0010 class TrackAndHCALLinker : public BlockElementLinkerBase {
0011 public:
0012   TrackAndHCALLinker(const edm::ParameterSet& conf)
0013       : BlockElementLinkerBase(conf),
0014         useKDTree_(conf.getParameter<bool>("useKDTree")),
0015         trajectoryLayerEntranceString_(conf.getParameter<std::string>("trajectoryLayerEntrance")),
0016         trajectoryLayerExitString_(conf.getParameter<std::string>("trajectoryLayerExit")),
0017         debug_(conf.getUntrackedParameter<bool>("debug", false)) {
0018     // convert TrajectoryLayers info from string to enum
0019     trajectoryLayerEntrance_ = reco::PFTrajectoryPoint::layerTypeByName(trajectoryLayerEntranceString_);
0020     trajectoryLayerExit_ = reco::PFTrajectoryPoint::layerTypeByName(trajectoryLayerExitString_);
0021     // make sure the requested setting is supported
0022     assert((trajectoryLayerEntrance_ == reco::PFTrajectoryPoint::HCALEntrance &&
0023             trajectoryLayerExit_ == reco::PFTrajectoryPoint::HCALExit) ||
0024            (trajectoryLayerEntrance_ == reco::PFTrajectoryPoint::HCALEntrance &&
0025             trajectoryLayerExit_ == reco::PFTrajectoryPoint::Unknown) ||
0026            (trajectoryLayerEntrance_ == reco::PFTrajectoryPoint::VFcalEntrance &&
0027             trajectoryLayerExit_ == reco::PFTrajectoryPoint::Unknown));
0028     // flag if exit layer should be checked or not
0029     checkExit_ = trajectoryLayerExit_ != reco::PFTrajectoryPoint::Unknown;
0030   }
0031
0032   bool linkPrefilter(const reco::PFBlockElement*, const reco::PFBlockElement*) const override;
0033
0034   double testLink(const reco::PFBlockElement*, const reco::PFBlockElement*) const override;
0035
0036 private:
0037   bool useKDTree_;
0038   std::string trajectoryLayerEntranceString_;
0039   std::string trajectoryLayerExitString_;
0040   reco::PFTrajectoryPoint::LayerType trajectoryLayerEntrance_;
0041   reco::PFTrajectoryPoint::LayerType trajectoryLayerExit_;
0042   bool debug_;
0043   bool checkExit_;
0044 };
0045
0046 DEFINE_EDM_PLUGIN(BlockElementLinkerFactory, TrackAndHCALLinker, "TrackAndHCALLinker");
0047
0048 bool TrackAndHCALLinker::linkPrefilter(const reco::PFBlockElement* elem1, const reco::PFBlockElement* elem2) const {
0049   bool result = false;
0050   // Track-HCAL KDTree multilinks are stored to track's elem
0051   switch (elem1->type()) {
0052     case reco::PFBlockElement::TRACK:
0053       result = (elem1->isMultilinksValide(elem2->type()) && !elem1->getMultilinks(elem2->type()).empty() &&
0054                 elem2->isMultilinksValide(elem1->type()));
0055       break;
0056     case reco::PFBlockElement::HCAL:
0057       result = (elem2->isMultilinksValide(elem1->type()) && !elem2->getMultilinks(elem1->type()).empty() &&
0058                 elem1->isMultilinksValide(elem2->type()));
0059     default:
0060       break;
0061   }
0062   return (useKDTree_ ? result : true);
0063 }
0064
0065 double TrackAndHCALLinker::testLink(const reco::PFBlockElement* elem1, const reco::PFBlockElement* elem2) const {
0066   const reco::PFBlockElementCluster* hcalelem(nullptr);
0067   const reco::PFBlockElementTrack* tkelem(nullptr);
0068   double dist(-1.0);
0069   if (elem1->type() < elem2->type()) {
0070     tkelem = static_cast<const reco::PFBlockElementTrack*>(elem1);
0071     hcalelem = static_cast<const reco::PFBlockElementCluster*>(elem2);
0072   } else {
0073     tkelem = static_cast<const reco::PFBlockElementTrack*>(elem2);
0074     hcalelem = static_cast<const reco::PFBlockElementCluster*>(elem1);
0075   }
0076   const reco::PFRecTrackRef& trackref = tkelem->trackRefPF();
0077   const reco::PFClusterRef& clusterref = hcalelem->clusterRef();
0078   const reco::PFCluster::REPPoint& hcalreppos = clusterref->positionREP();
0079   const reco::PFTrajectoryPoint& tkAtHCALEnt = trackref->extrapolatedPoint(trajectoryLayerEntrance_);
0080   const reco::PFCluster::REPPoint& tkreppos = tkAtHCALEnt.positionREP();
0081   // Check exit point
0082   double dHEta = 0.;
0083   double dHPhi = 0.;
0084   double dRHCALEx = 0.;
0085   if (checkExit_) {
0086     const reco::PFTrajectoryPoint& tkAtHCALEx = trackref->extrapolatedPoint(trajectoryLayerExit_);
0087     dHEta = (tkAtHCALEx.positionREP().Eta() - tkAtHCALEnt.positionREP().Eta());
0088     dHPhi = reco::deltaPhi(tkAtHCALEx.positionREP().Phi(), tkAtHCALEnt.positionREP().Phi());
0089     dRHCALEx = tkAtHCALEx.position().R();
0090   }
0091
0092   // Check if the linking has been done using the KDTree algo
0093   // Glowinski & Gouzevitch
0094   if (useKDTree_ && tkelem->isMultilinksValide(hcalelem->type())) {  //KDTree Algo
0095     const reco::PFMultilinksType& multilinks = tkelem->getMultilinks(hcalelem->type());
0096
0097     // Check if the link Track/Hcal exist
0098     reco::PFMultilinksType::const_iterator mlit = multilinks.begin();
0099     for (; mlit != multilinks.end(); ++mlit)
0100       if (mlit->clusterRef == clusterref)
0101         break;
0102
0103     // If the link exist, we fill dist and linktest.
0104
0105     if (mlit != multilinks.end()) {
0106       // when checkExit_ is false
0107       if (!checkExit_) {
0108         dist = LinkByRecHit::computeDist(hcalreppos.Eta(), hcalreppos.Phi(), tkreppos.Eta(), tkreppos.Phi());
0109       }
0110       // when checkExit_ is true
0111       else {
0112         //special case ! A looper  can exit the barrel inwards and hit the endcap
0113         //In this case calculate the distance based on the first crossing since
0114         //the looper will probably never make it to the endcap
0115         if (dRHCALEx < tkAtHCALEnt.position().R()) {
0116           dist = LinkByRecHit::computeDist(hcalreppos.Eta(), hcalreppos.Phi(), tkreppos.Eta(), tkreppos.Phi());
0117           edm::LogWarning("TrackHCALLinker ")
0118               << "Special case of linking with track hitting HCAL and looping back in the tracker ";
0119         } else {
0120           dist = LinkByRecHit::computeDist(
0121               hcalreppos.Eta(), hcalreppos.Phi(), tkreppos.Eta() + 0.1 * dHEta, tkreppos.Phi() + 0.1 * dHPhi);
0122         }
0123       }  // checkExit_
0124     }    // multilinks
0125
0126   } else {  // Old algorithm
0127     if (tkAtHCALEnt.isValid())
0128       dist = LinkByRecHit::testTrackAndClusterByRecHit(*trackref, *clusterref, false, debug_);
0129   }
0130   return dist;
0131 }