RecoTau/src/AntiElectronIDCut2.cc

0001 #include "RecoTauTag/RecoTau/interface/AntiElectronIDCut2.h"
0002
0003 #include "FWCore/Utilities/interface/Exception.h"
0004
0005 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidate.h"
0006 #include "DataFormats/ParticleFlowCandidate/interface/PFCandidateFwd.h"
0007 #include "DataFormats/TrackReco/interface/Track.h"
0008 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0009 #include <TMath.h>
0010
0011 AntiElectronIDCut2::AntiElectronIDCut2() {
0012   //default, keep all taus in cracks for HLT
0013   keepAllInEcalCrack_ = true;
0014   rejectAllInEcalCrack_ = false;
0015
0016   //Default Cuts to be applied
0017   Tau_applyCut_hcal3x3OverPLead_ = true;
0018   Tau_applyCut_leadPFChargedHadrEoP_ = true;
0019   Tau_applyCut_GammaEtaMom_ = true;
0020   Tau_applyCut_GammaPhiMom_ = false;
0021   Tau_applyCut_GammaEnFrac_ = false;
0022   Tau_applyCut_HLTSpecific_ = true;
0023
0024   verbosity_ = 0;
0025 }
0026
0027 AntiElectronIDCut2::~AntiElectronIDCut2() {}
0028
0029 double AntiElectronIDCut2::Discriminator(float TauPt,
0030                                          float TauEta,
0031                                          float TauLeadChargedPFCandPt,
0032                                          float TauLeadChargedPFCandEtaAtEcalEntrance,
0033                                          float TauLeadPFChargedHadrEoP,
0034                                          float TauHcal3x3OverPLead,
0035                                          float TauGammaEtaMom,
0036                                          float TauGammaPhiMom,
0037                                          float TauGammaEnFrac) {
0038   if (verbosity_) {
0039     std::cout << "Tau GammaEnFrac " << TauGammaEnFrac << std::endl;
0040     std::cout << "Tau GammaEtaMom " << TauGammaEtaMom << std::endl;
0041     std::cout << "Tau GammaPhiMom " << TauGammaPhiMom << std::endl;
0042     std::cout << "Tau Hcal3x3OverPLead " << TauHcal3x3OverPLead << std::endl;
0043     std::cout << "Tau LeadPFChargedHadrEoP " << TauLeadPFChargedHadrEoP << std::endl;
0044     std::cout << "Tau LeadChargedPFCandEtaAtEcalEntrance " << TauLeadChargedPFCandEtaAtEcalEntrance << std::endl;
0045     std::cout << "Tau LeadChargedPFCandPt " << TauLeadChargedPFCandPt << std::endl;
0046     std::cout << "Tau Eta " << TauEta << std::endl;
0047     std::cout << "Tau Pt " << TauPt << std::endl;
0048   }
0049
0050   // ensure tau has at least one charged object
0051   if (TauLeadChargedPFCandPt <= 0) {
0052     return 0.;
0053   } else {
0054     // Check if track goes to Ecal Crack
0055     if (isInEcalCrack(TauLeadChargedPFCandEtaAtEcalEntrance)) {
0056       if (keepAllInEcalCrack_)
0057         return 1.0;
0058       else if (rejectAllInEcalCrack_)
0059         return 0.;
0060     }
0061   }
0062
0063   bool decision = false;
0064   //Apply separate cuts for barrel and endcap
0065   if (TMath::Abs(TauEta) < 1.479) {
0066     //Apply cut for barrel
0067     if (Tau_applyCut_GammaEnFrac_ && TauGammaEnFrac > TauGammaEnFrac_barrel_max_)
0068       decision = true;
0069
0070     if (Tau_applyCut_GammaPhiMom_ && TauGammaPhiMom > TauGammaPhiMom_barrel_max_)
0071       decision = true;
0072
0073     if (Tau_applyCut_GammaEtaMom_ && TauGammaEtaMom > TauGammaEtaMom_barrel_max_)
0074       decision = true;
0075
0076     if (Tau_applyCut_hcal3x3OverPLead_ && TauHcal3x3OverPLead > TauHcal3x3OverPLead_barrel_max_)
0077       decision = true;
0078
0079     if (Tau_applyCut_leadPFChargedHadrEoP_ && (TauLeadPFChargedHadrEoP < TauLeadPFChargedHadrEoP_barrel_min_ ||
0080                                                TauLeadPFChargedHadrEoP > TauLeadPFChargedHadrEoP_barrel_max_))
0081       decision = true;
0082   } else {
0083     //Apply cut for endcap
0084     if (Tau_applyCut_GammaEnFrac_ && TauGammaEnFrac > TauGammaEnFrac_endcap_max_)
0085       decision = true;
0086
0087     if (Tau_applyCut_GammaPhiMom_ && TauGammaPhiMom > TauGammaPhiMom_endcap_max_)
0088       decision = true;
0089
0090     if (Tau_applyCut_GammaEtaMom_ && TauGammaEtaMom > TauGammaEtaMom_endcap_max_)
0091       decision = true;
0092
0093     //This cut is for both offline and HLT. For offline, use cut 0.99-1.01,
0094     //For HLT use cut 0.7-1.3
0095     if (Tau_applyCut_leadPFChargedHadrEoP_ && (TauLeadPFChargedHadrEoP < TauLeadPFChargedHadrEoP_endcap_min1_ ||
0096                                                TauLeadPFChargedHadrEoP > TauLeadPFChargedHadrEoP_endcap_max1_))
0097       decision = true;
0098
0099     //This cut is only for HLT. For HLT, use cut like 0.99-1.01 & H3x3/P>0.1
0100     //For offline, keep the values same as above : 0.99-1.01 & H3x3/P>0, otherwise it may select events in a wrong way.
0101     if (Tau_applyCut_HLTSpecific_ &&
0102         (TauLeadPFChargedHadrEoP < TauLeadPFChargedHadrEoP_endcap_min2_ ||
0103          TauLeadPFChargedHadrEoP > TauLeadPFChargedHadrEoP_endcap_max2_) &&
0104         TauHcal3x3OverPLead > TauHcal3x3OverPLead_endcap_max_)
0105       decision = true;
0106   }
0107
0108   return (decision ? 1. : 0.);
0109 }
0110
0111 double AntiElectronIDCut2::Discriminator(float TauPt,
0112                                          float TauEta,
0113                                          float TauLeadChargedPFCandPt,
0114                                          float TauLeadChargedPFCandEtaAtEcalEntrance,
0115                                          float TauLeadPFChargedHadrEoP,
0116                                          float TauHcal3x3OverPLead,
0117                                          const std::vector<float>& GammasdEta,
0118                                          const std::vector<float>& GammasdPhi,
0119                                          const std::vector<float>& GammasPt) {
0120   double sumPt = 0.;
0121   double dEta2 = 0.;
0122   double dPhi2 = 0.;
0123   for (unsigned int i = 0; i < GammasPt.size(); ++i) {
0124     double pt_i = GammasPt[i];
0125     double phi_i = GammasdPhi[i];
0126     if (GammasdPhi[i] > TMath::Pi())
0127       phi_i = GammasdPhi[i] - 2 * TMath::Pi();
0128     else if (GammasdPhi[i] < -TMath::Pi())
0129       phi_i = GammasdPhi[i] + 2 * TMath::Pi();
0130     double eta_i = GammasdEta[i];
0131     sumPt += pt_i;
0132     dEta2 += (pt_i * eta_i * eta_i);
0133     dPhi2 += (pt_i * phi_i * phi_i);
0134   }
0135
0136   float TauGammaEnFrac = sumPt / TauPt;
0137
0138   if (sumPt > 0.) {
0139     dEta2 /= sumPt;
0140     dPhi2 /= sumPt;
0141   }
0142
0143   float TauGammaEtaMom = TMath::Sqrt(dEta2) * TMath::Sqrt(TauGammaEnFrac) * TauPt;
0144   float TauGammaPhiMom = TMath::Sqrt(dPhi2) * TMath::Sqrt(TauGammaEnFrac) * TauPt;
0145
0146   return Discriminator(TauPt,
0147                        TauEta,
0148                        TauLeadChargedPFCandPt,
0149                        TauLeadChargedPFCandEtaAtEcalEntrance,
0150                        TauLeadPFChargedHadrEoP,
0151                        TauHcal3x3OverPLead,
0152                        TauGammaEtaMom,
0153                        TauGammaPhiMom,
0154                        TauGammaEnFrac);
0155 }
0156
0157 void AntiElectronIDCut2::SetBarrelCutValues(float TauLeadPFChargedHadrEoP_min,
0158                                             float TauLeadPFChargedHadrEoP_max,
0159                                             float TauHcal3x3OverPLead_max,
0160                                             float TauGammaEtaMom_max,
0161                                             float TauGammaPhiMom_max,
0162                                             float TauGammaEnFrac_max) {
0163   TauLeadPFChargedHadrEoP_barrel_min_ = TauLeadPFChargedHadrEoP_min;
0164   TauLeadPFChargedHadrEoP_barrel_max_ = TauLeadPFChargedHadrEoP_max;
0165   TauHcal3x3OverPLead_barrel_max_ = TauHcal3x3OverPLead_max;
0166   TauGammaEtaMom_barrel_max_ = TauGammaEtaMom_max;
0167   TauGammaPhiMom_barrel_max_ = TauGammaPhiMom_max;
0168   TauGammaEnFrac_barrel_max_ = TauGammaEnFrac_max;
0169 }
0170
0171 void AntiElectronIDCut2::SetEndcapCutValues(float TauLeadPFChargedHadrEoP_min_1,
0172                                             float TauLeadPFChargedHadrEoP_max_1,
0173                                             float TauLeadPFChargedHadrEoP_min_2,
0174                                             float TauLeadPFChargedHadrEoP_max_2,
0175                                             float TauHcal3x3OverPLead_max,
0176                                             float TauGammaEtaMom_max,
0177                                             float TauGammaPhiMom_max,
0178                                             float TauGammaEnFrac_max) {
0179   TauLeadPFChargedHadrEoP_endcap_min1_ = TauLeadPFChargedHadrEoP_min_1;
0180   TauLeadPFChargedHadrEoP_endcap_max1_ = TauLeadPFChargedHadrEoP_max_1;
0181   TauLeadPFChargedHadrEoP_endcap_min2_ = TauLeadPFChargedHadrEoP_min_2;
0182   TauLeadPFChargedHadrEoP_endcap_max2_ = TauLeadPFChargedHadrEoP_max_2;
0183   TauHcal3x3OverPLead_endcap_max_ = TauHcal3x3OverPLead_max;
0184   TauGammaEtaMom_endcap_max_ = TauGammaEtaMom_max;
0185   TauGammaPhiMom_endcap_max_ = TauGammaPhiMom_max;
0186   TauGammaEnFrac_endcap_max_ = TauGammaEnFrac_max;
0187 }
0188
0189 void AntiElectronIDCut2::ApplyCuts(bool applyCut_hcal3x3OverPLead,
0190                                    bool applyCut_leadPFChargedHadrEoP,
0191                                    bool applyCut_GammaEtaMom,
0192                                    bool applyCut_GammaPhiMom,
0193                                    bool applyCut_GammaEnFrac,
0194                                    bool applyCut_HLTSpecific) {
0195   Tau_applyCut_hcal3x3OverPLead_ = applyCut_hcal3x3OverPLead;
0196   Tau_applyCut_leadPFChargedHadrEoP_ = applyCut_leadPFChargedHadrEoP;
0197   Tau_applyCut_GammaEtaMom_ = applyCut_GammaEtaMom;
0198   Tau_applyCut_GammaPhiMom_ = applyCut_GammaPhiMom;
0199   Tau_applyCut_GammaEnFrac_ = applyCut_GammaEnFrac;
0200   Tau_applyCut_HLTSpecific_ = applyCut_HLTSpecific;
0201 }
0202
0203 bool AntiElectronIDCut2::isInEcalCrack(double eta) const {
0204   bool in_ecal_crack = false;
0205
0206   eta = fabs(eta);
0207   for (std::vector<pdouble>::const_iterator etaCrack = ecalCracks_.begin(); etaCrack != ecalCracks_.end(); ++etaCrack)
0208     if (eta >= etaCrack->first && eta < etaCrack->second)
0209       in_ecal_crack = true;
0210
0211   return in_ecal_crack;
0212 }