File indexing completed on 2024-09-07 04:37:23
0001 #ifndef PhysicsTools_PatUtils_interface_JetIDSelectionFunctor_h
0002 #define PhysicsTools_PatUtils_interface_JetIDSelectionFunctor_h
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016
0017
0018 #ifndef __GCCXML__
0019 #include "FWCore/Framework/interface/ConsumesCollector.h"
0020 #endif
0021 #include "DataFormats/PatCandidates/interface/Jet.h"
0022 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0023
0024 #include "PhysicsTools/SelectorUtils/interface/Selector.h"
0025 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0026
0027 #include <TMath.h>
0028 class JetIDSelectionFunctor : public Selector<pat::Jet> {
0029 public:
0030 enum Version_t { CRAFT08, PURE09, DQM09, N_VERSIONS };
0031 enum Quality_t { MINIMAL, LOOSE_AOD, LOOSE, TIGHT, N_QUALITY };
0032
0033 JetIDSelectionFunctor() {}
0034
0035 #ifndef __GCCXML__
0036 JetIDSelectionFunctor(edm::ParameterSet const& parameters, edm::ConsumesCollector& iC)
0037 : JetIDSelectionFunctor(parameters) {}
0038 #endif
0039
0040 JetIDSelectionFunctor(edm::ParameterSet const& parameters) {
0041 std::string versionStr = parameters.getParameter<std::string>("version");
0042 std::string qualityStr = parameters.getParameter<std::string>("quality");
0043 Quality_t quality = N_QUALITY;
0044
0045 if (qualityStr == "MINIMAL")
0046 quality = MINIMAL;
0047 else if (qualityStr == "LOOSE_AOD")
0048 quality = LOOSE_AOD;
0049 else if (qualityStr == "LOOSE")
0050 quality = LOOSE;
0051 else if (qualityStr == "TIGHT")
0052 quality = TIGHT;
0053 else
0054 throw cms::Exception("InvalidInput")
0055 << "Expect quality to be one of MINIMAL, LOOSE_AOD, LOOSE,TIGHT" << std::endl;
0056
0057 if (versionStr == "CRAFT08") {
0058 initialize(CRAFT08, quality);
0059 } else if (versionStr == "PURE09") {
0060 initialize(PURE09, quality);
0061 } else if (versionStr == "DQM09") {
0062 initialize(DQM09, quality);
0063 } else {
0064 throw cms::Exception("InvalidInput") << "Expect version to be one of CRAFT08, PURE09, DQM09" << std::endl;
0065 }
0066 }
0067
0068 JetIDSelectionFunctor(Version_t version, Quality_t quality) { initialize(version, quality); }
0069
0070 void initialize(Version_t version, Quality_t quality) {
0071 version_ = version;
0072 quality_ = quality;
0073
0074 push_back("MINIMAL_EMF");
0075
0076 push_back("LOOSE_AOD_fHPD");
0077 push_back("LOOSE_AOD_N90Hits");
0078 push_back("LOOSE_AOD_EMF");
0079
0080 push_back("LOOSE_fHPD");
0081 push_back("LOOSE_N90Hits");
0082 push_back("LOOSE_EMF");
0083
0084 push_back("TIGHT_fHPD");
0085 push_back("TIGHT_EMF");
0086
0087 push_back("LOOSE_nHit");
0088 push_back("LOOSE_als");
0089 push_back("LOOSE_fls");
0090 push_back("LOOSE_foot");
0091
0092 push_back("TIGHT_nHit");
0093 push_back("TIGHT_als");
0094 push_back("TIGHT_fls");
0095 push_back("TIGHT_foot");
0096 push_back("widths");
0097 push_back("EF_N90Hits");
0098 push_back("EF_EMF");
0099
0100 bool use_09_fwd_id = version_ != CRAFT08;
0101 bool use_dqm_09 = version_ == DQM09 && quality_ != LOOSE_AOD;
0102
0103
0104 set("MINIMAL_EMF");
0105 set("LOOSE_AOD_fHPD");
0106 set("LOOSE_AOD_N90Hits");
0107 set("LOOSE_AOD_EMF", !use_09_fwd_id);
0108 set("LOOSE_fHPD");
0109 set("LOOSE_N90Hits");
0110 set("LOOSE_EMF", !use_09_fwd_id);
0111 set("TIGHT_fHPD");
0112 set("TIGHT_EMF");
0113
0114 set("LOOSE_nHit", use_09_fwd_id);
0115 set("LOOSE_als", use_09_fwd_id);
0116 set("TIGHT_nHit", use_09_fwd_id);
0117 set("TIGHT_als", use_09_fwd_id);
0118 set("widths", use_09_fwd_id);
0119 set("EF_N90Hits", use_09_fwd_id);
0120 set("EF_EMF", use_09_fwd_id);
0121
0122 set("LOOSE_fls", use_dqm_09);
0123 set("LOOSE_foot", use_dqm_09);
0124 set("TIGHT_fls", use_dqm_09);
0125 set("TIGHT_foot", use_dqm_09);
0126
0127 index_MINIMAL_EMF_ = index_type(&bits_, "MINIMAL_EMF");
0128
0129 index_LOOSE_AOD_fHPD_ = index_type(&bits_, "LOOSE_AOD_fHPD");
0130 index_LOOSE_AOD_N90Hits_ = index_type(&bits_, "LOOSE_AOD_N90Hits");
0131 index_LOOSE_AOD_EMF_ = index_type(&bits_, "LOOSE_AOD_EMF");
0132
0133 index_LOOSE_fHPD_ = index_type(&bits_, "LOOSE_fHPD");
0134 index_LOOSE_N90Hits_ = index_type(&bits_, "LOOSE_N90Hits");
0135 index_LOOSE_EMF_ = index_type(&bits_, "LOOSE_EMF");
0136
0137 index_TIGHT_fHPD_ = index_type(&bits_, "TIGHT_fHPD");
0138 index_TIGHT_EMF_ = index_type(&bits_, "TIGHT_EMF");
0139
0140 index_LOOSE_nHit_ = index_type(&bits_, "LOOSE_nHit");
0141 index_LOOSE_als_ = index_type(&bits_, "LOOSE_als");
0142 index_LOOSE_fls_ = index_type(&bits_, "LOOSE_fls");
0143 index_LOOSE_foot_ = index_type(&bits_, "LOOSE_foot");
0144
0145 index_TIGHT_nHit_ = index_type(&bits_, "TIGHT_nHit");
0146 index_TIGHT_als_ = index_type(&bits_, "TIGHT_als");
0147 index_TIGHT_fls_ = index_type(&bits_, "TIGHT_fls");
0148 index_TIGHT_foot_ = index_type(&bits_, "TIGHT_foot");
0149 index_widths_ = index_type(&bits_, "widths");
0150 index_EF_N90Hits_ = index_type(&bits_, "EF_N90Hits");
0151 index_EF_EMF_ = index_type(&bits_, "EF_EMF");
0152
0153
0154 bool use_loose_aod = false;
0155 bool use_loose = false;
0156 bool use_tight = false;
0157 bool use_tight_09_fwd_id = false;
0158 bool use_loose_09_fwd_id = false;
0159
0160 if (quality_ == LOOSE) {
0161 use_loose = true;
0162 if (use_09_fwd_id)
0163 use_loose_09_fwd_id = true;
0164 }
0165 if (quality_ == LOOSE_AOD) {
0166 use_loose_aod = true;
0167 if (use_09_fwd_id)
0168 use_loose_09_fwd_id = true;
0169 }
0170 if (quality_ == TIGHT) {
0171 use_tight = true;
0172 if (use_09_fwd_id)
0173 use_tight_09_fwd_id = true;
0174 }
0175
0176 if (!use_loose_aod) {
0177 set("LOOSE_AOD_fHPD", false);
0178 set("LOOSE_AOD_N90Hits", false);
0179 set("LOOSE_AOD_EMF", false);
0180 }
0181
0182 if (!(use_loose || use_tight)) {
0183 set("LOOSE_N90Hits", false);
0184 set("LOOSE_fHPD", false);
0185 set("LOOSE_EMF", false);
0186 }
0187
0188 if (!use_tight) {
0189 set("TIGHT_fHPD", false);
0190 set("TIGHT_EMF", false);
0191 }
0192
0193 if (!use_loose_09_fwd_id) {
0194 set("LOOSE_nHit", false);
0195 set("LOOSE_als", false);
0196 if (use_dqm_09) {
0197 set("LOOSE_fls", false);
0198 set("LOOSE_foot", false);
0199 }
0200 }
0201
0202 if (!use_tight_09_fwd_id) {
0203 set("TIGHT_nHit", false);
0204 set("TIGHT_als", false);
0205 set("widths", false);
0206 set("EF_N90Hits", false);
0207 set("EF_EMF", false);
0208 if (use_dqm_09) {
0209 set("TIGHT_fls", false);
0210 set("TIGHT_foot", false);
0211 }
0212 }
0213
0214 retInternal_ = getBitTemplate();
0215 }
0216
0217
0218 unsigned int count_hits(const std::vector<CaloTowerPtr>& towers) {
0219 unsigned int nHit = 0;
0220 for (unsigned int iTower = 0; iTower < towers.size(); ++iTower) {
0221 const std::vector<DetId>& cellIDs = towers[iTower]->constituents();
0222 nHit += cellIDs.size();
0223 }
0224 return nHit;
0225 }
0226
0227
0228
0229
0230 bool operator()(const pat::Jet& jet, pat::strbitset& ret) override {
0231 if (!jet.isCaloJet() && !jet.isJPTJet()) {
0232 edm::LogWarning("NYI") << "Criteria for pat::Jet-s other than CaloJets and JPTJets are not yet implemented";
0233 return false;
0234 }
0235 if (version_ == CRAFT08)
0236 return craft08Cuts(jet.p4(), jet.emEnergyFraction(), jet.jetID(), ret);
0237 if (version_ == PURE09 || version_ == DQM09) {
0238 unsigned int nHit = count_hits(jet.getCaloConstituents());
0239 if (jet.currentJECLevel() == "Uncorrected") {
0240 return fwd09Cuts(
0241 jet.p4(), jet.emEnergyFraction(), jet.etaetaMoment(), jet.phiphiMoment(), nHit, jet.jetID(), ret);
0242 } else {
0243 return fwd09Cuts(jet.correctedP4("Uncorrected"),
0244 jet.emEnergyFraction(),
0245 jet.etaetaMoment(),
0246 jet.phiphiMoment(),
0247 nHit,
0248 jet.jetID(),
0249 ret);
0250 }
0251 }
0252 edm::LogWarning("BadInput | NYI") << "Requested version (" << version_ << ") is unknown";
0253 return false;
0254 }
0255
0256
0257 using Selector<pat::Jet>::operator();
0258
0259
0260
0261
0262
0263 bool operator()(reco::Candidate::LorentzVector const& correctedP4,
0264 double emEnergyFraction,
0265 reco::JetID const& jetID,
0266 pat::strbitset& ret) {
0267 if (version_ == CRAFT08)
0268 return craft08Cuts(correctedP4, emEnergyFraction, jetID, ret);
0269 edm::LogWarning("BadInput | NYI") << "Requested version (" << version_
0270 << ") is unknown or doesn't match the depreceated interface used";
0271 return false;
0272 }
0273
0274 virtual bool operator()(reco::Candidate::LorentzVector const& correctedP4,
0275 double emEnergyFraction,
0276 reco::JetID const& jetID) {
0277 retInternal_.set(false);
0278 operator()(correctedP4, emEnergyFraction, jetID, retInternal_);
0279 setIgnored(retInternal_);
0280 return (bool)retInternal_;
0281 }
0282
0283
0284
0285
0286
0287 bool operator()(reco::CaloJet const& jet, reco::JetID const& jetID, pat::strbitset& ret) {
0288 if (version_ == CRAFT08)
0289 return craft08Cuts(jet.p4(), jet.emEnergyFraction(), jetID, ret);
0290 if (version_ == PURE09 || version_ == DQM09) {
0291 unsigned int nHit = count_hits(jet.getCaloConstituents());
0292 return fwd09Cuts(jet.p4(), jet.emEnergyFraction(), jet.etaetaMoment(), jet.phiphiMoment(), nHit, jetID, ret);
0293 }
0294 edm::LogWarning("BadInput | NYI") << "Requested version (" << version_ << ") is unknown";
0295 return false;
0296 }
0297
0298
0299 virtual bool operator()(reco::CaloJet const& jet, reco::JetID const& jetID) {
0300 retInternal_.set(false);
0301 operator()(jet, jetID, retInternal_);
0302 setIgnored(retInternal_);
0303 return (bool)retInternal_;
0304 }
0305
0306
0307
0308
0309 bool craft08Cuts(reco::Candidate::LorentzVector const& correctedP4,
0310 double emEnergyFraction,
0311 reco::JetID const& jetID,
0312 pat::strbitset& ret) {
0313 ret.set(false);
0314
0315
0316 double abs_eta = TMath::Abs(correctedP4.eta());
0317 double corrPt = correctedP4.pt();
0318 double emf = emEnergyFraction;
0319
0320 if (ignoreCut(index_MINIMAL_EMF_) || abs_eta > 2.6 || emf > 0.01)
0321 passCut(ret, index_MINIMAL_EMF_);
0322
0323 if (quality_ == LOOSE_AOD) {
0324
0325 if (ignoreCut(index_LOOSE_AOD_fHPD_) || jetID.approximatefHPD < 0.98)
0326 passCut(ret, index_LOOSE_AOD_fHPD_);
0327
0328 if (ignoreCut(index_LOOSE_AOD_N90Hits_) || jetID.hitsInN90 > 1)
0329 passCut(ret, index_LOOSE_AOD_N90Hits_);
0330
0331
0332 bool emf_loose = true;
0333 if (abs_eta <= 2.6) {
0334 if (emEnergyFraction <= 0.01)
0335 emf_loose = false;
0336 } else {
0337 if (emEnergyFraction <= -0.9)
0338 emf_loose = false;
0339 if (corrPt > 80 && emEnergyFraction >= 1)
0340 emf_loose = false;
0341 }
0342 if (ignoreCut(index_LOOSE_AOD_EMF_) || emf_loose)
0343 passCut(ret, index_LOOSE_AOD_EMF_);
0344
0345 } else if (quality_ == LOOSE || quality_ == TIGHT) {
0346
0347 if (ignoreCut(index_LOOSE_fHPD_) || jetID.fHPD < 0.98)
0348 passCut(ret, index_LOOSE_fHPD_);
0349
0350 if (ignoreCut(index_LOOSE_N90Hits_) || jetID.n90Hits > 1)
0351 passCut(ret, index_LOOSE_N90Hits_);
0352
0353
0354 bool emf_loose = true;
0355 if (abs_eta <= 2.6) {
0356 if (emEnergyFraction <= 0.01)
0357 emf_loose = false;
0358 } else {
0359 if (emEnergyFraction <= -0.9)
0360 emf_loose = false;
0361 if (corrPt > 80 && emEnergyFraction >= 1)
0362 emf_loose = false;
0363 }
0364 if (ignoreCut(index_LOOSE_EMF_) || emf_loose)
0365 passCut(ret, index_LOOSE_EMF_);
0366
0367 if (quality_ == TIGHT) {
0368
0369 bool tight_fhpd = true;
0370 if (corrPt >= 25 && jetID.fHPD >= 0.95)
0371 tight_fhpd = false;
0372 if (ignoreCut(index_TIGHT_fHPD_) || tight_fhpd)
0373 passCut(ret, index_TIGHT_fHPD_);
0374
0375
0376 bool tight_emf = true;
0377 if (abs_eta >= 1 && corrPt >= 80 && emEnergyFraction >= 1)
0378 tight_emf = false;
0379 if (abs_eta >= 2.6) {
0380 if (emEnergyFraction <= -0.3)
0381 tight_emf = false;
0382 if (abs_eta < 3.25) {
0383 if (corrPt >= 50 && emEnergyFraction <= -0.2)
0384 tight_emf = false;
0385 if (corrPt >= 80 && emEnergyFraction <= -0.1)
0386 tight_emf = false;
0387 if (corrPt >= 340 && emEnergyFraction >= 0.95)
0388 tight_emf = false;
0389 } else {
0390 if (emEnergyFraction >= 0.9)
0391 tight_emf = false;
0392 if (corrPt >= 50 && emEnergyFraction <= -0.2)
0393 tight_emf = false;
0394 if (corrPt >= 50 && emEnergyFraction >= 0.8)
0395 tight_emf = false;
0396 if (corrPt >= 130 && emEnergyFraction <= -0.1)
0397 tight_emf = false;
0398 if (corrPt >= 130 && emEnergyFraction >= 0.7)
0399 tight_emf = false;
0400
0401 }
0402 }
0403 if (ignoreCut(index_TIGHT_EMF_) || tight_emf)
0404 passCut(ret, index_TIGHT_EMF_);
0405 }
0406 }
0407
0408 setIgnored(ret);
0409
0410 return (bool)ret;
0411 }
0412
0413
0414
0415
0416 bool fwd09Cuts(reco::Candidate::LorentzVector const& rawP4,
0417 double emEnergyFraction,
0418 double etaWidth,
0419 double phiWidth,
0420 unsigned int nHit,
0421 reco::JetID const& jetID,
0422 pat::strbitset& ret) {
0423 ret.set(false);
0424
0425
0426 double abs_eta = TMath::Abs(rawP4.eta());
0427 double rawPt = rawP4.pt();
0428 double emf = emEnergyFraction;
0429 double fhf = jetID.fLong + jetID.fShort;
0430 double lnpt = (rawPt > 0) ? TMath::Log(rawPt) : -10;
0431 double lnE = (rawP4.energy() > 0) ? TMath::Log(rawP4.energy()) : -10;
0432
0433 bool HB = abs_eta < 1.0;
0434 bool EF = 2.6 <= abs_eta && abs_eta < 3.4 && 0.1 <= fhf && fhf < 0.9;
0435 bool HBHE = abs_eta < 2.6 || (abs_eta < 3.4 && fhf < 0.1);
0436 bool HF = 3.4 <= abs_eta || (2.6 <= abs_eta && 0.9 <= fhf);
0437 bool HFa = HF && abs_eta < 3.8;
0438 bool HFb = HF && !HFa;
0439
0440
0441
0442
0443
0444 if ((!HBHE) || ignoreCut(index_MINIMAL_EMF_) || emf > 0.01)
0445 passCut(ret, index_MINIMAL_EMF_);
0446
0447
0448 if ((!HBHE) || ignoreCut(index_LOOSE_AOD_fHPD_) || jetID.approximatefHPD < 0.98)
0449 passCut(ret, index_LOOSE_AOD_fHPD_);
0450
0451 if ((!HBHE) || ignoreCut(index_LOOSE_AOD_N90Hits_) || jetID.hitsInN90 > 1)
0452 passCut(ret, index_LOOSE_AOD_N90Hits_);
0453
0454
0455 if ((!HBHE) || ignoreCut(index_LOOSE_fHPD_) || jetID.fHPD < 0.98)
0456 passCut(ret, index_LOOSE_fHPD_);
0457
0458 if ((!HBHE) || ignoreCut(index_LOOSE_N90Hits_) || jetID.n90Hits > 1)
0459 passCut(ret, index_LOOSE_N90Hits_);
0460
0461
0462 if ((!HBHE) || ignoreCut(index_TIGHT_fHPD_) || rawPt < 25 || jetID.fHPD < 0.95)
0463 passCut(ret, index_TIGHT_fHPD_);
0464
0465
0466 if ((!HBHE) || ignoreCut(index_TIGHT_EMF_) || HB || rawPt < 55 || emf < 1)
0467 passCut(ret, index_TIGHT_EMF_);
0468
0469
0470
0471 if ((!EF) || ignoreCut(index_EF_N90Hits_) || jetID.n90Hits > 1 + 1.5 * TMath::Max(0., lnpt - 1.5))
0472 passCut(ret, index_EF_N90Hits_);
0473
0474 if ((!EF) || ignoreCut(index_EF_EMF_) ||
0475 emf > TMath::Max(-0.9, -0.1 - 0.05 * TMath::Power(TMath::Max(0., 5 - lnpt), 2.)))
0476 passCut(ret, index_EF_EMF_);
0477
0478
0479
0480 if ((!(EF || HF)) || ignoreCut(index_TIGHT_fls_) ||
0481 (EF && jetID.fLS < TMath::Min(0.8, 0.1 + 0.016 * TMath::Power(TMath::Max(0., 6 - lnpt), 2.5))) ||
0482 (HFa && jetID.fLS < TMath::Min(0.6, 0.05 + 0.045 * TMath::Power(TMath::Max(0., 7.5 - lnE), 2.2))) ||
0483 (HFb && jetID.fLS < TMath::Min(0.1, 0.05 + 0.07 * TMath::Power(TMath::Max(0., 7.8 - lnE), 2.))))
0484 passCut(ret, index_TIGHT_fls_);
0485
0486 if ((!(EF || HF)) || ignoreCut(index_widths_) ||
0487 (1E-10 < etaWidth && etaWidth < 0.12 && 1E-10 < phiWidth && phiWidth < 0.12))
0488 passCut(ret, index_widths_);
0489
0490
0491
0492 if ((!HF) || ignoreCut(index_LOOSE_nHit_) || (HFa && nHit > 1 + 2.4 * (lnpt - 1.)) ||
0493 (HFb && nHit > 1 + 3. * (lnpt - 1.)))
0494 passCut(ret, index_LOOSE_nHit_);
0495
0496 if ((!HF) || ignoreCut(index_LOOSE_als_) ||
0497 (emf < 0.6 + 0.05 * TMath::Power(TMath::Max(0., 9 - lnE), 1.5) &&
0498 emf > -0.2 - 0.041 * TMath::Power(TMath::Max(0., 7.5 - lnE), 2.2)))
0499 passCut(ret, index_LOOSE_als_);
0500
0501 if ((!HF) || ignoreCut(index_LOOSE_fls_) ||
0502 (HFa && jetID.fLS < TMath::Min(0.9, 0.1 + 0.05 * TMath::Power(TMath::Max(0., 7.5 - lnE), 2.2))) ||
0503 (HFb && jetID.fLS < TMath::Min(0.6, 0.1 + 0.065 * TMath::Power(TMath::Max(0., 7.5 - lnE), 2.2))))
0504 passCut(ret, index_LOOSE_fls_);
0505
0506 if ((!HF) || ignoreCut(index_LOOSE_foot_) || jetID.fHFOOT < 0.9)
0507 passCut(ret, index_LOOSE_foot_);
0508
0509 if ((!HF) || ignoreCut(index_TIGHT_nHit_) || (HFa && nHit > 1 + 2.7 * (lnpt - 0.8)) ||
0510 (HFb && nHit > 1 + 3.5 * (lnpt - 0.8)))
0511 passCut(ret, index_TIGHT_nHit_);
0512
0513 if ((!HF) || ignoreCut(index_TIGHT_als_) ||
0514 (emf < 0.5 + 0.057 * TMath::Power(TMath::Max(0., 9 - lnE), 1.5) &&
0515 emf > TMath::Max(-0.6, -0.1 - 0.026 * TMath::Power(TMath::Max(0., 8 - lnE), 2.2))))
0516 passCut(ret, index_TIGHT_als_);
0517
0518 if ((!HF) || ignoreCut(index_TIGHT_foot_) || jetID.fLS < 0.5)
0519 passCut(ret, index_TIGHT_foot_);
0520
0521 setIgnored(ret);
0522
0523 return (bool)ret;
0524 }
0525
0526 private:
0527 Version_t version_;
0528 Quality_t quality_;
0529
0530 index_type index_MINIMAL_EMF_;
0531
0532 index_type index_LOOSE_AOD_fHPD_;
0533 index_type index_LOOSE_AOD_N90Hits_;
0534 index_type index_LOOSE_AOD_EMF_;
0535
0536 index_type index_LOOSE_fHPD_;
0537 index_type index_LOOSE_N90Hits_;
0538 index_type index_LOOSE_EMF_;
0539
0540 index_type index_TIGHT_fHPD_;
0541 index_type index_TIGHT_EMF_;
0542
0543 index_type index_LOOSE_nHit_;
0544 index_type index_LOOSE_als_;
0545 index_type index_LOOSE_fls_;
0546 index_type index_LOOSE_foot_;
0547
0548 index_type index_TIGHT_nHit_;
0549 index_type index_TIGHT_als_;
0550 index_type index_TIGHT_fls_;
0551 index_type index_TIGHT_foot_;
0552 index_type index_widths_;
0553 index_type index_EF_N90Hits_;
0554 index_type index_EF_EMF_;
0555 };
0556
0557 #endif