CMS 3D CMS Logo

List of all members | Public Member Functions | Static Public Member Functions | Private Types | Private Member Functions | Private Attributes | Static Private Attributes
DeepJetTableProducer< T > Class Template Reference
Inheritance diagram for DeepJetTableProducer< T >:
edm::stream::EDProducer<>

Public Member Functions

 DeepJetTableProducer (const edm::ParameterSet &)
 
 ~DeepJetTableProducer () override
 
- Public Member Functions inherited from edm::stream::EDProducer<>
 EDProducer ()=default
 
 EDProducer (const EDProducer &)=delete
 
bool hasAbilityToProduceInBeginLumis () const final
 
bool hasAbilityToProduceInBeginProcessBlocks () const final
 
bool hasAbilityToProduceInBeginRuns () const final
 
bool hasAbilityToProduceInEndLumis () const final
 
bool hasAbilityToProduceInEndProcessBlocks () const final
 
bool hasAbilityToProduceInEndRuns () const final
 
const EDProduceroperator= (const EDProducer &)=delete
 

Static Public Member Functions

static void fillDescriptions (edm::ConfigurationDescriptions &descriptions)
 

Private Types

typedef std::vector< reco::DeepFlavourTagInfoTagInfoCollection
 

Private Member Functions

void produce (edm::Event &, const edm::EventSetup &) override
 

Private Attributes

const std::string idx_nameDeepJet_
 
edm::EDGetTokenT< edm::View< T > > jet_token_
 
const unsigned int n_cpf_ = 25
 
const unsigned int n_npf_ = 25
 
const unsigned int n_sv_ = 12
 
const std::string nameDeepJet_
 
const edm::EDGetTokenT< TagInfoCollectiontag_info_src_
 

Static Private Attributes

static constexpr bool usePhysForLightAndUndefined = false
 

Additional Inherited Members

- Public Types inherited from edm::stream::EDProducer<>
using CacheTypes = CacheContexts< T... >
 
using GlobalCache = typename CacheTypes::GlobalCache
 
using HasAbility = AbilityChecker< T... >
 
using InputProcessBlockCache = typename CacheTypes::InputProcessBlockCache
 
using LuminosityBlockCache = typename CacheTypes::LuminosityBlockCache
 
using LuminosityBlockContext = LuminosityBlockContextT< LuminosityBlockCache, RunCache, GlobalCache >
 
using LuminosityBlockSummaryCache = typename CacheTypes::LuminosityBlockSummaryCache
 
using RunCache = typename CacheTypes::RunCache
 
using RunContext = RunContextT< RunCache, GlobalCache >
 
using RunSummaryCache = typename CacheTypes::RunSummaryCache
 

Detailed Description

template<typename T>
class DeepJetTableProducer< T >

Definition at line 29 of file DeepJetTableProducer.cc.

Member Typedef Documentation

◆ TagInfoCollection

template<typename T >
typedef std::vector<reco::DeepFlavourTagInfo> DeepJetTableProducer< T >::TagInfoCollection
private

Definition at line 47 of file DeepJetTableProducer.cc.

Constructor & Destructor Documentation

◆ DeepJetTableProducer()

template<typename T >
DeepJetTableProducer< T >::DeepJetTableProducer ( const edm::ParameterSet iConfig)
explicit

Definition at line 57 of file DeepJetTableProducer.cc.

References DeepJetTableProducer< T >::nameDeepJet_.

58  : nameDeepJet_(iConfig.getParameter<std::string>("nameDeepJet")),
59  idx_nameDeepJet_(iConfig.getParameter<std::string>("idx_nameDeepJet")),
60  n_cpf_(iConfig.getParameter<unsigned int>("n_cpf")),
61  n_npf_(iConfig.getParameter<unsigned int>("n_npf")),
62  n_sv_(iConfig.getParameter<unsigned int>("n_sv")),
63  jet_token_(consumes<edm::View<T>>(iConfig.getParameter<edm::InputTag>("jets"))),
64  tag_info_src_(consumes<TagInfoCollection>(iConfig.getParameter<edm::InputTag>("tagInfo_src"))) {
65  produces<nanoaod::FlatTable>(nameDeepJet_);
66 }
T getParameter(std::string const &) const
Definition: ParameterSet.h:307
const unsigned int n_npf_
const std::string nameDeepJet_
edm::EDGetTokenT< edm::View< T > > jet_token_
const unsigned int n_cpf_
const std::string idx_nameDeepJet_
const unsigned int n_sv_
const edm::EDGetTokenT< TagInfoCollection > tag_info_src_

