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