CMS 3D CMS Logo

/data/doxygen/doxygen-1.7.3/gen/CMSSW_4_2_8/src/RecoEcal/EgammaClusterProducers/src/EgammaSCCorrectionMaker.cc

Go to the documentation of this file.
00001 #include "RecoEcal/EgammaClusterProducers/interface/EgammaSCCorrectionMaker.h"
00002 #include "DataFormats/EgammaReco/interface/SuperCluster.h"
00003 #include "DataFormats/EgammaReco/interface/SuperClusterFwd.h"
00004 #include "DataFormats/EcalRecHit/interface/EcalRecHit.h"
00005 #include "DataFormats/EcalRecHit/interface/EcalRecHitCollections.h"
00006 
00007 #include "FWCore/MessageLogger/interface/MessageLogger.h"
00008 #include "FWCore/Utilities/interface/Exception.h"
00009 #include "DataFormats/Common/interface/Handle.h"
00010 
00011 #include "Geometry/Records/interface/CaloGeometryRecord.h"
00012 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
00013 #include "Geometry/CaloGeometry/interface/CaloCellGeometry.h"
00014 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
00015 #include "Geometry/CaloTopology/interface/EcalBarrelTopology.h"
00016 #include "Geometry/CaloTopology/interface/EcalEndcapTopology.h"
00017 #include "Geometry/CaloTopology/interface/EcalPreshowerTopology.h"
00018 
00019 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterFunctionBaseClass.h" 
00020 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterFunctionFactory.h" 
00021 
00022 #include <string>
00023 
00024 EgammaSCCorrectionMaker::EgammaSCCorrectionMaker(const edm::ParameterSet& ps)
00025 {
00026  
00027   // The verbosity level
00028   std::string debugString = ps.getParameter<std::string>("VerbosityLevel");
00029   if      (debugString == "DEBUG")   verbosity_ = EgammaSCEnergyCorrectionAlgo::pDEBUG;
00030   else if (debugString == "INFO")    verbosity_ = EgammaSCEnergyCorrectionAlgo::pINFO;
00031   else                               verbosity_ = EgammaSCEnergyCorrectionAlgo::pERROR;
00032 
00033   // the input producers
00034   rHInputProducer_ = ps.getParameter<edm::InputTag>("recHitProducer");
00035   sCInputProducer_ = ps.getParameter<edm::InputTag>("rawSuperClusterProducer");
00036   std::string sCAlgo_str = ps.getParameter<std::string>("superClusterAlgo");
00037 
00038   // determine which BasicCluster algo we are correcting for
00039   //And obtain forrection parameters form cfg file
00040   edm::ParameterSet fCorrPset;
00041   if (sCAlgo_str=="Hybrid") {
00042     sCAlgo_= reco::CaloCluster::hybrid;
00043     fCorrPset = ps.getParameter<edm::ParameterSet>("hyb_fCorrPset"); 
00044   } else if (sCAlgo_str=="Island") {
00045     sCAlgo_= reco::CaloCluster::island;
00046     fCorrPset = ps.getParameter<edm::ParameterSet>("isl_fCorrPset");
00047   } else if (sCAlgo_str=="DynamicHybrid") {
00048     sCAlgo_ = reco::CaloCluster::dynamicHybrid;
00049     fCorrPset = ps.getParameter<edm::ParameterSet>("dyn_fCorrPset"); 
00050   } else if (sCAlgo_str=="Multi5x5") {
00051     sCAlgo_ = reco::CaloCluster::multi5x5;
00052     fCorrPset = ps.getParameter<edm::ParameterSet>("fix_fCorrPset");
00053   } else {
00054     edm::LogError("EgammaSCCorrectionMakerError") 
00055       << "Error! SuperClusterAlgo in config file must be Hybrid or Island: " 
00056       << sCAlgo_str << "  Using Hybrid by default";
00057     sCAlgo_=reco::CaloCluster::hybrid;
00058   }
00059   
00060   // set correction algo parameters
00061   applyEnergyCorrection_ = ps.getParameter<bool>("applyEnergyCorrection");
00062   applyCrackCorrection_ = ps.getParameter<bool>("applyCrackCorrection");
00063 
00064   sigmaElectronicNoise_ =  ps.getParameter<double>("sigmaElectronicNoise");
00065 
00066   etThresh_ =  ps.getParameter<double>("etThresh");
00067 
00068   // set the producer parameters
00069   outputCollection_ = ps.getParameter<std::string>("corectedSuperClusterCollection");
00070   produces<reco::SuperClusterCollection>(outputCollection_);
00071 
00072   // instanciate the correction algo object
00073   energyCorrector_ = new EgammaSCEnergyCorrectionAlgo(sigmaElectronicNoise_, sCAlgo_, fCorrPset, verbosity_);
00074   
00075   // energy correction class
00076   if (applyEnergyCorrection_ )
00077     energyCorrectionFunction_ = EcalClusterFunctionFactory::get()->create("EcalClusterEnergyCorrection", ps);
00078   else
00079     energyCorrectionFunction_=0;
00080   if (applyCrackCorrection_ )
00081     crackCorrectionFunction_ = EcalClusterFunctionFactory::get()->create("EcalClusterCrackCorrection", ps);
00082   else
00083     crackCorrectionFunction_=0;
00084 
00085   
00086 }
00087 
00088 EgammaSCCorrectionMaker::~EgammaSCCorrectionMaker()
00089 {
00090   if (energyCorrectionFunction_) delete energyCorrectionFunction_;
00091   if (crackCorrectionFunction_) delete crackCorrectionFunction_;
00092   delete energyCorrector_;
00093 }
00094 
00095 void
00096 EgammaSCCorrectionMaker::produce(edm::Event& evt, const edm::EventSetup& es)
00097 {
00098   using namespace edm;
00099 
00100   // initialize energy correction class
00101   if(applyEnergyCorrection_) 
00102     energyCorrectionFunction_->init(es);
00103 
00104   // initialize energy correction class
00105   if(applyCrackCorrection_) 
00106     crackCorrectionFunction_->init(es);
00107   
00108 
00109   // get the collection geometry:
00110   edm::ESHandle<CaloGeometry> geoHandle;
00111   es.get<CaloGeometryRecord>().get(geoHandle);
00112   const CaloGeometry& geometry = *geoHandle;
00113   const CaloSubdetectorGeometry *geometry_p;
00114 
00115   std::string rHInputCollection = rHInputProducer_.instance();
00116   if(rHInputCollection == "EcalRecHitsEB") {
00117     geometry_p = geometry.getSubdetectorGeometry(DetId::Ecal, EcalBarrel);
00118   } else if(rHInputCollection == "EcalRecHitsEE") {
00119     geometry_p = geometry.getSubdetectorGeometry(DetId::Ecal, EcalEndcap);
00120   } else if(rHInputCollection == "EcalRecHitsPS") {
00121     geometry_p = geometry.getSubdetectorGeometry(DetId::Ecal, EcalPreshower);
00122   } else {
00123           std::string str = "\n\nSCCorrectionMaker encountered invalied ecalhitcollection type: " + rHInputCollection + ".\n\n";
00124           throw(std::runtime_error( str.c_str() ));
00125   }
00126   
00127   // Get raw SuperClusters from the event    
00128   Handle<reco::SuperClusterCollection> pRawSuperClusters;
00129   try { 
00130     evt.getByLabel(sCInputProducer_, pRawSuperClusters);
00131   } catch ( cms::Exception& ex ) {
00132     edm::LogError("EgammaSCCorrectionMakerError") 
00133       << "Error! can't get the rawSuperClusters " 
00134       << sCInputProducer_.label() ;
00135   }    
00136   
00137   // Get the RecHits from the event
00138   Handle<EcalRecHitCollection> pRecHits;
00139   try { 
00140     evt.getByLabel(rHInputProducer_, pRecHits);
00141   } catch ( cms::Exception& ex ) {
00142     edm::LogError("EgammaSCCorrectionMakerError") 
00143       << "Error! can't get the RecHits " 
00144       << rHInputProducer_.label();
00145   }    
00146   
00147   // Create a pointer to the RecHits and raw SuperClusters
00148   const EcalRecHitCollection *hitCollection = pRecHits.product();
00149   const reco::SuperClusterCollection *rawClusters = pRawSuperClusters.product();
00150    
00151   // Define a collection of corrected SuperClusters to put back into the event
00152   std::auto_ptr<reco::SuperClusterCollection> corrClusters(new reco::SuperClusterCollection);
00153   
00154   //  Loop over raw clusters and make corrected ones
00155   reco::SuperClusterCollection::const_iterator aClus;
00156   for(aClus = rawClusters->begin(); aClus != rawClusters->end(); aClus++)
00157     {
00158       reco::SuperCluster enecorrClus,crackcorrClus;
00159 
00160       if(applyEnergyCorrection_) 
00161         enecorrClus = energyCorrector_->applyCorrection(*aClus, *hitCollection, sCAlgo_, geometry_p, energyCorrectionFunction_);
00162       else
00163         enecorrClus=*aClus;
00164 
00165       if(applyCrackCorrection_)
00166         crackcorrClus=energyCorrector_->applyCrackCorrection(enecorrClus,crackCorrectionFunction_);
00167       else 
00168         crackcorrClus=enecorrClus;
00169 
00170       if(crackcorrClus.energy()*sin(crackcorrClus.position().theta())>etThresh_) {
00171         
00172         corrClusters->push_back(crackcorrClus);
00173       }
00174     }
00175   // Put collection of corrected SuperClusters into the event
00176   evt.put(corrClusters, outputCollection_);   
00177   
00178 }
00179