File indexing completed on 2024-04-09 02:22:30
0001 #include "FWCore/Framework/interface/Frameworkfwd.h"
0002 #include "FWCore/Framework/interface/stream/EDProducer.h"
0003
0004 #include "FWCore/Framework/interface/Event.h"
0005 #include "FWCore/Framework/interface/MakerMacros.h"
0006
0007 #include "FWCore/ParameterSet/interface/ParameterSet.h"
0008 #include "FWCore/Utilities/interface/StreamID.h"
0009
0010 #include "RecoBTag/FeatureTools/interface/TrackInfoBuilder.h"
0011
0012 using namespace btagbtvdeep;
0013
0014 #include "DataFormats/Math/interface/deltaR.h"
0015
0016 #include "CommonTools/Utils/interface/StringCutObjectSelector.h"
0017 #include "DataFormats/NanoAOD/interface/FlatTable.h"
0018 #include <string>
0019
0020
0021 #include "FWCore/Framework/interface/makeRefToBaseProdFrom.h"
0022 #include "DataFormats/BTauReco/interface/JetTag.h"
0023 #include "DataFormats/BTauReco/interface/DeepFlavourTagInfo.h"
0024
0025
0026 #include "DataFormats/PatCandidates/interface/Jet.h"
0027
0028 template <typename T>
0029 class DeepJetTableProducer : public edm::stream::EDProducer<> {
0030 public:
0031 explicit DeepJetTableProducer(const edm::ParameterSet&);
0032 ~DeepJetTableProducer() override;
0033
0034 static void fillDescriptions(edm::ConfigurationDescriptions& descriptions);
0035
0036 private:
0037 void produce(edm::Event&, const edm::EventSetup&) override;
0038
0039 const std::string nameDeepJet_;
0040 const std::string idx_nameDeepJet_;
0041 const unsigned int n_cpf_ = 25;
0042 const unsigned int n_npf_ = 25;
0043 const unsigned int n_sv_ = 12;
0044
0045 edm::EDGetTokenT<edm::View<T>> jet_token_;
0046
0047 typedef std::vector<reco::DeepFlavourTagInfo> TagInfoCollection;
0048 const edm::EDGetTokenT<TagInfoCollection> tag_info_src_;
0049
0050 constexpr static bool usePhysForLightAndUndefined = false;
0051 };
0052
0053
0054
0055
0056 template <typename T>
0057 DeepJetTableProducer<T>::DeepJetTableProducer(const edm::ParameterSet& iConfig)
0058 : nameDeepJet_(iConfig.getParameter<std::string>("nameDeepJet")),
0059 idx_nameDeepJet_(iConfig.getParameter<std::string>("idx_nameDeepJet")),
0060 n_cpf_(iConfig.getParameter<unsigned int>("n_cpf")),
0061 n_npf_(iConfig.getParameter<unsigned int>("n_npf")),
0062 n_sv_(iConfig.getParameter<unsigned int>("n_sv")),
0063 jet_token_(consumes<edm::View<T>>(iConfig.getParameter<edm::InputTag>("jets"))),
0064 tag_info_src_(consumes<TagInfoCollection>(iConfig.getParameter<edm::InputTag>("tagInfo_src"))) {
0065 produces<nanoaod::FlatTable>(nameDeepJet_);
0066 }
0067
0068 template <typename T>
0069 DeepJetTableProducer<T>::~DeepJetTableProducer() {}
0070
0071 template <typename T>
0072 void DeepJetTableProducer<T>::produce(edm::Event& iEvent, const edm::EventSetup& iSetup) {
0073
0074
0075
0076
0077
0078 auto jets = iEvent.getHandle(jet_token_);
0079
0080 edm::Handle<TagInfoCollection> tag_infos;
0081 iEvent.getByToken(tag_info_src_, tag_infos);
0082
0083 unsigned nJets = jets->size();
0084
0085 std::vector<int> jet_N_CPFCands(nJets);
0086 std::vector<int> jet_N_NPFCands(nJets);
0087 std::vector<int> jet_N_PVs(nJets);
0088 std::vector<int> jet_N_SVs(nJets);
0089
0090
0091 std::vector<std::vector<float>> Cpfcan_BtagPf_trackEtaRel_nCpf(n_cpf_, std::vector<float>(nJets));
0092 std::vector<std::vector<float>> Cpfcan_BtagPf_trackPtRel_nCpf(n_cpf_, std::vector<float>(nJets));
0093 std::vector<std::vector<float>> Cpfcan_BtagPf_trackPPar_nCpf(n_cpf_, std::vector<float>(nJets));
0094 std::vector<std::vector<float>> Cpfcan_BtagPf_trackDeltaR_nCpf(n_cpf_, std::vector<float>(nJets));
0095 std::vector<std::vector<float>> Cpfcan_BtagPf_trackPParRatio_nCpf(n_cpf_, std::vector<float>(nJets));
0096 std::vector<std::vector<float>> Cpfcan_BtagPf_trackSip2dVal_nCpf(n_cpf_, std::vector<float>(nJets));
0097 std::vector<std::vector<float>> Cpfcan_BtagPf_trackSip2dSig_nCpf(n_cpf_, std::vector<float>(nJets));
0098 std::vector<std::vector<float>> Cpfcan_BtagPf_trackSip3dVal_nCpf(n_cpf_, std::vector<float>(nJets));
0099 std::vector<std::vector<float>> Cpfcan_BtagPf_trackSip3dSig_nCpf(n_cpf_, std::vector<float>(nJets));
0100 std::vector<std::vector<float>> Cpfcan_BtagPf_trackJetDistVal_nCpf(n_cpf_, std::vector<float>(nJets));
0101 std::vector<std::vector<float>> Cpfcan_ptrel_nCpf(n_cpf_, std::vector<float>(nJets));
0102 std::vector<std::vector<float>> Cpfcan_drminsv_nCpf(n_cpf_, std::vector<float>(nJets));
0103 std::vector<std::vector<int>> Cpfcan_VTX_ass_nCpf(n_cpf_, std::vector<int>(nJets));
0104 std::vector<std::vector<float>> Cpfcan_puppiw_nCpf(n_cpf_, std::vector<float>(nJets));
0105 std::vector<std::vector<float>> Cpfcan_chi2_nCpf(n_cpf_, std::vector<float>(nJets));
0106 std::vector<std::vector<int>> Cpfcan_quality_nCpf(n_cpf_, std::vector<int>(nJets));
0107
0108
0109 std::vector<std::vector<float>> Npfcan_ptrel_nNpf(n_npf_, std::vector<float>(nJets));
0110 std::vector<std::vector<float>> Npfcan_deltaR_nNpf(n_npf_, std::vector<float>(nJets));
0111 std::vector<std::vector<int>> Npfcan_isGamma_nNpf(n_npf_, std::vector<int>(nJets));
0112 std::vector<std::vector<float>> Npfcan_HadFrac_nNpf(n_npf_, std::vector<float>(nJets));
0113 std::vector<std::vector<float>> Npfcan_drminsv_nNpf(n_npf_, std::vector<float>(nJets));
0114 std::vector<std::vector<float>> Npfcan_puppiw_nNpf(n_npf_, std::vector<float>(nJets));
0115
0116
0117
0118
0119
0120
0121
0122 std::vector<std::vector<float>> sv_mass_nSV(n_sv_, std::vector<float>(nJets));
0123 std::vector<std::vector<float>> sv_pt_nSV(n_sv_, std::vector<float>(nJets));
0124 std::vector<std::vector<int>> sv_ntracks_nSV(n_sv_, std::vector<int>(nJets));
0125 std::vector<std::vector<float>> sv_chi2_nSV(n_sv_, std::vector<float>(nJets));
0126 std::vector<std::vector<float>> sv_normchi2_nSV(n_sv_, std::vector<float>(nJets));
0127 std::vector<std::vector<float>> sv_dxy_nSV(n_sv_, std::vector<float>(nJets));
0128 std::vector<std::vector<float>> sv_dxysig_nSV(n_sv_, std::vector<float>(nJets));
0129 std::vector<std::vector<float>> sv_d3d_nSV(n_sv_, std::vector<float>(nJets));
0130 std::vector<std::vector<float>> sv_d3dsig_nSV(n_sv_, std::vector<float>(nJets));
0131 std::vector<std::vector<float>> sv_costhetasvpv_nSV(n_sv_, std::vector<float>(nJets));
0132
0133
0134
0135
0136
0137 std::vector<std::vector<float>> sv_deltaR_nSV(n_sv_, std::vector<float>(nJets));
0138 std::vector<std::vector<float>> sv_enratio_nSV(n_sv_, std::vector<float>(nJets));
0139
0140 if (!tag_infos->empty()) {
0141 for (unsigned i_jet = 0; i_jet < nJets; ++i_jet) {
0142
0143
0144 const auto& taginfo = (*tag_infos)[i_jet];
0145 const auto& features = taginfo.features();
0146
0147
0148
0149 jet_N_CPFCands[i_jet] = features.c_pf_features.size();
0150 jet_N_NPFCands[i_jet] = features.n_pf_features.size();
0151 jet_N_SVs[i_jet] = features.sv_features.size();
0152 jet_N_PVs[i_jet] = features.npv;
0153
0154
0155 auto max_c_pf_n = std::min(features.c_pf_features.size(), (std::size_t)n_cpf_);
0156 for (std::size_t c_pf_n = 0; c_pf_n < max_c_pf_n; c_pf_n++) {
0157 const auto& c_pf_features = features.c_pf_features.at(c_pf_n);
0158 Cpfcan_BtagPf_trackEtaRel_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackEtaRel;
0159 Cpfcan_BtagPf_trackPtRel_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackPtRel;
0160 Cpfcan_BtagPf_trackPPar_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackPPar;
0161 Cpfcan_BtagPf_trackDeltaR_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackDeltaR;
0162 Cpfcan_BtagPf_trackPParRatio_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackPParRatio;
0163 Cpfcan_BtagPf_trackSip2dVal_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackSip2dVal;
0164 Cpfcan_BtagPf_trackSip2dSig_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackSip2dSig;
0165 Cpfcan_BtagPf_trackSip3dVal_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackSip3dVal;
0166 Cpfcan_BtagPf_trackSip3dSig_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackSip3dSig;
0167 Cpfcan_BtagPf_trackJetDistVal_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackJetDistVal;
0168 Cpfcan_ptrel_nCpf[c_pf_n][i_jet] = c_pf_features.ptrel;
0169 Cpfcan_drminsv_nCpf[c_pf_n][i_jet] = c_pf_features.drminsv;
0170 Cpfcan_VTX_ass_nCpf[c_pf_n][i_jet] = c_pf_features.vtx_ass;
0171 Cpfcan_puppiw_nCpf[c_pf_n][i_jet] = c_pf_features.puppiw;
0172 Cpfcan_chi2_nCpf[c_pf_n][i_jet] = c_pf_features.chi2;
0173 Cpfcan_quality_nCpf[c_pf_n][i_jet] = c_pf_features.quality;
0174 }
0175
0176
0177 auto max_n_pf_n = std::min(features.n_pf_features.size(), (std::size_t)n_npf_);
0178 for (std::size_t n_pf_n = 0; n_pf_n < max_n_pf_n; n_pf_n++) {
0179 const auto& n_pf_features = features.n_pf_features.at(n_pf_n);
0180 Npfcan_ptrel_nNpf[n_pf_n][i_jet] = n_pf_features.ptrel;
0181 Npfcan_deltaR_nNpf[n_pf_n][i_jet] = n_pf_features.deltaR;
0182 Npfcan_isGamma_nNpf[n_pf_n][i_jet] = n_pf_features.isGamma;
0183 Npfcan_HadFrac_nNpf[n_pf_n][i_jet] = n_pf_features.hadFrac;
0184 Npfcan_drminsv_nNpf[n_pf_n][i_jet] = n_pf_features.drminsv;
0185 Npfcan_puppiw_nNpf[n_pf_n][i_jet] = n_pf_features.puppiw;
0186
0187
0188
0189
0190
0191 }
0192
0193
0194 auto max_sv_n = std::min(features.sv_features.size(), (std::size_t)n_sv_);
0195 for (std::size_t sv_n = 0; sv_n < max_sv_n; sv_n++) {
0196 const auto& sv_features = features.sv_features.at(sv_n);
0197 sv_pt_nSV[sv_n][i_jet] = sv_features.pt;
0198 sv_deltaR_nSV[sv_n][i_jet] = sv_features.deltaR;
0199 sv_mass_nSV[sv_n][i_jet] = sv_features.mass;
0200 sv_ntracks_nSV[sv_n][i_jet] = sv_features.ntracks;
0201 sv_chi2_nSV[sv_n][i_jet] = sv_features.chi2;
0202 sv_normchi2_nSV[sv_n][i_jet] = sv_features.normchi2;
0203 sv_dxy_nSV[sv_n][i_jet] = sv_features.dxy;
0204 sv_dxysig_nSV[sv_n][i_jet] = sv_features.dxysig;
0205 sv_d3d_nSV[sv_n][i_jet] = sv_features.d3d;
0206 sv_d3dsig_nSV[sv_n][i_jet] = sv_features.d3dsig;
0207 sv_costhetasvpv_nSV[sv_n][i_jet] = sv_features.costhetasvpv;
0208 sv_enratio_nSV[sv_n][i_jet] = sv_features.enratio;
0209
0210
0211
0212
0213
0214 }
0215 }
0216 }
0217
0218
0219 auto djTable = std::make_unique<nanoaod::FlatTable>(jet_N_CPFCands.size(), nameDeepJet_, false, true);
0220
0221
0222 djTable->addColumn<int>("DeepJet_nCpfcand", jet_N_CPFCands, "Number of charged PF candidates in the jet");
0223 djTable->addColumn<int>("DeepJet_nNpfcand", jet_N_NPFCands, "Number of neutral PF candidates in the jet");
0224 djTable->addColumn<int>("DeepJet_nsv", jet_N_SVs, "Number of secondary vertices in the jet");
0225 djTable->addColumn<int>("DeepJet_npv", jet_N_PVs, "Number of primary vertices");
0226
0227
0228 for (unsigned int p = 0; p < n_cpf_; p++) {
0229 auto s = std::to_string(p);
0230
0231 djTable->addColumn<float>("DeepJet_Cpfcan_puppiw_" + s,
0232 Cpfcan_puppiw_nCpf[p],
0233 "charged candidate PUPPI weight of the " + s + ". cpf",
0234 10);
0235 djTable->addColumn<int>(
0236 "DeepJet_Cpfcan_VTX_ass_" + s,
0237 Cpfcan_VTX_ass_nCpf[p],
0238 "integer flag that indicates whether the track was used in the primary vertex fit for the " + s + ". cpf",
0239 10);
0240 djTable->addColumn<float>("DeepJet_Cpfcan_drminsv_" + s,
0241 Cpfcan_drminsv_nCpf[p],
0242 "track pseudoangular distance from the closest secondary vertex of the " + s + ". cpf",
0243 10);
0244 djTable->addColumn<float>("DeepJet_Cpfcan_ptrel_" + s,
0245 Cpfcan_ptrel_nCpf[p],
0246 "fraction of the jet momentum carried by the track for the " + s + ". cpf",
0247 10);
0248 djTable->addColumn<int>(
0249 "DeepJet_Cpfcan_quality_" + s,
0250 Cpfcan_quality_nCpf[p],
0251 "integer flag which indicates the quality of the fitted track, based on number of detector hits used for the "
0252 "reconstruction as well as the overall chi2 of the charged track fit for the " +
0253 s + ". cpf",
0254 10);
0255 djTable->addColumn<float>(
0256 "DeepJet_Cpfcan_chi2_" + s, Cpfcan_chi2_nCpf[p], "chi2 of the charged track fit for the " + s + ". cpf", 10);
0257
0258 djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackDeltaR_" + s,
0259 Cpfcan_BtagPf_trackDeltaR_nCpf[p],
0260 "track pseudoangular distance from the jet axis for the " + s + ". cpf",
0261 10);
0262 djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackEtaRel_" + s,
0263 Cpfcan_BtagPf_trackEtaRel_nCpf[p],
0264 "track pseudorapidity, relative to the jet axis for the " + s + ". cpf",
0265 10);
0266 djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackJetDistVal_" + s,
0267 Cpfcan_BtagPf_trackJetDistVal_nCpf[p],
0268 "minimum track approach distance to jet axis for the " + s + ". cpf",
0269 10);
0270 djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackPPar_" + s,
0271 Cpfcan_BtagPf_trackPPar_nCpf[p],
0272 "dot product of the jet and track momentum for the " + s + ". cpf",
0273 10);
0274 djTable->addColumn<float>(
0275 "DeepJet_Cpfcan_BtagPf_trackPParRatio_" + s,
0276 Cpfcan_BtagPf_trackPParRatio_nCpf[p],
0277 "dot product of the jet and track momentum divided by the magnitude of the jet momentum for the " + s + ". cpf",
0278 10);
0279 djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackPtRel_" + s,
0280 Cpfcan_BtagPf_trackPtRel_nCpf[p],
0281 "track transverse momentum, relative to the jet axis for the " + s + ". cpf",
0282 10);
0283 djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackSip2dSig_" + s,
0284 Cpfcan_BtagPf_trackSip2dSig_nCpf[p],
0285 "track 2D signed impact parameter significance for the " + s + ". cpf",
0286 10);
0287 djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackSip3dSig_" + s,
0288 Cpfcan_BtagPf_trackSip3dSig_nCpf[p],
0289 "track 3D signed impact parameter significance for the " + s + ". cpf",
0290 10);
0291 djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackSip2dVal_" + s,
0292 Cpfcan_BtagPf_trackSip2dVal_nCpf[p],
0293 "track 2D signed impact parameter for the " + s + ". cpf",
0294 10);
0295 djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackSip3dVal_" + s,
0296 Cpfcan_BtagPf_trackSip3dVal_nCpf[p],
0297 "track 3D signed impact parameter for the " + s + ". cpf",
0298 10);
0299 }
0300
0301
0302 for (unsigned int p = 0; p < n_npf_; p++) {
0303 auto s = std::to_string(p);
0304
0305 djTable->addColumn<float>("DeepJet_Npfcan_puppiw_" + s,
0306 Npfcan_puppiw_nNpf[p],
0307 "neutral candidate PUPPI weight for the " + s + ". npf",
0308 10);
0309 djTable->addColumn<float>(
0310 "DeepJet_Npfcan_deltaR_" + s,
0311 Npfcan_deltaR_nNpf[p],
0312 "pseudoangular distance between the neutral candidate and the jet axis for the " + s + ". npf",
0313 10);
0314 djTable->addColumn<float>(
0315 "DeepJet_Npfcan_drminsv_" + s,
0316 Npfcan_drminsv_nNpf[p],
0317 "pseudoangular distance between the neutral candidate and the closest secondary vertex for the " + s + ". npf",
0318 10);
0319 djTable->addColumn<float>(
0320 "DeepJet_Npfcan_HadFrac_" + s,
0321 Npfcan_HadFrac_nNpf[p],
0322 "fraction of the neutral candidate energy deposited in the hadronic calorimeter for the " + s + ". npf",
0323 10);
0324 djTable->addColumn<float>("DeepJet_Npfcan_ptrel_" + s,
0325 Npfcan_ptrel_nNpf[p],
0326 "fraction of the jet momentum carried by the neutral candidate for the " + s + ". npf",
0327 10);
0328 djTable->addColumn<int>("DeepJet_Npfcan_isGamma_" + s,
0329 Npfcan_isGamma_nNpf[p],
0330 "integer flag indicating whether the neutral candidate is a photon for the " + s + ". npf",
0331 10);
0332
0333
0334
0335
0336
0337
0338
0339
0340
0341
0342
0343 }
0344
0345
0346 for (unsigned int p = 0; p < n_sv_; p++) {
0347 auto s = std::to_string(p);
0348
0349 djTable->addColumn<float>("DeepJet_sv_mass_" + s, sv_mass_nSV[p], "SV mass of the " + s + ". SV", 10);
0350 djTable->addColumn<float>("DeepJet_sv_pt_" + s, sv_pt_nSV[p], "SV pt of the " + s + ". SV", 10);
0351 djTable->addColumn<float>(
0352 "DeepJet_sv_ntracks_" + s, sv_ntracks_nSV[p], "Number of tracks asociated to the " + s + ". SV", 10);
0353 djTable->addColumn<float>("DeepJet_sv_chi2_" + s, sv_chi2_nSV[p], "chi2 of the " + s + ". SV", 10);
0354 djTable->addColumn<float>("DeepJet_sv_normchi2_" + s, sv_normchi2_nSV[p], "chi2/dof of the " + s + ". SV", 10);
0355 djTable->addColumn<float>(
0356 "DeepJet_sv_dxy_" + s, sv_dxy_nSV[p], "2D impact parameter (flight distance) value of the " + s + ". SV", 10);
0357 djTable->addColumn<float>("DeepJet_sv_dxysig_" + s,
0358 sv_dxysig_nSV[p],
0359 "2D impact parameter (flight distance) significance of the " + s + ". SV",
0360 10);
0361 djTable->addColumn<float>(
0362 "DeepJet_sv_d3d_" + s, sv_d3d_nSV[p], "3D impact parameter (flight distance) value of the " + s + ". SV", 10);
0363 djTable->addColumn<float>("DeepJet_sv_d3dsig_" + s,
0364 sv_d3dsig_nSV[p],
0365 "3D impact parameter (flight distance) significance of the " + s + ". SV",
0366 10);
0367 djTable->addColumn<float>("DeepJet_sv_costhetasvpv_" + s,
0368 sv_costhetasvpv_nSV[p],
0369 "cosine of the angle between the " + s +
0370 ". SV flight direction and the direction of the " + s + ". SV momentum",
0371 10);
0372
0373
0374
0375
0376
0377
0378
0379
0380
0381
0382
0383 djTable->addColumn<float>("DeepJet_sv_deltaR_" + s,
0384 sv_deltaR_nSV[p],
0385 "pseudoangular distance between jet axis and the " + s + ". SV direction",
0386 10);
0387 djTable->addColumn<float>(
0388 "DeepJet_sv_enratio_" + s, sv_enratio_nSV[p], "ratio of the " + s + ". SV energy ratio to the jet energy", 10);
0389 }
0390
0391 iEvent.put(std::move(djTable), nameDeepJet_);
0392 }
0393
0394 template <typename T>
0395 void DeepJetTableProducer<T>::fillDescriptions(edm::ConfigurationDescriptions& descriptions) {
0396 edm::ParameterSetDescription desc;
0397 desc.add<std::string>("nameDeepJet", "Jet");
0398 desc.add<std::string>("idx_nameDeepJet", "djIdx");
0399
0400 desc.add<unsigned int>("n_cpf", 3);
0401 desc.add<unsigned int>("n_npf", 3);
0402 desc.add<unsigned int>("n_sv", 4);
0403 desc.add<edm::InputTag>("jets", edm::InputTag("slimmedJetsPuppi"));
0404 desc.add<edm::InputTag>("tagInfo_src", edm::InputTag("pfDeepFlavourTagInfosPuppiWithDeepInfo"));
0405 descriptions.addWithDefaultLabel(desc);
0406 }
0407
0408 typedef DeepJetTableProducer<pat::Jet> PatJetDeepJetTableProducer;
0409
0410 DEFINE_FWK_MODULE(PatJetDeepJetTableProducer);