File indexing completed on 2025-01-31 02:19:53
0001 #include <limits>
0002 #include <memory>
0003
0004 #include <random>
0005
0006 #include "CondFormats/DataRecord/interface/BTauGenericMVAJetTagComputerRcd.h"
0007 #include "CondFormats/DataRecord/interface/GBRWrapperRcd.h"
0008 #include "RecoBTau/JetTagComputer/interface/JetTagComputerRecord.h"
0009 #include "DataFormats/BTauReco/interface/SoftLeptonTagInfo.h"
0010 #include "RecoBTag/SoftLepton/interface/LeptonSelector.h"
0011 #include "RecoBTag/SoftLepton/interface/ElectronTagger.h"
0012 #include "DataFormats/BTauReco/interface/CandSoftLeptonTagInfo.h"
0013 #include <iostream>
0014
0015 ElectronTagger::Tokens::Tokens(const edm::ParameterSet& cfg, edm::ESConsumesCollector&& cc) {
0016 if (cfg.getParameter<bool>("useCondDB")) {
0017 gbrForest_ = cc.consumes(edm::ESInputTag{"", cfg.getParameter<std::string>("gbrForestLabel")});
0018 }
0019 }
0020
0021 ElectronTagger::ElectronTagger(const edm::ParameterSet& cfg, Tokens tokens)
0022 : m_selector(cfg),
0023 m_weightFile(cfg.getParameter<edm::FileInPath>("weightFile")),
0024 m_useGBRForest(cfg.getParameter<bool>("useGBRForest")),
0025 m_useAdaBoost(cfg.getParameter<bool>("useAdaBoost")),
0026 m_tokens{tokens} {
0027 uses("seTagInfos");
0028 mvaID = std::make_unique<TMVAEvaluator>();
0029 }
0030
0031 void ElectronTagger::initialize(const JetTagComputerRecord& record) {
0032
0033 std::vector<std::string> variables({"sip3d", "sip2d", "ptRel", "deltaR", "ratio", "mva_e_pi"});
0034 std::vector<std::string> spectators;
0035
0036 if (m_tokens.gbrForest_.isInitialized()) {
0037 mvaID->initializeGBRForest(&record.get(m_tokens.gbrForest_), variables, spectators, m_useAdaBoost);
0038 } else
0039 mvaID->initialize(
0040 "Color:Silent:Error", "BDT", m_weightFile.fullPath(), variables, spectators, m_useGBRForest, m_useAdaBoost);
0041 }
0042
0043
0044 float ElectronTagger::discriminator(const TagInfoHelper& tagInfo) const {
0045
0046 float bestTag = -std::numeric_limits<float>::infinity();
0047 const reco::CandSoftLeptonTagInfo& info = tagInfo.get<reco::CandSoftLeptonTagInfo>();
0048
0049 std::mt19937_64 random;
0050 std::uniform_real_distribution<float> dist(0.f, 1.f);
0051
0052
0053 for (unsigned int i = 0; i < info.leptons(); i++) {
0054 const reco::SoftLeptonProperties& properties = info.properties(i);
0055 if (m_selector(properties)) {
0056 int theSeed = 1 + round(10000.0 * std::abs(properties.deltaR));
0057 random.seed(theSeed);
0058 float rndm = dist(random);
0059
0060 float sip3dsig = (m_selector.isNegative() && rndm < 0.5) ? -properties.sip3dsig : properties.sip3dsig;
0061 float sip2dsig = (m_selector.isNegative() && rndm < 0.5) ? -properties.sip2dsig : properties.sip2dsig;
0062
0063 std::map<std::string, float> inputs;
0064 inputs["sip3d"] = sip3dsig;
0065 inputs["sip2d"] = sip2dsig;
0066 inputs["ptRel"] = properties.ptRel;
0067 inputs["deltaR"] = properties.deltaR;
0068 inputs["ratio"] = properties.ratio;
0069 inputs["mva_e_pi"] = properties.elec_mva;
0070
0071 float tag = mvaID->evaluate(inputs);
0072
0073 tag = (tag + 1.0) / 2.0;
0074 if (tag > bestTag)
0075 bestTag = tag;
0076 }
0077 }
0078 return bestTag;
0079 }
0080
0081 void ElectronTagger::fillPSetDescription(edm::ParameterSetDescription& desc) {
0082 btag::LeptonSelector::fillPSetDescription(desc);
0083 desc.add<bool>("useCondDB", false);
0084 desc.add<std::string>("gbrForestLabel", "");
0085 desc.add<edm::FileInPath>("weightFile", edm::FileInPath());
0086 desc.add<bool>("useGBRForest", false);
0087 desc.add<bool>("useAdaBoost", false);
0088 }