◆ ~DeepJetTableProducer()

template<typename T >
DeepJetTableProducer< T >::~DeepJetTableProducer ( )
override

Definition at line 69 of file DeepJetTableProducer.cc.

69 {}

Member Function Documentation

◆ fillDescriptions()

template<typename T >
void DeepJetTableProducer< T >::fillDescriptions ( edm::ConfigurationDescriptions descriptions)
static

Definition at line 395 of file DeepJetTableProducer.cc.

References edm::ConfigurationDescriptions::addWithDefaultLabel(), submitPVResolutionJobs::desc, ProducerED_cfi::InputTag, and AlCaHLTBitMon_QueryRunRegistry::string.

395  {
397  desc.add<std::string>("nameDeepJet", "Jet");
398  desc.add<std::string>("idx_nameDeepJet", "djIdx");
399 
400  desc.add<unsigned int>("n_cpf", 3);
401  desc.add<unsigned int>("n_npf", 3);
402  desc.add<unsigned int>("n_sv", 4);
403  desc.add<edm::InputTag>("jets", edm::InputTag("slimmedJetsPuppi"));
404  desc.add<edm::InputTag>("tagInfo_src", edm::InputTag("pfDeepFlavourTagInfosPuppiWithDeepInfo"));
405  descriptions.addWithDefaultLabel(desc);
406 }
void addWithDefaultLabel(ParameterSetDescription const &psetDescription)

◆ produce()

template<typename T >
void DeepJetTableProducer< T >::produce ( edm::Event iEvent,
const edm::EventSetup iSetup 
)
overrideprivate

Definition at line 72 of file DeepJetTableProducer.cc.

References lowptgsfeleseed::features(), iEvent, PDWG_EXODelayedJetMET_cff::jets, SiStripPI::min, eostools::move(), nJets, AlCaHLTBitMon_ParallelJobs::p, alignCSCRings::s, and to_string().

