ATLAS Offline Software
DiTauIDVarCalculator.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2024 CERN for the benefit of the ATLAS collaboration
3 */
4 
6 
7 // Core include(s):
8 #include "AthLinks/ElementLink.h"
11 
12 // EDM include(s):
15 
16 #include "xAODJet/Jet.h"
17 #include "xAODJet/JetContainer.h"
19 
20 #include "xAODTau/DiTauJet.h"
24 
26 
27 // fastjet
28 #include "fastjet/PseudoJet.hh"
29 #include "fastjet/JetDefinition.hh"
30 #include "fastjet/AreaDefinition.hh"
31 #include "fastjet/ClusterSequenceArea.hh"
32 #include "fastjet/tools/Filter.hh"
33 #include "fastjet/tools/MassDropTagger.hh"
34 
35 using namespace DiTauRecTools;
36 using namespace fastjet;
37 
38 using TrackParticleLinks_t = std::vector<ElementLink<xAOD::TrackParticleContainer>>;
40 
41 //=================================PUBLIC-PART==================================
42 //______________________________________________________________________________
43 DiTauIDVarCalculator::DiTauIDVarCalculator( const std::string& name )
44  : AsgTool(name)
45 {
46  declareProperty( "DefaultValue", m_dDefault = -1234);
47 }
48 
49 //______________________________________________________________________________
50 DiTauIDVarCalculator::~DiTauIDVarCalculator( )
51 = default;
52 
53 //______________________________________________________________________________
55 {
56  ATH_MSG_INFO( "Initializing DiTauIDVarCalculator" );
57 
58  return StatusCode::SUCCESS;
59 }
60 
62 // Wrapper functions //
64 
66  return execute(xDiTau);
67 }
68 
70 {
71 
72  ATH_MSG_DEBUG("Calculate DiTau ID variables");
73 
74  static const SG::Decorator< int > n_subjetsDec("n_subjets");
75  n_subjetsDec(xDiTau) = n_subjets(xDiTau);
76  ATH_CHECK( decorNtracks(xDiTau) );
77 
78  static const SG::Decorator< float > ditau_ptDec("ditau_pt");
79  static const SG::Decorator< float > f_core_leadDec("f_core_lead");
80  static const SG::Decorator< float > f_core_sublDec("f_core_subl");
81  static const SG::Decorator< float > f_subjet_leadDec("f_subjet_lead");
82  static const SG::Decorator< float > f_subjet_sublDec("f_subjet_subl");
83  static const SG::Decorator< float > f_subjetsDec("f_subjets");
84  static const SG::Decorator< float > f_track_leadDec("f_track_lead");
85  static const SG::Decorator< float > f_track_sublDec("f_track_subl");
86  static const SG::Decorator< float > R_max_leadDec("R_max_lead");
87  static const SG::Decorator< float > R_max_sublDec("R_max_subl");
88  static const SG::Decorator< int > n_trackDec("n_track");
89  static const SG::Decorator< int > n_tracks_leadDec("n_tracks_lead");
90  static const SG::Decorator< int > n_tracks_sublDec("n_tracks_subl");
91  static const SG::Decorator< int > n_isotrackDec("n_isotrack");
92  static const SG::Decorator< float > R_trackDec("R_track");
93  static const SG::Decorator< float > R_track_coreDec("R_track_core");
94  static const SG::Decorator< float > R_track_allDec("R_track_all");
95  static const SG::Decorator< float > R_isotrackDec("R_isotrack");
96  static const SG::Decorator< float > R_core_leadDec("R_core_lead");
97  static const SG::Decorator< float > R_core_sublDec("R_core_subl");
98  static const SG::Decorator< float > R_tracks_leadDec("R_tracks_lead");
99  static const SG::Decorator< float > R_tracks_sublDec("R_tracks_subl");
100  static const SG::Decorator< float > M_trackDec("m_track");
101  static const SG::Decorator< float > M_track_coreDec("m_track_core");
102  static const SG::Decorator< float > M_core_leadDec("m_core_lead");
103  static const SG::Decorator< float > M_core_sublDec("m_core_subl");
104  static const SG::Decorator< float > M_track_allDec("m_track_all");
105  static const SG::Decorator< float > M_tracks_leadDec("m_tracks_lead");
106  static const SG::Decorator< float > M_tracks_sublDec("m_tracks_subl");
107  static const SG::Decorator< float > E_frac_sublDec("E_frac_subl");
108  static const SG::Decorator< float > E_frac_subsublDec("E_frac_subsubl");
109  static const SG::Decorator< float > R_subjets_sublDec("R_subjets_subl");
110  static const SG::Decorator< float > R_subjets_subsublDec("R_subjets_subsubl");
111  static const SG::Decorator< float > d0_leadtrack_leadDec("d0_leadtrack_lead");
112  static const SG::Decorator< float > d0_leadtrack_sublDec("d0_leadtrack_subl");
113  static const SG::Decorator< float > f_isotracksDec("f_isotracks");
114 
115  ditau_ptDec(xDiTau) = ditau_pt(xDiTau);
116  f_core_leadDec(xDiTau) = f_core(xDiTau, 0);
117  f_core_sublDec(xDiTau) = f_core(xDiTau, 1);
118  f_subjet_leadDec(xDiTau) = f_subjet(xDiTau, 0);
119  f_subjet_sublDec(xDiTau) = f_subjet(xDiTau, 1);
120  f_subjetsDec(xDiTau) = f_subjets(xDiTau);
121  f_track_leadDec(xDiTau) = f_track(xDiTau, 0);
122  f_track_sublDec(xDiTau) = f_track(xDiTau, 1);
123  R_max_leadDec(xDiTau) = R_max(xDiTau, 0);
124  R_max_sublDec(xDiTau) = R_max(xDiTau, 1);
125  n_trackDec(xDiTau) = n_track(xDiTau);
126  n_tracks_leadDec(xDiTau) = n_tracks(xDiTau, 0);
127  n_tracks_sublDec(xDiTau) = n_tracks(xDiTau, 1);
128  n_isotrackDec(xDiTau) = n_isotrack(xDiTau);
129  R_trackDec(xDiTau) = R_track(xDiTau);
130  R_track_coreDec(xDiTau) = R_track_core(xDiTau);
131  R_track_allDec(xDiTau) = R_track_all(xDiTau);
132  R_isotrackDec(xDiTau) = R_isotrack(xDiTau);
133  R_core_leadDec(xDiTau) = R_core(xDiTau, 0);
134  R_core_sublDec(xDiTau) = R_core(xDiTau, 1);
135  R_tracks_leadDec(xDiTau) = R_tracks(xDiTau, 0);
136  R_tracks_sublDec(xDiTau) = R_tracks(xDiTau, 1);
137  M_trackDec(xDiTau) = mass_track(xDiTau);
138  M_track_coreDec(xDiTau) = mass_track_core(xDiTau);
139  M_core_leadDec(xDiTau) = mass_core(xDiTau, 0);
140  M_core_sublDec(xDiTau) = mass_core(xDiTau, 1);
141  M_track_allDec(xDiTau) = mass_track_all(xDiTau);
142  M_tracks_leadDec(xDiTau) = mass_tracks(xDiTau, 0);
143  M_tracks_sublDec(xDiTau) = mass_tracks(xDiTau, 1);
144  E_frac_sublDec(xDiTau) = E_frac(xDiTau,1);
145  E_frac_subsublDec(xDiTau) = E_frac(xDiTau, 2);
146  R_subjets_sublDec(xDiTau) = R_subjets(xDiTau, 1);
147  R_subjets_subsublDec(xDiTau) = R_subjets(xDiTau, 2);
148  d0_leadtrack_leadDec(xDiTau) = d0_leadtrack(xDiTau, 0);
149  d0_leadtrack_sublDec(xDiTau) = d0_leadtrack(xDiTau, 1);
150  f_isotracksDec(xDiTau) = f_isotracks(xDiTau);
151 
152  return StatusCode::SUCCESS;
153 }
154 
155 //=================================PRIVATE-PART=================================
156 //______________________________________________________________________________
157 
159 {
160  int nSubjet = 0;
161  while (xDiTau.subjetPt(nSubjet) > 0. )
162  {
163  nSubjet++;
164  }
165 
166  return nSubjet;
167 }
168 
169 
171 {
172  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
173  if (n_subjetsAcc(xDiTau) < 2 ) {
174  return m_dDefault;
175  }
176 
177  return xDiTau.subjetPt(0)+xDiTau.subjetPt(1);
178 }
179 
180 
181 //______________________________________________________________________________;
182 float DiTauIDVarCalculator::f_core(const xAOD::DiTauJet& xDiTau, int iSubjet) const
183 {
184  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
185  if (iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
186  return m_dDefault;
187  }
188 
189  return xDiTau.fCore(iSubjet);
190 }
191 
192 
193 //______________________________________________________________________________;
194 float DiTauIDVarCalculator::f_subjet(const xAOD::DiTauJet& xDiTau, int iSubjet) const
195 {
196  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
197  if (iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
198  return m_dDefault;
199  }
200 
201  return xDiTau.subjetPt(iSubjet) / xDiTau.pt();
202 }
203 
204 
205 //______________________________________________________________________________;
207 {
208  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
209  if (n_subjetsAcc(xDiTau) < 2 ) {
210  return m_dDefault;
211  }
212 
213  return (xDiTau.subjetPt(0) + xDiTau.subjetPt(1))/ xDiTau.pt();
214 }
215 
216 
217 //______________________________________________________________________________;
218 float DiTauIDVarCalculator::f_track(const xAOD::DiTauJet& xDiTau, int iSubjet) const
219 {
220  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
221  if (iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
222  return m_dDefault;
223  }
224 
225  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
226  if (!trackLinksAcc.isAvailable(xDiTau) )
227  {
228  ATH_MSG_WARNING("Link not available");
229  }
230 
231  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
232 
233  TLorentzVector tlvSubjet;
234  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(iSubjet),
235  xDiTau.subjetEta(iSubjet),
236  xDiTau.subjetPhi(iSubjet),
237  xDiTau.subjetE(iSubjet) );
238 
239  TLorentzVector tlvTrack;
240  TLorentzVector tlvLeadTrack;
241  tlvLeadTrack.SetPtEtaPhiE( 0,0,0, 0);
242 
243  for (const auto &xTrack: xTracks)
244  {
245  if (!xTrack)
246  {
247  ATH_MSG_ERROR("Could not read Track");
248  continue;
249  }
250  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
251  (*xTrack)->eta(),
252  (*xTrack)->phi(),
253  (*xTrack)->e() );
254 
255  if ( tlvSubjet.DeltaR(tlvTrack) < 0.2 )
256  {
257  if (tlvLeadTrack.Pt() < tlvTrack.Pt())
258  {
259  tlvLeadTrack = tlvTrack;
260  }
261  }
262  }
263 
264  return tlvLeadTrack.Pt() / tlvSubjet.Pt();
265 }
266 
267 
268 //______________________________________________________________________________;
269 float DiTauIDVarCalculator::R_max(const xAOD::DiTauJet& xDiTau, int iSubjet) const
270 {
271  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
272  if (iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
273  return m_dDefault;
274  }
275 
276  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
277 
278  TLorentzVector tlvSubjet;
279  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(iSubjet),
280  xDiTau.subjetEta(iSubjet),
281  xDiTau.subjetPhi(iSubjet),
282  xDiTau.subjetE(iSubjet) );
283 
284  TLorentzVector tlvTrack;
285  TLorentzVector tlvRmaxTrack;
286  double Rmax = 0;
287  static const SG::ConstAccessor<float> R_subjetAcc("R_subjet");
288  for (const auto &xTrack: xTracks)
289  {
290  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
291  (*xTrack)->eta(),
292  (*xTrack)->phi(),
293  (*xTrack)->e() );
294 
295  if ( tlvSubjet.DeltaR(tlvTrack) < R_subjetAcc(xDiTau) )
296  {
297  if (tlvTrack.DeltaR(tlvSubjet) > Rmax)
298  {
299  Rmax = tlvTrack.DeltaR(tlvSubjet);
300  }
301  }
302  }
303 
304  return Rmax;
305 }
306 
307 
308 //______________________________________________________________________________;
310 {
311  return xDiTau.nTracks();
312 }
313 
314 //______________________________________________________________________________;
315 int DiTauIDVarCalculator::n_tracks(const xAOD::DiTauJet& xDiTau, int iSubjet) const
316 {
317  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
318  if (iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
319  return m_dDefault;
320  }
321 
322  static const SG::ConstAccessor<std::vector<int> > n_tracksAcc("n_tracks");
323  if (!n_tracksAcc.isAvailable(xDiTau))
324  {
325  ATH_MSG_DEBUG("n_tracks decoration not available. Try with track links.");
326 
327  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
328  if (!trackLinksAcc.isAvailable(xDiTau) )
329  {
330  ATH_MSG_WARNING("Track links not available. Return 0.");
331  return (int)m_dDefault;
332  }
333 
334  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
335 
336  TLorentzVector tlvSubjet;
337  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(iSubjet),
338  xDiTau.subjetEta(iSubjet),
339  xDiTau.subjetPhi(iSubjet),
340  xDiTau.subjetE(iSubjet) );
341 
342  TLorentzVector tlvTrack;
343  int nTracks = 0;
344  for (const auto &xTrack: xTracks)
345  {
346  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
347  (*xTrack)->eta(),
348  (*xTrack)->phi(),
349  (*xTrack)->e() );
350  if ( tlvSubjet.DeltaR(tlvTrack) < 0.2 ) nTracks++;
351  }
352 
353  return nTracks;
354  }
355 
356  return n_tracksAcc(xDiTau).at(iSubjet);
357 
358 }
359 
360 //______________________________________________________________________________;
362 {
363  return xDiTau.nIsoTracks();
364 }
365 
366 //______________________________________________________________________________;
367 float DiTauIDVarCalculator::R_tracks(const xAOD::DiTauJet& xDiTau, int iSubjet) const
368 {
369  double R_sum = 0;
370  double pt = 0;
371 
372  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
373  if (!trackLinksAcc.isAvailable(xDiTau) )
374  {
375  ATH_MSG_WARNING("Link not available");
376  }
377 
378  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
379  if (iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
380  return m_dDefault;
381  }
382 
383  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
384 
385  TLorentzVector tlvSubjet;
386  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(iSubjet),
387  xDiTau.subjetEta(iSubjet),
388  xDiTau.subjetPhi(iSubjet),
389  xDiTau.subjetE(iSubjet) );
390 
391  TLorentzVector tlvTrack;
392 
393  for (const auto& xTrack: xTracks)
394  {
395  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
396  (*xTrack)->eta(),
397  (*xTrack)->phi(),
398  (*xTrack)->e() );
399 
400  if ( tlvSubjet.DeltaR(tlvTrack) < 0.2 )
401  {
402  //ATH_MSG_DEBUG("smaller");
403  R_sum += tlvSubjet.DeltaR(tlvTrack)*tlvTrack.Pt();
404  pt += tlvTrack.Pt();
405  }
406  }
407 
408  if (pt == 0)
409  {
410  return m_dDefault;
411  }
412 
413  return R_sum / pt;
414 }
415 
416 //______________________________________________________________________________;
417 float DiTauIDVarCalculator::R_core(const xAOD::DiTauJet& xDiTau, int iSubjet) const
418 {
419  double R_sum = 0;
420  double pt = 0;
421 
422  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
423  if (!trackLinksAcc.isAvailable(xDiTau) )
424  {
425  ATH_MSG_WARNING("Link not available");
426  }
427  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
428  if (iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
429  return m_dDefault;
430  }
431 
432  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
433 
434  TLorentzVector tlvSubjet;
435  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(iSubjet),
436  xDiTau.subjetEta(iSubjet),
437  xDiTau.subjetPhi(iSubjet),
438  xDiTau.subjetE(iSubjet) );
439 
440  TLorentzVector tlvTrack;
441 
442  static const SG::ConstAccessor<float> R_coreAcc("R_core");
443  for (const auto& xTrack: xTracks)
444  {
445  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
446  (*xTrack)->eta(),
447  (*xTrack)->phi(),
448  (*xTrack)->e() );
449 
450  if ( tlvSubjet.DeltaR(tlvTrack) < R_coreAcc(xDiTau) )
451  {
452  R_sum += tlvSubjet.DeltaR(tlvTrack)*tlvTrack.Pt();
453  pt += tlvTrack.Pt();
454  }
455  }
456 
457  if (pt == 0)
458  {
459  return m_dDefault;
460  }
461 
462  return R_sum / pt;
463 }
464 
465 //______________________________________________________________________________;
467 {
468  double R_sum = 0;
469  double pt = 0;
470 
471  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
472  if (!trackLinksAcc.isAvailable(xDiTau) )
473  {
474  ATH_MSG_WARNING("Link not available");
475  }
476  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
477  if (n_subjetsAcc(xDiTau) < 2) {
478  return m_dDefault;
479  }
480 
481 
482  for (int i = 0; i<=1; i++)
483  {
484 
485  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
486 
487  TLorentzVector tlvSubjet;
488  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(i),
489  xDiTau.subjetEta(i),
490  xDiTau.subjetPhi(i),
491  xDiTau.subjetE(i) );
492 
493  TLorentzVector tlvTrack;
494 
495  static const SG::ConstAccessor<float> R_coreAcc("R_core");
496  for (const auto& xTrack: xTracks)
497  {
498  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
499  (*xTrack)->eta(),
500  (*xTrack)->phi(),
501  (*xTrack)->e() );
502  if ( tlvSubjet.DeltaR(tlvTrack) < R_coreAcc(xDiTau) )
503  {
504  //ATH_MSG_DEBUG("smaller");
505  R_sum += tlvSubjet.DeltaR(tlvTrack)*tlvTrack.Pt();
506  pt += tlvTrack.Pt();
507  }
508  }
509  }
510  if (pt == 0)
511  {
512  return m_dDefault;
513  }
514 
515  return R_sum / pt;
516 }
517 
518 //______________________________________________________________________________;
520 {
521  double R_sum = 0;
522  double pt = 0;
523 
524  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
525  if (!trackLinksAcc.isAvailable(xDiTau) )
526  {
527  ATH_MSG_WARNING("Link not available");
528  }
529  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
530  if (n_subjetsAcc(xDiTau) < 2) {
531  return m_dDefault;
532  }
533 
534  for (int i = 0; i<=1; i++)
535  {
536 
537  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
538 
539  TLorentzVector tlvSubjet;
540  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(i),
541  xDiTau.subjetEta(i),
542  xDiTau.subjetPhi(i),
543  xDiTau.subjetE(i) );
544 
545  TLorentzVector tlvTrack;
546 
547  for (const auto& xTrack: xTracks)
548  {
549  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
550  (*xTrack)->eta(),
551  (*xTrack)->phi(),
552  (*xTrack)->e() );
553 
554  if (tlvSubjet.DeltaR(tlvTrack) < 0.2)
555  {
556  R_sum += tlvSubjet.DeltaR(tlvTrack)*tlvTrack.Pt();
557  pt += tlvTrack.Pt();
558  }
559  }
560  }
561  if (pt == 0)
562  {
563  return m_dDefault;
564  }
565 
566  return R_sum / pt;
567 }
568 //______________________________________________________________________________;
570 {
571  double R_sum = 0;
572  double pt = 0;
573 
574  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
575  if (!trackLinksAcc.isAvailable(xDiTau) )
576  {
577  ATH_MSG_WARNING("Link not available");
578  }
579 
580  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
581  for (int i = 0; i<n_subjetsAcc(xDiTau); i++)
582  {
583 
584  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
585 
586  TLorentzVector tlvSubjet;
587  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(i),
588  xDiTau.subjetEta(i),
589  xDiTau.subjetPhi(i),
590  xDiTau.subjetE(i) );
591 
592  TLorentzVector tlvTrack;
593 
594  for (const auto& xTrack: xTracks)
595  {
596  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
597  (*xTrack)->eta(),
598  (*xTrack)->phi(),
599  (*xTrack)->e() );
600 
601  if (tlvSubjet.DeltaR(tlvTrack) <= 0.2)
602  {
603  R_sum += tlvSubjet.DeltaR(tlvTrack)*tlvTrack.Pt();
604  pt += tlvTrack.Pt();
605  }
606  }
607  }
608 
609  if (pt == 0)
610  {
611  return m_dDefault;
612  }
613 
614  return R_sum / pt;
615 }
616 
617 //______________________________________________________________________________;
619 {
620  double R_sum = 0;
621  double pt = 0;
622 
623  static const SG::ConstAccessor<TrackParticleLinks_t> isoTrackLinksAcc("isoTrackLinks");
624  if (!isoTrackLinksAcc.isAvailable(xDiTau) )
625  {
626  ATH_MSG_WARNING("Link not available");
627  }
628 
629  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
630  if (n_subjetsAcc(xDiTau) < 2) {
631  return m_dDefault;
632  }
633 
634  for (int i = 0; i<=1; i++)
635  {
636 
637  TrackParticleLinks_t xIsoTracks = xDiTau.isoTrackLinks();
638 
639  TLorentzVector tlvSubjet;
640  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(i),
641  xDiTau.subjetEta(i),
642  xDiTau.subjetPhi(i),
643  xDiTau.subjetE(i) );
644 
645  TLorentzVector tlvIsoTrack;
646 
647  for (const auto& xIsoTrack: xIsoTracks)
648  {
649  tlvIsoTrack.SetPtEtaPhiE( (*xIsoTrack)->pt(),
650  (*xIsoTrack)->eta(),
651  (*xIsoTrack)->phi(),
652  (*xIsoTrack)->e() );
653 
654  if (tlvSubjet.DeltaR(tlvIsoTrack) < 0.4)
655  {
656  R_sum += tlvSubjet.DeltaR(tlvIsoTrack)*tlvIsoTrack.Pt();
657  pt += tlvIsoTrack.Pt();
658  }
659  }
660  }
661 
662  if (pt == 0)
663  {
664  return m_dDefault;
665  }
666 
667  return R_sum / pt;
668 }
669 
670 //______________________________________________________________________________;
672 {
673 
674  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
675  if (!trackLinksAcc.isAvailable(xDiTau) )
676  {
677  ATH_MSG_WARNING("Link not available");
678  }
679 
680  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
681  if (n_subjetsAcc(xDiTau) < 2) {
682  return m_dDefault;
683  }
684 
685  TLorentzVector tlvallTracks;
686 
687  static const SG::ConstAccessor<float> R_coreAcc("R_core");
688  for (int i = 0; i<=1; i++)
689  {
690 
691  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
692 
693  TLorentzVector tlvSubjet;
694  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(i),
695  xDiTau.subjetEta(i),
696  xDiTau.subjetPhi(i),
697  xDiTau.subjetE(i) );
698 
699  TLorentzVector tlvTrack;
700 
701  for (const auto& xTrack: xTracks)
702  {
703  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
704  (*xTrack)->eta(),
705  (*xTrack)->phi(),
706  (*xTrack)->e() );
707  if ( tlvSubjet.DeltaR(tlvTrack) < R_coreAcc(xDiTau) )
708  {
709  //ATH_MSG_DEBUG("smaller");
710  tlvallTracks += tlvTrack;
711  }
712  }
713  }
714  if (tlvallTracks.M() < 0)
715  {
716  return m_dDefault;
717  }
718 
719  return tlvallTracks.M();
720 }
721 
722 //______________________________________________________________________________;
723 float DiTauIDVarCalculator::mass_core(const xAOD::DiTauJet& xDiTau, int iSubjet) const
724 {
725 
726  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
727  if (!trackLinksAcc.isAvailable(xDiTau) )
728  {
729  ATH_MSG_WARNING("Link not available");
730  }
731 
732  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
733  if ( iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
734  return m_dDefault;
735  }
736 
737  TLorentzVector tlvallTracks;
738 
739 
740  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
741 
742  TLorentzVector tlvSubjet;
743  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(iSubjet),
744  xDiTau.subjetEta(iSubjet),
745  xDiTau.subjetPhi(iSubjet),
746  xDiTau.subjetE(iSubjet) );
747 
748  TLorentzVector tlvTrack;
749 
750  static const SG::ConstAccessor<float> R_coreAcc("R_core");
751  for (const auto& xTrack: xTracks)
752  {
753  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
754  (*xTrack)->eta(),
755  (*xTrack)->phi(),
756  (*xTrack)->e() );
757  if ( tlvSubjet.DeltaR(tlvTrack) < R_coreAcc(xDiTau) )
758  {
759  //ATH_MSG_DEBUG("smaller");
760  tlvallTracks += tlvTrack;
761  }
762  }
763 
764  if (tlvallTracks.M() < 0)
765  {
766  return m_dDefault;
767  }
768 
769  return tlvallTracks.M();
770 }
771 
772 //______________________________________________________________________________;
773 float DiTauIDVarCalculator::mass_tracks(const xAOD::DiTauJet& xDiTau, int iSubjet) const
774 {
775 
776  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
777  if (!trackLinksAcc.isAvailable(xDiTau) )
778  {
779  ATH_MSG_WARNING("Link not available");
780  }
781 
782  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
783  if ( iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
784  return m_dDefault;
785  }
786 
787  TLorentzVector tlvallTracks;
788 
789  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
790 
791  TLorentzVector tlvSubjet;
792  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(iSubjet),
793  xDiTau.subjetEta(iSubjet),
794  xDiTau.subjetPhi(iSubjet),
795  xDiTau.subjetE(iSubjet) );
796 
797  TLorentzVector tlvTrack;
798 
799  for (const auto& xTrack: xTracks)
800  {
801  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
802  (*xTrack)->eta(),
803  (*xTrack)->phi(),
804  (*xTrack)->e() );
805  if ( tlvSubjet.DeltaR(tlvTrack) < 0.2 )
806  {
807  tlvallTracks += tlvTrack;
808  }
809  }
810 
811  if (tlvallTracks.M() < 0)
812  {
813  return m_dDefault;
814  }
815 
816  return tlvallTracks.M();
817 }
818 //______________________________________________________________________________;
820 {
821 
822  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
823  if (!trackLinksAcc.isAvailable(xDiTau) )
824  {
825  ATH_MSG_WARNING("Link not available");
826  }
827 
828  TLorentzVector tlvallTracks;
829 
830  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
831 
832  TLorentzVector tlvTrack;
833 
834  for (const auto& xTrack: xTracks)
835  {
836  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
837  (*xTrack)->eta(),
838  (*xTrack)->phi(),
839  (*xTrack)->e() );
840 
841  tlvallTracks += tlvTrack;
842  }
843 
844  if (tlvallTracks.M() < 0)
845  {
846  return m_dDefault;
847  }
848  return tlvallTracks.M();
849 }
850 //______________________________________________________________________________;
852 {
853 
854  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
855  if (!trackLinksAcc.isAvailable(xDiTau) )
856  {
857  ATH_MSG_WARNING("Link not available");
858  }
859 
860  TLorentzVector tlvallTracks;
861 
862  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
863 
864  TLorentzVector tlvTrack;
865 
866  for (const auto& xTrack: xTracks)
867  {
868  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
869  (*xTrack)->eta(),
870  (*xTrack)->phi(),
871  (*xTrack)->e() );
872 
873  tlvallTracks += tlvTrack;
874  }
875 
876 
877  TrackParticleLinks_t xIsoTracks = xDiTau.isoTrackLinks();
878 
879  TLorentzVector tlvIsoTrack;
880 
881  for (const auto& xIsoTrack: xIsoTracks)
882  {
883  tlvIsoTrack.SetPtEtaPhiE( (*xIsoTrack)->pt(),
884  (*xIsoTrack)->eta(),
885  (*xIsoTrack)->phi(),
886  (*xIsoTrack)->e() );
887 
888  tlvallTracks += tlvIsoTrack;
889  }
890 
891  if (tlvallTracks.M() < 0)
892  {
893  return m_dDefault;
894  }
895 
896  return tlvallTracks.M();
897 }
898 
899 //______________________________________________________________________________;
900 float DiTauIDVarCalculator::E_frac(const xAOD::DiTauJet& xDiTau, int iSubjet) const
901 {
902  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
903  if ( iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
904  return m_dDefault;
905  }
906 
907  return xDiTau.subjetE(iSubjet) / xDiTau.subjetE(0);
908 }
909 
910 //______________________________________________________________________________;
911 float DiTauIDVarCalculator::R_subjets(const xAOD::DiTauJet& xDiTau, int iSubjet) const
912 {
913 
914  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
915  if ( iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
916  return m_dDefault;
917  }
918 
919  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
920  if (!trackLinksAcc.isAvailable(xDiTau) )
921  {
922  ATH_MSG_WARNING("Track links not available");
923  }
924 
925  TLorentzVector tlvLeadSubjet;
926  tlvLeadSubjet.SetPtEtaPhiE( xDiTau.subjetPt(0),
927  xDiTau.subjetEta(0),
928  xDiTau.subjetPhi(0),
929  xDiTau.subjetE(0) );
930 
931  TLorentzVector tlvSubjet;
932  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(iSubjet),
933  xDiTau.subjetEta(iSubjet),
934  xDiTau.subjetPhi(iSubjet),
935  xDiTau.subjetE(iSubjet) );
936  return tlvLeadSubjet.DeltaR(tlvSubjet);
937 }
938 
939 //______________________________________________________________________________;
940 float DiTauIDVarCalculator::d0_leadtrack(const xAOD::DiTauJet& xDiTau, int iSubjet) const
941 {
942  double pt_leadtrk = 0;
943  double d0 = m_dDefault;
944  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
945  if (!trackLinksAcc.isAvailable(xDiTau) )
946  {
947  ATH_MSG_WARNING("Track links not available");
948  }
949 
950  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
951  if ( iSubjet < 0 || iSubjet >= n_subjetsAcc(xDiTau)) {
952  return m_dDefault;
953  }
954 
955  TLorentzVector tlvSubjet;
956  tlvSubjet.SetPtEtaPhiE( xDiTau.subjetPt(iSubjet),
957  xDiTau.subjetEta(iSubjet),
958  xDiTau.subjetPhi(iSubjet),
959  xDiTau.subjetE(iSubjet) );
960 
961  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
962 
963  TLorentzVector tlvTrack;
964 
965  static const SG::ConstAccessor<float> R_coreAcc("R_core");
966  for (auto &xTrack: xTracks)
967  {
968  tlvTrack.SetPtEtaPhiE( (*xTrack)->pt(),
969  (*xTrack)->eta(),
970  (*xTrack)->phi(),
971  (*xTrack)->e() );
972 
973  if (tlvTrack.DeltaR(tlvSubjet) < R_coreAcc(xDiTau))
974  {
975  if (tlvTrack.Pt() > pt_leadtrk)
976  {
977  pt_leadtrk = tlvTrack.Pt();
978  d0 = (*xTrack)->d0();
979  }
980  }
981  }
982  return d0;
983 }
984 
985 //______________________________________________________________________________;
987 {
988  double iso_pt = 0;
989  static const SG::ConstAccessor<TrackParticleLinks_t> isoTrackLinksAcc("isoTrackLinks");
990  if (!isoTrackLinksAcc.isAvailable(xDiTau) )
991  {
992  ATH_MSG_WARNING("Track links not available");
993  }
994 
995  TrackParticleLinks_t xIsoTracks = xDiTau.isoTrackLinks();
996 
997  TLorentzVector tlvIsoTrack;
998 
999  for (const auto& xIsoTrack: xIsoTracks)
1000  {
1001  tlvIsoTrack.SetPtEtaPhiE( (*xIsoTrack)->pt(),
1002  (*xIsoTrack)->eta(),
1003  (*xIsoTrack)->phi(),
1004  (*xIsoTrack)->e() );
1005 
1006  iso_pt += tlvIsoTrack.Pt();
1007  }
1008 
1009  return iso_pt / xDiTau.pt();
1010 }
1011 
1012 //______________________________________________________________________________;
1014 {
1015  static const SG::ConstAccessor<TrackParticleLinks_t> trackLinksAcc("trackLinks");
1016  if (!trackLinksAcc.isAvailable(xDiTau) )
1017  {
1018  Warning("decorNtracks()", "Track links not available.");
1019  return StatusCode::FAILURE;
1020  }
1021 
1022  static const SG::ConstAccessor<int> n_subjetsAcc("n_subjets");
1023  int nSubjets = n_subjetsAcc(xDiTau);
1024 
1025  static const SG::ConstAccessor<float> R_subjetAcc("R_subjet");
1026  float Rsubjet = R_subjetAcc(xDiTau);
1027  std::vector<int> nTracks(nSubjets, 0);
1028 
1029  TrackParticleLinks_t xTracks = xDiTau.trackLinks();
1030  for (const auto &xTrack: xTracks)
1031  {
1032  double dRmin = 1111;
1033  double itrmin = -1;
1034 
1035  for (int i=0; i<nSubjets; ++i)
1036  {
1037  TLorentzVector tlvSubjet = TLorentzVector();
1038  tlvSubjet.SetPtEtaPhiE(xDiTau.subjetPt(i),
1039  xDiTau.subjetEta(i),
1040  xDiTau.subjetPhi(i),
1041  xDiTau.subjetE(i));
1042  double dR = tlvSubjet.DeltaR((*xTrack)->p4());
1043 
1044 
1045  if ((dR < Rsubjet) && (dR < dRmin))
1046  {
1047  dRmin = dR;
1048  itrmin = i;
1049  }
1050  } // loop over subjets
1051  if (itrmin > -1) nTracks[itrmin]++;
1052  } // loop over tracks
1053 
1054  static const SG::Decorator< std::vector<int> > n_tracksDec("n_tracks");
1055  n_tracksDec(xDiTau) = nTracks;
1056 
1057  return StatusCode::SUCCESS;
1058 }
DiTauRecTools::DiTauIDVarCalculator::R_max
float R_max(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:269
DiTauRecTools::DiTauIDVarCalculator::ditau_pt
float ditau_pt(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:170
xAOD::DiTauJet_v1::pt
virtual double pt() const
The transverse momentum ( ) of the particle.
Jet.h
DiTauRecTools::DiTauIDVarCalculator::f_core
float f_core(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:182
xAOD::DiTauJet_v1::fCore
float fCore(unsigned int numSubjet) const
Definition: DiTauJet_v1.cxx:168
DiTauRecTools::DiTauIDVarCalculator::R_track
float R_track(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:519
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
fastjet
Definition: FastJetLinkBase.h:22
DiTauRecTools::DiTauIDVarCalculator::f_isotracks
float f_isotracks(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:986
initialize
void initialize()
Definition: run_EoverP.cxx:894
DiTauRecTools::DiTauIDVarCalculator::mass_track_all
float mass_track_all(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:851
xAOD::DiTauJet_v1::subjetPhi
float subjetPhi(unsigned int numSubjet) const
Definition: DiTauJet_v1.cxx:112
DiTauRecTools::DiTauIDVarCalculator::n_track
static int n_track(const xAOD::DiTauJet &xDiTau)
Definition: DiTauIDVarCalculator.cxx:309
test_pyathena.pt
pt
Definition: test_pyathena.py:11
DiTauRecTools::DiTauIDVarCalculator::R_tracks
float R_tracks(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:367
DiTauRecTools::DiTauIDVarCalculator::d0_leadtrack
float d0_leadtrack(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:940
DiTauJetContainer.h
SG::ConstAccessor< int >
JetConstituentVector.h
This file defines helper classes to deal with jet constituents.
DiTauJetAuxContainer.h
DiTauRecTools::DiTauIDVarCalculator::R_subjets
float R_subjets(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:911
DiTauIDVarCalculator.h
DiTauRecTools::DiTauIDVarCalculator::n_isotrack
static int n_isotrack(const xAOD::DiTauJet &xDiTau)
Definition: DiTauIDVarCalculator.cxx:361
DiTauRecTools::DiTauIDVarCalculator::R_isotrack
float R_isotrack(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:618
xAOD::DiTauJet_v1::subjetE
float subjetE(unsigned int numSubjet) const
Definition: DiTauJet_v1.cxx:122
DiTauRecTools::DiTauIDVarCalculator::R_core
float R_core(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:417
DiTauRecTools::DiTauIDVarCalculator::R_track_all
float R_track_all(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:569
DiTauRecTools::DiTauIDVarCalculator::mass_track
float mass_track(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:819
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
DiTauRecTools::DiTauIDVarCalculator::mass_core
float mass_core(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:723
SG::Decorator< int >
lumiFormat.i
int i
Definition: lumiFormat.py:85
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
ATH_MSG_DEBUG
#define ATH_MSG_DEBUG(x)
Definition: AthMsgStreamMacros.h:29
TrackParticleLinks_t
std::vector< ElementLink< xAOD::TrackParticleContainer > > TrackParticleLinks_t
Definition: DiTauIDVarCalculator.cxx:38
TRT::Track::d0
@ d0
Definition: InnerDetector/InDetCalibEvent/TRT_CalibData/TRT_CalibData/TrackInfo.h:62
DiTauRecTools::DiTauIDVarCalculator::mass_track_core
float mass_track_core(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:671
DiTauRecTools::DiTauIDVarCalculator::E_frac
float E_frac(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:900
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
DiTauRecTools::DiTauIDVarCalculator::m_dDefault
float m_dDefault
Definition: DiTauIDVarCalculator.h:94
DiTauRecTools::DiTauIDVarCalculator::f_subjet
float f_subjet(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:194
TauJetContainer.h
DiTauJet.h
xAOD::DiTauJet_v1::subjetEta
float subjetEta(unsigned int numSubjet) const
Definition: DiTauJet_v1.cxx:102
xAOD::DiTauJet_v1::nIsoTracks
size_t nIsoTracks() const
Definition: DiTauJet_v1.cxx:290
xAOD::DiTauJet_v1::nTracks
size_t nTracks() const
Definition: DiTauJet_v1.cxx:258
name
std::string name
Definition: Control/AthContainers/Root/debug.cxx:221
DiTauRecTools
Implementation of boosted di-tau ID.
Definition: DiTauDiscriminantTool.h:36
DiTauRecTools::DiTauIDVarCalculator::R_track_core
float R_track_core(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:466
DiTauRecTools::DiTauIDVarCalculator::calculateIDVariables
virtual StatusCode calculateIDVariables(const xAOD::DiTauJet &xDiTau)
Definition: DiTauIDVarCalculator.cxx:65
DiTauRecTools::DiTauIDVarCalculator::f_track
float f_track(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:218
EventInfo.h
TrackParticle.h
JetContainer.h
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
DiTauRecTools::DiTauIDVarCalculator::decorNtracks
static StatusCode decorNtracks(const xAOD::DiTauJet &xDiTau)
Definition: DiTauIDVarCalculator.cxx:1013
DiTauRecTools::DiTauIDVarCalculator::f_subjets
float f_subjets(const xAOD::DiTauJet &xDiTau) const
Definition: DiTauIDVarCalculator.cxx:206
xAOD::DiTauJet_v1
Definition: DiTauJet_v1.h:31
xAOD::DiTauJet_v1::isoTrackLinks
const TrackParticleLinks_t & isoTrackLinks() const
SG::ConstAccessor::isAvailable
bool isAvailable(const ELT &e) const
Test to see if this variable exists in the store.
xAOD::DiTauJet_v1::subjetPt
float subjetPt(unsigned int numSubjet) const
Definition: DiTauJet_v1.cxx:92
ConstAccessor.h
Helper class to provide constant type-safe access to aux data.
DiTauRecTools::DiTauIDVarCalculator::n_tracks
int n_tracks(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:315
xAOD::DiTauJet_v1::trackLinks
const TrackParticleLinks_t & trackLinks() const
DiTauRecTools::DiTauIDVarCalculator::n_subjets
static float n_subjets(const xAOD::DiTauJet &xDiTau)
Definition: DiTauIDVarCalculator.cxx:158
Decorator.h
Helper class to provide type-safe access to aux data.
TrackParticleContainer.h
DiTauRecTools::DiTauIDVarCalculator::execute
virtual StatusCode execute(const xAOD::DiTauJet &xDiTau) override
Declare the interface that the class provides.
Definition: DiTauIDVarCalculator.cxx:69
DiTauRecTools::DiTauIDVarCalculator::mass_tracks
float mass_tracks(const xAOD::DiTauJet &xDiTau, int iSubjet) const
Definition: DiTauIDVarCalculator.cxx:773