EgammaTools/src/HGCalEgammaIDHelper.cc

0001 #include "RecoEgamma/EgammaTools/interface/HGCalEgammaIDHelper.h"
0002
0003 #include <iostream>
0004
0005 HGCalEgammaIDHelper::HGCalEgammaIDHelper(const edm::ParameterSet& iConfig, edm::ConsumesCollector&& iC)
0006     : eeRecHitInputTag_(iConfig.getParameter<edm::InputTag>("EERecHits")),
0007       fhRecHitInputTag_(iConfig.getParameter<edm::InputTag>("FHRecHits")),
0008       bhRecHitInputTag_(iConfig.getParameter<edm::InputTag>("BHRecHits")),
0009       hitMapInputTag_(iConfig.getParameter<edm::InputTag>("hitMapTag")),
0010       dEdXWeights_(iConfig.getParameter<std::vector<double>>("dEdXWeights")) {
0011   isoHelper_.setDeltaR(iConfig.getParameter<double>("isoDeltaR"));
0012   isoHelper_.setNRings(iConfig.getParameter<unsigned int>("isoNRings"));
0013   isoHelper_.setMinDeltaR(iConfig.getParameter<double>("isoDeltaRmin"));
0014
0015   recHitsEE_ = iC.consumes<HGCRecHitCollection>(eeRecHitInputTag_);
0016   recHitsFH_ = iC.consumes<HGCRecHitCollection>(fhRecHitInputTag_);
0017   recHitsBH_ = iC.consumes<HGCRecHitCollection>(bhRecHitInputTag_);
0018   hitMap_ = iC.consumes<std::unordered_map<DetId, const unsigned int>>(hitMapInputTag_);
0019   caloGeometry_ = iC.esConsumes();
0020   pcaHelper_.setdEdXWeights(dEdXWeights_);
0021   debug_ = iConfig.getUntrackedParameter<bool>("debug", false);
0022 }
0023
0024 void HGCalEgammaIDHelper::eventInit(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0025   auto recHitHandleEE = iEvent.getHandle(recHitsEE_);
0026   auto recHitHandleFH = iEvent.getHandle(recHitsFH_);
0027   auto recHitHandleBH = iEvent.getHandle(recHitsBH_);
0028
0029   recHitTools_.setGeometry(iSetup.getData(caloGeometry_));
0030   pcaHelper_.setRecHitTools(&recHitTools_);
0031   isoHelper_.setRecHitTools(&recHitTools_);
0032   pcaHelper_.setHitMap(&iEvent.get(hitMap_));
0033   pcaHelper_.setRecHits(recHitHandleEE, recHitHandleFH, recHitHandleBH);
0034   isoHelper_.setRecHits(recHitHandleEE, recHitHandleFH, recHitHandleBH);
0035 }
0036
0037 void HGCalEgammaIDHelper::computeHGCAL(const reco::Photon& thePhoton, float radius) {
0038   if (thePhoton.isEB()) {
0039     if (debug_)
0040       std::cout << "The photon is in the barrel" << std::endl;
0041     pcaHelper_.clear();
0042     return;
0043   }
0044
0045   pcaHelper_.storeRecHits(*thePhoton.superCluster()->seed());
0046   if (debug_)
0047     std::cout << " Stored the hits belonging to the photon superCluster seed " << std::endl;
0048
0049   // initial computation, no radius cut, but halo hits not taken
0050   if (debug_)
0051     std::cout << " Calling PCA initial computation" << std::endl;
0052   pcaHelper_.pcaInitialComputation();
0053   // first computation within cylinder, halo hits included
0054   pcaHelper_.computePCA(radius);
0055   // second computation within cylinder, halo hits included
0056   pcaHelper_.computePCA(radius);
0057   pcaHelper_.computeShowerWidth(radius);
0058
0059   // isolation
0060   isoHelper_.produceHGCalIso(thePhoton.superCluster()->seed());
0061 }
0062
0063 void HGCalEgammaIDHelper::computeHGCAL(const reco::GsfElectron& theElectron, float radius) {
0064   if (theElectron.isEB()) {
0065     if (debug_)
0066       std::cout << "The electron is in the barrel" << std::endl;
0067     pcaHelper_.clear();
0068     return;
0069   }
0070
0071   pcaHelper_.storeRecHits(*theElectron.electronCluster());
0072   if (debug_)
0073     std::cout << " Stored the hits belonging to the electronCluster " << std::endl;
0074
0075   // initial computation, no radius cut, but halo hits not taken
0076   if (debug_)
0077     std::cout << " Calling PCA initial computation" << std::endl;
0078   pcaHelper_.pcaInitialComputation();
0079   // first computation within cylinder, halo hits included
0080   pcaHelper_.computePCA(radius);
0081   // second computation within cylinder, halo hits included
0082   pcaHelper_.computePCA(radius);
0083   pcaHelper_.computeShowerWidth(radius);
0084   isoHelper_.produceHGCalIso(theElectron.electronCluster());
0085 }