ATLAS Offline Software
Reconstruction/VKalVrt/VrtSecInclusive/VrtSecInclusive/details/Utilities.h
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2022 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 #ifndef _VrtSecInclusive_VrtSecInclusive_Utilities_H
6 #define _VrtSecInclusive_VrtSecInclusive_Utilities_H
7 
8 #include "VrtSecInclusive/IntersectionPos.h"
9 #include <numeric>
10 
11 #include <TString.h>
12 
13 namespace VKalVrtAthena {
14 
15  using namespace std;
16 
17  //____________________________________________________________________________________________________
18  template<class TrackT> void VrtSecInclusive::getIntersection( TrackT *trk, vector<IntersectionPos*>& layers, const Trk::Perigee* per) {
19  //
20  // get intersection point of track with various surfaces
21  //
22 
23  //--------------------
24  // main loop
25  for( auto *layer : layers ) {
26 
27  ATH_MSG_VERBOSE( " >> getIntersection(): attempt to loop " << layer->name() );
28  setIntersection( trk, layer, per );
29 
30  }
31 
32 
33  ATH_MSG_VERBOSE( " >> getIntersection(): End of getIntersection" );
34  }
35 
36 
37  //____________________________________________________________________________________________________
38  template<class TrackT>
39  void
40  VrtSecInclusive::setIntersection(TrackT* trk,
41  IntersectionPos* layer,
42  const Trk::Perigee* per)
43  {
44  const EventContext& ctx = Gaudi::Hive::currentContext();
45  std::unique_ptr<Trk::TrackParameters> Output;
46 
47  if( layer->bec() == IntersectionPos::barrel ) {
48 
49  IntersectionPos_barrel *barrel = dynamic_cast<IntersectionPos_barrel*>( layer );
50 
51  ATH_MSG_VERBOSE( " >>> setIntersection: name = " << barrel->name() << ", Radius = " << barrel->radius() );
52 
53  Amg::Transform3D trnsf;
54  trnsf.setIdentity();
55  Trk::CylinderSurface surfBorder(trnsf, barrel->radius(), 300000.); //
56  Output = m_extrapolator->extrapolateDirectly(ctx, *per,surfBorder,Trk::anyDirection,true,Trk::pion);
57 
58  barrel->x()->push_back( Output? Output->position().x() : -10000 );
59  barrel->y()->push_back( Output? Output->position().y() : -10000 );
60  barrel->z()->push_back( Output? Output->position().z() : -10000 );
61 
62 
63  } else if( layer->bec() == IntersectionPos::endcap ) {
64 
65  IntersectionPos_endcap *endcap = dynamic_cast<IntersectionPos_endcap*>( layer );
66 
67  ATH_MSG_VERBOSE( " >>> setIntersection: name = " << endcap->name() <<
68  ", zpos = " << endcap->zpos() <<
69  ", Rmin = " << endcap->rmin() <<
70  ", Rmax = " << endcap->rmax() );
71 
72  Amg::Transform3D trnsf;
73  trnsf.setIdentity();
74  trnsf.translate( Amg::Vector3D(0.,0.,endcap->zpos()) );
75  Trk::DiscSurface surfBorder(trnsf, endcap->rmin(), endcap->rmax() );
76  Output = m_extrapolator->extrapolateDirectly(ctx, *per, surfBorder, Trk::anyDirection, true, Trk::pion);
77 
78  endcap->x()->push_back( Output? Output->position().x() : -10000 );
79  endcap->y()->push_back( Output? Output->position().y() : -10000 );
80  endcap->z()->push_back( Output? Output->position().z() : -10000 );
81 
82  }
83 
84 
85  if( Output ) {
86  trk->template auxdecor<float>( Form("intersection_%s_x", layer->name().c_str()) ) = Output->position().x();
87  trk->template auxdecor<float>( Form("intersection_%s_y", layer->name().c_str()) ) = Output->position().y();
88  trk->template auxdecor<float>( Form("intersection_%s_z", layer->name().c_str()) ) = Output->position().z();
89  }
90  }
91 
92  //____________________________________________________________________________________________________
93  template<class LeptonFlavor>
94  void genSequence( const LeptonFlavor*, std::vector<unsigned>& ) {}
95 
96  template<> void genSequence( const xAOD::Muon*, std::vector<unsigned>& trackTypes );
97  template<> void genSequence( const xAOD::Electron* electron, std::vector<unsigned>& trackTypes );
98 
99  //____________________________________________________________________________________________________
100  template<class LeptonFlavor>
101  const xAOD::TrackParticle* getLeptonTrackParticle( const LeptonFlavor*, const unsigned& ) { return nullptr; }
102 
103  template<> const xAOD::TrackParticle* getLeptonTrackParticle( const xAOD::Muon* muon, const unsigned& trackType );
104  template<> const xAOD::TrackParticle* getLeptonTrackParticle( const xAOD::Electron* electron, const unsigned& trackType );
105 
106  //____________________________________________________________________________________________________
107  template<class LeptonFlavor>
108  StatusCode VrtSecInclusive::augmentDVimpactParametersToLeptons( const std::string& containerName )
109  {
110 
111  const xAOD::VertexContainer *secondaryVertexContainer( nullptr );
112  ATH_CHECK( evtStore()->retrieve( secondaryVertexContainer, "VrtSecInclusive_" + m_jp.secondaryVerticesContainerName ) );
113 
114  using LeptonContainer = DataVector<LeptonFlavor>;
115 
116  const LeptonContainer *leptonContainer( nullptr );
117  ATH_CHECK( evtStore()->retrieve( leptonContainer, containerName ) );
118 
119  if( !m_decor_d0_wrtSVs ) m_decor_d0_wrtSVs = std::make_unique<IPDecoratorType>( "d0_wrtSVs" + m_jp.augVerString );
120  if( !m_decor_z0_wrtSVs ) m_decor_z0_wrtSVs = std::make_unique<IPDecoratorType>( "z0_wrtSVs" + m_jp.augVerString );
121  if( !m_decor_pt_wrtSVs ) m_decor_pt_wrtSVs = std::make_unique<IPDecoratorType>( "pt_wrtSVs" + m_jp.augVerString );
122  if( !m_decor_eta_wrtSVs ) m_decor_eta_wrtSVs = std::make_unique<IPDecoratorType>( "eta_wrtSVs" + m_jp.augVerString );
123  if( !m_decor_phi_wrtSVs ) m_decor_phi_wrtSVs = std::make_unique<IPDecoratorType>( "phi_wrtSVs" + m_jp.augVerString );
124  if( !m_decor_d0err_wrtSVs ) m_decor_d0err_wrtSVs = std::make_unique<IPDecoratorType>( "d0err_wrtSVs" + m_jp.augVerString );
125  if( !m_decor_z0err_wrtSVs ) m_decor_z0err_wrtSVs = std::make_unique<IPDecoratorType>( "z0err_wrtSVs" + m_jp.augVerString );
126 
127  // Grouping decorators
128  std::vector< IPDecoratorType* > decor_ipWrtSVs { m_decor_d0_wrtSVs.get(), m_decor_z0_wrtSVs.get(), m_decor_pt_wrtSVs.get(), m_decor_eta_wrtSVs.get(), m_decor_phi_wrtSVs.get(), m_decor_d0err_wrtSVs.get(), m_decor_z0err_wrtSVs.get() };
129  enum { k_ip_d0, k_ip_z0, k_ip_pt, k_ip_eta, k_ip_phi, k_ip_d0err, k_ip_z0err };
130 
131  if( !m_decor_svLink ) m_decor_svLink = std::make_unique< VertexELType >( "svLinks" + m_jp.augVerString );
132 
133  // Loop over leptons
134  for( const auto& lepton : *leptonContainer ) {
135 
136  std::vector< std::vector< std::vector<float> > > ip_wrtSVs( decor_ipWrtSVs.size() ); // triple nest of { ip parameters, tracks, DVs }
137 
138  bool linkFlag { false };
139 
140  std::vector<unsigned> trackTypes;
141  genSequence<LeptonFlavor>( lepton, trackTypes );
142 
143  // Loop over lepton types
144  for( auto& trackType : trackTypes ) {
145 
146  std::vector< std::vector<float> > ip_wrtSV( decor_ipWrtSVs.size() ); // nest of { tracks, DVs }
147 
148  const auto* trk = getLeptonTrackParticle<LeptonFlavor>( lepton, trackType );
149 
150  if( !trk ) continue;
151 
152  std::map< const xAOD::Vertex*, std::vector<double> > distanceMap;
153 
154  std::vector<ElementLink< xAOD::VertexContainer > > links;
155 
156  // Loop over vertices
157  for( const auto vtx : *secondaryVertexContainer ) {
158 
159  std::vector<double> impactParameters;
160  std::vector<double> impactParErrors;
161 
162  m_fitSvc->VKalGetImpact( trk, vtx->position(), static_cast<int>( lepton->charge() ), impactParameters, impactParErrors );
163 
164  enum { k_d0, k_z0, k_theta, k_phi, k_qOverP }; // for the impact parameter
165  enum { k_d0d0, k_d0z0, k_z0z0 }; // for the par errors
166 
167  const auto& theta = impactParameters.at( k_theta );
168  const auto& phi = impactParameters.at( k_phi );
169  const auto p = fabs( 1.0 / impactParameters.at(k_qOverP) );
170  const auto pt = fabs( p * sin( theta ) );
171  const auto eta = -log( tan(theta/2.) );
172 
173  // filling the parameters to the corresponding container
174  ip_wrtSV.at( k_ip_d0 ) .emplace_back( impactParameters.at(k_d0) );
175  ip_wrtSV.at( k_ip_z0 ) .emplace_back( impactParameters.at(k_z0) );
176  ip_wrtSV.at( k_ip_pt ) .emplace_back( pt );
177  ip_wrtSV.at( k_ip_eta ) .emplace_back( eta );
178  ip_wrtSV.at( k_ip_phi ) .emplace_back( phi );
179  ip_wrtSV.at( k_ip_d0err ) .emplace_back( impactParErrors.at(k_d0d0) );
180  ip_wrtSV.at( k_ip_z0err ) .emplace_back( impactParErrors.at(k_z0z0) );
181 
182  if( !linkFlag ) {
183 
184  ElementLink<xAOD::VertexContainer> link_SV( *( dynamic_cast<const xAOD::VertexContainer*>( vtx->container() ) ), static_cast<size_t>( vtx->index() ) );
185  links.emplace_back( link_SV );
186 
187  }
188 
189  } // end of vertex loop
190 
191  // The linking to the vertices need to be done only once
192  if( !linkFlag ) {
193  ( *m_decor_svLink )( *lepton ) = links;
194  linkFlag = true;
195  }
196 
197  for( size_t ipar = 0; ipar < ip_wrtSVs.size(); ipar++ ) ip_wrtSVs.at( ipar ).emplace_back( ip_wrtSV.at( ipar ) );
198 
199  } // end of track type loop
200 
201  // decoration
202  for( size_t ipar = 0; ipar < decor_ipWrtSVs.size(); ipar++ ) ( *( decor_ipWrtSVs.at( ipar ) ) )( *lepton ) = ip_wrtSVs.at( ipar );
203 
204  } // end of lepton container loop
205 
206  return StatusCode::SUCCESS;
207  }
208 
209 
210 }
211 
212 
213 #endif
python.PyKernel.retrieve
def retrieve(aClass, aKey=None)
Definition: PyKernel.py:110
VKalVrtAthena::VrtSecInclusive::setIntersection
void setIntersection(Track *trk, IntersectionPos *bec, const Trk::Perigee *per)
Trk::anyDirection
@ anyDirection
Definition: PropDirection.h:22
xAOD::muon
@ muon
Definition: TrackingPrimitives.h:195
VKalVrtAthena::IntersectionPos_barrel
Definition: IntersectionPos.h:34
VKalVrtAthena::getLeptonTrackParticle
const xAOD::TrackParticle * getLeptonTrackParticle(const xAOD::Muon *muon, const unsigned &trackType)
Definition: Reconstruction/VKalVrt/VrtSecInclusive/src/Utilities.cxx:101
VKalVrtAthena::IntersectionPos_endcap::rmax
const double & rmax()
Definition: IntersectionPos.h:73
python.PerfMonSerializer.p
def p
Definition: PerfMonSerializer.py:743
phi
Scalar phi() const
phi method
Definition: AmgMatrixBasePlugin.h:64
VKalVrtAthena::IntersectionPos
Definition: IntersectionPos.h:18
eta
Scalar eta() const
pseudorapidity method
Definition: AmgMatrixBasePlugin.h:79
Trk::ParametersT
Dummy class used to allow special convertors to be called for surfaces owned by a detector element.
Definition: EMErrorDetail.h:25
theta
Scalar theta() const
theta method
Definition: AmgMatrixBasePlugin.h:71
module_driven_slicing.layers
layers
Definition: module_driven_slicing.py:114
test_pyathena.pt
pt
Definition: test_pyathena.py:11
VKalVrtAthena::IntersectionPos_endcap::name
const std::string & name()
Definition: IntersectionPos.h:67
Trk::DiscSurface
Definition: DiscSurface.h:54
ATH_MSG_VERBOSE
#define ATH_MSG_VERBOSE(x)
Definition: AthMsgStreamMacros.h:28
xAOD::Muon_v1
Class describing a Muon.
Definition: Muon_v1.h:38
VKalVrtAthena
Definition: AANT_Tools.cxx:24
VKalVrtAthena::IntersectionPos_endcap
Definition: IntersectionPos.h:58
EL::StatusCode
::StatusCode StatusCode
StatusCode definition for legacy code.
Definition: PhysicsAnalysis/D3PDTools/EventLoop/EventLoop/StatusCode.h:22
DMTest::links
links
Definition: CLinks_v1.cxx:22
VKalVrtAthena::genSequence
void genSequence(const xAOD::Muon *, std::vector< unsigned > &trackTypes)
Definition: Reconstruction/VKalVrt/VrtSecInclusive/src/Utilities.cxx:84
Trk::CylinderSurface
Definition: CylinderSurface.h:55
TRT::Hit::layer
@ layer
Definition: HitInfo.h:79
Trk::pion
@ pion
Definition: ParticleHypothesis.h:29
Amg::Transform3D
Eigen::Affine3d Transform3D
Definition: GeoPrimitives.h:46
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
drawFromPickle.tan
tan
Definition: drawFromPickle.py:36
DataVector
Derived DataVector<T>.
Definition: DataVector.h:581
VKalVrtAthena::IntersectionPos_endcap::x
std::vector< double > * x()
Definition: IntersectionPos.h:68
VKalVrtAthena::IntersectionPos_endcap::zpos
const double & zpos()
Definition: IntersectionPos.h:71
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition: GeoPrimitives.h:47
xAOD::Electron_v1
Definition: Electron_v1.h:34
VKalVrtAthena::IntersectionPos::endcap
@ endcap
Definition: IntersectionPos.h:29
skel.Output
Output
Definition: skel.ABtoEVGEN.py:112
python.CaloCondTools.log
log
Definition: CaloCondTools.py:20
xAOD::EgammaParameters::electron
@ electron
Definition: EgammaEnums.h:18
VKalVrtAthena::IntersectionPos::barrel
@ barrel
Definition: IntersectionPos.h:29
VKalVrtAthena::IntersectionPos_endcap::y
std::vector< double > * y()
Definition: IntersectionPos.h:69
VKalVrtAthena::IntersectionPos_endcap::rmin
const double & rmin()
Definition: IntersectionPos.h:72
xAOD::TrackParticle_v1
Class describing a TrackParticle.
Definition: TrackParticle_v1.h:43
drawFromPickle.sin
sin
Definition: drawFromPickle.py:36
VKalVrtAthena::IntersectionPos_endcap::z
std::vector< double > * z()
Definition: IntersectionPos.h:70
IntersectionPos.h
VKalVrtAthena::VrtSecInclusive::augmentDVimpactParametersToLeptons
StatusCode augmentDVimpactParametersToLeptons(const std::string &containerName)
Definition: Reconstruction/VKalVrt/VrtSecInclusive/VrtSecInclusive/details/Utilities.h:108
VKalVrtAthena::VrtSecInclusive::getIntersection
void getIntersection(Track *trk, std::vector< IntersectionPos * > &layers, const Trk::Perigee *per)
Generated on Thu Jul 4 2024 21:32:19 for ATLAS Offline Software by   doxygen 1.8.18