CMS 3D CMS Logo

EgammaEcalIsolation.cc

Go to the documentation of this file.
00001 //*****************************************************************************
00002 // File:      EgammaEcalIsolation.cc
00003 // ----------------------------------------------------------------------------
00004 // OrigAuth:  Gilles De Lentdecker
00005 // Institute: IIHE-ULB
00006 //=============================================================================
00007 //*****************************************************************************
00008 
00009 //C++ includes
00010 #include <vector>
00011 #include <functional>
00012 
00013 //ROOT includes
00014 #include <Math/VectorUtil.h>
00015 
00016 //CMSSW includes
00017 #include "RecoEgamma/EgammaIsolationAlgos/interface/EgammaEcalIsolation.h"
00018 #include "DataFormats/RecoCandidate/interface/RecoCandidate.h"
00019 
00020 using namespace ROOT::Math::VectorUtil ;
00021 
00022 
00023 EgammaEcalIsolation::EgammaEcalIsolation (double extRadius,
00024                                           double etLow,
00025                                           const reco::BasicClusterCollection* basicClusterCollection,
00026                                           const reco::SuperClusterCollection* superClusterCollection):
00027   etMin(etLow),
00028   conesize(extRadius),
00029   basicClusterCollection_(basicClusterCollection),
00030   superClusterCollection_(superClusterCollection)
00031 {
00032 }
00033 
00034 EgammaEcalIsolation::~EgammaEcalIsolation(){}
00035 
00036 double EgammaEcalIsolation::getEcalEtSum(const reco::Candidate* candidate){
00037 
00038   
00039   double ecalIsol=0.;
00040   reco::SuperClusterRef sc = candidate->get<reco::SuperClusterRef>();
00041   math::XYZVector position(sc.get()->position().x(),
00042                    sc.get()->position().y(),
00043                    sc.get()->position().z());
00044   
00045   // match the photon hybrid supercluster with those with Algo==0 (island)
00046   double delta1=1000.;
00047   double deltacur=1000.;
00048   const reco::SuperCluster *matchedsupercluster=0;
00049   bool MATCHEDSC = false;
00050   
00051   for(reco::SuperClusterCollection::const_iterator scItr = superClusterCollection_->begin(); scItr != superClusterCollection_->end(); ++scItr){
00052     
00053     const reco::SuperCluster *supercluster = &(*scItr);
00054  
00055     math::XYZVector currentPosition(supercluster->position().x(),
00056                      supercluster->position().y(),
00057                      supercluster->position().z());
00058  
00059     
00060     if(supercluster->seed()->algo() == 0){
00061       deltacur = DeltaR(currentPosition, position); 
00062       
00063       if (deltacur < delta1) {
00064         delta1=deltacur;
00065         matchedsupercluster = supercluster;
00066         MATCHEDSC = true;
00067       }
00068     }
00069   }
00070 
00071   const reco::BasicCluster *cluster= 0;
00072   
00073   //loop over basic clusters
00074   for(reco::BasicClusterCollection::const_iterator cItr = basicClusterCollection_->begin(); cItr != basicClusterCollection_->end(); ++cItr){
00075  
00076     cluster = &(*cItr);
00077     double ebc_bcchi2 = cluster->chi2();
00078     int   ebc_bcalgo = cluster->algo();
00079     double ebc_bce    = cluster->energy();
00080     double ebc_bceta  = cluster->eta();
00081     double ebc_bcet   = ebc_bce*sin(2*atan(exp(ebc_bceta)));
00082     double newDelta = 0.;
00083 
00084 
00085     if (ebc_bcet > etMin && ebc_bcalgo == 0) {
00086       if (ebc_bcchi2 < 30.) {
00087         
00088         if(MATCHEDSC){
00089           bool inSuperCluster = false;
00090 
00091           reco::basicCluster_iterator theEclust = matchedsupercluster->clustersBegin();
00092           // loop over the basic clusters of the matched supercluster
00093           for(;theEclust != matchedsupercluster->clustersEnd();
00094               theEclust++) {
00095             if ((**theEclust) ==  (*cluster) ) inSuperCluster = true;
00096           }
00097           if (!inSuperCluster) {
00098 
00099             math::XYZVector basicClusterPosition(cluster->position().x(),
00100                                                   cluster->position().y(),
00101                                                   cluster->position().z());
00102             newDelta=DeltaR(basicClusterPosition,position);
00103             if(newDelta < conesize) {
00104               ecalIsol+=ebc_bcet;
00105             }
00106           }
00107         }
00108       } // matches ebc_bcchi2
00109     } // matches ebc_bcet && ebc_bcalgo
00110 
00111   }
00112   
00113   //  std::cout << "Will return ecalIsol = " << ecalIsol << std::endl; 
00114   return ecalIsol;
00115   
00116 }

Generated on Tue Jun 9 17:43:25 2009 for CMSSW by  doxygen 1.5.4