Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-05-18 09:15:20

0001 
0002 #ifndef RecoParticleFlow_PFClusterProducer_PFRecHitCaloNavigatorWithTime_h
0003 #define RecoParticleFlow_PFClusterProducer_PFRecHitCaloNavigatorWithTime_h
0004 
0005 #include <memory>
0006 
0007 #include "FWCore/Framework/interface/ConsumesCollector.h"
0008 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0009 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0010 #include "RecoParticleFlow/PFClusterProducer/interface/PFRecHitNavigatorBase.h"
0011 
0012 #include "RecoCaloTools/Navigation/interface/CaloNavigator.h"
0013 #include "DataFormats/EcalDetId/interface/EBDetId.h"
0014 #include "DataFormats/EcalDetId/interface/EEDetId.h"
0015 #include "DataFormats/EcalDetId/interface/ESDetId.h"
0016 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0017 
0018 #include "Geometry/CaloTopology/interface/EcalEndcapTopology.h"
0019 #include "Geometry/CaloTopology/interface/EcalBarrelTopology.h"
0020 #include "Geometry/CaloTopology/interface/EcalPreshowerTopology.h"
0021 #include "Geometry/CaloTopology/interface/HcalTopology.h"
0022 
0023 #include "Geometry/CaloTopology/interface/CaloTowerTopology.h"
0024 #include "DataFormats/CaloTowers/interface/CaloTowerDetId.h"
0025 
0026 #include "RecoParticleFlow/PFClusterProducer/interface/CaloRecHitResolutionProvider.h"
0027 
0028 template <typename D, typename T, bool ownsTopo = true>
0029 class PFRecHitCaloNavigatorWithTime : public PFRecHitNavigatorBase {
0030 public:
0031   PFRecHitCaloNavigatorWithTime(const edm::ParameterSet& iConfig, edm::ConsumesCollector& cc) {
0032     sigmaCut2_ = pow(iConfig.getParameter<double>("sigmaCut"), 2);
0033     const edm::ParameterSet& timeResConf = iConfig.getParameterSet("timeResolutionCalc");
0034     _timeResolutionCalc = std::make_unique<CaloRecHitResolutionProvider>(timeResConf);
0035   }
0036 
0037   ~PFRecHitCaloNavigatorWithTime() override {
0038     if (!ownsTopo) {
0039       topology_.release();
0040     }
0041   }
0042 
0043   void associateNeighbours(reco::PFRecHit& hit,
0044                            std::unique_ptr<reco::PFRecHitCollection>& hits,
0045                            edm::RefProd<reco::PFRecHitCollection>& refProd) override {
0046     DetId detid(hit.detId());
0047 
0048     CaloNavigator<D> navigator(detid, topology_.get());
0049 
0050     DetId N(0);
0051     DetId E(0);
0052     DetId S(0);
0053     DetId W(0);
0054     DetId NW(0);
0055     DetId NE(0);
0056     DetId SW(0);
0057     DetId SE(0);
0058 
0059     N = navigator.north();
0060     associateNeighbour(N, hit, hits, refProd, 0, 1);
0061 
0062     if (N != DetId(0)) {
0063       NE = navigator.east();
0064     } else {
0065       navigator.home();
0066       E = navigator.east();
0067       NE = navigator.north();
0068     }
0069     associateNeighbour(NE, hit, hits, refProd, 1, 1);
0070     navigator.home();
0071 
0072     S = navigator.south();
0073     associateNeighbour(S, hit, hits, refProd, 0, -1);
0074 
0075     if (S != DetId(0)) {
0076       SW = navigator.west();
0077     } else {
0078       navigator.home();
0079       W = navigator.west();
0080       SW = navigator.south();
0081     }
0082     associateNeighbour(SW, hit, hits, refProd, -1, -1);
0083     navigator.home();
0084 
0085     E = navigator.east();
0086     associateNeighbour(E, hit, hits, refProd, 1, 0);
0087 
0088     if (E != DetId(0)) {
0089       SE = navigator.south();
0090     } else {
0091       navigator.home();
0092       S = navigator.south();
0093       SE = navigator.east();
0094     }
0095     associateNeighbour(SE, hit, hits, refProd, 1, -1);
0096     navigator.home();
0097 
0098     W = navigator.west();
0099     associateNeighbour(W, hit, hits, refProd, -1, 0);
0100 
0101     if (W != DetId(0)) {
0102       NW = navigator.north();
0103     } else {
0104       navigator.home();
0105       N = navigator.north();
0106       NW = navigator.west();
0107     }
0108     associateNeighbour(NW, hit, hits, refProd, -1, 1);
0109   }
0110 
0111 protected:
0112   double sigmaCut2_;
0113   std::unique_ptr<const T> topology_;
0114   std::unique_ptr<CaloRecHitResolutionProvider> _timeResolutionCalc;
0115 
0116   void associateNeighbour(const DetId& id,
0117                           reco::PFRecHit& hit,
0118                           std::unique_ptr<reco::PFRecHitCollection>& hits,
0119                           edm::RefProd<reco::PFRecHitCollection>& refProd,
0120                           short eta,
0121                           short phi) {
0122     double sigma2 = 10000.0;
0123 
0124     auto found_hit = std::lower_bound(hits->begin(), hits->end(), id, [](const reco::PFRecHit& a, DetId b) {
0125       if (DetId(a.detId()).det() == DetId::Hcal || b.det() == DetId::Hcal)
0126         return (HcalDetId)(a.detId()) < (HcalDetId)(b);
0127 
0128       else
0129         return a.detId() < b.rawId();
0130     });
0131 
0132     if (found_hit != hits->end() &&
0133         ((found_hit->detId() == id.rawId()) ||
0134          ((id.det() == DetId::Hcal) && ((HcalDetId)(found_hit->detId()) == (HcalDetId)(id))))) {
0135       sigma2 = _timeResolutionCalc->timeResolution2(hit.energy()) +
0136                _timeResolutionCalc->timeResolution2(found_hit->energy());
0137       const double deltaTime = hit.time() - found_hit->time();
0138       if (deltaTime * deltaTime < sigmaCut2_ * sigma2) {
0139         hit.addNeighbour(eta, phi, 0, found_hit - hits->begin());
0140       }
0141     }
0142   }
0143 };
0144 
0145 #endif