File indexing completed on 2023-03-17 11:28:00
0001
0002
0003
0004
0005
0006
0007
0008
0009
0010
0011
0012
0013
0014
0015
0016 #include <memory>
0017 #include <chrono>
0018 #include <ctime>
0019
0020
0021 #include "FWCore/Framework/interface/Event.h"
0022 #include "FWCore/Framework/interface/MakerMacros.h"
0023
0024 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0025 #include "FWCore/Utilities/interface/InputTag.h"
0026
0027
0028 #include "DataFormats/HepMCCandidate/interface/GenParticle.h"
0029
0030
0031 #include "DataFormats/METReco/interface/GenMETCollection.h"
0032 #include "DataFormats/METReco/interface/GenMET.h"
0033
0034
0035 #include "FWCore/ServiceRegistry/interface/Service.h"
0036 #include "CommonTools/UtilAlgos/interface/TFileService.h"
0037
0038
0039 #include "DQMOffline/Trigger/interface/FunctionDefs.h"
0040
0041
0042 #include "Validation/HLTrigger/interface/HLTGenValHistCollPath.h"
0043
0044
0045 #include "DQMServices/Core/interface/DQMStore.h"
0046 #include "DQMServices/Core/interface/DQMEDAnalyzer.h"
0047
0048
0049 #include "HLTrigger/HLTcore/interface/HLTConfigProvider.h"
0050 #include "DataFormats/HLTReco/interface/TriggerEvent.h"
0051
0052
0053 #include "Validation/HLTrigger/interface/HLTGenValObject.h"
0054
0055 class HLTGenValSource : public DQMEDAnalyzer {
0056 public:
0057 explicit HLTGenValSource(const edm::ParameterSet&);
0058 ~HLTGenValSource() override = default;
0059 HLTGenValSource(const HLTGenValSource&) = delete;
0060 HLTGenValSource& operator=(const HLTGenValSource&) = delete;
0061
0062 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0063
0064 private:
0065 void analyze(const edm::Event&, const edm::EventSetup&) override;
0066 void bookHistograms(DQMStore::IBooker&, edm::Run const& run, edm::EventSetup const& c) override;
0067 void dqmBeginRun(const edm::Run&, const edm::EventSetup&) override;
0068
0069
0070 std::vector<HLTGenValObject> getObjectCollection(const edm::Event&);
0071 std::vector<HLTGenValObject> getGenParticles(const edm::Event&);
0072 reco::GenParticle getLastCopyPreFSR(reco::GenParticle part);
0073 reco::GenParticle getLastCopy(reco::GenParticle part);
0074
0075
0076
0077
0078 const edm::EDGetTokenT<reco::GenParticleCollection> genParticleToken_;
0079 const edm::EDGetTokenT<reco::GenMETCollection> genMETToken_;
0080 const edm::EDGetTokenT<reco::GenJetCollection> ak4genJetToken_;
0081 const edm::EDGetTokenT<reco::GenJetCollection> ak8genJetToken_;
0082 const edm::EDGetTokenT<trigger::TriggerEvent> trigEventToken_;
0083
0084
0085 std::string objType_;
0086 std::string dirName_;
0087 std::vector<edm::ParameterSet> histConfigs_;
0088 std::vector<edm::ParameterSet> histConfigs2D_;
0089 std::vector<edm::ParameterSet> binnings_;
0090 std::string hltProcessName_;
0091
0092
0093
0094 std::string infoString_ = "{";
0095
0096
0097 std::vector<HLTGenValHistCollPath> collectionPath_;
0098
0099
0100 HLTConfigProvider hltConfig_;
0101
0102
0103 std::vector<std::string> hltPathsToCheck_;
0104 std::vector<std::string> hltPaths;
0105 std::vector<std::string> hltPathSpecificCuts;
0106 double dR2limit_;
0107 bool doOnlyLastFilter_;
0108 };
0109
0110 HLTGenValSource::HLTGenValSource(const edm::ParameterSet& iConfig)
0111 : genParticleToken_(consumes<reco::GenParticleCollection>(
0112 iConfig.getParameterSet("inputCollections").getParameter<edm::InputTag>("genParticles"))),
0113 genMETToken_(consumes<reco::GenMETCollection>(
0114 iConfig.getParameterSet("inputCollections").getParameter<edm::InputTag>("genMET"))),
0115 ak4genJetToken_(consumes<reco::GenJetCollection>(
0116 iConfig.getParameterSet("inputCollections").getParameter<edm::InputTag>("ak4GenJets"))),
0117 ak8genJetToken_(consumes<reco::GenJetCollection>(
0118 iConfig.getParameterSet("inputCollections").getParameter<edm::InputTag>("ak8GenJets"))),
0119 trigEventToken_(consumes<trigger::TriggerEvent>(
0120 iConfig.getParameterSet("inputCollections").getParameter<edm::InputTag>("TrigEvent"))) {
0121
0122 histConfigs_ = iConfig.getParameterSetVector("histConfigs");
0123 histConfigs2D_ = iConfig.getParameterSetVector("histConfigs2D");
0124 binnings_ = iConfig.getParameterSetVector("binnings");
0125
0126
0127 dirName_ = iConfig.getParameter<std::string>("dqmDirName");
0128 objType_ = iConfig.getParameter<std::string>("objType");
0129 dR2limit_ = iConfig.getParameter<double>("dR2limit");
0130 doOnlyLastFilter_ = iConfig.getParameter<bool>("doOnlyLastFilter");
0131 hltProcessName_ = iConfig.getParameter<std::string>("hltProcessName");
0132 hltPathsToCheck_ = iConfig.getParameter<std::vector<std::string>>("hltPathsToCheck");
0133 }
0134
0135 void HLTGenValSource::dqmBeginRun(const edm::Run& iRun, const edm::EventSetup& iSetup) {
0136
0137 auto t = std::time(nullptr);
0138 auto tm = *std::localtime(&t);
0139
0140
0141 std::ostringstream timeStringStream;
0142 timeStringStream << std::put_time(&tm, "%d-%m-%Y %H-%M-%S");
0143 auto timeString = timeStringStream.str();
0144 infoString_ += "\"date & time\":\"" + timeString + "\",";
0145
0146
0147 std::string cmsswVersion = std::getenv("CMSSW_VERSION");
0148 infoString_ += std::string("\"CMSSW release\":\"") + cmsswVersion + "\",";
0149
0150
0151 bool changedConfig;
0152 if (!hltConfig_.init(iRun, iSetup, hltProcessName_, changedConfig)) {
0153 edm::LogError("HLTGenValSource") << "Initialization of HLTConfigProvider failed!";
0154 return;
0155 }
0156
0157
0158 infoString_ += std::string("\"global tag\":\"") + hltConfig_.globalTag() + "\",";
0159
0160
0161 infoString_ += std::string("\"HLT ConfDB table\":\"") + hltConfig_.tableName() + "\",";
0162
0163
0164 std::vector<std::string> notFoundPaths;
0165 for (auto const& pathToCheck : hltPathsToCheck_) {
0166
0167
0168 std::string cleanedPathToCheck;
0169 std::string pathSpecificCuts = "";
0170 if (pathToCheck.find(':') != std::string::npos) {
0171
0172 std::stringstream hltPathToCheckInputStream(pathToCheck);
0173 std::string hltPathToCheckInputSegment;
0174 std::vector<std::string> hltPathToCheckInputSeglist;
0175 while (std::getline(hltPathToCheckInputStream, hltPathToCheckInputSegment, ':')) {
0176 hltPathToCheckInputSeglist.push_back(hltPathToCheckInputSegment);
0177 }
0178
0179
0180 if (hltPathToCheckInputSeglist.size() != 2)
0181 throw cms::Exception("InputError")
0182 << "Path string can not be properly split into path and cuts: please use exactly one colon!.\n";
0183
0184
0185 cleanedPathToCheck = hltPathToCheckInputSeglist.at(0);
0186
0187
0188 pathSpecificCuts = hltPathToCheckInputSeglist.at(1);
0189
0190 } else {
0191 cleanedPathToCheck = pathToCheck;
0192 }
0193
0194 bool pathfound = false;
0195 for (auto const& pathFromConfig : hltConfig_.triggerNames()) {
0196 if (pathFromConfig.find(cleanedPathToCheck) != std::string::npos) {
0197 hltPaths.push_back(pathFromConfig);
0198
0199
0200 int count = std::count(hltPaths.begin(), hltPaths.end(), pathFromConfig);
0201 if (count > 1) {
0202 pathSpecificCuts += std::string(",autotag=v") + std::to_string(count);
0203 }
0204 hltPathSpecificCuts.push_back(pathSpecificCuts);
0205 pathfound = true;
0206 }
0207 }
0208 if (!pathfound)
0209 notFoundPaths.push_back(cleanedPathToCheck);
0210 }
0211 if (!notFoundPaths.empty()) {
0212
0213 std::string notFoundPathsMessage = "";
0214 for (const auto& path : notFoundPaths)
0215 notFoundPathsMessage += "- " + path + "\n";
0216 edm::LogError("HLTGenValSource") << "The following paths could not be found and will not be used: \n"
0217 << notFoundPathsMessage << std::endl;
0218 }
0219
0220
0221
0222
0223 edm::ParameterSet pathCollConfig;
0224 pathCollConfig.addParameter<std::string>("objType", objType_);
0225 pathCollConfig.addParameter<double>("dR2limit", dR2limit_);
0226 pathCollConfig.addParameter<bool>("doOnlyLastFilter", doOnlyLastFilter_);
0227 pathCollConfig.addParameter<std::string>("hltProcessName", hltProcessName_);
0228
0229
0230 for (const auto& path : hltPaths) {
0231 edm::ParameterSet pathCollConfigStep = pathCollConfig;
0232 pathCollConfigStep.addParameter<std::string>("triggerPath", path);
0233 collectionPath_.emplace_back(HLTGenValHistCollPath(pathCollConfigStep, hltConfig_));
0234 }
0235 }
0236
0237
0238 void HLTGenValSource::analyze(const edm::Event& iEvent, const edm::EventSetup& iSetup) {
0239
0240 const std::vector<HLTGenValObject> objects = getObjectCollection(iEvent);
0241
0242
0243 edm::Handle<trigger::TriggerEvent> triggerEvent;
0244 iEvent.getByToken(trigEventToken_, triggerEvent);
0245
0246
0247 for (const auto& object : objects) {
0248 for (auto& collection_path : collectionPath_) {
0249 collection_path.fillHists(object, triggerEvent);
0250 }
0251 }
0252 }
0253
0254
0255 void HLTGenValSource::bookHistograms(DQMStore::IBooker& iBooker, const edm::Run& run, const edm::EventSetup& setup) {
0256 iBooker.setCurrentFolder(dirName_);
0257
0258 if (infoString_.back() == ',')
0259 infoString_.pop_back();
0260 infoString_ += "}";
0261 iBooker.bookString("HLTGenValInfo", infoString_);
0262
0263
0264 for (long unsigned int i = 0; i < collectionPath_.size(); i++) {
0265 std::vector<edm::ParameterSet> histConfigs = histConfigs_;
0266 for (auto& histConfig : histConfigs) {
0267 histConfig.addParameter<std::string>("pathSpecificCuts", hltPathSpecificCuts.at(i));
0268 histConfig.addParameter<std::vector<edm::ParameterSet>>("binnings",
0269 binnings_);
0270 }
0271
0272 collectionPath_.at(i).bookHists(iBooker, histConfigs, histConfigs2D_);
0273 }
0274 }
0275
0276
0277 void HLTGenValSource::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0278 edm::ParameterSetDescription desc;
0279
0280
0281 desc.add<std::string>(
0282 "objType");
0283 desc.add<std::vector<std::string>>(
0284 "hltPathsToCheck");
0285 desc.add<std::string>("dqmDirName", "HLTGenVal");
0286 desc.add<std::string>("hltProcessName", "HLT");
0287 desc.add<double>("dR2limit", 0.1);
0288 desc.add<bool>("doOnlyLastFilter", false);
0289
0290
0291 edm::ParameterSetDescription inputCollections;
0292 inputCollections.add<edm::InputTag>("genParticles", edm::InputTag("genParticles"));
0293 inputCollections.add<edm::InputTag>("genMET", edm::InputTag("genMetTrue"));
0294 inputCollections.add<edm::InputTag>("ak4GenJets", edm::InputTag("ak4GenJets"));
0295 inputCollections.add<edm::InputTag>("ak8GenJets", edm::InputTag("ak8GenJets"));
0296 inputCollections.add<edm::InputTag>("TrigEvent", edm::InputTag("hltTriggerSummaryAOD"));
0297 desc.add<edm::ParameterSetDescription>("inputCollections", inputCollections);
0298
0299
0300
0301
0302
0303 edm::ParameterSetDescription histConfig;
0304 histConfig.add<std::string>("vsVar");
0305 histConfig.add<std::vector<double>>("binLowEdges");
0306 histConfig.addVPSet(
0307 "rangeCuts", VarRangeCut<HLTGenValObject>::makePSetDescription(), std::vector<edm::ParameterSet>());
0308
0309
0310 std::vector<edm::ParameterSet> histConfigDefaults;
0311
0312 edm::ParameterSet histConfigDefault0;
0313 histConfigDefault0.addParameter<std::string>("vsVar", "pt");
0314 std::vector<double> defaultPtBinning{0, 5, 10, 12.5, 15, 17.5, 20, 22.5, 25, 30, 35, 40,
0315 45, 50, 60, 80, 100, 150, 200, 250, 300, 350, 400};
0316 histConfigDefault0.addParameter<std::vector<double>>("binLowEdges", defaultPtBinning);
0317 histConfigDefaults.push_back(histConfigDefault0);
0318
0319 edm::ParameterSet histConfigDefault1;
0320 histConfigDefault1.addParameter<std::string>("vsVar", "eta");
0321 std::vector<double> defaultetaBinning{-10, -8, -6, -5, -4, -3, -2, -1, 0, 1, 2, 3, 4, 5, 6, 8, 10};
0322 histConfigDefault1.addParameter<std::vector<double>>("binLowEdges", defaultetaBinning);
0323 histConfigDefaults.push_back(histConfigDefault1);
0324
0325 desc.addVPSet("histConfigs", histConfig, histConfigDefaults);
0326
0327
0328 edm::ParameterSetDescription histConfig2D;
0329 histConfig2D.add<std::string>("vsVarX");
0330 histConfig2D.add<std::string>("vsVarY");
0331 histConfig2D.add<std::vector<double>>("binLowEdgesX");
0332 histConfig2D.add<std::vector<double>>("binLowEdgesY");
0333
0334
0335 std::vector<edm::ParameterSet> histConfigDefaults2D;
0336
0337 edm::ParameterSet histConfigDefault2D0;
0338 histConfigDefault2D0.addParameter<std::string>("vsVarX", "pt");
0339 histConfigDefault2D0.addParameter<std::string>("vsVarY", "eta");
0340 histConfigDefault2D0.addParameter<std::vector<double>>("binLowEdgesX", defaultPtBinning);
0341 histConfigDefault2D0.addParameter<std::vector<double>>("binLowEdgesY", defaultetaBinning);
0342 histConfigDefaults2D.push_back(histConfigDefault2D0);
0343
0344 desc.addVPSet("histConfigs2D", histConfig2D, histConfigDefaults2D);
0345
0346
0347
0348 edm::ParameterSetDescription binningConfig;
0349 binningConfig.add<std::string>("name");
0350 binningConfig.add<std::string>("vsVar");
0351 binningConfig.add<std::vector<double>>("binLowEdges");
0352
0353
0354 std::vector<edm::ParameterSet> binningConfigDefaults;
0355
0356 desc.addVPSet("binnings", binningConfig, binningConfigDefaults);
0357
0358 descriptions.addDefault(desc);
0359 }
0360
0361
0362 std::vector<HLTGenValObject> HLTGenValSource::getObjectCollection(const edm::Event& iEvent) {
0363 std::vector<HLTGenValObject> objects;
0364
0365
0366 std::vector<std::string> implementedGenParticles = {"ele", "pho", "mu", "tau"};
0367 if (std::find(implementedGenParticles.begin(), implementedGenParticles.end(), objType_) !=
0368 implementedGenParticles.end()) {
0369 objects = getGenParticles(iEvent);
0370 } else if (objType_ == "AK4jet") {
0371 const auto& genJets = iEvent.getHandle(ak4genJetToken_);
0372 for (size_t i = 0; i < genJets->size(); i++) {
0373 const reco::GenJet p = (*genJets)[i];
0374 objects.emplace_back(p);
0375 }
0376 } else if (objType_ == "AK8jet") {
0377 const auto& genJets = iEvent.getHandle(ak8genJetToken_);
0378 for (size_t i = 0; i < genJets->size(); i++) {
0379 const reco::GenJet p = (*genJets)[i];
0380 objects.emplace_back(p);
0381 }
0382 } else if (objType_ == "AK4HT") {
0383 const auto& genJets = iEvent.getHandle(ak4genJetToken_);
0384 if (!genJets->empty()) {
0385 double HTsum = 0.;
0386 for (const auto& genJet : *genJets) {
0387 if (genJet.pt() > 30 && std::abs(genJet.eta()) < 2.5)
0388 HTsum += genJet.pt();
0389 }
0390 if (HTsum > 0)
0391 objects.emplace_back(reco::Candidate::PolarLorentzVector(HTsum, 0, 0, 0));
0392 }
0393 } else if (objType_ == "AK8HT") {
0394 const auto& genJets = iEvent.getHandle(ak8genJetToken_);
0395 if (!genJets->empty()) {
0396 double HTsum = 0.;
0397 for (const auto& genJet : *genJets) {
0398 if (genJet.pt() > 200 && std::abs(genJet.eta()) < 2.5)
0399 HTsum += genJet.pt();
0400 }
0401 if (HTsum > 0)
0402 objects.emplace_back(reco::Candidate::PolarLorentzVector(HTsum, 0, 0, 0));
0403 }
0404 } else if (objType_ == "MET") {
0405 const auto& genMET = iEvent.getHandle(genMETToken_);
0406 if (!genMET->empty()) {
0407 auto genMETpt = (*genMET)[0].pt();
0408 objects.emplace_back(reco::Candidate::PolarLorentzVector(genMETpt, 0, 0, 0));
0409 }
0410 } else
0411 throw cms::Exception("InputError") << "Generator-level validation is not available for type " << objType_ << ".\n"
0412 << "Please check for a potential spelling error.\n";
0413
0414 return objects;
0415 }
0416
0417
0418 std::vector<HLTGenValObject> HLTGenValSource::getGenParticles(const edm::Event& iEvent) {
0419 std::vector<HLTGenValObject> objects;
0420
0421 const auto& genParticles = iEvent.getHandle(genParticleToken_);
0422
0423
0424 int pdgID = -1;
0425 if (objType_ == "ele")
0426 pdgID = 11;
0427 else if (objType_ == "pho")
0428 pdgID = 22;
0429 else if (objType_ == "mu")
0430 pdgID = 13;
0431 else if (objType_ == "tau")
0432 pdgID = 15;
0433
0434
0435 for (size_t i = 0; i < genParticles->size(); ++i) {
0436 const reco::GenParticle p = (*genParticles)[i];
0437
0438
0439 if (std::abs(p.pdgId()) != pdgID)
0440 continue;
0441
0442
0443 if (p.isHardProcess()) {
0444
0445 if ((objType_ == "ele") || (objType_ == "pho"))
0446 objects.emplace_back(getLastCopyPreFSR(p));
0447 else if ((objType_ == "mu") || (objType_ == "tau"))
0448 objects.emplace_back(getLastCopy(p));
0449 }
0450 }
0451
0452 return objects;
0453 }
0454
0455
0456 reco::GenParticle HLTGenValSource::getLastCopyPreFSR(reco::GenParticle part) {
0457 const auto& daughters = part.daughterRefVector();
0458 if (daughters.size() == 1 && daughters.at(0)->pdgId() == part.pdgId())
0459 return getLastCopyPreFSR(*daughters.at(0).get());
0460 else
0461 return part;
0462 }
0463
0464
0465 reco::GenParticle HLTGenValSource::getLastCopy(reco::GenParticle part) {
0466 for (const auto& daughter : part.daughterRefVector()) {
0467 if (daughter->pdgId() == part.pdgId())
0468 return getLastCopy(*daughter.get());
0469 }
0470 return part;
0471 }
0472
0473
0474 DEFINE_FWK_MODULE(HLTGenValSource);