CMS 3D CMS Logo
LowPtGsfElectronFeatures.cc
Go to the documentation of this file.
1 #include "RecoEgamma/EgammaElectronProducers/interface/LowPtGsfElectronFeatures.h"
2 #include "CommonTools/BaseParticlePropagator/interface/BaseParticlePropagator.h"
3 #include "DataFormats/GsfTrackReco/interface/GsfTrack.h"
4 #include "DataFormats/ParticleFlowReco/interface/PFCluster.h"
5 #include "DataFormats/ParticleFlowReco/interface/PFClusterFwd.h"
6 #include "DataFormats/PatCandidates/interface/Electron.h"
7 #include "DataFormats/TrackReco/interface/Track.h"
8 #include "DataFormats/TrackReco/interface/TrackFwd.h"
9 #include "TVector3.h"
10 #include <algorithm>
11 #include <iostream>
12 
13 namespace lowptgsfeleseed {
14 
16  //
17  std::vector<float> features(const reco::PreId& ecal,
18  const reco::PreId& hcal,
19  double rho,
20  const reco::BeamSpot& spot,
21  noZS::EcalClusterLazyTools& tools) {
22  float trk_pt_ = -1.;
23  float trk_eta_ = -1.;
24  float trk_phi_ = -1.;
25  float trk_p_ = -1.;
26  float trk_nhits_ = -1.;
27  float trk_high_quality_ = -1.;
28  float trk_chi2red_ = -1.;
29  float rho_ = -1.;
30  float ktf_ecal_cluster_e_ = -1.;
31  float ktf_ecal_cluster_deta_ = -42.;
32  float ktf_ecal_cluster_dphi_ = -42.;
33  float ktf_ecal_cluster_e3x3_ = -1.;
34  float ktf_ecal_cluster_e5x5_ = -1.;
35  float ktf_ecal_cluster_covEtaEta_ = -42.;
36  float ktf_ecal_cluster_covEtaPhi_ = -42.;
37  float ktf_ecal_cluster_covPhiPhi_ = -42.;
38  float ktf_ecal_cluster_r9_ = -0.1;
39  float ktf_ecal_cluster_circularity_ = -0.1;
40  float ktf_hcal_cluster_e_ = -1.;
41  float ktf_hcal_cluster_deta_ = -42.;
42  float ktf_hcal_cluster_dphi_ = -42.;
43  float preid_gsf_dpt_ = -1.;
44  float preid_trk_gsf_chiratio_ = -1.;
45  float preid_gsf_chi2red_ = -1.;
46  float trk_dxy_sig_ = -1.; // must be last (not used by unbiased model)
47 
48  // Tracks
49  const auto& trk = ecal.trackRef(); // reco::TrackRef
50  if (trk.isNonnull()) {
51  trk_pt_ = trk->pt();
52  trk_eta_ = trk->eta();
53  trk_phi_ = trk->phi();
54  trk_p_ = trk->p();
55  trk_nhits_ = static_cast<float>(trk->found());
56  trk_high_quality_ = static_cast<float>(trk->quality(reco::TrackBase::qualityByName("highPurity")));
57  trk_chi2red_ = trk->normalizedChi2();
58  if (trk->dxy(spot) > 0.) {
59  trk_dxy_sig_ = trk->dxyError() / trk->dxy(spot); //@@ to be consistent with the training based on 94X MC
60  }
61  ktf_ecal_cluster_dphi_ *= trk->charge(); //@@ to be consistent with the training based on 94X MC
62  }
63 
64  // Rho
65  rho_ = static_cast<float>(rho);
66 
67  // ECAL clusters
68  const auto& ecal_clu = ecal.clusterRef(); // reco::PFClusterRef
69  if (ecal_clu.isNonnull()) {
70  ktf_ecal_cluster_e_ = ecal_clu->energy();
71  ktf_ecal_cluster_deta_ = ecal.geomMatching()[0];
72  ktf_ecal_cluster_dphi_ = ecal.geomMatching()[1];
73  ktf_ecal_cluster_e3x3_ = tools.e3x3(*ecal_clu);
74  ktf_ecal_cluster_e5x5_ = tools.e5x5(*ecal_clu);
75  const auto& covs = tools.localCovariances(*ecal_clu);
76  ktf_ecal_cluster_covEtaEta_ = covs[0];
77  ktf_ecal_cluster_covEtaPhi_ = covs[1];
78  ktf_ecal_cluster_covPhiPhi_ = covs[2];
79  if (ktf_ecal_cluster_e_ > 0.) {
80  ktf_ecal_cluster_r9_ = ktf_ecal_cluster_e3x3_ / ktf_ecal_cluster_e_;
81  }
82  if (ktf_ecal_cluster_e5x5_ > 0.) {
83  ktf_ecal_cluster_circularity_ = 1. - tools.e1x5(*ecal_clu) / ktf_ecal_cluster_e5x5_;
84  } else {
85  ktf_ecal_cluster_circularity_ = -0.1;
86  }
87  }
88 
89  // HCAL clusters
90  const auto& hcal_clu = hcal.clusterRef(); // reco::PFClusterRef
91  if (hcal_clu.isNonnull()) {
92  ktf_hcal_cluster_e_ = hcal_clu->energy();
93  ktf_hcal_cluster_deta_ = hcal.geomMatching()[0];
94  ktf_hcal_cluster_dphi_ = hcal.geomMatching()[1];
95  }
96 
97  // PreId
98  preid_gsf_dpt_ = ecal.dpt();
99  preid_trk_gsf_chiratio_ = ecal.chi2Ratio();
100  preid_gsf_chi2red_ = ecal.gsfChi2();
101 
102  // Set contents of vector
103  std::vector<float> output = {trk_pt_,
104  trk_eta_,
105  trk_phi_,
106  trk_p_,
107  trk_nhits_,
108  trk_high_quality_,
109  trk_chi2red_,
110  rho_,
111  ktf_ecal_cluster_e_,
112  ktf_ecal_cluster_deta_,
113  ktf_ecal_cluster_dphi_,
114  ktf_ecal_cluster_e3x3_,
115  ktf_ecal_cluster_e5x5_,
116  ktf_ecal_cluster_covEtaEta_,
117  ktf_ecal_cluster_covEtaPhi_,
118  ktf_ecal_cluster_covPhiPhi_,
119  ktf_ecal_cluster_r9_,
120  ktf_ecal_cluster_circularity_,
121  ktf_hcal_cluster_e_,
122  ktf_hcal_cluster_deta_,
123  ktf_hcal_cluster_dphi_,
124  preid_gsf_dpt_,
125  preid_trk_gsf_chiratio_,
126  preid_gsf_chi2red_,
127  trk_dxy_sig_};
128  return output;
129  };
130 
131 } // namespace lowptgsfeleseed
132 
133 namespace lowptgsfeleid {
134 
135  std::vector<float> features_V1(
136  reco::GsfElectron const& ele, float rho, float unbiased, float field_z, const reco::Track* trk) {
137  float eid_rho = -999.;
138  float eid_sc_eta = -999.;
139  float eid_shape_full5x5_r9 = -999.;
140  float eid_sc_etaWidth = -999.;
141  float eid_sc_phiWidth = -999.;
142  float eid_shape_full5x5_HoverE = -999.;
143  float eid_trk_nhits = -999.;
144  float eid_trk_chi2red = -999.;
145  float eid_gsf_chi2red = -999.;
146  float eid_brem_frac = -999.;
147  float eid_gsf_nhits = -999.;
148  float eid_match_SC_EoverP = -999.;
149  float eid_match_eclu_EoverP = -999.;
150  float eid_match_SC_dEta = -999.;
151  float eid_match_SC_dPhi = -999.;
152  float eid_match_seed_dEta = -999.;
153  float eid_sc_E = -999.;
154  float eid_trk_p = -999.;
155  float gsf_mode_p = -999.;
156  float core_shFracHits = -999.;
157  float gsf_bdtout1 = -999.;
158  float gsf_dr = -999.;
159  float trk_dr = -999.;
160  float sc_Nclus = -999.;
161  float sc_clus1_nxtal = -999.;
162  float sc_clus1_dphi = -999.;
163  float sc_clus2_dphi = -999.;
164  float sc_clus1_deta = -999.;
165  float sc_clus2_deta = -999.;
166  float sc_clus1_E = -999.;
167  float sc_clus2_E = -999.;
168  float sc_clus1_E_ov_p = -999.;
169  float sc_clus2_E_ov_p = -999.;
170 
171  // KF tracks
172  if (trk != nullptr || (ele.core().isNonnull() && ele.closestCtfTrackRef().isNonnull())) {
173  const reco::Track* tk = trk ? trk : &*ele.closestCtfTrackRef();
174  eid_trk_p = tk->p();
175  eid_trk_nhits = tk->found();
176  eid_trk_chi2red = tk->normalizedChi2();
177  trk_dr = reco::deltaR(*tk, ele);
178  }
179 
180  // GSF tracks
181  if (ele.core().isNonnull()) {
182  reco::GsfTrackRef gsf = ele.gsfTrack();
183  if (gsf.isNonnull()) {
184  gsf_mode_p = gsf->pMode();
185  eid_gsf_nhits = (float)gsf->found();
186  eid_gsf_chi2red = gsf->normalizedChi2();
187  TVector3 gsfTV3(0, 0, 0);
188  gsfTV3.SetPtEtaPhi(gsf->ptMode(), gsf->etaMode(), gsf->phiMode());
189  TVector3 eleTV3(0, 0, 0);
190  eleTV3.SetPtEtaPhi(ele.pt(), ele.eta(), ele.phi());
191  gsf_dr = eleTV3.DeltaR(gsfTV3);
192  }
193  }
194 
195  // Super clusters
196  if (ele.core().isNonnull()) {
197  reco::SuperClusterRef sc = ele.superCluster();
198  if (sc.isNonnull()) {
199  eid_sc_E = sc->energy();
200  eid_sc_eta = sc->eta();
201  eid_sc_etaWidth = sc->etaWidth();
202  eid_sc_phiWidth = sc->phiWidth();
203  sc_Nclus = sc->clustersSize();
204  }
205  }
206 
207  // Track-cluster matching
208  eid_match_seed_dEta = ele.deltaEtaSeedClusterTrackAtCalo();
209  eid_match_eclu_EoverP = (1. / ele.ecalEnergy()) - (1. / ele.p());
210  eid_match_SC_EoverP = ele.eSuperClusterOverP();
211  eid_match_SC_dEta = ele.deltaEtaSuperClusterTrackAtVtx();
212  eid_match_SC_dPhi = ele.deltaPhiSuperClusterTrackAtVtx();
213 
214  // Shower shape vars
215  eid_shape_full5x5_HoverE = ele.full5x5_hcalOverEcal();
216  eid_shape_full5x5_r9 = ele.full5x5_r9();
217 
218  // Misc
219  eid_rho = rho;
220 
221  eid_brem_frac = ele.fbrem();
222  core_shFracHits = ele.shFracInnerHits();
223 
224  // Unbiased BDT from ElectronSeed
225  gsf_bdtout1 = unbiased;
226 
227  // Clusters
228  if (ele.core().isNonnull()) {
229  reco::GsfTrackRef gsf = ele.gsfTrack();
230  if (gsf.isNonnull()) {
231  reco::SuperClusterRef sc = ele.superCluster();
232  if (sc.isNonnull()) {
233  // Propagate electron track to ECAL surface
234  double mass2 = 0.000511 * 0.000511;
235  float p2 = pow(gsf->p(), 2);
236  float energy = sqrt(mass2 + p2);
237  XYZTLorentzVector mom = XYZTLorentzVector(gsf->px(), gsf->py(), gsf->pz(), energy);
238  XYZTLorentzVector pos = XYZTLorentzVector(gsf->vx(), gsf->vy(), gsf->vz(), 0.);
239  BaseParticlePropagator propagator(RawParticle(mom, pos, gsf->charge()), 0, 0, field_z);
240  propagator.propagateToEcalEntrance(true); // true only first half loop , false more than one loop
241  bool reach_ECAL = propagator.getSuccess(); // 0 does not reach ECAL, 1 yes barrel, 2 yes endcaps
242  // ECAL entry point for track
243  GlobalPoint ecal_pos(propagator.particle().x(), propagator.particle().y(), propagator.particle().z());
244 
245  // Track-cluster matching for most energetic clusters
246  sc_clus1_nxtal = -999;
247  sc_clus1_dphi = -999.;
248  sc_clus2_dphi = -999.;
249  sc_clus1_deta = -999.;
250  sc_clus2_deta = -999.;
251  sc_clus1_E = -999.;
252  sc_clus2_E = -999.;
253  sc_clus1_E_ov_p = -999.;
254  sc_clus2_E_ov_p = -999.;
255  trackClusterMatching(*sc,
256  *gsf,
257  reach_ECAL,
258  ecal_pos,
259  sc_clus1_nxtal,
260  sc_clus1_dphi,
261  sc_clus2_dphi,
262  sc_clus1_deta,
263  sc_clus2_deta,
264  sc_clus1_E,
265  sc_clus2_E,
266  sc_clus1_E_ov_p,
267  sc_clus2_E_ov_p);
268  sc_clus1_nxtal = (int)sc_clus1_nxtal;
269 
270  } // sc.isNonnull()
271  } // gsf.isNonnull()
272  } // clusters
273 
274  // Out-of-range
275  eid_rho = std::clamp(eid_rho, 0.f, 100.f);
276  eid_sc_eta = std::clamp(eid_sc_eta, -5.f, 5.f);
277  eid_shape_full5x5_r9 = std::clamp(eid_shape_full5x5_r9, 0.f, 2.f);
278  eid_sc_etaWidth = std::clamp(eid_sc_etaWidth, 0.f, 3.14f);
279  eid_sc_phiWidth = std::clamp(eid_sc_phiWidth, 0.f, 3.14f);
280  eid_shape_full5x5_HoverE = std::clamp(eid_shape_full5x5_HoverE, 0.f, 50.f);
281  eid_trk_nhits = std::clamp(eid_trk_nhits, -1.f, 50.f);
282  eid_trk_chi2red = std::clamp(eid_trk_chi2red, -1.f, 50.f);
283  eid_gsf_chi2red = std::clamp(eid_gsf_chi2red, -1.f, 100.f);
284  if (eid_brem_frac < 0.)
285  eid_brem_frac = -1.; //
286  if (eid_brem_frac > 1.)
287  eid_brem_frac = 1.; //
288  eid_gsf_nhits = std::clamp(eid_gsf_nhits, -1.f, 50.f);
289  eid_match_SC_EoverP = std::clamp(eid_match_SC_EoverP, 0.f, 100.f);
290  eid_match_eclu_EoverP = std::clamp(eid_match_eclu_EoverP, -1.f, 1.f);
291  eid_match_SC_dEta = std::clamp(eid_match_SC_dEta, -10.f, 10.f);
292  eid_match_SC_dPhi = std::clamp(eid_match_SC_dPhi, -3.14f, 3.14f);
293  eid_match_seed_dEta = std::clamp(eid_match_seed_dEta, -10.f, 10.f);
294  eid_sc_E = std::clamp(eid_sc_E, 0.f, 1000.f);
295  eid_trk_p = std::clamp(eid_trk_p, -1.f, 1000.f);
296  gsf_mode_p = std::clamp(gsf_mode_p, 0.f, 1000.f);
297  core_shFracHits = std::clamp(core_shFracHits, 0.f, 1.f);
298  gsf_bdtout1 = std::clamp(gsf_bdtout1, -20.f, 20.f);
299  if (gsf_dr < 0.)
300  gsf_dr = 5.; //
301  if (gsf_dr > 5.)
302  gsf_dr = 5.; //
303  if (trk_dr < 0.)
304  trk_dr = 5.; //
305  if (trk_dr > 5.)
306  trk_dr = 5.; //
307  sc_Nclus = std::clamp(sc_Nclus, 0.f, 20.f);
308  sc_clus1_nxtal = std::clamp(sc_clus1_nxtal, 0.f, 100.f);
309  sc_clus1_dphi = std::clamp(sc_clus1_dphi, -3.14f, 3.14f);
310  sc_clus2_dphi = std::clamp(sc_clus2_dphi, -3.14f, 3.14f);
311  sc_clus1_deta = std::clamp(sc_clus1_deta, -5.f, 5.f);
312  sc_clus2_deta = std::clamp(sc_clus2_deta, -5.f, 5.f);
313  sc_clus1_E = std::clamp(sc_clus1_E, 0.f, 1000.f);
314  sc_clus2_E = std::clamp(sc_clus2_E, 0.f, 1000.f);
315  if (sc_clus1_E_ov_p < 0.)
316  sc_clus1_E_ov_p = -1.; //
317  if (sc_clus2_E_ov_p < 0.)
318  sc_clus2_E_ov_p = -1.; //
319 
320  // Set contents of vector
321  std::vector<float> output = {eid_rho,
322  eid_sc_eta,
323  eid_shape_full5x5_r9,
324  eid_sc_etaWidth,
325  eid_sc_phiWidth,
326  eid_shape_full5x5_HoverE,
327  eid_trk_nhits,
328  eid_trk_chi2red,
329  eid_gsf_chi2red,
330  eid_brem_frac,
331  eid_gsf_nhits,
332  eid_match_SC_EoverP,
333  eid_match_eclu_EoverP,
334  eid_match_SC_dEta,
335  eid_match_SC_dPhi,
336  eid_match_seed_dEta,
337  eid_sc_E,
338  eid_trk_p,
339  gsf_mode_p,
340  core_shFracHits,
341  gsf_bdtout1,
342  gsf_dr,
343  trk_dr,
344  sc_Nclus,
345  sc_clus1_nxtal,
346  sc_clus1_dphi,
347  sc_clus2_dphi,
348  sc_clus1_deta,
349  sc_clus2_deta,
350  sc_clus1_E,
351  sc_clus2_E,
352  sc_clus1_E_ov_p,
353  sc_clus2_E_ov_p};
354  return output;
355  }
356 
358  // feature list for original models (2019Aug07 and earlier)
359  std::vector<float> features_V0(reco::GsfElectron const& ele, float rho, float unbiased) {
360  float eid_rho = -999.;
361  float eid_sc_eta = -999.;
362  float eid_shape_full5x5_r9 = -999.;
363  float eid_sc_etaWidth = -999.;
364  float eid_sc_phiWidth = -999.;
365  float eid_shape_full5x5_HoverE = -999.;
366  float eid_trk_nhits = -999.;
367  float eid_trk_chi2red = -999.;
368  float eid_gsf_chi2red = -999.;
369  float eid_brem_frac = -999.;
370  float eid_gsf_nhits = -999.;
371  float eid_match_SC_EoverP = -999.;
372  float eid_match_eclu_EoverP = -999.;
373  float eid_match_SC_dEta = -999.;
374  float eid_match_SC_dPhi = -999.;
375  float eid_match_seed_dEta = -999.;
376  float eid_sc_E = -999.;
377  float eid_trk_p = -999.;
378  float gsf_mode_p = -999.;
379  float core_shFracHits = -999.;
380  float gsf_bdtout1 = -999.;
381  float gsf_dr = -999.;
382  float trk_dr = -999.;
383  float sc_Nclus = -999.;
384  float sc_clus1_nxtal = -999.;
385  float sc_clus1_dphi = -999.;
386  float sc_clus2_dphi = -999.;
387  float sc_clus1_deta = -999.;
388  float sc_clus2_deta = -999.;
389  float sc_clus1_E = -999.;
390  float sc_clus2_E = -999.;
391  float sc_clus1_E_ov_p = -999.;
392  float sc_clus2_E_ov_p = -999.;
393 
394  // KF tracks
395  if (ele.core().isNonnull()) {
396  const auto& trk = ele.closestCtfTrackRef(); // reco::TrackRef
397  if (trk.isNonnull()) {
398  eid_trk_p = (float)trk->p();
399  eid_trk_nhits = (float)trk->found();
400  eid_trk_chi2red = (float)trk->normalizedChi2();
401  TVector3 trkTV3(0, 0, 0);
402  trkTV3.SetPtEtaPhi(trk->pt(), trk->eta(), trk->phi());
403  TVector3 eleTV3(0, 0, 0);
404  eleTV3.SetPtEtaPhi(ele.pt(), ele.eta(), ele.phi());
405  trk_dr = eleTV3.DeltaR(trkTV3);
406  }
407  }
408 
409  // GSF tracks
410  if (ele.core().isNonnull()) {
411  const auto& gsf = ele.core()->gsfTrack(); // reco::GsfTrackRef
412  if (gsf.isNonnull()) {
413  gsf_mode_p = gsf->pMode();
414  eid_gsf_nhits = (float)gsf->found();
415  eid_gsf_chi2red = gsf->normalizedChi2();
416  TVector3 gsfTV3(0, 0, 0);
417  gsfTV3.SetPtEtaPhi(gsf->ptMode(), gsf->etaMode(), gsf->phiMode());
418  TVector3 eleTV3(0, 0, 0);
419  eleTV3.SetPtEtaPhi(ele.pt(), ele.eta(), ele.phi());
420  gsf_dr = eleTV3.DeltaR(gsfTV3);
421  }
422  }
423 
424  // Super clusters
425  if (ele.core().isNonnull()) {
426  const auto& sc = ele.core()->superCluster(); // reco::SuperClusterRef
427  if (sc.isNonnull()) {
428  eid_sc_E = sc->energy();
429  eid_sc_eta = sc->eta();
430  eid_sc_etaWidth = sc->etaWidth();
431  eid_sc_phiWidth = sc->phiWidth();
432  sc_Nclus = (float)sc->clustersSize();
433  }
434  }
435 
436  // Track-cluster matching
437  eid_match_seed_dEta = ele.deltaEtaSeedClusterTrackAtCalo();
438  eid_match_eclu_EoverP = (1. / ele.ecalEnergy()) - (1. / ele.p());
439  eid_match_SC_EoverP = ele.eSuperClusterOverP();
440  eid_match_SC_dEta = ele.deltaEtaSuperClusterTrackAtVtx();
441  eid_match_SC_dPhi = ele.deltaPhiSuperClusterTrackAtVtx();
442 
443  // Shower shape vars
444  eid_shape_full5x5_HoverE = ele.full5x5_hcalOverEcal();
445  eid_shape_full5x5_r9 = ele.full5x5_r9();
446 
447  // Misc
448  eid_rho = rho;
449 
450  eid_brem_frac = ele.fbrem();
451  core_shFracHits = (float)ele.shFracInnerHits();
452 
453  // Unbiased BDT from ElectronSeed
454  gsf_bdtout1 = unbiased;
455 
456  // Clusters
457  if (ele.core().isNonnull()) {
458  const auto& gsf = ele.core()->gsfTrack(); // reco::GsfTrackRef
459  if (gsf.isNonnull()) {
460  const auto& sc = ele.core()->superCluster(); // reco::SuperClusterRef
461  if (sc.isNonnull()) {
462  // Propagate electron track to ECAL surface
463  double mass2 = 0.000511 * 0.000511;
464  float p2 = pow(gsf->p(), 2);
465  float energy = sqrt(mass2 + p2);
466  math::XYZTLorentzVector mom = math::XYZTLorentzVector(gsf->px(), gsf->py(), gsf->pz(), energy);
467  math::XYZTLorentzVector pos = math::XYZTLorentzVector(gsf->vx(), gsf->vy(), gsf->vz(), 0.);
468  float field_z = 3.8;
469  BaseParticlePropagator mypart(RawParticle(mom, pos, gsf->charge()), 0, 0, field_z);
470  mypart.propagateToEcalEntrance(true); // true only first half loop , false more than one loop
471  bool reach_ECAL = mypart.getSuccess(); // 0 does not reach ECAL, 1 yes barrel, 2 yes endcaps
472 
473  // ECAL entry point for track
474  GlobalPoint ecal_pos(
475  mypart.particle().vertex().x(), mypart.particle().vertex().y(), mypart.particle().vertex().z());
476 
477  // Track-cluster matching for most energetic clusters
478  sc_clus1_nxtal = -999.;
479  sc_clus1_dphi = -999.;
480  sc_clus2_dphi = -999.;
481  sc_clus1_deta = -999.;
482  sc_clus2_deta = -999.;
483  sc_clus1_E = -999.;
484  sc_clus2_E = -999.;
485  sc_clus1_E_ov_p = -999.;
486  sc_clus2_E_ov_p = -999.;
487  trackClusterMatching(*sc,
488  *gsf,
489  reach_ECAL,
490  ecal_pos,
491  sc_clus1_nxtal,
492  sc_clus1_dphi,
493  sc_clus2_dphi,
494  sc_clus1_deta,
495  sc_clus2_deta,
496  sc_clus1_E,
497  sc_clus2_E,
498  sc_clus1_E_ov_p,
499  sc_clus2_E_ov_p);
500 
501  } // sc.isNonnull()
502  } // gsf.isNonnull()
503  } // clusters
504 
505  // Out-of-range
506  eid_sc_eta = std::clamp(eid_sc_eta, -5.f, 5.f);
507  eid_shape_full5x5_r9 = std::clamp(eid_shape_full5x5_r9, 0.f, 2.f);
508  eid_sc_etaWidth = std::clamp(eid_sc_etaWidth, 0.f, 3.14f);
509  eid_sc_phiWidth = std::clamp(eid_sc_phiWidth, 0.f, 3.14f);
510  eid_shape_full5x5_HoverE = std::clamp(eid_shape_full5x5_HoverE, 0.f, 50.f);
511  eid_trk_nhits = std::clamp(eid_trk_nhits, -1.f, 50.f);
512  eid_trk_chi2red = std::clamp(eid_trk_chi2red, -1.f, 50.f);
513  eid_gsf_chi2red = std::clamp(eid_gsf_chi2red, -1.f, 100.f);
514  if (eid_brem_frac < 0.)
515  eid_brem_frac = -1.; //
516  if (eid_brem_frac > 1.)
517  eid_brem_frac = 1.; //
518  eid_gsf_nhits = std::clamp(eid_gsf_nhits, -1.f, 50.f);
519  eid_match_SC_EoverP = std::clamp(eid_match_SC_EoverP, 0.f, 100.f);
520  eid_match_eclu_EoverP = std::clamp(eid_match_eclu_EoverP, -1.f, 1.f);
521  eid_match_SC_dEta = std::clamp(eid_match_SC_dEta, -10.f, 10.f);
522  eid_match_SC_dPhi = std::clamp(eid_match_SC_dPhi, -3.14f, 3.14f);
523  eid_match_seed_dEta = std::clamp(eid_match_seed_dEta, -10.f, 10.f);
524  eid_sc_E = std::clamp(eid_sc_E, 0.f, 1000.f);
525  eid_trk_p = std::clamp(eid_trk_p, -1.f, 1000.f);
526  gsf_mode_p = std::clamp(gsf_mode_p, 0.f, 1000.f);
527  core_shFracHits = std::clamp(core_shFracHits, 0.f, 1.f);
528  gsf_bdtout1 = std::clamp(gsf_bdtout1, -20.f, 20.f);
529  if (gsf_dr < 0.)
530  gsf_dr = 5.; //
531  if (gsf_dr > 5.)
532  gsf_dr = 5.; //
533  if (trk_dr < 0.)
534  trk_dr = 5.; //
535  if (trk_dr > 5.)
536  trk_dr = 5.; //
537  sc_Nclus = std::clamp(sc_Nclus, 0.f, 20.f);
538  sc_clus1_nxtal = std::clamp(sc_clus1_nxtal, 0.f, 100.f);
539  if (sc_clus1_dphi < -3.14)
540  sc_clus1_dphi = -5.; //
541  if (sc_clus1_dphi > 3.14)
542  sc_clus1_dphi = 5.; //
543  if (sc_clus2_dphi < -3.14)
544  sc_clus2_dphi = -5.; //
545  if (sc_clus2_dphi > 3.14)
546  sc_clus2_dphi = 5.; //
547  sc_clus1_deta = std::clamp(sc_clus1_deta, -5.f, 5.f);
548  sc_clus2_deta = std::clamp(sc_clus2_deta, -5.f, 5.f);
549  sc_clus1_E = std::clamp(sc_clus1_E, 0.f, 1000.f);
550  sc_clus2_E = std::clamp(sc_clus2_E, 0.f, 1000.f);
551  if (sc_clus1_E_ov_p < 0.)
552  sc_clus1_E_ov_p = -1.; //
553  if (sc_clus2_E_ov_p < 0.)
554  sc_clus2_E_ov_p = -1.; //
555 
556  // Set contents of vector
557  std::vector<float> output = {eid_rho,
558  eid_sc_eta,
559  eid_shape_full5x5_r9,
560  eid_sc_etaWidth,
561  eid_sc_phiWidth,
562  eid_shape_full5x5_HoverE,
563  eid_trk_nhits,
564  eid_trk_chi2red,
565  eid_gsf_chi2red,
566  eid_brem_frac,
567  eid_gsf_nhits,
568  eid_match_SC_EoverP,
569  eid_match_eclu_EoverP,
570  eid_match_SC_dEta,
571  eid_match_SC_dPhi,
572  eid_match_seed_dEta,
573  eid_sc_E,
574  eid_trk_p,
575  gsf_mode_p,
576  core_shFracHits,
577  gsf_bdtout1,
578  gsf_dr,
579  trk_dr,
580  sc_Nclus,
581  sc_clus1_nxtal,
582  sc_clus1_dphi,
583  sc_clus2_dphi,
584  sc_clus1_deta,
585  sc_clus2_deta,
586  sc_clus1_E,
587  sc_clus2_E,
588  sc_clus1_E_ov_p,
589  sc_clus2_E_ov_p};
590  return output;
591  }
592 
594  // Find most energetic clusters
595  void findEnergeticClusters(
596  reco::SuperCluster const& sc, int& clusNum, float& maxEne1, float& maxEne2, int& i1, int& i2) {
597  if (sc.clustersSize() > 0 && sc.clustersBegin() != sc.clustersEnd()) {
598  for (auto const& cluster : sc.clusters()) {
599  if (cluster->energy() > maxEne1) {
600  maxEne1 = cluster->energy();
601  i1 = clusNum;
602  }
603  clusNum++;
604  }
605  if (sc.clustersSize() > 1) {
606  clusNum = 0;
607  for (auto const& cluster : sc.clusters()) {
608  if (clusNum != i1) {
609  if (cluster->energy() > maxEne2) {
610  maxEne2 = cluster->energy();
611  i2 = clusNum;
612  }
613  }
614  clusNum++;
615  }
616  }
617  } // loop over clusters
618  }
619 
621  // Track-cluster matching for most energetic clusters
622  void trackClusterMatching(reco::SuperCluster const& sc,
623  reco::GsfTrack const& gsf,
624  bool const& reach_ECAL,
625  GlobalPoint const& ecal_pos,
626  float& sc_clus1_nxtal,
627  float& sc_clus1_dphi,
628  float& sc_clus2_dphi,
629  float& sc_clus1_deta,
630  float& sc_clus2_deta,
631  float& sc_clus1_E,
632  float& sc_clus2_E,
633  float& sc_clus1_E_ov_p,
634  float& sc_clus2_E_ov_p) {
635  // Iterate through ECAL clusters and sort in energy
636  int clusNum = 0;
637  float maxEne1 = -1;
638  float maxEne2 = -1;
639  int i1 = -1;
640  int i2 = -1;
641  findEnergeticClusters(sc, clusNum, maxEne1, maxEne2, i1, i2);
642 
643  // track-clusters match
644  clusNum = 0;
645  if (sc.clustersSize() > 0 && sc.clustersBegin() != sc.clustersEnd()) {
646  for (auto const& cluster : sc.clusters()) {
647  float deta = ecal_pos.eta() - cluster->eta();
648  float dphi = reco::deltaPhi(ecal_pos.phi(), cluster->phi());
649  if (clusNum == i1) {
650  sc_clus1_E = cluster->energy();
651  if (gsf.pMode() > 0)
652  sc_clus1_E_ov_p = cluster->energy() / gsf.pMode();
653  sc_clus1_nxtal = (float)cluster->size();
654  if (reach_ECAL > 0) {
655  sc_clus1_deta = deta;
656  sc_clus1_dphi = dphi;
657  }
658  } else if (clusNum == i2) {
659  sc_clus2_E = cluster->energy();
660  if (gsf.pMode() > 0)
661  sc_clus2_E_ov_p = cluster->energy() / gsf.pMode();
662  if (reach_ECAL > 0) {
663  sc_clus2_deta = deta;
664  sc_clus2_dphi = dphi;
665  }
666  }
667  clusNum++;
668  }
669  }
670  }
671 
672 } // namespace lowptgsfeleid
reco::deltaPhi
constexpr double deltaPhi(double phi1, double phi2)
Definition: deltaPhi.h:26
lowptgsfeleid::findEnergeticClusters
void findEnergeticClusters(reco::SuperCluster const &, int &, float &, float &, int &, int &)
Definition: LowPtGsfElectronFeatures.cc:595
reco::LeafCandidate::pt
double pt() const final
transverse momentum
Definition: LeafCandidate.h:146
SiStripOfflineCRack_cfg.p2
p2
Definition: SiStripOfflineCRack_cfg.py:140
TrackCandidateProducer_cfi.propagator
propagator
Definition: TrackCandidateProducer_cfi.py:17
PV3DBase::phi
Geom::Phi< T > phi() const
Definition: PV3DBase.h:66
reco::TrackBase::p
double p() const
momentum vector magnitude
Definition: TrackBase.h:631
PV3DBase::eta
T eta() const
Definition: PV3DBase.h:73
reco::SuperCluster::clustersBegin
CaloCluster_iterator clustersBegin() const
fist iterator over BasicCluster constituents
Definition: SuperCluster.h:88
testProducerWithPsetDescEmpty_cfi.i2
i2
Definition: testProducerWithPsetDescEmpty_cfi.py:46
edm::Ref::isNonnull
bool isNonnull() const
Checks for non-null.
Definition: Ref.h:238
convertSQLitetoXML_cfg.output
output
Definition: convertSQLitetoXML_cfg.py:72
lowptgsfeleid::trackClusterMatching
void trackClusterMatching(reco::SuperCluster const &, reco::GsfTrack const &, bool const &, GlobalPoint const &, float &, float &, float &, float &, float &, float &, float &, float &, float &)
Definition: LowPtGsfElectronFeatures.cc:622
hcalRecHitTable_cff.energy
energy
Definition: hcalRecHitTable_cff.py:13
reco::SuperCluster::clustersEnd
CaloCluster_iterator clustersEnd() const
last iterator over BasicCluster constituents
Definition: SuperCluster.h:91
reco::GsfElectron::eSuperClusterOverP
float eSuperClusterOverP() const
Definition: GsfElectron.h:221
reco::Track::found
unsigned short found() const
Number of valid hits on track.
Definition: Track.h:142
math::XYZTLorentzVector
XYZTLorentzVectorD XYZTLorentzVector
Lorentz vector with cylindrical internal representation using pseudorapidity.
Definition: LorentzVector.h:29
reco::GsfElectron::gsfTrack
GsfTrackRef gsfTrack() const override
reference to a GsfTrack
Definition: GsfElectron.h:156
reco::LeafCandidate::p
double p() const final
magnitude of momentum vector
Definition: LeafCandidate.h:123
nano_mu_digi_cff.float
float
Definition: nano_mu_digi_cff.py:13
reco::GsfElectron::shFracInnerHits
float shFracInnerHits() const
Definition: GsfElectron.h:179
lowptgsfeleseed::features
std::vector< float > features(const reco::PreId &ecal, const reco::PreId &hcal, double rho, const reco::BeamSpot &spot, noZS::EcalClusterLazyTools &ecalTools)
Definition: LowPtGsfElectronFeatures.cc:17
mathSSE::sqrt
T sqrt(T t)
Definition: SSEVec.h:19
tools
Definition: tools.py:1
BaseParticlePropagator::particle
RawParticle const & particle() const
The particle being propagated.
Definition: BaseParticlePropagator.h:164
f
double f[11][100]
Definition: MuScleFitUtils.cc:78
lowptgsfeleid::features_V0
std::vector< float > features_V0(reco::GsfElectron const &ele, float rho, float unbiased)
Definition: LowPtGsfElectronFeatures.cc:359
BaseParticlePropagator::propagateToEcalEntrance
bool propagateToEcalEntrance(bool first=true)
Definition: BaseParticlePropagator.cc:450
reco::GsfElectron::core
virtual GsfElectronCoreRef core() const
Definition: GsfElectron.cc:8
reco::deltaR
constexpr auto deltaR(const T1 &t1, const T2 &t2) -> decltype(t1.eta())
Definition: deltaR.h:30
reco::GsfElectron::deltaPhiSuperClusterTrackAtVtx
float deltaPhiSuperClusterTrackAtVtx() const
Definition: GsfElectron.h:228
lowPtElectrons_cff.unbiased
unbiased
Definition: lowPtElectrons_cff.py:75
lowptgsfeleid::features_V1
std::vector< float > features_V1(reco::GsfElectron const &ele, float rho, float unbiased, float field_z, const reco::Track *trk=nullptr)
Definition: LowPtGsfElectronFeatures.cc:135
reco::TrackBase::qualityByName
static TrackQuality qualityByName(const std::string &name)
Definition: TrackBase.cc:126
reco::GsfElectron::deltaEtaSuperClusterTrackAtVtx
float deltaEtaSuperClusterTrackAtVtx() const
Definition: GsfElectron.h:225
BaseParticlePropagator::getSuccess
int getSuccess() const
Has propagation been performed and was barrel or endcap reached ?
Definition: BaseParticlePropagator.h:296
Geom::Phi::phi
T1 phi() const
Definition: Phi.h:78
reco::GsfElectron::ecalEnergy
float ecalEnergy() const
Definition: GsfElectron.h:896
reco::TrackBase::normalizedChi2
double normalizedChi2() const
chi-squared divided by n.d.o.f. (or chi-squared * 1e6 if n.d.o.f. is zero)
Definition: TrackBase.h:593
reco::SuperCluster::clusters
const CaloClusterPtrVector & clusters() const
const access to the cluster list itself
Definition: SuperCluster.h:82
testProducerWithPsetDescEmpty_cfi.i1
i1
Definition: testProducerWithPsetDescEmpty_cfi.py:45
reco::GsfElectron::closestCtfTrackRef
virtual TrackRef closestCtfTrackRef() const
Definition: GsfElectron.h:180
reco::GsfElectron::fbrem
float fbrem() const
Definition: GsfElectron.h:809
reco::GsfElectron::full5x5_hcalOverEcal
float full5x5_hcalOverEcal(int depth=0) const
Definition: GsfElectron.h:477
reco::GsfElectron::full5x5_r9
float full5x5_r9() const
Definition: GsfElectron.h:476
output
Definition: output.py:1
reco::GsfElectron::deltaEtaSeedClusterTrackAtCalo
float deltaEtaSeedClusterTrackAtCalo() const
Definition: GsfElectron.h:226
reco::SuperCluster::clustersSize
size_t clustersSize() const
number of BasicCluster constituents
Definition: SuperCluster.h:100
reco::LeafCandidate::phi
double phi() const final
momentum azimuthal angle
Definition: LeafCandidate.h:148
funct::pow
Power< A, B >::type pow(const A &a, const B &b)
Definition: Power.h:29
reco::GsfElectron::superCluster
SuperClusterRef superCluster() const override
reference to a SuperCluster
Definition: GsfElectron.h:155
XYZTLorentzVector
math::XYZTLorentzVector XYZTLorentzVector
Definition: RawParticle.h:25
RawParticle::vertex
const XYZTLorentzVector & vertex() const
the vertex fourvector
Definition: RawParticle.h:320
reco::LeafCandidate::eta
double eta() const final
momentum pseudorapidity
Definition: LeafCandidate.h:152
Generated for CMSSW Reference Manual by doxygen 1.8.14