72  {
73  // elements in all these collections must have the same order!
74 
75  // only necessary to explicitly check correct matching of jets
76  // std::vector<int> jetIdx_dj;
77 
78  auto jets = iEvent.getHandle(jet_token_);
79 
81  iEvent.getByToken(tag_info_src_, tag_infos);
82 
83  unsigned nJets = jets->size();
84 
85  std::vector<int> jet_N_CPFCands(nJets);
86  std::vector<int> jet_N_NPFCands(nJets);
87  std::vector<int> jet_N_PVs(nJets);
88  std::vector<int> jet_N_SVs(nJets);
89 
90  // should default to 0 if less than nCpf cpf with information
91  std::vector<std::vector<float>> Cpfcan_BtagPf_trackEtaRel_nCpf(n_cpf_, std::vector<float>(nJets));
92  std::vector<std::vector<float>> Cpfcan_BtagPf_trackPtRel_nCpf(n_cpf_, std::vector<float>(nJets));
93  std::vector<std::vector<float>> Cpfcan_BtagPf_trackPPar_nCpf(n_cpf_, std::vector<float>(nJets));
94  std::vector<std::vector<float>> Cpfcan_BtagPf_trackDeltaR_nCpf(n_cpf_, std::vector<float>(nJets));
95  std::vector<std::vector<float>> Cpfcan_BtagPf_trackPParRatio_nCpf(n_cpf_, std::vector<float>(nJets));
96  std::vector<std::vector<float>> Cpfcan_BtagPf_trackSip2dVal_nCpf(n_cpf_, std::vector<float>(nJets));
97  std::vector<std::vector<float>> Cpfcan_BtagPf_trackSip2dSig_nCpf(n_cpf_, std::vector<float>(nJets));
98  std::vector<std::vector<float>> Cpfcan_BtagPf_trackSip3dVal_nCpf(n_cpf_, std::vector<float>(nJets));
99  std::vector<std::vector<float>> Cpfcan_BtagPf_trackSip3dSig_nCpf(n_cpf_, std::vector<float>(nJets));
100  std::vector<std::vector<float>> Cpfcan_BtagPf_trackJetDistVal_nCpf(n_cpf_, std::vector<float>(nJets));
101  std::vector<std::vector<float>> Cpfcan_ptrel_nCpf(n_cpf_, std::vector<float>(nJets));
102  std::vector<std::vector<float>> Cpfcan_drminsv_nCpf(n_cpf_, std::vector<float>(nJets));
103  std::vector<std::vector<int>> Cpfcan_VTX_ass_nCpf(n_cpf_, std::vector<int>(nJets));
104  std::vector<std::vector<float>> Cpfcan_puppiw_nCpf(n_cpf_, std::vector<float>(nJets));
105  std::vector<std::vector<float>> Cpfcan_chi2_nCpf(n_cpf_, std::vector<float>(nJets));
106  std::vector<std::vector<int>> Cpfcan_quality_nCpf(n_cpf_, std::vector<int>(nJets));
107 
108  // should default to 0 if less than nNpf npf with information
109  std::vector<std::vector<float>> Npfcan_ptrel_nNpf(n_npf_, std::vector<float>(nJets));
110  std::vector<std::vector<float>> Npfcan_deltaR_nNpf(n_npf_, std::vector<float>(nJets));
111  std::vector<std::vector<int>> Npfcan_isGamma_nNpf(n_npf_, std::vector<int>(nJets));
112  std::vector<std::vector<float>> Npfcan_HadFrac_nNpf(n_npf_, std::vector<float>(nJets));
113  std::vector<std::vector<float>> Npfcan_drminsv_nNpf(n_npf_, std::vector<float>(nJets));
114  std::vector<std::vector<float>> Npfcan_puppiw_nNpf(n_npf_, std::vector<float>(nJets));
115  /*
116  // only after dataformat updated as well
117  std::vector<std::vector<float>> Npfcan_etarel_nNpf(n_npf_, std::vector<float>(nJets));
118  std::vector<std::vector<float>> Npfcan_phirel_nNpf(n_npf_, std::vector<float>(nJets));
119  */
120 
121  // should default to 0 if less than nSv SVs with information
122  std::vector<std::vector<float>> sv_mass_nSV(n_sv_, std::vector<float>(nJets));
123  std::vector<std::vector<float>> sv_pt_nSV(n_sv_, std::vector<float>(nJets));
124  std::vector<std::vector<int>> sv_ntracks_nSV(n_sv_, std::vector<int>(nJets));
125  std::vector<std::vector<float>> sv_chi2_nSV(n_sv_, std::vector<float>(nJets));
126  std::vector<std::vector<float>> sv_normchi2_nSV(n_sv_, std::vector<float>(nJets));
127  std::vector<std::vector<float>> sv_dxy_nSV(n_sv_, std::vector<float>(nJets));
128  std::vector<std::vector<float>> sv_dxysig_nSV(n_sv_, std::vector<float>(nJets));
129  std::vector<std::vector<float>> sv_d3d_nSV(n_sv_, std::vector<float>(nJets));
130  std::vector<std::vector<float>> sv_d3dsig_nSV(n_sv_, std::vector<float>(nJets));
131  std::vector<std::vector<float>> sv_costhetasvpv_nSV(n_sv_, std::vector<float>(nJets));
132  /*
133  // only after dataformat updated as well
134  std::vector<std::vector<float>> sv_etarel_nSV(n_sv_, std::vector<float>(nJets));
135  std::vector<std::vector<float>> sv_phirel_nSV(n_sv_, std::vector<float>(nJets));
136  */
137  std::vector<std::vector<float>> sv_deltaR_nSV(n_sv_, std::vector<float>(nJets));
138  std::vector<std::vector<float>> sv_enratio_nSV(n_sv_, std::vector<float>(nJets));
139 
140  if (!tag_infos->empty()) {
141  for (unsigned i_jet = 0; i_jet < nJets; ++i_jet) {
142  // jet loop reads tag info instead of constituent info
143 
144  const auto& taginfo = (*tag_infos)[i_jet];
145  const auto& features = taginfo.features();
146 
147  // jet.pt and jet.eta as well as other jet variables (ShallowTagInfo) already included (via DeepCSV)
148  // number of elements in different collections
149  jet_N_CPFCands[i_jet] = features.c_pf_features.size();
150  jet_N_NPFCands[i_jet] = features.n_pf_features.size();
151  jet_N_SVs[i_jet] = features.sv_features.size();
152  jet_N_PVs[i_jet] = features.npv;
153 
154  // c_pf candidates
155  auto max_c_pf_n = std::min(features.c_pf_features.size(), (std::size_t)n_cpf_);
156  for (std::size_t c_pf_n = 0; c_pf_n < max_c_pf_n; c_pf_n++) {
157  const auto& c_pf_features = features.c_pf_features.at(c_pf_n);
158  Cpfcan_BtagPf_trackEtaRel_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackEtaRel;
159  Cpfcan_BtagPf_trackPtRel_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackPtRel;
160  Cpfcan_BtagPf_trackPPar_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackPPar;
161  Cpfcan_BtagPf_trackDeltaR_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackDeltaR;
162  Cpfcan_BtagPf_trackPParRatio_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackPParRatio;
163  Cpfcan_BtagPf_trackSip2dVal_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackSip2dVal;
164  Cpfcan_BtagPf_trackSip2dSig_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackSip2dSig;
165  Cpfcan_BtagPf_trackSip3dVal_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackSip3dVal;
166  Cpfcan_BtagPf_trackSip3dSig_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackSip3dSig;
167  Cpfcan_BtagPf_trackJetDistVal_nCpf[c_pf_n][i_jet] = c_pf_features.btagPf_trackJetDistVal;
168  Cpfcan_ptrel_nCpf[c_pf_n][i_jet] = c_pf_features.ptrel;
169  Cpfcan_drminsv_nCpf[c_pf_n][i_jet] = c_pf_features.drminsv;
170  Cpfcan_VTX_ass_nCpf[c_pf_n][i_jet] = c_pf_features.vtx_ass;
171  Cpfcan_puppiw_nCpf[c_pf_n][i_jet] = c_pf_features.puppiw;
172  Cpfcan_chi2_nCpf[c_pf_n][i_jet] = c_pf_features.chi2;
173  Cpfcan_quality_nCpf[c_pf_n][i_jet] = c_pf_features.quality;
174  }
175 
176  // n_pf candidates
177  auto max_n_pf_n = std::min(features.n_pf_features.size(), (std::size_t)n_npf_);
178  for (std::size_t n_pf_n = 0; n_pf_n < max_n_pf_n; n_pf_n++) {
179  const auto& n_pf_features = features.n_pf_features.at(n_pf_n);
180  Npfcan_ptrel_nNpf[n_pf_n][i_jet] = n_pf_features.ptrel;
181  Npfcan_deltaR_nNpf[n_pf_n][i_jet] = n_pf_features.deltaR;
182  Npfcan_isGamma_nNpf[n_pf_n][i_jet] = n_pf_features.isGamma;
183  Npfcan_HadFrac_nNpf[n_pf_n][i_jet] = n_pf_features.hadFrac;
184  Npfcan_drminsv_nNpf[n_pf_n][i_jet] = n_pf_features.drminsv;
185  Npfcan_puppiw_nNpf[n_pf_n][i_jet] = n_pf_features.puppiw;
186  /*
187  // only after dataformat updated as well
188  Npfcan_etarel_nNpf[n_pf_n][i_jet] = n_pf_features.etarel;
189  Npfcan_phirel_nNpf[n_pf_n][i_jet] = n_pf_features.phirel;
190  */
191  }
192 
193  // sv candidates
194  auto max_sv_n = std::min(features.sv_features.size(), (std::size_t)n_sv_);
195  for (std::size_t sv_n = 0; sv_n < max_sv_n; sv_n++) {
196  const auto& sv_features = features.sv_features.at(sv_n);
197  sv_pt_nSV[sv_n][i_jet] = sv_features.pt;
198  sv_deltaR_nSV[sv_n][i_jet] = sv_features.deltaR;
199  sv_mass_nSV[sv_n][i_jet] = sv_features.mass;
200  sv_ntracks_nSV[sv_n][i_jet] = sv_features.ntracks;
201  sv_chi2_nSV[sv_n][i_jet] = sv_features.chi2;
202  sv_normchi2_nSV[sv_n][i_jet] = sv_features.normchi2;
203  sv_dxy_nSV[sv_n][i_jet] = sv_features.dxy;
204  sv_dxysig_nSV[sv_n][i_jet] = sv_features.dxysig;
205  sv_d3d_nSV[sv_n][i_jet] = sv_features.d3d;
206  sv_d3dsig_nSV[sv_n][i_jet] = sv_features.d3dsig;
207  sv_costhetasvpv_nSV[sv_n][i_jet] = sv_features.costhetasvpv;
208  sv_enratio_nSV[sv_n][i_jet] = sv_features.enratio;
209  /*
210  // only after dataformat updated as well
211  sv_etarel_nSV[sv_n][i_jet] = sv_features.etarel;
212  sv_phirel_nSV[sv_n][i_jet] = sv_features.phirel;
213  */
214  }
215  }
216  }
217 
218  // DeepJetInputs table
219  auto djTable = std::make_unique<nanoaod::FlatTable>(jet_N_CPFCands.size(), nameDeepJet_, false, true);
220  //djTable->addColumn<int>("DeepJet_jetIdx", jetIdx_dj, "Index of the parent jet", );
221 
222  djTable->addColumn<int>("DeepJet_nCpfcand", jet_N_CPFCands, "Number of charged PF candidates in the jet");
223  djTable->addColumn<int>("DeepJet_nNpfcand", jet_N_NPFCands, "Number of neutral PF candidates in the jet");
224  djTable->addColumn<int>("DeepJet_nsv", jet_N_SVs, "Number of secondary vertices in the jet");
225  djTable->addColumn<int>("DeepJet_npv", jet_N_PVs, "Number of primary vertices");
226 
227  // ============================================================== Cpfs ===================================================================
228  for (unsigned int p = 0; p < n_cpf_; p++) {
229  auto s = std::to_string(p);
230 
231  djTable->addColumn<float>("DeepJet_Cpfcan_puppiw_" + s,
232  Cpfcan_puppiw_nCpf[p],
233  "charged candidate PUPPI weight of the " + s + ". cpf",
234  10);
235  djTable->addColumn<int>(
236  "DeepJet_Cpfcan_VTX_ass_" + s,
237  Cpfcan_VTX_ass_nCpf[p],
238  "integer flag that indicates whether the track was used in the primary vertex fit for the " + s + ". cpf",
239  10);
240  djTable->addColumn<float>("DeepJet_Cpfcan_drminsv_" + s,
241  Cpfcan_drminsv_nCpf[p],
242  "track pseudoangular distance from the closest secondary vertex of the " + s + ". cpf",
243  10);
244  djTable->addColumn<float>("DeepJet_Cpfcan_ptrel_" + s,
245  Cpfcan_ptrel_nCpf[p],
246  "fraction of the jet momentum carried by the track for the " + s + ". cpf",
247  10);
248  djTable->addColumn<int>(
249  "DeepJet_Cpfcan_quality_" + s,
250  Cpfcan_quality_nCpf[p],
251  "integer flag which indicates the quality of the fitted track, based on number of detector hits used for the "
252  "reconstruction as well as the overall chi2 of the charged track fit for the " +
253  s + ". cpf",
254  10);
255  djTable->addColumn<float>(
256  "DeepJet_Cpfcan_chi2_" + s, Cpfcan_chi2_nCpf[p], "chi2 of the charged track fit for the " + s + ". cpf", 10);
257 
258  djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackDeltaR_" + s,
259  Cpfcan_BtagPf_trackDeltaR_nCpf[p],
260  "track pseudoangular distance from the jet axis for the " + s + ". cpf",
261  10);
262  djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackEtaRel_" + s,
263  Cpfcan_BtagPf_trackEtaRel_nCpf[p],
264  "track pseudorapidity, relative to the jet axis for the " + s + ". cpf",
265  10);
266  djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackJetDistVal_" + s,
267  Cpfcan_BtagPf_trackJetDistVal_nCpf[p],
268  "minimum track approach distance to jet axis for the " + s + ". cpf",
269  10);
270  djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackPPar_" + s,
271  Cpfcan_BtagPf_trackPPar_nCpf[p],
272  "dot product of the jet and track momentum for the " + s + ". cpf",
273  10);
274  djTable->addColumn<float>(
275  "DeepJet_Cpfcan_BtagPf_trackPParRatio_" + s,
276  Cpfcan_BtagPf_trackPParRatio_nCpf[p],
277  "dot product of the jet and track momentum divided by the magnitude of the jet momentum for the " + s + ". cpf",
278  10);
279  djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackPtRel_" + s,
280  Cpfcan_BtagPf_trackPtRel_nCpf[p],
281  "track transverse momentum, relative to the jet axis for the " + s + ". cpf",
282  10);
283  djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackSip2dSig_" + s,
284  Cpfcan_BtagPf_trackSip2dSig_nCpf[p],
285  "track 2D signed impact parameter significance for the " + s + ". cpf",
286  10);
287  djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackSip3dSig_" + s,
288  Cpfcan_BtagPf_trackSip3dSig_nCpf[p],
289  "track 3D signed impact parameter significance for the " + s + ". cpf",
290  10);
291  djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackSip2dVal_" + s,
292  Cpfcan_BtagPf_trackSip2dVal_nCpf[p],
293  "track 2D signed impact parameter for the " + s + ". cpf",
294  10);
295  djTable->addColumn<float>("DeepJet_Cpfcan_BtagPf_trackSip3dVal_" + s,
296  Cpfcan_BtagPf_trackSip3dVal_nCpf[p],
297  "track 3D signed impact parameter for the " + s + ". cpf",
298  10);
299  }
300 
301  // ============================================================== Npfs ===================================================================
302  for (unsigned int p = 0; p < n_npf_; p++) {
303  auto s = std::to_string(p);
304 
305  djTable->addColumn<float>("DeepJet_Npfcan_puppiw_" + s,
306  Npfcan_puppiw_nNpf[p],
307  "neutral candidate PUPPI weight for the " + s + ". npf",
308  10);
309  djTable->addColumn<float>(
310  "DeepJet_Npfcan_deltaR_" + s,
311  Npfcan_deltaR_nNpf[p],
312  "pseudoangular distance between the neutral candidate and the jet axis for the " + s + ". npf",
313  10);
314  djTable->addColumn<float>(
315  "DeepJet_Npfcan_drminsv_" + s,
316  Npfcan_drminsv_nNpf[p],
317  "pseudoangular distance between the neutral candidate and the closest secondary vertex for the " + s + ". npf",
318  10);
319  djTable->addColumn<float>(
320  "DeepJet_Npfcan_HadFrac_" + s,
321  Npfcan_HadFrac_nNpf[p],
322  "fraction of the neutral candidate energy deposited in the hadronic calorimeter for the " + s + ". npf",
323  10);
324  djTable->addColumn<float>("DeepJet_Npfcan_ptrel_" + s,
325  Npfcan_ptrel_nNpf[p],
326  "fraction of the jet momentum carried by the neutral candidate for the " + s + ". npf",
327  10);
328  djTable->addColumn<int>("DeepJet_Npfcan_isGamma_" + s,
329  Npfcan_isGamma_nNpf[p],
330  "integer flag indicating whether the neutral candidate is a photon for the " + s + ". npf",
331  10);
332  /*
333  // only after dataformat updated as well
334  djTable->addColumn<float>("DeepJetExtra_Npfcan_etarel_" + s,
335  Npfcan_etarel_nNpf[p],
336  "pseudorapidity relative to parent jet for the " + s + ". npf",
337  10);
338  djTable->addColumn<float>("DeepJetExtra_Npfcan_phirel_" + s,
339  Npfcan_phirel_nNpf[p],
340  "DeltaPhi(npf, jet) for the " + s + ". npf",
341  10);
342  */
343  }
344 
345  // ============================================================== SVs ===================================================================
346  for (unsigned int p = 0; p < n_sv_; p++) {
347  auto s = std::to_string(p);
348 
349  djTable->addColumn<float>("DeepJet_sv_mass_" + s, sv_mass_nSV[p], "SV mass of the " + s + ". SV", 10);
350  djTable->addColumn<float>("DeepJet_sv_pt_" + s, sv_pt_nSV[p], "SV pt of the " + s + ". SV", 10);
351  djTable->addColumn<float>(
352  "DeepJet_sv_ntracks_" + s, sv_ntracks_nSV[p], "Number of tracks asociated to the " + s + ". SV", 10);
353  djTable->addColumn<float>("DeepJet_sv_chi2_" + s, sv_chi2_nSV[p], "chi2 of the " + s + ". SV", 10);
354  djTable->addColumn<float>("DeepJet_sv_normchi2_" + s, sv_normchi2_nSV[p], "chi2/dof of the " + s + ". SV", 10);
355  djTable->addColumn<float>(
356  "DeepJet_sv_dxy_" + s, sv_dxy_nSV[p], "2D impact parameter (flight distance) value of the " + s + ". SV", 10);
357  djTable->addColumn<float>("DeepJet_sv_dxysig_" + s,
358  sv_dxysig_nSV[p],
359  "2D impact parameter (flight distance) significance of the " + s + ". SV",
360  10);
361  djTable->addColumn<float>(
362  "DeepJet_sv_d3d_" + s, sv_d3d_nSV[p], "3D impact parameter (flight distance) value of the " + s + ". SV", 10);
363  djTable->addColumn<float>("DeepJet_sv_d3dsig_" + s,
364  sv_d3dsig_nSV[p],
365  "3D impact parameter (flight distance) significance of the " + s + ". SV",
366  10);
367  djTable->addColumn<float>("DeepJet_sv_costhetasvpv_" + s,
368  sv_costhetasvpv_nSV[p],
369  "cosine of the angle between the " + s +
370  ". SV flight direction and the direction of the " + s + ". SV momentum",
371  10);
372  /*
373  // only after dataformat updated as well
374  djTable->addColumn<float>("DeepJetExtra_sv_etarel_" + s,
375  sv_etarel_nSV[p],
376  "pseudorapidity relative to parent jet for the " + s + ". SV",
377  10);
378  djTable->addColumn<float>("DeepJetExtra_sv_phirel_" + s,
379  sv_phirel_nSV[p],
380  "DeltaPhi(sv, jet) for the " + s + ". SV",
381  10);
382  */
383  djTable->addColumn<float>("DeepJet_sv_deltaR_" + s,
384  sv_deltaR_nSV[p],
385  "pseudoangular distance between jet axis and the " + s + ". SV direction",
386  10);
387  djTable->addColumn<float>(
388  "DeepJet_sv_enratio_" + s, sv_enratio_nSV[p], "ratio of the " + s + ". SV energy ratio to the jet energy", 10);
389  }
390 
391  iEvent.put(std::move(djTable), nameDeepJet_);
392 }
const unsigned int n_npf_
const std::string nameDeepJet_
edm::EDGetTokenT< edm::View< T > > jet_token_
static constexpr int nJets
static std::string to_string(const XMLCh *ch)
const unsigned int n_cpf_
int iEvent
Definition: GenABIO.cc:224
std::vector< float > features(const reco::PreId &ecal, const reco::PreId &hcal, double rho, const reco::BeamSpot &spot, noZS::EcalClusterLazyTools &ecalTools)
const unsigned int n_sv_
const edm::EDGetTokenT< TagInfoCollection > tag_info_src_
def move(src, dest)
Definition: eostools.py:511

