Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2021-02-14 14:25:27

0001 // -*- C++ -*-
0002 //
0003 // Package:    METAlgorithms
0004 // Class:      CaloSpecificAlgo
0005 //
0006 // Original Author:  R. Cavanaugh (taken from F.Ratnikov, UMd)
0007 //         Created:  June 6, 2006
0008 //
0009 //
0010 #include "RecoMET/METAlgorithms/interface/CaloSpecificAlgo.h"
0011 #include "DataFormats/HcalDetId/interface/HcalDetId.h"
0012 #include "DataFormats/EcalDetId/interface/EcalSubdetector.h"
0013 #include "DataFormats/CaloTowers/interface/CaloTower.h"
0014 
0015 #include <iostream>
0016 
0017 using namespace reco;
0018 
0019 //____________________________________________________________________________||
0020 // This algorithm adds calorimeter specific global event information to
0021 // the MET object which may be useful/needed for MET Data Quality Monitoring
0022 // and MET cleaning.
0023 //____________________________________________________________________________||
0024 
0025 //____________________________________________________________________________||
0026 reco::CaloMET CaloSpecificAlgo::addInfo(edm::Handle<edm::View<Candidate> > towers,
0027                                         const CommonMETData &umet,
0028                                         bool noHF,
0029                                         double globalThreshold) {
0030   SpecificCaloMETData specific;
0031   CommonMETData met = umet;
0032   double totalEt = 0.0;
0033   double totalEm = 0.0;
0034 
0035   double sumEtInpHF = 0.0;
0036   double sumEtInmHF = 0.0;
0037   double MExInpHF = 0.0;
0038   double MEyInpHF = 0.0;
0039   double MExInmHF = 0.0;
0040   double MEyInmHF = 0.0;
0041 
0042   for (edm::View<Candidate>::const_iterator towerCand = towers->begin(); towerCand != towers->end(); ++towerCand) {
0043     const CaloTower *calotower = dynamic_cast<const CaloTower *>(&(*towerCand));
0044     if (!calotower)
0045       continue;
0046     if (calotower->et() < globalThreshold)
0047       continue;
0048     update_totalEt_totalEm(totalEt, totalEm, calotower, noHF);
0049     update_MaxTowerEm_MaxTowerHad(specific.MaxEtInEmTowers, specific.MaxEtInHadTowers, calotower, noHF);
0050     update_EmEtInEB_EmEtInEE(specific.EmEtInEB, specific.EmEtInEE, calotower);
0051     update_HadEtInHB_HadEtInHE_HadEtInHO_HadEtInHF_EmEtInHF(specific.HadEtInHB,
0052                                                             specific.HadEtInHE,
0053                                                             specific.HadEtInHO,
0054                                                             specific.HadEtInHF,
0055                                                             specific.EmEtInHF,
0056                                                             calotower,
0057                                                             noHF);
0058     update_sumEtInpHF_MExInpHF_MEyInpHF_sumEtInmHF_MExInmHF_MEyInmHF(
0059         sumEtInpHF, MExInpHF, MEyInpHF, sumEtInmHF, MExInmHF, MEyInmHF, calotower);
0060   }
0061 
0062   double totalHad = totalEt - totalEm;
0063 
0064   if (noHF)
0065     remove_HF_from_MET(met, sumEtInpHF, MExInpHF, MEyInpHF, sumEtInmHF, MExInmHF, MEyInmHF);
0066 
0067   if (!noHF)
0068     add_MET_in_HF(specific, sumEtInpHF, MExInpHF, MEyInpHF, sumEtInmHF, MExInmHF, MEyInmHF);
0069 
0070   specific.EtFractionHadronic = (totalEt == 0.0) ? 0.0 : totalHad / totalEt;
0071   specific.EtFractionEm = (totalEt == 0.0) ? 0.0 : totalEm / totalEt;
0072 
0073   const LorentzVector p4(met.mex, met.mey, 0.0, met.met);
0074   const Point vtx(0.0, 0.0, 0.0);
0075   CaloMET caloMET(specific, met.sumet, p4, vtx);
0076   return caloMET;
0077 }
0078 
0079 //____________________________________________________________________________||
0080 void CaloSpecificAlgo::update_totalEt_totalEm(double &totalEt, double &totalEm, const CaloTower *calotower, bool noHF) {
0081   if (noHF) {
0082     DetId detIdHcal = find_DetId_of_HCAL_cell_in_constituent_of(calotower);
0083     if (!detIdHcal.null()) {
0084       HcalSubdetector subdet = HcalDetId(detIdHcal).subdet();
0085       if (subdet == HcalForward)
0086         return;
0087     }
0088   }
0089 
0090   totalEt += calotower->et();
0091   totalEm += calotower->emEt();
0092 }
0093 
0094 //____________________________________________________________________________||
0095 void CaloSpecificAlgo::update_MaxTowerEm_MaxTowerHad(float &MaxTowerEm,
0096                                                      float &MaxTowerHad,
0097                                                      const CaloTower *calotower,
0098                                                      bool noHF) {
0099   DetId detIdHcal = find_DetId_of_HCAL_cell_in_constituent_of(calotower);
0100   DetId detIdEcal = find_DetId_of_ECAL_cell_in_constituent_of(calotower);
0101 
0102   if (!detIdHcal.null()) {
0103     HcalSubdetector subdet = HcalDetId(detIdHcal).subdet();
0104     if (subdet == HcalBarrel || subdet == HcalOuter || subdet == HcalEndcap || (!noHF && subdet == HcalForward)) {
0105       if (calotower->hadEt() > MaxTowerHad)
0106         MaxTowerHad = calotower->hadEt();
0107       if (calotower->emEt() > MaxTowerEm)
0108         MaxTowerEm = calotower->emEt();
0109     }
0110   }
0111 
0112   if (!detIdEcal.null()) {
0113     if (calotower->emEt() > MaxTowerEm)
0114       MaxTowerEm = calotower->emEt();
0115   }
0116 }
0117 
0118 //____________________________________________________________________________||
0119 void CaloSpecificAlgo::update_EmEtInEB_EmEtInEE(float &EmEtInEB, float &EmEtInEE, const CaloTower *calotower) {
0120   DetId detIdEcal = find_DetId_of_ECAL_cell_in_constituent_of(calotower);
0121   if (detIdEcal.null())
0122     return;
0123 
0124   EcalSubdetector subdet = EcalSubdetector(detIdEcal.subdetId());
0125   if (subdet == EcalBarrel) {
0126     EmEtInEB += calotower->emEt();
0127   } else if (subdet == EcalEndcap) {
0128     EmEtInEE += calotower->emEt();
0129   }
0130 }
0131 
0132 //____________________________________________________________________________||
0133 void CaloSpecificAlgo::update_HadEtInHB_HadEtInHE_HadEtInHO_HadEtInHF_EmEtInHF(float &HadEtInHB,
0134                                                                                float &HadEtInHE,
0135                                                                                float &HadEtInHO,
0136                                                                                float &HadEtInHF,
0137                                                                                float &EmEtInHF,
0138                                                                                const CaloTower *calotower,
0139                                                                                bool noHF) {
0140   DetId detIdHcal = find_DetId_of_HCAL_cell_in_constituent_of(calotower);
0141   if (detIdHcal.null())
0142     return;
0143 
0144   HcalSubdetector subdet = HcalDetId(detIdHcal).subdet();
0145   if (subdet == HcalBarrel || subdet == HcalOuter) {
0146     HadEtInHB += calotower->hadEt();
0147     HadEtInHO += calotower->outerEt();
0148   }
0149 
0150   if (subdet == HcalEndcap) {
0151     HadEtInHE += calotower->hadEt();
0152   }
0153 
0154   if (subdet == HcalForward && !noHF) {
0155     HadEtInHF += calotower->hadEt();
0156     EmEtInHF += calotower->emEt();
0157   }
0158 }
0159 
0160 //____________________________________________________________________________||
0161 void CaloSpecificAlgo::update_sumEtInpHF_MExInpHF_MEyInpHF_sumEtInmHF_MExInmHF_MEyInmHF(double &sumEtInpHF,
0162                                                                                         double &MExInpHF,
0163                                                                                         double &MEyInpHF,
0164                                                                                         double &sumEtInmHF,
0165                                                                                         double &MExInmHF,
0166                                                                                         double &MEyInmHF,
0167                                                                                         const CaloTower *calotower) {
0168   DetId detIdHcal = find_DetId_of_HCAL_cell_in_constituent_of(calotower);
0169   if (detIdHcal.null())
0170     return;
0171 
0172   HcalSubdetector subdet = HcalDetId(detIdHcal).subdet();
0173   if (!(subdet == HcalForward))
0174     return;
0175 
0176   if (calotower->eta() >= 0) {
0177     sumEtInpHF += calotower->et();
0178     MExInpHF -= (calotower->et() * cos(calotower->phi()));
0179     MEyInpHF -= (calotower->et() * sin(calotower->phi()));
0180   } else {
0181     sumEtInmHF += calotower->et();
0182     MExInmHF -= (calotower->et() * cos(calotower->phi()));
0183     MEyInmHF -= (calotower->et() * sin(calotower->phi()));
0184   }
0185 }
0186 
0187 //____________________________________________________________________________||
0188 void CaloSpecificAlgo::remove_HF_from_MET(CommonMETData &met,
0189                                           double sumEtInpHF,
0190                                           double MExInpHF,
0191                                           double MEyInpHF,
0192                                           double sumEtInmHF,
0193                                           double MExInmHF,
0194                                           double MEyInmHF) {
0195   met.mex -= (MExInmHF + MExInpHF);
0196   met.mey -= (MEyInmHF + MEyInpHF);
0197   met.sumet -= (sumEtInpHF + sumEtInmHF);
0198   met.met = sqrt(met.mex * met.mex + met.mey * met.mey);
0199 }
0200 
0201 //____________________________________________________________________________||
0202 void CaloSpecificAlgo::add_MET_in_HF(SpecificCaloMETData &specific,
0203                                      double sumEtInpHF,
0204                                      double MExInpHF,
0205                                      double MEyInpHF,
0206                                      double sumEtInmHF,
0207                                      double MExInmHF,
0208                                      double MEyInmHF) {
0209   LorentzVector METpHF(MExInpHF, MEyInpHF, 0, sqrt(MExInpHF * MExInpHF + MEyInpHF * MEyInpHF));
0210   LorentzVector METmHF(MExInmHF, MEyInmHF, 0, sqrt(MExInmHF * MExInmHF + MEyInmHF * MEyInmHF));
0211   specific.CaloMETInpHF = METpHF.pt();
0212   specific.CaloMETInmHF = METmHF.pt();
0213   specific.CaloMETPhiInpHF = METpHF.Phi();
0214   specific.CaloMETPhiInmHF = METmHF.Phi();
0215   specific.CaloSETInpHF = sumEtInpHF;
0216   specific.CaloSETInmHF = sumEtInmHF;
0217 }
0218 
0219 //____________________________________________________________________________||
0220 DetId CaloSpecificAlgo::find_DetId_of_HCAL_cell_in_constituent_of(const CaloTower *calotower) {
0221   DetId ret;
0222   for (int cell = calotower->constituentsSize() - 1; cell >= 0; --cell) {
0223     DetId id = calotower->constituent(cell);
0224     if (id.det() == DetId::Hcal) {
0225       ret = id;
0226       break;
0227     }
0228   }
0229   return ret;
0230 }
0231 
0232 //____________________________________________________________________________||
0233 DetId CaloSpecificAlgo::find_DetId_of_ECAL_cell_in_constituent_of(const CaloTower *calotower) {
0234   DetId ret;
0235   for (int cell = calotower->constituentsSize() - 1; cell >= 0; --cell) {
0236     DetId id = calotower->constituent(cell);
0237     if (id.det() == DetId::Ecal) {
0238       ret = id;
0239       break;
0240     }
0241   }
0242   return ret;
0243 }
0244 
0245 //____________________________________________________________________________||