HiEgammaAlgos/plugins/HiEgammaSCCorrectionMaker.cc

0001 #include "RecoHI/HiEgammaAlgos/plugins/HiEgammaSCCorrectionMaker.h"
0002 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
0003 #include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
0004
0005 #include "FWCore/PluginManager/interface/ModuleDef.h"
0006 #include "FWCore/Framework/interface/MakerMacros.h"
0007
0008 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0009 #include "FWCore/Utilities/interface/Exception.h"
0010 #include "DataFormats/Common/interface/Handle.h"
0011
0012 #include "Geometry/Records/interface/CaloGeometryRecord.h"
0013 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0014 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
0015 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
0016 #include "Geometry/CaloTopology/interface/EcalBarrelTopology.h"
0017 #include "Geometry/CaloTopology/interface/EcalEndcapTopology.h"
0018 #include "Geometry/CaloTopology/interface/EcalPreshowerTopology.h"
0019 #include "Geometry/CaloTopology/interface/CaloTopology.h"
0020 #include "Geometry/Records/interface/CaloTopologyRecord.h"
0021
0022 #include <string>
0023
0024 HiEgammaSCCorrectionMaker::HiEgammaSCCorrectionMaker(const edm::ParameterSet& ps)
0025     : geoToken_(esConsumes()), topologyToken_(esConsumes()) {
0026   // The verbosity level
0027   std::string debugString = ps.getParameter<std::string>("VerbosityLevel");
0028   if (debugString == "DEBUG")
0029     verbosity_ = HiEgammaSCEnergyCorrectionAlgo::pDEBUG;
0030   else if (debugString == "INFO")
0031     verbosity_ = HiEgammaSCEnergyCorrectionAlgo::pINFO;
0032   else
0033     verbosity_ = HiEgammaSCEnergyCorrectionAlgo::pERROR;
0034
0035   // the input producers
0036   rHInputProducerTag_ = ps.getParameter<edm::InputTag>("recHitProducer");
0037   sCInputProducerTag_ = ps.getParameter<edm::InputTag>("rawSuperClusterProducer");
0038   rHInputProducer_ = consumes<EcalRecHitCollection>(rHInputProducerTag_);
0039   sCInputProducer_ = consumes<reco::SuperClusterCollection>(sCInputProducerTag_);
0040   std::string sCAlgo_str = ps.getParameter<std::string>("superClusterAlgo");
0041
0042   // determine which BasicCluster algo we are correcting for
0043   //And obtain forrection parameters form cfg file
0044   edm::ParameterSet fCorrPset;
0045   if (sCAlgo_str == "Hybrid") {
0046     sCAlgo_ = reco::CaloCluster::hybrid;
0047     fCorrPset = ps.getParameter<edm::ParameterSet>("hyb_fCorrPset");
0048   } else if (sCAlgo_str == "Island") {
0049     sCAlgo_ = reco::CaloCluster::island;
0050     fCorrPset = ps.getParameter<edm::ParameterSet>("isl_fCorrPset");
0051   } else if (sCAlgo_str == "DynamicHybrid") {
0052     sCAlgo_ = reco::CaloCluster::dynamicHybrid;
0053     fCorrPset = ps.getParameter<edm::ParameterSet>("dyn_fCorrPset");
0054   } else if (sCAlgo_str == "Multi5x5") {
0055     sCAlgo_ = reco::CaloCluster::multi5x5;
0056     fCorrPset = ps.getParameter<edm::ParameterSet>("fix_fCorrPset");
0057   } else {
0058     edm::LogError("HiEgammaSCCorrectionMakerError")
0059         << "Error! SuperClusterAlgo in config file must be Hybrid or Island: " << sCAlgo_str
0060         << "  Using Hybrid by default";
0061     sCAlgo_ = reco::CaloCluster::hybrid;
0062   }
0063
0064   // set correction algo parameters
0065   applyEnergyCorrection_ = ps.getParameter<bool>("applyEnergyCorrection");
0066   sigmaElectronicNoise_ = ps.getParameter<double>("sigmaElectronicNoise");
0067
0068   etThresh_ = ps.getParameter<double>("etThresh");
0069
0070   // set the producer parameters
0071   outputCollection_ = ps.getParameter<std::string>("corectedSuperClusterCollection");
0072   produces<reco::SuperClusterCollection>(outputCollection_);
0073
0074   // instanciate the correction algo object
0075   energyCorrector_ = std::make_unique<HiEgammaSCEnergyCorrectionAlgo>(sigmaElectronicNoise_, fCorrPset, verbosity_);
0076 }
0077
0078 HiEgammaSCCorrectionMaker::~HiEgammaSCCorrectionMaker() = default;
0079
0080 void HiEgammaSCCorrectionMaker::produce(edm::Event& evt, const edm::EventSetup& es) {
0081   using namespace edm;
0082
0083   // get the collection geometry:
0084   const CaloGeometry& geometry = es.getData(geoToken_);
0085   const CaloSubdetectorGeometry* geometry_p;
0086
0087   const CaloTopology& topology = es.getData(topologyToken_);
0088
0089   std::string rHInputCollection = rHInputProducerTag_.instance();
0090   if (rHInputCollection == "EcalRecHitsEB") {
0091     geometry_p = geometry.getSubdetectorGeometry(DetId::Ecal, EcalBarrel);
0092   } else if (rHInputCollection == "EcalRecHitsEE") {
0093     geometry_p = geometry.getSubdetectorGeometry(DetId::Ecal, EcalEndcap);
0094   } else if (rHInputCollection == "EcalRecHitsPS") {
0095     geometry_p = geometry.getSubdetectorGeometry(DetId::Ecal, EcalPreshower);
0096   } else {
0097     std::string str =
0098         "\n\nSCCorrectionMaker encountered invalied ecalhitcollection type: " + rHInputCollection + ".\n\n";
0099     throw(std::runtime_error(str.c_str()));
0100   }
0101
0102   // Get raw SuperClusters from the event
0103   Handle<reco::SuperClusterCollection> pRawSuperClusters;
0104   try {
0105     evt.getByToken(sCInputProducer_, pRawSuperClusters);
0106   } catch (cms::Exception& ex) {
0107     edm::LogError("HiEgammaSCCorrectionMakerError")
0108         << "Error! can't get the rawSuperClusters " << sCInputProducerTag_.label();
0109   }
0110
0111   // Get the RecHits from the event
0112   Handle<EcalRecHitCollection> pRecHits;
0113   try {
0114     evt.getByToken(rHInputProducer_, pRecHits);
0115   } catch (cms::Exception& ex) {
0116     edm::LogError("HiEgammaSCCorrectionMakerError") << "Error! can't get the RecHits " << rHInputProducerTag_.label();
0117   }
0118
0119   // Create a pointer to the RecHits and raw SuperClusters
0120   const EcalRecHitCollection* hitCollection = pRecHits.product();
0121   const reco::SuperClusterCollection* rawClusters = pRawSuperClusters.product();
0122
0123   // Define a collection of corrected SuperClusters to put back into the event
0124   auto corrClusters = std::make_unique<reco::SuperClusterCollection>();
0125
0126   //  Loop over raw clusters and make corrected ones
0127   reco::SuperClusterCollection::const_iterator aClus;
0128   for (aClus = rawClusters->begin(); aClus != rawClusters->end(); aClus++) {
0129     reco::SuperCluster newClus;
0130     if (applyEnergyCorrection_)
0131       newClus = energyCorrector_->applyCorrection(*aClus, *hitCollection, sCAlgo_, *geometry_p, topology);
0132     else
0133       newClus = *aClus;
0134
0135     if (newClus.energy() * sin(newClus.position().theta()) > etThresh_) {
0136       corrClusters->push_back(newClus);
0137     }
0138   }
0139   // Put collection of corrected SuperClusters into the event
0140   evt.put(std::move(corrClusters), outputCollection_);
0141 }
0142
0143 DEFINE_FWK_MODULE(HiEgammaSCCorrectionMaker);