Member Data Documentation

◆ idx_nameDeepJet_

template<typename T >
const std::string DeepJetTableProducer< T >::idx_nameDeepJet_
private

Definition at line 40 of file DeepJetTableProducer.cc.

◆ jet_token_

template<typename T >
edm::EDGetTokenT<edm::View<T> > DeepJetTableProducer< T >::jet_token_
private

Definition at line 45 of file DeepJetTableProducer.cc.

◆ n_cpf_

template<typename T >
const unsigned int DeepJetTableProducer< T >::n_cpf_ = 25
private

Definition at line 41 of file DeepJetTableProducer.cc.

◆ n_npf_

template<typename T >
const unsigned int DeepJetTableProducer< T >::n_npf_ = 25
private

Definition at line 42 of file DeepJetTableProducer.cc.

◆ n_sv_

template<typename T >
const unsigned int DeepJetTableProducer< T >::n_sv_ = 12
private

Definition at line 43 of file DeepJetTableProducer.cc.

◆ nameDeepJet_

template<typename T >
const std::string DeepJetTableProducer< T >::nameDeepJet_
private

◆ tag_info_src_

template<typename T >
const edm::EDGetTokenT<TagInfoCollection> DeepJetTableProducer< T >::tag_info_src_
private

Definition at line 48 of file DeepJetTableProducer.cc.

◆ usePhysForLightAndUndefined

template<typename T >
constexpr bool DeepJetTableProducer< T >::usePhysForLightAndUndefined = false
staticprivate

Definition at line 50 of file DeepJetTableProducer.cc.