CMS 3D CMS Logo

/data/refman/pasoursint/CMSSW_5_3_3/src/RecoEcal/EgammaCoreTools/src/EcalClusterPUCleaningTools.cc

Go to the documentation of this file.
00001 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterPUCleaningTools.h"
00002 
00003 #include "Geometry/CaloGeometry/interface/CaloGeometry.h"
00004 #include "Geometry/CaloTopology/interface/CaloTopology.h"
00005 #include "Geometry/Records/interface/CaloGeometryRecord.h"
00006 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
00007 
00008 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
00009 
00010 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterTools.h"
00011 
00012 #include "RecoEcal/EgammaCoreTools/interface/SuperClusterShapeAlgo.h"
00013 
00014 EcalClusterPUCleaningTools::EcalClusterPUCleaningTools( const edm::Event &ev, const edm::EventSetup &es, edm::InputTag redEBRecHits, edm::InputTag redEERecHits )
00015 {
00016   getGeometry( es );
00017   getEBRecHits( ev, redEBRecHits );
00018   getEERecHits( ev, redEERecHits );
00019 }
00020 
00021 
00022 
00023 EcalClusterPUCleaningTools::~EcalClusterPUCleaningTools()
00024 {
00025 }
00026 
00027 
00028 void EcalClusterPUCleaningTools::getGeometry( const edm::EventSetup &es )
00029 {
00030   edm::ESHandle<CaloGeometry> pGeometry;
00031   es.get<CaloGeometryRecord>().get(pGeometry);
00032   geometry_ = pGeometry.product();
00033 }
00034 
00035 void EcalClusterPUCleaningTools::getEBRecHits( const edm::Event &ev, edm::InputTag redEBRecHits )
00036 {
00037   edm::Handle< EcalRecHitCollection > pEBRecHits;
00038   ev.getByLabel( redEBRecHits, pEBRecHits );
00039   ebRecHits_ = pEBRecHits.product();
00040 }
00041 
00042 
00043 
00044 void EcalClusterPUCleaningTools::getEERecHits( const edm::Event &ev, edm::InputTag redEERecHits )
00045 {
00046   edm::Handle< EcalRecHitCollection > pEERecHits;
00047   ev.getByLabel( redEERecHits, pEERecHits );
00048   eeRecHits_ = pEERecHits.product();
00049 }
00050 
00051 
00052 reco::SuperCluster EcalClusterPUCleaningTools::CleanedSuperCluster(float xi, const reco::SuperCluster &scluster, const edm::Event &ev ){
00053   
00054   //std::cout << "\nEcalClusterPUCleaningTools::CleanedSuperCluster called, this will give you back a cleaned supercluster" << std::endl;
00055 
00056   // seed basic cluster of initial SC: this will remain in the cleaned SC, by construction
00057   reco::CaloClusterPtr seed = scluster.seed();
00058 
00059   float seedBCEnergy      = (scluster.seed())->energy(); // this should be replaced by the 5x5 around the seed; a good approx of how E_seed is defined 
00060   float eSeed             = 0.35;                        // standard eSeed in EB ; see CMS IN-2010/008
00061   int   numBcRemoved      = 0;  
00062 
00063   double ClusterE = 0; //Sum of cluster energies for supercluster.
00064   //Holders for position of this supercluster.
00065   double posX=0;
00066   double posY=0;
00067   double posZ=0;
00068 
00069   reco::CaloClusterPtrVector thissc;
00070 
00071   // looping on basic clusters within the Supercluster
00072   for(reco::CaloCluster_iterator bcIt = scluster.clustersBegin(); bcIt!=scluster.clustersEnd(); bcIt++)
00073     {
00074       // E_seed is an Et  selection on 5x1 dominoes (around which BC's are built), see CMS IN-2010/008
00075       // here such selection is implemented on the full BC around it 
00076       if( (*bcIt)->energy() > sqrt( eSeed*eSeed + xi*xi*seedBCEnergy*seedBCEnergy/cosh((*bcIt)->eta())/cosh((*bcIt)->eta())  ) ) 
00077         {
00078           ;
00079         }// the sum only of the BC's that pass the Esee selection
00080       else{
00081         numBcRemoved++;
00082         continue;
00083       }// count how many basic cluster get removed by the cleaning
00084       
00085       // if BC passes dynamic selection, include it in the 'cleaned' supercluster
00086       thissc.push_back( (*bcIt) );
00087       // cumulate energy and position of the cleaned supercluster
00088       ClusterE += (*bcIt)->energy();
00089       posX += (*bcIt)->energy() * (*bcIt)->position().X();
00090       posY += (*bcIt)->energy() * (*bcIt)->position().Y();
00091       posZ += (*bcIt)->energy() * (*bcIt)->position().Z();
00092       
00093     }// loop on basic clusters of the original SC
00094 
00095   posX /= ClusterE;
00096   posY /= ClusterE;
00097   posZ /= ClusterE;
00098 
00099   // make temporary 'cleaned' supercluster in oder to compute the covariances 
00100   double Epreshower=scluster.preshowerEnergy();
00101   double phiWidth=0.; 
00102   double etaWidth=0.; 
00103   reco::SuperCluster suCltmp(ClusterE, math::XYZPoint(posX, posY, posZ), seed, thissc, Epreshower, phiWidth, etaWidth);
00104   
00105 
00106   // construct cluster shape to compute ieta and iphi covariances of the SC 
00107   const CaloSubdetectorGeometry *geometry_p=0;
00108   if (seed->seed().det() == DetId::Ecal && seed->seed().subdetId() == EcalBarrel){
00109     geometry_p = geometry_->getSubdetectorGeometry(DetId::Ecal, EcalBarrel);
00110     SuperClusterShapeAlgo  SCShape(   ebRecHits_ , geometry_p); 
00111     SCShape.Calculate_Covariances( suCltmp );
00112     phiWidth= SCShape.phiWidth(); 
00113     etaWidth= SCShape.etaWidth(); 
00114   }
00115   else if (seed->seed().det() == DetId::Ecal && seed->seed().subdetId() == EcalEndcap){
00116     geometry_p = geometry_->getSubdetectorGeometry(DetId::Ecal, EcalEndcap);
00117     SuperClusterShapeAlgo  SCShape(   eeRecHits_ , geometry_p); 
00118     SCShape.Calculate_Covariances( suCltmp );
00119     phiWidth= SCShape.phiWidth(); 
00120     etaWidth= SCShape.etaWidth(); 
00121   }
00122   else {
00123     std::cout << "The seed crystal of this SC is neither in EB nor in EE. This is a problem. Bailing out " << std::endl;
00124     assert(-1);
00125   }
00126   
00127   // return the cleaned supercluster SCluster, with covariances updated 
00128   return reco::SuperCluster(ClusterE, math::XYZPoint(posX, posY, posZ), seed, thissc, Epreshower, phiWidth, etaWidth);
00129 
00130 }