Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-09-09 02:18:52

0001 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0002 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterPUCleaningTools.h"
0003 
0004 #include "Geometry/CaloTopology/interface/CaloTopology.h"
0005 #include "Geometry/CaloGeometry/interface/CaloSubdetectorGeometry.h"
0006 
0007 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0008 
0009 #include "RecoEcal/EgammaCoreTools/interface/EcalClusterTools.h"
0010 
0011 #include "RecoEcal/EgammaCoreTools/interface/SuperClusterShapeAlgo.h"
0012 
0013 EcalClusterPUCleaningTools::EcalClusterPUCleaningTools(edm::ConsumesCollector &cc,
0014                                                        const edm::InputTag &redEBRecHits,
0015                                                        const edm::InputTag &redEERecHits)
0016     : pEBRecHitsToken_(cc.consumes<EcalRecHitCollection>(redEBRecHits)),
0017       pEERecHitsToken_(cc.consumes<EcalRecHitCollection>(redEERecHits)),
0018       geometryToken_(cc.esConsumes()) {}
0019 
0020 EcalClusterPUCleaningTools::~EcalClusterPUCleaningTools() {}
0021 
0022 void EcalClusterPUCleaningTools::getEBRecHits(const edm::Event &ev) {
0023   edm::Handle<EcalRecHitCollection> pEBRecHits;
0024   ev.getByToken(pEBRecHitsToken_, pEBRecHits);
0025   ebRecHits_ = pEBRecHits.product();
0026 }
0027 
0028 void EcalClusterPUCleaningTools::getEERecHits(const edm::Event &ev) {
0029   edm::Handle<EcalRecHitCollection> pEERecHits;
0030   ev.getByToken(pEERecHitsToken_, pEERecHits);
0031   eeRecHits_ = pEERecHits.product();
0032 }
0033 
0034 reco::SuperCluster EcalClusterPUCleaningTools::CleanedSuperCluster(float xi,
0035                                                                    const reco::SuperCluster &scluster,
0036                                                                    const edm::Event &ev,
0037                                                                    const edm::EventSetup &es) {
0038   // get the geometry
0039   const auto &geometry = es.getData(geometryToken_);
0040 
0041   // get the RecHits
0042   getEBRecHits(ev);
0043   getEERecHits(ev);
0044 
0045   // seed basic cluster of initial SC: this will remain in the cleaned SC, by construction
0046   const reco::CaloClusterPtr &seed = scluster.seed();
0047 
0048   // this should be replaced by the 5x5 around the seed; a good approx of how E_seed is defined
0049   float seedBCEnergy = (scluster.seed())->energy();
0050   float eSeed = 0.35;  // standard eSeed in EB ; see CMS IN-2010/008
0051   int numBcRemoved = 0;
0052 
0053   double ClusterE = 0;  //Sum of cluster energies for supercluster.
0054   //Holders for position of this supercluster.
0055   double posX = 0;
0056   double posY = 0;
0057   double posZ = 0;
0058 
0059   reco::CaloClusterPtrVector thissc;
0060 
0061   // looping on basic clusters within the Supercluster
0062   for (reco::CaloCluster_iterator bcIt = scluster.clustersBegin(); bcIt != scluster.clustersEnd(); bcIt++) {
0063     // E_seed is an Et  selection on 5x1 dominoes (around which BC's are built), see CMS IN-2010/008
0064     // here such selection is implemented on the full BC around it
0065     if ((*bcIt)->energy() >
0066         sqrt(eSeed * eSeed + xi * xi * seedBCEnergy * seedBCEnergy / cosh((*bcIt)->eta()) / cosh((*bcIt)->eta()))) {
0067       ;
0068     }  // the sum only of the BC's that pass the Esee selection
0069     else {
0070       numBcRemoved++;
0071       continue;
0072     }  // count how many basic cluster get removed by the cleaning
0073 
0074     // if BC passes dynamic selection, include it in the 'cleaned' supercluster
0075     thissc.push_back((*bcIt));
0076     // cumulate energy and position of the cleaned supercluster
0077     ClusterE += (*bcIt)->energy();
0078     posX += (*bcIt)->energy() * (*bcIt)->position().X();
0079     posY += (*bcIt)->energy() * (*bcIt)->position().Y();
0080     posZ += (*bcIt)->energy() * (*bcIt)->position().Z();
0081 
0082   }  // loop on basic clusters of the original SC
0083 
0084   posX /= ClusterE;
0085   posY /= ClusterE;
0086   posZ /= ClusterE;
0087 
0088   // make temporary 'cleaned' supercluster in oder to compute the covariances
0089   double Epreshower = scluster.preshowerEnergy();
0090   double phiWidth = 0.;
0091   double etaWidth = 0.;
0092   reco::SuperCluster suCltmp(ClusterE, math::XYZPoint(posX, posY, posZ), seed, thissc, Epreshower, phiWidth, etaWidth);
0093 
0094   // construct cluster shape to compute ieta and iphi covariances of the SC
0095   const CaloSubdetectorGeometry *geometry_p = nullptr;
0096   if (seed->seed().det() == DetId::Ecal && seed->seed().subdetId() == EcalBarrel) {
0097     geometry_p = geometry.getSubdetectorGeometry(DetId::Ecal, EcalBarrel);
0098     SuperClusterShapeAlgo SCShape(ebRecHits_, geometry_p);
0099     SCShape.Calculate_Covariances(suCltmp);
0100     phiWidth = SCShape.phiWidth();
0101     etaWidth = SCShape.etaWidth();
0102   } else if (seed->seed().det() == DetId::Ecal && seed->seed().subdetId() == EcalEndcap) {
0103     geometry_p = geometry.getSubdetectorGeometry(DetId::Ecal, EcalEndcap);
0104     SuperClusterShapeAlgo SCShape(eeRecHits_, geometry_p);
0105     SCShape.Calculate_Covariances(suCltmp);
0106     phiWidth = SCShape.phiWidth();
0107     etaWidth = SCShape.etaWidth();
0108   } else {
0109     edm::LogError("SeedError")
0110         << "The seed crystal of this SC is neither in EB nor in EE. This is a problem. Bailing out";
0111     assert(-1);
0112   }
0113 
0114   // return the cleaned supercluster SCluster, with covariances updated
0115   return reco::SuperCluster(ClusterE, math::XYZPoint(posX, posY, posZ), seed, thissc, Epreshower, phiWidth, etaWidth);
0116 }