Back to home page

Project CMSSW displayed by LXR

 
 

    


File indexing completed on 2022-09-21 02:12:37

0001 /** \class HLTPrescaleProvider
0002  *
0003  * See header file for documentation
0004  *
0005  *
0006  *  \author Martin Grunewald
0007  *
0008  */
0009 
0010 #include "HLTrigger/HLTcore/interface/HLTPrescaleProvider.h"
0011 #include "FWCore/MessageLogger/interface/MessageLogger.h"
0012 #include "FWCore/Utilities/interface/Exception.h"
0013 #include "FWCore/ParameterSet/interface/ParameterSetDescription.h"
0014 
0015 #include <cassert>
0016 #include <sstream>
0017 
0018 static const bool useL1EventSetup(true);
0019 static const bool useL1GtTriggerMenuLite(false);
0020 static const unsigned char countMax(2);
0021 
0022 bool HLTPrescaleProvider::init(const edm::Run& iRun,
0023                                const edm::EventSetup& iSetup,
0024                                const std::string& processName,
0025                                bool& changed) {
0026   inited_ = true;
0027 
0028   count_[0] = 0;
0029   count_[1] = 0;
0030   count_[2] = 0;
0031   count_[3] = 0;
0032   count_[4] = 0;
0033 
0034   const bool result(hltConfigProvider_.init(iRun, iSetup, processName, changed));
0035 
0036   const unsigned int l1tType(hltConfigProvider_.l1tType());
0037   if (l1tType == 1) {
0038     checkL1GtUtils();
0039     /// L1 GTA V3: https://twiki.cern.ch/twiki/bin/view/CMSPublic/SWGuideL1TriggerL1GtUtils#Version_3
0040     l1GtUtils_->getL1GtRunCache(iRun, iSetup, useL1EventSetup, useL1GtTriggerMenuLite);
0041   } else if (l1tType == 2) {
0042     checkL1TGlobalUtil();
0043     l1tGlobalUtil_->retrieveL1Setup(iSetup);
0044   } else {
0045     edm::LogError("HLTPrescaleProvider") << " Unknown L1T Type " << l1tType << " - prescales will not be avaiable!";
0046   }
0047 
0048   return result;
0049 }
0050 
0051 L1GtUtils const& HLTPrescaleProvider::l1GtUtils() const {
0052   checkL1GtUtils();
0053   return *l1GtUtils_;
0054 }
0055 
0056 l1t::L1TGlobalUtil const& HLTPrescaleProvider::l1tGlobalUtil() const {
0057   checkL1TGlobalUtil();
0058   return *l1tGlobalUtil_;
0059 }
0060 
0061 int HLTPrescaleProvider::prescaleSet(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0062   if (!inited_) {
0063     throw cms::Exception("LogicError") << "HLTPrescaleProvider::prescaleSet,\n"
0064                                           "HLTPrescaleProvider::init was not called at beginRun\n";
0065   }
0066   const unsigned int l1tType(hltConfigProvider_.l1tType());
0067   if (l1tType == 1) {
0068     checkL1GtUtils();
0069 
0070     // return hltPrescaleTable_.set();
0071     l1GtUtils_->getL1GtRunCache(iEvent, iSetup, useL1EventSetup, useL1GtTriggerMenuLite);
0072     int errorTech(0);
0073     const int psfsiTech(l1GtUtils_->prescaleFactorSetIndex(iEvent, L1GtUtils::TechnicalTrigger, errorTech));
0074     int errorPhys(0);
0075     const int psfsiPhys(l1GtUtils_->prescaleFactorSetIndex(iEvent, L1GtUtils::AlgorithmTrigger, errorPhys));
0076     assert(psfsiTech == psfsiPhys);
0077     if ((errorTech == 0) && (errorPhys == 0) && (psfsiTech >= 0) && (psfsiPhys >= 0) && (psfsiTech == psfsiPhys)) {
0078       return psfsiPhys;
0079     } else {
0080       /// error - notify user!
0081       if (count_[0] < countMax) {
0082         count_[0] += 1;
0083         edm::LogError("HLTPrescaleProvider")
0084             << " Using processName '" << hltConfigProvider_.processName() << "':"
0085             << " Error in determining HLT prescale set index from L1 data using L1GtUtils:"
0086             << " Tech/Phys error = " << errorTech << "/" << errorPhys << " Tech/Phys psfsi = " << psfsiTech << "/"
0087             << psfsiPhys;
0088       }
0089       return -1;
0090     }
0091   } else if (l1tType == 2) {
0092     checkL1TGlobalUtil();
0093     l1tGlobalUtil_->retrieveL1Event(iEvent, iSetup);
0094     return static_cast<int>(l1tGlobalUtil_->prescaleColumn());
0095   } else {
0096     if (count_[0] < countMax) {
0097       count_[0] += 1;
0098       edm::LogError("HLTPrescaleProvider")
0099           << " Using processName '" << hltConfigProvider_.processName() << "':"
0100           << " Unknown L1T Type " << l1tType << " - can not determine prescale set index!";
0101     }
0102     return -1;
0103   }
0104 }
0105 
0106 template <>
0107 FractionalPrescale HLTPrescaleProvider::convertL1PS(double val) const {
0108   int numer = static_cast<int>(val * kL1PrescaleDenominator_ + 0.5);
0109   static constexpr double kL1RoundingEpsilon = 0.001;
0110   if (std::abs(numer - val * kL1PrescaleDenominator_) > kL1RoundingEpsilon) {
0111     edm::LogWarning("ValueError") << " Error, L1 prescale val " << val
0112                                   << " does not appear to precisely expressable as int / " << kL1PrescaleDenominator_
0113                                   << ", using a FractionalPrescale is a loss of precision";
0114   }
0115 
0116   return {numer, kL1PrescaleDenominator_};
0117 }
0118 
0119 double HLTPrescaleProvider::getL1PrescaleValue(const edm::Event& iEvent,
0120                                                const edm::EventSetup& iSetup,
0121                                                const std::string& trigger) {
0122   // get L1T prescale - works only for those hlt trigger paths with
0123   // exactly one L1GT seed module which has exactly one L1T name as seed
0124 
0125   double result = -1;
0126 
0127   const unsigned int l1tType(hltConfigProvider_.l1tType());
0128   if (l1tType == 1) {
0129     checkL1GtUtils();
0130     const unsigned int nL1GTSeedModules(hltConfigProvider_.hltL1GTSeeds(trigger).size());
0131     if (nL1GTSeedModules == 0) {
0132       // no L1 seed module on path hence no L1 seed hence formally no L1 prescale
0133       result = 1;
0134     } else if (nL1GTSeedModules == 1) {
0135       l1GtUtils_->getL1GtRunCache(iEvent, iSetup, useL1EventSetup, useL1GtTriggerMenuLite);
0136       const std::string l1tname(hltConfigProvider_.hltL1GTSeeds(trigger).at(0).second);
0137       int l1error(0);
0138       result = l1GtUtils_->prescaleFactor(iEvent, l1tname, l1error);
0139       if (l1error != 0) {
0140         if (count_[1] < countMax) {
0141           count_[1] += 1;
0142           edm::LogError("HLTPrescaleProvider")
0143               << " Error in determining L1T prescale for HLT path: '" << trigger << "' with L1T seed: '" << l1tname
0144               << "' using L1GtUtils: error code = " << l1error << "." << std::endl
0145               << " Note: for this method ('prescaleValues'), only a single L1T name (and not a bit number) is allowed "
0146                  "as seed!"
0147               << std::endl
0148               << " For seeds being complex logical expressions, try the new method 'prescaleValuesInDetail'."
0149               << std::endl;
0150         }
0151         result = -1;
0152       }
0153     } else {
0154       /// error - can't handle properly multiple L1GTSeed modules
0155       if (count_[2] < countMax) {
0156         count_[2] += 1;
0157         std::string dump("'" + hltConfigProvider_.hltL1GTSeeds(trigger).at(0).second + "'");
0158         for (unsigned int i = 1; i != nL1GTSeedModules; ++i) {
0159           dump += " * '" + hltConfigProvider_.hltL1GTSeeds(trigger).at(i).second + "'";
0160         }
0161         edm::LogError("HLTPrescaleProvider")
0162             << " Error in determining L1T prescale for HLT path: '" << trigger << "' has multiple L1GTSeed modules, "
0163             << nL1GTSeedModules << ", with L1 seeds: " << dump
0164             << ". (Note: at most one L1GTSeed module is allowed for a proper determination of the L1T prescale!)";
0165       }
0166       result = -1;
0167     }
0168   } else if (l1tType == 2) {
0169     checkL1TGlobalUtil();
0170     const unsigned int nL1TSeedModules(hltConfigProvider_.hltL1TSeeds(trigger).size());
0171     if (nL1TSeedModules == 0) {
0172       // no L1 seed module on path hence no L1 seed hence formally no L1 prescale
0173       result = 1;
0174     } else if (nL1TSeedModules == 1) {
0175       //    l1tGlobalUtil_->retrieveL1Event(iEvent,iSetup);
0176       const std::string l1tname(hltConfigProvider_.hltL1TSeeds(trigger).at(0));
0177       bool l1error(!l1tGlobalUtil_->getPrescaleByName(l1tname, result));
0178       if (l1error) {
0179         if (count_[1] < countMax) {
0180           count_[1] += 1;
0181           edm::LogError("HLTPrescaleProvider")
0182               << " Error in determining L1T prescale for HLT path: '" << trigger << "' with L1T seed: '" << l1tname
0183               << "' using L1TGlobalUtil: error cond = " << l1error << "." << std::endl
0184               << " Note: for this method ('prescaleValues'), only a single L1T name (and not a bit number) is allowed "
0185                  "as seed!"
0186               << std::endl
0187               << " For seeds being complex logical expressions, try the new method 'prescaleValuesInDetail'."
0188               << std::endl;
0189         }
0190         result = -1;
0191       }
0192     } else {
0193       /// error - can't handle properly multiple L1TSeed modules
0194       if (count_[2] < countMax) {
0195         count_[2] += 1;
0196         std::string dump("'" + hltConfigProvider_.hltL1TSeeds(trigger).at(0) + "'");
0197         for (unsigned int i = 1; i != nL1TSeedModules; ++i) {
0198           dump += " * '" + hltConfigProvider_.hltL1TSeeds(trigger).at(i) + "'";
0199         }
0200         edm::LogError("HLTPrescaleProvider")
0201             << " Error in determining L1T prescale for HLT path: '" << trigger << "' has multiple L1TSeed modules, "
0202             << nL1TSeedModules << ", with L1T seeds: " << dump
0203             << ". (Note: at most one L1TSeed module is allowed for a proper determination of the L1T prescale!)";
0204       }
0205       result = -1;
0206     }
0207   } else {
0208     if (count_[1] < countMax) {
0209       count_[1] += 1;
0210       edm::LogError("HLTPrescaleProvider") << " Unknown L1T Type " << l1tType << " - can not determine L1T prescale! ";
0211     }
0212     result = -1;
0213   }
0214 
0215   return result;
0216 }
0217 
0218 std::vector<std::pair<std::string, double>> HLTPrescaleProvider::getL1PrescaleValueInDetail(
0219     const edm::Event& iEvent, const edm::EventSetup& iSetup, const std::string& trigger) {
0220   std::vector<std::pair<std::string, double>> result;
0221 
0222   const unsigned int l1tType(hltConfigProvider_.l1tType());
0223   if (l1tType == 1) {
0224     checkL1GtUtils();
0225 
0226     const unsigned int nL1GTSeedModules(hltConfigProvider_.hltL1GTSeeds(trigger).size());
0227     if (nL1GTSeedModules == 0) {
0228       // no L1 seed module on path hence no L1 seed hence formally no L1 prescale
0229       result.clear();
0230     } else if (nL1GTSeedModules == 1) {
0231       l1GtUtils_->getL1GtRunCache(iEvent, iSetup, useL1EventSetup, useL1GtTriggerMenuLite);
0232       const std::string l1tname(hltConfigProvider_.hltL1GTSeeds(trigger).at(0).second);
0233       L1GtUtils::LogicalExpressionL1Results l1Logical(l1tname, *l1GtUtils_);
0234       l1Logical.logicalExpressionRunUpdate(iEvent.getRun(), iSetup, l1tname);
0235       const std::vector<std::pair<std::string, int>>& errorCodes(l1Logical.errorCodes(iEvent));
0236       auto resultInt = l1Logical.prescaleFactors();
0237       result.clear();
0238       for (const auto& entry : resultInt) {
0239         result.push_back(entry);
0240       }
0241 
0242       int l1error(l1Logical.isValid() ? 0 : 1);
0243       for (auto const& errorCode : errorCodes) {
0244         l1error += std::abs(errorCode.second);
0245       }
0246       if (l1error != 0) {
0247         if (count_[3] < countMax) {
0248           count_[3] += 1;
0249           std::ostringstream message;
0250           message << " Error in determining L1T prescales for HLT path: '" << trigger << "' with complex L1T seed: '"
0251                   << l1tname << "' using L1GtUtils: " << std::endl
0252                   << " isValid=" << l1Logical.isValid() << " l1tname/error/prescale " << errorCodes.size() << std::endl;
0253           for (unsigned int i = 0; i < errorCodes.size(); ++i) {
0254             message << " " << i << ":" << errorCodes[i].first << "/" << errorCodes[i].second << "/" << result[i].second;
0255           }
0256           message << ".";
0257           edm::LogError("HLTPrescaleProvider") << message.str();
0258         }
0259         result.clear();
0260       }
0261     } else {
0262       /// error - can't handle properly multiple L1GTSeed modules
0263       if (count_[4] < countMax) {
0264         count_[4] += 1;
0265         std::string dump("'" + hltConfigProvider_.hltL1GTSeeds(trigger).at(0).second + "'");
0266         for (unsigned int i = 1; i != nL1GTSeedModules; ++i) {
0267           dump += " * '" + hltConfigProvider_.hltL1GTSeeds(trigger).at(i).second + "'";
0268         }
0269         edm::LogError("HLTPrescaleProvider")
0270             << " Error in determining L1T prescale for HLT path: '" << trigger << "' has multiple L1GTSeed modules, "
0271             << nL1GTSeedModules << ", with L1 seeds: " << dump
0272             << ". (Note: at most one L1GTSeed module is allowed for a proper determination of the L1T prescale!)";
0273       }
0274       result.clear();
0275     }
0276   } else if (l1tType == 2) {
0277     checkL1TGlobalUtil();
0278     const unsigned int nL1TSeedModules(hltConfigProvider_.hltL1TSeeds(trigger).size());
0279     if (nL1TSeedModules == 0) {
0280       // no L1 seed module on path hence no L1 seed hence formally no L1 prescale
0281       result.clear();
0282     } else if (nL1TSeedModules == 1) {
0283       //    l1tGlobalUtil_->retrieveL1Event(iEvent,iSetup);
0284       std::string l1tname(hltConfigProvider_.hltL1TSeeds(trigger).at(0));
0285       GlobalLogicParser l1tGlobalLogicParser = GlobalLogicParser(l1tname);
0286       const std::vector<GlobalLogicParser::OperandToken> l1tSeeds = l1tGlobalLogicParser.expressionSeedsOperandList();
0287       int l1error(0);
0288       double l1tPrescale(-1);
0289       for (auto const& i : l1tSeeds) {
0290         const string& l1tSeed = i.tokenName;
0291         if (!l1tGlobalUtil_->getPrescaleByName(l1tSeed, l1tPrescale)) {
0292           l1error += 1;
0293         }
0294         result.push_back(std::pair<std::string, double>(l1tSeed, l1tPrescale));
0295       }
0296       if (l1error != 0) {
0297         if (count_[3] < countMax) {
0298           count_[3] += 1;
0299           string l1name = l1tname;
0300           std::ostringstream message;
0301           message << " Error in determining L1T prescales for HLT path: '" << trigger << "' with complex L1T seed: '"
0302                   << l1tname << "' using L1TGlobalUtil: " << std::endl
0303                   << " isValid=" << l1tGlobalLogicParser.checkLogicalExpression(l1name) << " l1tname/error/prescale "
0304                   << l1tSeeds.size() << std::endl;
0305           for (unsigned int i = 0; i < l1tSeeds.size(); ++i) {
0306             const string& l1tSeed = l1tSeeds[i].tokenName;
0307             message << " " << i << ":" << l1tSeed << "/" << l1tGlobalUtil_->getPrescaleByName(l1tSeed, l1tPrescale)
0308                     << "/" << result[i].second;
0309           }
0310           message << ".";
0311           edm::LogError("HLTPrescaleProvider") << message.str();
0312         }
0313         result.clear();
0314       }
0315     } else {
0316       /// error - can't handle properly multiple L1TSeed modules
0317       if (count_[4] < countMax) {
0318         count_[4] += 1;
0319         std::string dump("'" + hltConfigProvider_.hltL1TSeeds(trigger).at(0) + "'");
0320         for (unsigned int i = 1; i != nL1TSeedModules; ++i) {
0321           dump += " * '" + hltConfigProvider_.hltL1TSeeds(trigger).at(i) + "'";
0322         }
0323         edm::LogError("HLTPrescaleProvider")
0324             << " Error in determining L1T prescale for HLT path: '" << trigger << "' has multiple L1TSeed modules, "
0325             << nL1TSeedModules << ", with L1T seeds: " << dump
0326             << ". (Note: at most one L1TSeed module is allowed for a proper determination of the L1T prescale!)";
0327       }
0328       result.clear();
0329     }
0330   } else {
0331     if (count_[3] < countMax) {
0332       count_[3] += 1;
0333       edm::LogError("HLTPrescaleProvider") << " Unknown L1T Type " << l1tType << " - can not determine L1T prescale! ";
0334     }
0335     result.clear();
0336   }
0337 
0338   return result;
0339 }
0340 
0341 bool HLTPrescaleProvider::rejectedByHLTPrescaler(const edm::TriggerResults& triggerResults, unsigned int i) const {
0342   return hltConfigProvider_.moduleType(hltConfigProvider_.moduleLabel(i, triggerResults.index(i))) == "HLTPrescaler";
0343 }
0344 
0345 void HLTPrescaleProvider::checkL1GtUtils() const {
0346   if (!l1GtUtils_) {
0347     throw cms::Exception("Configuration") << "HLTPrescaleProvider::checkL1GtUtils(),\n"
0348                                              "Attempt to use L1GtUtils object when none was constructed.\n"
0349                                              "Possibly the proper era is not configured or\n"
0350                                              "the module configuration does not use the era properly\n"
0351                                              "or input is from mixed eras";
0352   }
0353 }
0354 
0355 void HLTPrescaleProvider::checkL1TGlobalUtil() const {
0356   if (!l1tGlobalUtil_) {
0357     throw cms::Exception("Configuration") << "HLTPrescaleProvider:::checkL1TGlobalUtil(),\n"
0358                                              "Attempt to use L1TGlobalUtil object when none was constructed.\n"
0359                                              "Possibly the proper era is not configured or\n"
0360                                              "the module configuration does not use the era properly\n"
0361                                              "or input is from mixed eras";
0362   }
0363 }
0364 
0365 void HLTPrescaleProvider::fillPSetDescription(edm::ParameterSetDescription& desc,
0366                                               unsigned int stageL1Trigger,
0367                                               edm::InputTag const& l1tAlgBlkInputTag,
0368                                               edm::InputTag const& l1tExtBlkInputTag,
0369                                               bool readPrescalesFromFile) {
0370   desc.add<unsigned int>("stageL1Trigger", stageL1Trigger);
0371   L1GtUtils::fillDescription(desc);
0372   l1t::L1TGlobalUtil::fillDescription(desc, l1tAlgBlkInputTag, l1tExtBlkInputTag, readPrescalesFromFile);
0373 }