ATLAS Offline Software
FastCaloSimCaloTransportation.cxx
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2023 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 /* Athena includes */
7 #include "GaudiKernel/ToolHandle.h"
8 #include "GaudiKernel/IPartPropSvc.h"
9 
10 #include "CxxUtils/inline_hints.h"
11 /* Header include */
13 
14 /* ISF includes */
17 
18 /* Tracking includes */
20 
21 /* Geometry primitives */
23 
24 /* Particle data */
25 #include "HepPDT/ParticleDataTable.hh"
26 
27 /* Transport steps will be return as G4FieldTracks*/
28 #include "G4FieldTrack.hh"
30 
31 
32 FastCaloSimCaloTransportation::FastCaloSimCaloTransportation(const std::string& t, const std::string& n, const IInterface* p)
33  : base_class(t,n,p)
34 {
35 
36 }
37 
39 {
40  ATH_MSG_INFO( "Initializing FastCaloSimCaloTransportation" );
41 
42  // Get TimedExtrapolator
43  if(!m_extrapolator.empty()){
44  ATH_CHECK(m_extrapolator.retrieve());
45  ATH_MSG_INFO("Extrapolator retrieved "<< m_extrapolator);
46  }
47 
48  return StatusCode::SUCCESS;
49 
50 }
51 
52 
54  ATH_MSG_INFO( "Finalizing FastCaloSimCaloTransportation" );
55  return StatusCode::SUCCESS;
56 }
57 
58 
59 std::vector<G4FieldTrack> FastCaloSimCaloTransportation::transport(const TFCSTruthState* truth, bool forceNeutral) const{
60  // Start calo extrapolation
61  ATH_MSG_DEBUG ("[ fastCaloSim transport ] processing particle "<<truth->pdgid() );
62 
63  auto hitVector = std::make_unique<std::vector<Trk::HitInfo>>();
64 
65  int pdgId = truth->pdgid();
66  double charge = HepPDT::ParticleID(pdgId).charge();
67  if (forceNeutral) charge = 0.;
68 
69  // particle Hypothesis for the extrapolation
71 
72  ATH_MSG_DEBUG ("particle hypothesis "<< pHypothesis);
73 
74  // geantinos not handled by PdgToParticleHypothesis - fix there
75  if (pdgId == 999) pHypothesis = Trk::geantino;
76 
77  Amg::Vector3D pos = Amg::Vector3D(truth->vertex().X(), truth->vertex().Y(), truth->vertex().Z());
78 
79  Amg::Vector3D mom(truth->X(), truth->Y(), truth->Z());
80 
81  ATH_MSG_DEBUG( "[ fastCaloSim transport ] x from position eta="<<pos.eta()<<" phi="<<pos.phi()<<" d="<<pos.mag()<<" pT="<<mom.perp() );
82 
83  // input parameters : curvilinear parameters
85 
86  double freepath = -1.;
87  double tDec = freepath > 0. ? freepath : -1.;
88  int decayProc = 0;
89 
90  Trk::TimeLimit timeLim(tDec, 0., decayProc);
91  Trk::PathLimit pathLim(-1., 0);
93 
94  // first extrapolation to reach the ID boundary
95  ATH_MSG_DEBUG( "[ fastCaloSim transport ] before calo entrance ");
96 
97  // get CaloEntrance if not done already
98  if(!m_caloEntrance.get()){
99  m_caloEntrance.set(m_extrapolator->trackingGeometry()->trackingVolume(m_caloEntranceName));
100  if(!m_caloEntrance.get())
101  ATH_MSG_WARNING("CaloEntrance not found");
102  else
103  ATH_MSG_DEBUG("CaloEntrance found");
104  }
105 
106  ATH_MSG_DEBUG( "[ fastCaloSim transport ] after calo entrance ");
107 
108  std::unique_ptr<const Trk::TrackParameters> caloEntry = nullptr;
109 
110  if(m_caloEntrance.get() && m_caloEntrance.get()->inside(pos, 0.001) && !m_extrapolator->trackingGeometry()->atVolumeBoundary(pos,m_caloEntrance.get(), 0.001)){
111  std::vector<Trk::HitInfo>* dummyHitVector = nullptr;
112  if (charge == 0){
113  caloEntry =
114  m_extrapolator->transportNeutralsWithPathLimit(inputPar,
115  pathLim,
116  timeLim,
118  pHypothesis,
119  dummyHitVector,
120  nextGeoID,
121  m_caloEntrance.get());
122  }else{
123  caloEntry = m_extrapolator->extrapolateWithPathLimit(inputPar,
124  pathLim,
125  timeLim,
127  pHypothesis,
128  dummyHitVector,
129  nextGeoID,
130  m_caloEntrance.get());
131  }
132  } else
133  caloEntry = inputPar.uniqueClone();
134 
135  ATH_MSG_DEBUG( "[ fastCaloSim transport ] after calo caloEntry ");
136 
137  if(caloEntry){
138  std::unique_ptr<const Trk::TrackParameters> eParameters = nullptr;
139 
140  // save Calo entry hit (fallback info)
141  hitVector->push_back(Trk::HitInfo(caloEntry->uniqueClone(), timeLim.time, nextGeoID, 0.));
142 
144  "[ fastCaloSim transport ] starting Calo transport from position eta="
145  << caloEntry->position().eta() << " phi=" << caloEntry->position().phi()
146  << " d=" << caloEntry->position().mag());
147 
148  std::vector<Trk::HitInfo>* rawHitVector = hitVector.get();
149  if (charge == 0){
150  eParameters =
151  m_extrapolator->transportNeutralsWithPathLimit(*caloEntry,
152  pathLim,
153  timeLim,
155  pHypothesis,
156  rawHitVector,
157  nextGeoID);
158  }else{
159  eParameters =
160  m_extrapolator->extrapolateWithPathLimit(*caloEntry,
161  pathLim,
162  timeLim,
164  pHypothesis,
165  rawHitVector,
166  nextGeoID);
167  }
168  // save Calo exit hit (fallback info)
169  if(eParameters) hitVector->push_back(Trk::HitInfo(std::move(eParameters), timeLim.time, nextGeoID, 0.));
170  } //if caloEntry
171 
172  //used to identify ID-Calo boundary in tracking tools
173  int IDCaloBoundary = 3000;
174 
175  if(msgLvl(MSG::DEBUG)){
177  while (it < hitVector->end()){
178  int sample=(*it).detID;
179  Amg::Vector3D hitPos = (*it).trackParms->position();
180  ATH_MSG_DEBUG(" HIT: layer="<<sample<<" sample="<<sample-IDCaloBoundary<<" eta="<<hitPos.eta()<<" phi="<<hitPos.phi()<<" d="<<hitPos.mag());
181  ++it;
182  }
183  }
184 
186  while(it2 < hitVector->end()){
187  int sample=(*it2).detID;
188  Amg::Vector3D hitPos = (*it2).trackParms->position();
189  ATH_MSG_DEBUG(" HIT: layer="<<sample<<" sample="<<sample-IDCaloBoundary<<" eta="<<hitPos.eta()<<" phi="<<hitPos.phi()<<" r="<<hitPos.perp()<<" z="<<hitPos[Amg::z]);
190  ++it2;
191  }
192 
193  // Extrapolation may fail for very low pT charged particles. Enforce charge 0 to prevent this
194  if (!forceNeutral && hitVector->empty()){
195  ATH_MSG_DEBUG("forcing neutral charge in FastCaloSimCaloTransportation::caloHits");
196  transport(truth, true);
197  }
198  // Don't expect this ever to happen. Nevertheless, error handling should be improved.
199  // This may require changes in periphery (adjustments after setting function type to StatusCode)
200  else if(hitVector->empty()) ATH_MSG_ERROR("Empty hitVector even after forcing neutral charge. This may cause a segfault soon.");
201 
202 
203  std::vector<G4FieldTrack> caloSteps = convertToFieldTrack(*hitVector);
204 
205  return caloSteps;
206 
207 }
208 
209 
210 std::vector<G4FieldTrack> FastCaloSimCaloTransportation::convertToFieldTrack(const std::vector<Trk::HitInfo>& vec) const{
211 
212  std::vector<G4FieldTrack> caloSteps;
213  for (auto& step : vec){
214  G4FieldTrack track = G4FieldTrack(' ');
215  track.SetPosition(Amg::EigenToHep3Vector(step.trackParms->position()));
216  track.SetMomentum(Amg::EigenToHep3Vector(step.trackParms->momentum()));
217  track.SetChargeAndMoments(step.trackParms->charge());
218  caloSteps.push_back(track);
219  }
220 
221  return caloSteps;
222 
223 }
xAOD::iterator
JetConstituentVector::iterator iterator
Definition: JetConstituentVector.cxx:68
FastCaloSimCaloTransportation::finalize
virtual StatusCode finalize() override final
Definition: FastCaloSimCaloTransportation.cxx:53
inline_hints.h
FastCaloSimCaloTransportation::m_extrapolator
PublicToolHandle< Trk::ITimedExtrapolator > m_extrapolator
Definition: FastCaloSimCaloTransportation.h:38
ATH_MSG_INFO
#define ATH_MSG_INFO(x)
Definition: AthMsgStreamMacros.h:31
Trk::PdgToParticleHypothesis::convert
Trk::ParticleHypothesis convert(int pdg, bool &stable, bool &exiting, double charge=1.) const
Converter method : PDG -> Particle Hyptothesis.
FastCaloSimCaloTransportation::initialize
virtual StatusCode initialize() override final
Definition: FastCaloSimCaloTransportation.cxx:38
FastCaloSimCaloTransportation::transport
virtual std::vector< G4FieldTrack > transport(const TFCSTruthState *truth, bool forceNeutral=false) const override final
Definition: FastCaloSimCaloTransportation.cxx:59
FastCaloSimCaloTransportation::m_caloEntrance
CxxUtils::CachedPointer< const Trk::TrackingVolume > m_caloEntrance
Definition: FastCaloSimCaloTransportation.h:40
Trk::ParametersBase::position
const Amg::Vector3D & position() const
Access method for the position.
Trk::ParametersBase::uniqueClone
std::unique_ptr< ParametersBase< DIM, T > > uniqueClone() const
clone method for polymorphic deep copy returning unique_ptr; it is not overriden, but uses the existi...
Definition: ParametersBase.h:97
hitVector
std::vector< FPGATrackSimHit > hitVector
Definition: FPGATrackSimCluster.h:22
skel.it
it
Definition: skel.GENtoEVGEN.py:396
Amg::EigenToHep3Vector
CLHEP::Hep3Vector EigenToHep3Vector(const Amg::Vector3D &eigenvector)
Converts an Eigen-based Amg::Vector3D into a CLHEP-based CLHEP::Hep3Vector.
Definition: CLHEPtoEigenConverter.h:147
Trk::alongMomentum
@ alongMomentum
Definition: PropDirection.h:20
vec
std::vector< size_t > vec
Definition: CombinationsGeneratorTest.cxx:12
read_hist_ntuple.t
t
Definition: read_hist_ntuple.py:5
Trk::Calo
@ Calo
Definition: GeometrySignature.h:28
TFCSTruthState::vertex
const TLorentzVector & vertex() const
Definition: TFCSTruthState.h:28
CxxUtils::CachedPointer::set
void set(pointer_t elt) const
Set the element, assuming it is currently null.
mergePhysValFiles.end
end
Definition: DataQuality/DataQualityUtils/scripts/mergePhysValFiles.py:93
Trk::ParticleHypothesis
ParticleHypothesis
Definition: ParticleHypothesis.h:25
Amg::z
@ z
Definition: GeoPrimitives.h:36
Trk::GeometrySignature
GeometrySignature
Definition: GeometrySignature.h:24
FastCaloSimCaloTransportation::m_caloEntranceName
StringProperty m_caloEntranceName
Definition: FastCaloSimCaloTransportation.h:41
python.utils.AtlRunQueryDQUtils.p
p
Definition: AtlRunQueryDQUtils.py:210
ParticleGun_EoverP_Config.mom
mom
Definition: ParticleGun_EoverP_Config.py:63
ATH_MSG_ERROR
#define ATH_MSG_ERROR(x)
Definition: AthMsgStreamMacros.h:33
FullCPAlgorithmsTest_eljob.sample
sample
Definition: FullCPAlgorithmsTest_eljob.py:113
beamspotman.n
n
Definition: beamspotman.py:731
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
AthAlgTool.h
Trk::geantino
@ geantino
Definition: ParticleHypothesis.h:26
FastCaloSimCaloTransportation::convertToFieldTrack
std::vector< G4FieldTrack > convertToFieldTrack(const std::vector< Trk::HitInfo > &vec) const
Definition: FastCaloSimCaloTransportation.cxx:210
jobOptions.ParticleID
ParticleID
Definition: jobOptions.decayer.py:85
ATH_CHECK
#define ATH_CHECK
Definition: AthCheckMacros.h:40
Trk::CurvilinearParametersT
Definition: CurvilinearParametersT.h:48
Trk::TimeLimit
Definition: HelperStructs.h:58
CLHEPtoEigenConverter.h
charge
double charge(const T &p)
Definition: AtlasPID.h:756
Trk::PathLimit
Definition: HelperStructs.h:34
TFCSTruthState::pdgid
int pdgid() const
Definition: TFCSTruthState.h:25
Amg::Vector3D
Eigen::Matrix< double, 3, 1 > Vector3D
Definition: GeoPrimitives.h:47
python.LumiBlobConversion.pos
pos
Definition: LumiBlobConversion.py:18
TFCSTruthState.h
ATH_MSG_WARNING
#define ATH_MSG_WARNING(x)
Definition: AthMsgStreamMacros.h:32
GeoPrimitivesHelpers.h
DEBUG
#define DEBUG
Definition: page_access.h:11
FastCaloSimCaloTransportation.h
Trk::TimeLimit::time
float time
Definition: HelperStructs.h:60
CxxUtils::CachedPointer::get
pointer_t get() const
Return the current value of the element.
LArCellBinning.step
step
Definition: LArCellBinning.py:158
FastCaloSimCaloTransportation::FastCaloSimCaloTransportation
FastCaloSimCaloTransportation(const std::string &t, const std::string &n, const IInterface *p)
Definition: FastCaloSimCaloTransportation.cxx:32
TrackingGeometry.h
Trk::HitInfo
Definition: HelperStructs.h:12
FastCaloSim_CaloCell_ID.h
xAOD::track
@ track
Definition: TrackingPrimitives.h:512
TFCSTruthState
Definition: TFCSTruthState.h:13
FastCaloSimCaloTransportation::m_pdgToParticleHypothesis
Trk::PdgToParticleHypothesis m_pdgToParticleHypothesis
Definition: FastCaloSimCaloTransportation.h:42