Back to home page

Project CMSSW displayed by LXR
 
 

    

File indexing completed on 2025-02-05 23:51:42 

0001 #include "CutBasedElectronID.h"
0002 
0003 #include "DataFormats/EgammaReco/interface/BasicCluster.h"
0004 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
0005 #include "DataFormats/BeamSpot/interface/BeamSpot.h"
0006 
0007 #include <algorithm>
0008 
0009 CutBasedElectronID::CutBasedElectronID(const edm::ParameterSet& conf, edm::ConsumesCollector& iC) {
0010   verticesCollection_ = iC.consumes<std::vector<reco::Vertex>>(conf.getParameter<edm::InputTag>("verticesCollection"));
0011 }
0012 
0013 void CutBasedElectronID::setup(const edm::ParameterSet& conf) {
0014   // Get all the parameters
0015   //baseSetup(conf);
0016 
0017   type_ = conf.getParameter<std::string>("electronIDType");
0018   quality_ = conf.getParameter<std::string>("electronQuality");
0019   version_ = conf.getParameter<std::string>("electronVersion");
0020 
0021   if (type_ == "classbased" and (version_ == "V06")) {
0022     newCategories_ = conf.getParameter<bool>("additionalCategories");
0023   }
0024 
0025   if (type_ == "classbased" and (version_ == "V03" or version_ == "V04" or version_ == "V05" or version_.empty())) {
0026     wantBinning_ = conf.getParameter<bool>("etBinning");
0027     newCategories_ = conf.getParameter<bool>("additionalCategories");
0028   }
0029 
0030   if (type_ == "robust" || type_ == "classbased") {
0031     std::string stringCut = type_ + quality_ + "EleIDCuts" + version_;
0032     cuts_ = conf.getParameter<edm::ParameterSet>(stringCut);
0033   } else {
0034     throw cms::Exception("Configuration")
0035         << "Invalid electronType parameter in CutBasedElectronID: must be robust or classbased\n";
0036   }
0037 }
0038 
0039 void CutBasedElectronID::fillPSetDescription(edm::ParameterSetDescription& desc) {
0040   // General parameters
0041   desc.add<std::string>("electronIDType", "classbased")->setComment("Type of electron ID: classbased or robust");
0042   desc.add<std::string>("electronQuality", "loose")->setComment("Quality of electron ID: loose or tight");
0043   desc.add<std::string>("electronVersion", "V06")->setComment("Version of electron ID");
0044 
0045   // Optional parameters
0046   desc.addOptional<bool>("additionalCategories", false)->setComment("Flag to enable additional categories");
0047   desc.addOptional<bool>("etBinning", false)->setComment("Flag to enable ET binning");
0048 
0049   // Parameters for classbased ID
0050   edm::ParameterSetDescription classbasedCuts;
0051   classbasedCuts.add<std::vector<double>>("hOverE", std::vector<double>())->setComment("H/E cut values");
0052   classbasedCuts.add<std::vector<double>>("sigmaEtaEta", std::vector<double>())->setComment("SigmaEtaEta cut values");
0053   classbasedCuts.add<std::vector<double>>("deltaPhiIn", std::vector<double>())->setComment("DeltaPhiIn cut values");
0054   classbasedCuts.add<std::vector<double>>("deltaEtaIn", std::vector<double>())->setComment("DeltaEtaIn cut values");
0055   classbasedCuts.add<std::vector<double>>("eSeedOverPin", std::vector<double>())->setComment("ESeedOverPin cut values");
0056   classbasedCuts.add<std::vector<double>>("cutiso_sum", std::vector<double>())->setComment("Isolation sum cut values");
0057   classbasedCuts.add<std::vector<double>>("cutiso_sumoet", std::vector<double>())
0058       ->setComment("Isolation sum over ET cut values");
0059   classbasedCuts.add<std::vector<double>>("cutfmishits", std::vector<double>())->setComment("Missing hits cut values");
0060   classbasedCuts.add<std::vector<double>>("cutdcotdist", std::vector<double>())->setComment("DCOT distance cut values");
0061   classbasedCuts.add<std::vector<double>>("cutip_gsf", std::vector<double>())
0062       ->setComment("Impact parameter cut values");
0063   desc.add<edm::ParameterSetDescription>("classbasedLooseEleIDCutsV06", classbasedCuts)
0064       ->setComment("Classbased loose electron ID cuts for version V06");
0065   desc.add<edm::ParameterSetDescription>("classbasedTightEleIDCutsV06", classbasedCuts)
0066       ->setComment("Classbased tight electron ID cuts for version V06");
0067 
0068   // Parameters for robust ID
0069   edm::ParameterSetDescription robustCuts;
0070   robustCuts.add<std::vector<double>>("barrel", std::vector<double>())->setComment("Cut values for barrel region");
0071   robustCuts.add<std::vector<double>>("endcap", std::vector<double>())->setComment("Cut values for endcap region");
0072   desc.add<edm::ParameterSetDescription>("robustLooseEleIDCuts", robustCuts)
0073       ->setComment("Robust loose electron ID cuts");
0074   desc.add<edm::ParameterSetDescription>("robustTightEleIDCuts", robustCuts)
0075       ->setComment("Robust tight electron ID cuts");
0076 
0077   // Add the vertices collection parameter
0078   desc.add<edm::InputTag>("verticesCollection", edm::InputTag("offlinePrimaryVertices"))
0079       ->setComment("Input tag for the vertices collection");
0080 }
0081 
0082 double CutBasedElectronID::result(const reco::GsfElectron* electron, const edm::Event& e, const edm::EventSetup& es) {
0083   if (type_ == "classbased")
0084     return cicSelection(electron, e, es);
0085   else if (type_ == "robust")
0086     return robustSelection(electron, e, es);
0087 
0088   return 0;
0089 }
0090 
0091 int CutBasedElectronID::classify(const reco::GsfElectron* electron) {
0092   double eta = fabs(electron->superCluster()->eta());
0093   double eOverP = electron->eSuperClusterOverP();
0094   double fBrem = electron->fbrem();
0095 
0096   int cat = -1;
0097   if (version_ == "V00" || version_ == "V01") {
0098     if ((electron->isEB() && fBrem < 0.06) || (electron->isEE() && fBrem < 0.1))
0099       cat = 1;
0100     else if (eOverP < 1.2 && eOverP > 0.8)
0101       cat = 0;
0102     else
0103       cat = 2;
0104 
0105     return cat;
0106 
0107   } else if (version_ == "V02") {
0108     if (electron->isEB()) {  // BARREL
0109       if (fBrem < 0.12)
0110         cat = 1;
0111       else if (eOverP < 1.2 && eOverP > 0.9)
0112         cat = 0;
0113       else
0114         cat = 2;
0115     } else {  // ENDCAP
0116       if (fBrem < 0.2)
0117         cat = 1;
0118       else if (eOverP < 1.22 && eOverP > 0.82)
0119         cat = 0;
0120       else
0121         cat = 2;
0122     }
0123 
0124     return cat;
0125 
0126   } else {
0127     if (electron->isEB()) {
0128       if ((fBrem >= 0.12) and (eOverP > 0.9) and (eOverP < 1.2))
0129         cat = 0;
0130       else if (((eta > .445 and eta < .45) or (eta > .79 and eta < .81) or (eta > 1.137 and eta < 1.157) or
0131                 (eta > 1.47285 and eta < 1.4744)) and
0132                newCategories_)
0133         cat = 6;
0134       else if (electron->trackerDrivenSeed() and !electron->ecalDrivenSeed() and newCategories_)
0135         cat = 8;
0136       else if (fBrem < 0.12)
0137         cat = 1;
0138       else
0139         cat = 2;
0140     } else {
0141       if ((fBrem >= 0.2) and (eOverP > 0.82) and (eOverP < 1.22))
0142         cat = 3;
0143       else if (eta > 1.5 and eta < 1.58 and newCategories_)
0144         cat = 7;
0145       else if (electron->trackerDrivenSeed() and !electron->ecalDrivenSeed() and newCategories_)
0146         cat = 8;
0147       else if (fBrem < 0.2)
0148         cat = 4;
0149       else
0150         cat = 5;
0151     }
0152 
0153     return cat;
0154   }
0155 
0156   return -1;
0157 }
0158 
0159 double CutBasedElectronID::cicSelection(const reco::GsfElectron* electron,
0160                                         const edm::Event& e,
0161                                         const edm::EventSetup& es) {
0162   double scTheta = (2 * atan(exp(-electron->superCluster()->eta())));
0163   double scEt = electron->superCluster()->energy() * sin(scTheta);
0164 
0165   double eta = fabs(electron->superCluster()->eta());
0166 
0167   double eOverP = electron->eSuperClusterOverP();
0168   double eSeedOverPin = electron->eSeedClusterOverP();
0169   double fBrem = electron->fbrem();
0170   double hOverE = electron->hadronicOverEm();
0171   double sigmaee = electron->sigmaIetaIeta();  //sqrt(vLocCov[0]);
0172   double deltaPhiIn = electron->deltaPhiSuperClusterTrackAtVtx();
0173   double deltaEtaIn = electron->deltaEtaSuperClusterTrackAtVtx();
0174 
0175   double ip = 0;
0176   int mishits = electron->gsfTrack()->missingInnerHits();
0177   double tkIso = electron->dr03TkSumPt();
0178   double ecalIso = electron->dr04EcalRecHitSumEt();
0179   double hcalIso = electron->dr04HcalTowerSumEt();
0180 
0181   if (version_ == "V00") {
0182     sigmaee = electron->sigmaEtaEta();  //sqrt(vCov[0]);
0183     if (electron->isEE())
0184       sigmaee = sigmaee - 0.02 * (fabs(eta) - 2.3);  //correct sigmaetaeta dependence on eta in endcap
0185   }
0186 
0187   if (version_ != "V01" or version_ != "V00") {
0188     edm::Handle<reco::VertexCollection> vtxH;
0189     e.getByToken(verticesCollection_, vtxH);
0190     if (!vtxH->empty()) {
0191       reco::VertexRef vtx(vtxH, 0);
0192       ip = fabs(electron->gsfTrack()->dxy(math::XYZPoint(vtx->x(), vtx->y(), vtx->z())));
0193     } else
0194       ip = fabs(electron->gsfTrack()->dxy());
0195 
0196     if (electron->isEB()) {
0197       sigmaee = electron->sigmaIetaIeta();  //sqrt(vCov[0]);
0198     }
0199   }
0200 
0201   std::vector<double> cut;
0202 
0203   int cat = classify(electron);
0204   int eb;
0205 
0206   if (electron->isEB())
0207     eb = 0;
0208   else
0209     eb = 1;
0210 
0211   // LOOSE and TIGHT Selections
0212   if (type_ == "classbased" && (version_ == "V01" || version_ == "V00")) {
0213     if ((eOverP < 0.8) && (fBrem < 0.2))
0214       return 0.;
0215 
0216     cut = cuts_.getParameter<std::vector<double>>("hOverE");
0217     if (hOverE > cut[cat + 4 * eb])
0218       return 0.;
0219 
0220     cut = cuts_.getParameter<std::vector<double>>("sigmaEtaEta");
0221     if (sigmaee > cut[cat + 4 * eb])
0222       return 0.;
0223 
0224     cut = cuts_.getParameter<std::vector<double>>("deltaPhiIn");
0225     if (eOverP < 1.5) {
0226       if (fabs(deltaPhiIn) > cut[cat + 4 * eb])
0227         return 0.;
0228     } else {
0229       if (fabs(deltaPhiIn) > cut[3 + 4 * eb])
0230         return 0.;
0231     }
0232 
0233     cut = cuts_.getParameter<std::vector<double>>("deltaEtaIn");
0234     if (fabs(deltaEtaIn) > cut[cat + 4 * eb])
0235       return 0.;
0236 
0237     cut = cuts_.getParameter<std::vector<double>>("eSeedOverPin");
0238     if (eSeedOverPin < cut[cat + 4 * eb])
0239       return 0.;
0240 
0241     if (quality_ == "tight")
0242       if (eOverP < 0.9 * (1 - fBrem))
0243         return 0.;
0244 
0245     return 1.;
0246   }
0247 
0248   if (type_ == "classbased" and version_ == "V02") {
0249     double result = 0.;
0250 
0251     int bin = 0;
0252 
0253     if (scEt < 20.)
0254       bin = 2;
0255     else if (scEt > 30.)
0256       bin = 0;
0257     else
0258       bin = 1;
0259 
0260     if (fBrem > 0)
0261       eSeedOverPin = eSeedOverPin + fBrem;
0262 
0263     if (bin != 2) {
0264       tkIso = tkIso * pow(40. / scEt, 2);
0265       ecalIso = ecalIso * pow(40. / scEt, 2);
0266       hcalIso = hcalIso * pow(40. / scEt, 2);
0267     }
0268 
0269     std::vector<double> cutTk = cuts_.getParameter<std::vector<double>>("cutisotk");
0270     std::vector<double> cutEcal = cuts_.getParameter<std::vector<double>>("cutisoecal");
0271     std::vector<double> cutHcal = cuts_.getParameter<std::vector<double>>("cutisohcal");
0272     if ((tkIso > cutTk[cat + 3 * eb + bin * 6]) || (ecalIso > cutEcal[cat + 3 * eb + bin * 6]) ||
0273         (hcalIso > cutHcal[cat + 3 * eb + bin * 6]))
0274       result = 0.;
0275     else
0276       result = 2.;
0277 
0278     if (fBrem > -2) {
0279       std::vector<double> cuthoe = cuts_.getParameter<std::vector<double>>("cuthoe");
0280       std::vector<double> cutsee = cuts_.getParameter<std::vector<double>>("cutsee");
0281       std::vector<double> cutdphi = cuts_.getParameter<std::vector<double>>("cutdphiin");
0282       std::vector<double> cutdeta = cuts_.getParameter<std::vector<double>>("cutdetain");
0283       std::vector<double> cuteopin = cuts_.getParameter<std::vector<double>>("cuteseedopcor");
0284       std::vector<double> cutet = cuts_.getParameter<std::vector<double>>("cutet");
0285       std::vector<double> cutip = cuts_.getParameter<std::vector<double>>("cutip");
0286       std::vector<double> cutmishits = cuts_.getParameter<std::vector<double>>("cutmishits");
0287       if ((hOverE < cuthoe[cat + 3 * eb + bin * 6]) and (sigmaee < cutsee[cat + 3 * eb + bin * 6]) and
0288           (fabs(deltaPhiIn) < cutdphi[cat + 3 * eb + bin * 6]) and
0289           (fabs(deltaEtaIn) < cutdeta[cat + 3 * eb + bin * 6]) and (eSeedOverPin > cuteopin[cat + 3 * eb + bin * 6]) and
0290           (ip < cutip[cat + 3 * eb + bin * 6]) and (mishits < cutmishits[cat + 3 * eb + bin * 6]))
0291         result = result + 1.;
0292     }
0293     return result;
0294   }
0295 
0296   if (version_ == "V03" or version_ == "V04" or version_ == "V05") {
0297     double result = 0.;
0298 
0299     int bin = 0;
0300 
0301     if (wantBinning_) {
0302       if (scEt < 20.)
0303         bin = 2;
0304       else if (scEt > 30.)
0305         bin = 0;
0306       else
0307         bin = 1;
0308     }
0309 
0310     if (fBrem > 0)
0311       eSeedOverPin = eSeedOverPin + fBrem;
0312 
0313     float iso_sum = tkIso + ecalIso + hcalIso;
0314     float iso_sum_corrected = iso_sum * pow(40. / scEt, 2);
0315 
0316     std::vector<double> cutIsoSum = cuts_.getParameter<std::vector<double>>("cutiso_sum");
0317     std::vector<double> cutIsoSumCorr = cuts_.getParameter<std::vector<double>>("cutiso_sumoet");
0318     if ((iso_sum < cutIsoSum[cat + bin * 9]) and (iso_sum_corrected < cutIsoSumCorr[cat + bin * 9]))
0319       result += 2.;
0320 
0321     if (fBrem > -2) {
0322       std::vector<double> cuthoe = cuts_.getParameter<std::vector<double>>("cuthoe");
0323       std::vector<double> cutsee = cuts_.getParameter<std::vector<double>>("cutsee");
0324       std::vector<double> cutdphi = cuts_.getParameter<std::vector<double>>("cutdphiin");
0325       std::vector<double> cutdeta = cuts_.getParameter<std::vector<double>>("cutdetain");
0326       std::vector<double> cuteopin = cuts_.getParameter<std::vector<double>>("cuteseedopcor");
0327       std::vector<double> cutet = cuts_.getParameter<std::vector<double>>("cutet");
0328 
0329       if ((hOverE < cuthoe[cat + bin * 9]) and (sigmaee < cutsee[cat + bin * 9]) and
0330           (fabs(deltaPhiIn) < cutdphi[cat + bin * 9]) and (fabs(deltaEtaIn) < cutdeta[cat + bin * 9]) and
0331           (eSeedOverPin > cuteopin[cat + bin * 9]) and (scEt > cutet[cat + bin * 9]))
0332         result += 1.;
0333     }
0334 
0335     std::vector<double> cutip = cuts_.getParameter<std::vector<double>>("cutip_gsf");
0336     if (ip < cutip[cat + bin * 9])
0337       result += 8;
0338 
0339     std::vector<double> cutmishits = cuts_.getParameter<std::vector<double>>("cutfmishits");
0340     std::vector<double> cutdcotdist = cuts_.getParameter<std::vector<double>>("cutdcotdist");
0341 
0342     float dist = (electron->convDist() == -9999. ? 9999 : electron->convDist());
0343     float dcot = (electron->convDcot() == -9999. ? 9999 : electron->convDcot());
0344 
0345     float dcotdistcomb =
0346         ((0.04 - std::max(fabs(dist), fabs(dcot))) > 0 ? (0.04 - std::max(fabs(dist), fabs(dcot))) : 0);
0347 
0348     if ((mishits < cutmishits[cat + bin * 9]) and (dcotdistcomb < cutdcotdist[cat + bin * 9]))
0349       result += 4;
0350 
0351     return result;
0352   }
0353 
0354   if (type_ == "classbased" && (version_ == "V06" || version_.empty())) {
0355     std::vector<double> cutIsoSum = cuts_.getParameter<std::vector<double>>("cutiso_sum");
0356     std::vector<double> cutIsoSumCorr = cuts_.getParameter<std::vector<double>>("cutiso_sumoet");
0357     std::vector<double> cuthoe = cuts_.getParameter<std::vector<double>>("cuthoe");
0358     std::vector<double> cutsee = cuts_.getParameter<std::vector<double>>("cutsee");
0359     std::vector<double> cutdphi = cuts_.getParameter<std::vector<double>>("cutdphiin");
0360     std::vector<double> cutdeta = cuts_.getParameter<std::vector<double>>("cutdetain");
0361     std::vector<double> cuteopin = cuts_.getParameter<std::vector<double>>("cuteseedopcor");
0362     std::vector<double> cutmishits = cuts_.getParameter<std::vector<double>>("cutfmishits");
0363     std::vector<double> cutdcotdist = cuts_.getParameter<std::vector<double>>("cutdcotdist");
0364     std::vector<double> cutip = cuts_.getParameter<std::vector<double>>("cutip_gsf");
0365     std::vector<double> cutIsoSumCorrl = cuts_.getParameter<std::vector<double>>("cutiso_sumoetl");
0366     std::vector<double> cuthoel = cuts_.getParameter<std::vector<double>>("cuthoel");
0367     std::vector<double> cutseel = cuts_.getParameter<std::vector<double>>("cutseel");
0368     std::vector<double> cutdphil = cuts_.getParameter<std::vector<double>>("cutdphiinl");
0369     std::vector<double> cutdetal = cuts_.getParameter<std::vector<double>>("cutdetainl");
0370     std::vector<double> cutipl = cuts_.getParameter<std::vector<double>>("cutip_gsfl");
0371 
0372     int result = 0;
0373 
0374     const int ncuts = 10;
0375     std::vector<bool> cut_results(ncuts, false);
0376 
0377     float iso_sum = tkIso + ecalIso + hcalIso;
0378     float scEta = electron->superCluster()->eta();
0379     if (fabs(scEta) > 1.5)
0380       iso_sum += (fabs(scEta) - 1.5) * 1.09;
0381 
0382     float iso_sumoet = iso_sum * (40. / scEt);
0383 
0384     float eseedopincor = eSeedOverPin + fBrem;
0385     if (fBrem < 0)
0386       eseedopincor = eSeedOverPin;
0387 
0388     float dist = (electron->convDist() == -9999. ? 9999 : electron->convDist());
0389     float dcot = (electron->convDcot() == -9999. ? 9999 : electron->convDcot());
0390 
0391     float dcotdistcomb =
0392         ((0.04 - std::max(fabs(dist), fabs(dcot))) > 0 ? (0.04 - std::max(fabs(dist), fabs(dcot))) : 0);
0393 
0394     for (int cut = 0; cut < ncuts; cut++) {
0395       switch (cut) {
0396         case 0:
0397           cut_results[cut] = compute_cut(fabs(deltaEtaIn), scEt, cutdetal[cat], cutdeta[cat]);
0398           break;
0399         case 1:
0400           cut_results[cut] = compute_cut(fabs(deltaPhiIn), scEt, cutdphil[cat], cutdphi[cat]);
0401           break;
0402         case 2:
0403           cut_results[cut] = (eseedopincor > cuteopin[cat]);
0404           break;
0405         case 3:
0406           cut_results[cut] = compute_cut(hOverE, scEt, cuthoel[cat], cuthoe[cat]);
0407           break;
0408         case 4:
0409           cut_results[cut] = compute_cut(sigmaee, scEt, cutseel[cat], cutsee[cat]);
0410           break;
0411         case 5:
0412           cut_results[cut] = compute_cut(iso_sumoet, scEt, cutIsoSumCorrl[cat], cutIsoSumCorr[cat]);
0413           break;
0414         case 6:
0415           cut_results[cut] = (iso_sum < cutIsoSum[cat]);
0416           break;
0417         case 7:
0418           cut_results[cut] = compute_cut(fabs(ip), scEt, cutipl[cat], cutip[cat]);
0419           break;
0420         case 8:
0421           cut_results[cut] = (mishits < cutmishits[cat]);
0422           break;
0423         case 9:
0424           cut_results[cut] = (dcotdistcomb < cutdcotdist[cat]);
0425           break;
0426       }
0427     }
0428 
0429     // ID part
0430     if (cut_results[0] && cut_results[1] && cut_results[2] && cut_results[3] && cut_results[4])
0431       result = result + 1;
0432 
0433     // ISO part
0434     if (cut_results[5] && cut_results[6])
0435       result = result + 2;
0436 
0437     // IP part
0438     if (cut_results[7])
0439       result = result + 8;
0440 
0441     // Conversion part
0442     if (cut_results[8] && cut_results[9])
0443       result = result + 4;
0444 
0445     return result;
0446   }
0447 
0448   return -1.;
0449 }
0450 
0451 bool CutBasedElectronID::compute_cut(double x, double et, double cut_min, double cut_max, bool gtn) {
0452   float et_min = 10;
0453   float et_max = 40;
0454 
0455   bool accept = false;
0456   float cut = cut_max;  //  the cut at et=40 GeV
0457 
0458   if (et < et_max) {
0459     cut = cut_min + (1 / et_min - 1 / et) * (cut_max - cut_min) / (1 / et_min - 1 / et_max);
0460   }
0461 
0462   if (et < et_min) {
0463     cut = cut_min;
0464   }
0465 
0466   if (gtn) {  // useful for e/p cut which is gt
0467     accept = (x >= cut);
0468   } else {
0469     accept = (x <= cut);
0470   }
0471 
0472   //std::cout << x << " " << cut_min << " " << cut << " " << cut_max << " " << et << " " << accept << std::endl;
0473   return accept;
0474 }
0475 
0476 double CutBasedElectronID::robustSelection(const reco::GsfElectron* electron,
0477                                            const edm::Event& e,
0478                                            const edm::EventSetup& es) {
0479   double scTheta = (2 * atan(exp(-electron->superCluster()->eta())));
0480   double scEt = electron->superCluster()->energy() * sin(scTheta);
0481   double eta = electron->p4().Eta();
0482   double eOverP = electron->eSuperClusterOverP();
0483   double hOverE = electron->hadronicOverEm();
0484   double sigmaee = electron->sigmaIetaIeta();
0485   double e25Max = electron->e2x5Max();
0486   double e15 = electron->e1x5();
0487   double e55 = electron->e5x5();
0488   double e25Maxoe55 = e25Max / e55;
0489   double e15oe55 = e15 / e55;
0490   double deltaPhiIn = electron->deltaPhiSuperClusterTrackAtVtx();
0491   double deltaEtaIn = electron->deltaEtaSuperClusterTrackAtVtx();
0492 
0493   double ip = 0;
0494   int mishits = electron->gsfTrack()->missingInnerHits();
0495   double tkIso = electron->dr03TkSumPt();
0496   double ecalIso = electron->dr04EcalRecHitSumEt();
0497   double ecalIsoPed = (electron->isEB()) ? std::max(0., ecalIso - 1.) : ecalIso;
0498   double hcalIso = electron->dr04HcalTowerSumEt();
0499   double hcalIso1 = electron->dr04HcalTowerSumEt(1);
0500   double hcalIso2 = electron->dr04HcalTowerSumEt(2);
0501 
0502   if (version_ == "V00") {
0503     sigmaee = electron->sigmaEtaEta();
0504     if (electron->isEE())
0505       sigmaee = sigmaee - 0.02 * (fabs(eta) - 2.3);  //correct sigmaetaeta dependence on eta in endcap
0506   }
0507 
0508   if (version_ == "V03" or version_ == "V04") {
0509     edm::Handle<reco::BeamSpot> pBeamSpot;
0510     // uses the same name for the vertex collection to avoid adding more new names
0511     e.getByToken(verticesCollection_, pBeamSpot);
0512     if (pBeamSpot.isValid()) {
0513       const reco::BeamSpot* bspot = pBeamSpot.product();
0514       const math::XYZPoint& bspotPosition = bspot->position();
0515       ip = fabs(electron->gsfTrack()->dxy(bspotPosition));
0516     } else
0517       ip = fabs(electron->gsfTrack()->dxy());
0518   }
0519 
0520   if (version_ == "V04" or version_ == "V05") {
0521     ecalIso = electron->dr03EcalRecHitSumEt();
0522     ecalIsoPed = (electron->isEB()) ? std::max(0., ecalIso - 1.) : ecalIso;
0523     hcalIso = electron->dr03HcalTowerSumEt();
0524     hcalIso1 = electron->dr03HcalTowerSumEt(1);
0525     hcalIso2 = electron->dr03HcalTowerSumEt(2);
0526   }
0527 
0528   if (version_ == "V05") {
0529     edm::Handle<reco::VertexCollection> vtxH;
0530     e.getByToken(verticesCollection_, vtxH);
0531     if (!vtxH->empty()) {
0532       reco::VertexRef vtx(vtxH, 0);
0533       ip = fabs(electron->gsfTrack()->dxy(math::XYZPoint(vtx->x(), vtx->y(), vtx->z())));
0534     } else
0535       ip = fabs(electron->gsfTrack()->dxy());
0536   }
0537 
0538   // .....................................................................................
0539   std::vector<double> cut;
0540   // ROBUST Selection
0541   if (type_ == "robust") {
0542     double result = 0;
0543 
0544     // hoe, sigmaEtaEta, dPhiIn, dEtaIn
0545     if (electron->isEB())
0546       cut = cuts_.getParameter<std::vector<double>>("barrel");
0547     else
0548       cut = cuts_.getParameter<std::vector<double>>("endcap");
0549     // check isolations: if only isolation passes result = 2
0550     if (quality_ == "highenergy") {
0551       if ((tkIso > cut[6] || hcalIso2 > cut[11]) ||
0552           (electron->isEB() && ((ecalIso + hcalIso1) > cut[7] + cut[8] * scEt)) ||
0553           (electron->isEE() && (scEt >= 50.) && ((ecalIso + hcalIso1) > cut[7] + cut[8] * (scEt - 50))) ||
0554           (electron->isEE() && (scEt < 50.) && ((ecalIso + hcalIso1) > cut[9] + cut[10] * (scEt - 50))))
0555         result = 0;
0556       else
0557         result = 2;
0558     } else {
0559       if ((tkIso > cut[6]) || (ecalIso > cut[7]) || (hcalIso > cut[8]) || (hcalIso1 > cut[9]) || (hcalIso2 > cut[10]) ||
0560           (tkIso / electron->p4().Pt() > cut[11]) || (ecalIso / electron->p4().Pt() > cut[12]) ||
0561           (hcalIso / electron->p4().Pt() > cut[13]) || ((tkIso + ecalIso + hcalIso) > cut[14]) ||
0562           (((tkIso + ecalIso + hcalIso) / electron->p4().Pt()) > cut[15]) ||
0563           ((tkIso + ecalIsoPed + hcalIso) > cut[16]) ||
0564           (((tkIso + ecalIsoPed + hcalIso) / electron->p4().Pt()) > cut[17]))
0565         result = 0.;
0566       else
0567         result = 2.;
0568     }
0569 
0570     if ((hOverE < cut[0]) && (sigmaee < cut[1]) && (fabs(deltaPhiIn) < cut[2]) && (fabs(deltaEtaIn) < cut[3]) &&
0571         (e25Maxoe55 > cut[4] && e15oe55 > cut[5]) && (sigmaee >= cut[18]) && (eOverP > cut[19] && eOverP < cut[20])) {
0572       result = result + 1;
0573     }
0574 
0575     if (ip > cut[21])
0576       return result;
0577     if (mishits > cut[22])  // expected missing hits
0578       return result;
0579     // positive cut[23] means to demand a valid hit in 1st layer PXB
0580     if (cut[23] > 0 &&
0581         !electron->gsfTrack()->hitPattern().hasValidHitInPixelLayer(PixelSubdetector::SubDetector::PixelBarrel, 1))
0582       return result;
0583 
0584     // cut[24]: Dist cut[25]: dcot
0585     float dist = fabs(electron->convDist());
0586     float dcot = fabs(electron->convDcot());
0587     bool isConversion = (cut[24] > 99. || cut[25] > 99.) ? false : (dist < cut[24] && dcot < cut[25]);
0588     if (isConversion)
0589       return result;
0590 
0591     result += 4;
0592 
0593     return result;
0594   }
0595 
0596   return -1.;
0597 }
[ Source navigation ] [ Diff markup ] [ Identifier search ] [ general search ]

This page was automatically generated by the 2.2.1 LXR engine.

The LXR team
Valid CSS 2.1! Valid HTML 4.01!