ATLAS Offline Software
ExtrapolationEngineTest.icc
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2025 CERN for the benefit of the ATLAS collaboration
3 */
4 
5 // Root includes
6 #include "TTree.h"
7 #include "TString.h"
8 // Amg
9 #include "GeoPrimitives/GeoPrimitivesToStringConverter.h"
10 // Trk includes
11 #include "TrkParameters/TrackParameters.h"
12 #include "TrkNeutralParameters/NeutralParameters.h"
13 #include "TrkValInterfaces/IPositionMomentumWriter.h"
14 #include "TrkEventPrimitives/ParticleHypothesis.h"
15 #include "TrkEventPrimitives/PdgToParticleHypothesis.h"
16 #include "TrkExUtils/MaterialInteraction.h"
17 #include "TrkGeometry/Layer.h"
18 #include "InDetIdentifier/PixelID.h"
19 #include "InDetIdentifier/SCT_ID.h"
20 #include "HGTD_Identifier/HGTD_ID.h"
21 #include "AtlasDetDescr/AtlasDetectorID.h"
22 #include "Identifier/Identifier.h"
23 #include "Identifier/IdentifierHash.h"
24 #include "InDetReadoutGeometry/SiDetectorElement.h"
25 
26 template < class T, class P > StatusCode Trk::ExtrapolationEngineTest::runTestT() {
27  // ----------------- creation of the surfaces & track parameters -------------
28  for (size_t it = 0; it < TrkExUnitTestBase::m_numTests; ++it) {
29  // verbose output
30  ATH_MSG_DEBUG("===> starting test " << it << " <<===");
31  // create the curvilinear parameters
32  double eta = (m_scanMode) ? m_currentEta : m_etaMin + (m_etaMax - m_etaMin) *
33  Trk::TrkExUnitTestBase::m_flatDist->shoot();
34  double theta = 2. * atan(exp(-eta));
35  double phi = (m_scanMode) ? m_currentPhi : m_phiMin + (m_phiMax - m_phiMin) *
36  Trk::TrkExUnitTestBase::m_flatDist->shoot();
37  double pt = m_ptMin + (m_ptMax - m_ptMin) * Trk::TrkExUnitTestBase::m_flatDist->shoot();
38  double p = pt / sin(theta);
39  double q = m_splitCharge ? m_charge *
40  -1. : (m_parametersMode ? (Trk::TrkExUnitTestBase::m_flatDist->shoot() > 0.5 ? 1. : -1) : 1.); // charge
41  // or
42  // neutral
43 
44  // initializa the validation variables
45  m_endSuccessful = 0;
46  m_endPositionX = 0.;
47  m_endPositionY = 0.;
48  m_endPositionZ = 0.;
49  m_endPositionR = 0.;
50  m_endPhi = 0.;
51  m_endEta = 0.;
52  m_endTheta = 0.;
53  m_endP = 0.;
54  m_endPt = 0.;
55 
56  // material collection
57  m_materialThicknessInX0 = 0.;
58  m_materialThicknessInL0 = 0.;
59  m_materialThicknessZARho = 0.;
60  m_materialEmulatedIonizationLoss = 0.;
61  m_materialThicknessInX0Sensitive = 0.;
62  m_materialThicknessInX0Passive = 0.;
63  m_materialThicknessInX0Boundary = 0.;
64 
65  m_materialThicknessInX0Bwd = 0.;
66  m_materialThicknessInL0Bwd = 0.;
67 
68  m_materialThicknessInX0Cylinder = 0.;
69  m_materialThicknessInX0Disc = 0.;
70  m_materialThicknessInX0Plane = 0.;
71 
72  if (m_collectMaterial) {
73  m_materialThicknessInX0Accumulated->clear();
74  m_materialThicknessInX0Steps->clear();
75  m_materialThicknessInX0Steps->clear();
76  m_materialPositionX->clear();
77  m_materialPositionY->clear();
78  m_materialPositionZ->clear();
79  m_materialPositionR->clear();
80  m_materialPositionP->clear();
81  m_materialPositionPt->clear();
82  m_materialScaling->clear();
83  m_stepDirection->clear();
84  }
85 
86  if (m_collectSensitive) {
87  m_sensitiveLayerIndex->clear();
88  m_sensitiveSurfaceType->clear();
89  m_sensitiveCenterPosX->clear();
90  m_sensitiveCenterPosY->clear();
91  m_sensitiveCenterPosZ->clear();
92  m_sensitiveCenterPosR->clear();
93  m_sensitiveCenterPosPhi->clear();
94  m_sensitiveLocalPosX->clear();
95  m_sensitiveLocalPosY->clear();
96  m_sensitiveLocalPosR->clear();
97  m_sensitiveLocalPosPhi->clear();
98  m_sensitiveDetector->clear();
99  m_sensitiveIsInnermost->clear();
100  m_sensitiveIsNextToInnermost->clear();
101  m_sensitiveBarrelEndcap->clear();
102  m_sensitiveLayerDisc->clear();
103  m_sensitiveEtaModule->clear();
104  m_sensitivePhiModule->clear();
105  m_sensitiveSide->clear();
106  m_sensitiveNearBondGap->clear();
107  m_sensitiveisInside->clear();
108  m_sensitiveisInsideBound->clear();
109  m_materialtInX0AccumulatedUpTo->clear();
110  }
111 
112  for (size_t ip = 0; ip < m_parameterNames.size(); ++ip) {
113  // clear
114  m_pPositionX[ip]->clear();
115  m_pPositionY[ip]->clear();
116  m_pPositionZ[ip]->clear();
117  m_pPositionR[ip]->clear();
118  m_pPhi[ip]->clear();
119  m_pTheta[ip]->clear();
120  m_pEta[ip]->clear();
121  m_pP[ip]->clear();
122  m_pPt[ip]->clear();
123  }
124 
125  Amg::Vector3D momentum(p* sin(theta)* cos(phi), p* sin(theta)* sin(phi), p* cos(theta));
126  // create the start parameters
127  double d0 =
128  m_smearProductionVertex ? (m_smearFlatOriginT ? (m_d0Min + (m_d0Max - m_d0Min) *
129  Trk::TrkExUnitTestBase::m_flatDist->shoot()) : Trk::
130  TrkExUnitTestBase::m_gaussDist->shoot() * m_sigmaOriginT) : 0.;
131  double z0 =
132  m_smearProductionVertex ? (m_smearFlatOriginZ ? (m_z0Min + (m_z0Max - m_z0Min) *
133  Trk::TrkExUnitTestBase::m_flatDist->shoot()) : Trk::
134  TrkExUnitTestBase::m_gaussDist->shoot() * m_sigmaOriginZ) : 0.;
135  // create the start parameters
136  P startParameters(d0, z0, phi, theta, q / p, Trk::PerigeeSurface(Amg::Vector3D(0., 0., 0.)));
137 
138  m_startPositionX = startParameters.position().x();
139  m_startPositionY = startParameters.position().y();
140  m_startPositionZ = startParameters.position().z();
141  m_startPositionR = startParameters.position().perp();
142  m_startPhi = phi;
143  m_startEta = eta;
144  m_startTheta = theta;
145  m_startPt = pt;
146  m_startP = p;
147  m_charge = q;
148 
149  // the step surfaces for stepwise extrapolation afterwards
150  std::vector < const Trk::Surface * > stepSurfaces;
151 
152  // position / momentum writer
153  if (!m_posmomWriter.empty() && m_writeTTree) {
154  // masses and pdg
155  Trk::PdgToParticleHypothesis pdgFromHypothesis;
156  // now convert
157  double mass = Trk::ParticleMasses::mass[m_particleHypothesis];
158  int pdg = pdgFromHypothesis.convert(static_cast<Trk::ParticleHypothesis>(m_particleHypothesis.value()), m_charge);
159  m_posmomWriter->initializeTrack(startParameters.position(),
160  startParameters.momentum(),
161  mass, pdg);
162  }
163 
164  // setup the extrapolation how you'd like it
165  Trk::ExtrapolationCell < T > ecc(startParameters);
166  ecc.setParticleHypothesis(static_cast<Trk::ParticleHypothesis>(m_particleHypothesis.value()));
167  ecc.addConfigurationMode(Trk::ExtrapolationMode::StopAtBoundary);
168  ecc.addConfigurationMode(Trk::ExtrapolationMode::FATRAS);
169  if (m_collectSensitive) ecc.addConfigurationMode(Trk::ExtrapolationMode::CollectSensitive);
170  if (m_collectPassive) ecc.addConfigurationMode(Trk::ExtrapolationMode::CollectPassive);
171  if (m_collectBoundary) ecc.addConfigurationMode(Trk::ExtrapolationMode::CollectBoundary);
172  if (m_collectMaterial) ecc.addConfigurationMode(Trk::ExtrapolationMode::CollectMaterial);
173  // force a stop in the extrapoaltion mode
174  if (m_pathLimit > 0.) {
175  ecc.pathLimit = m_pathLimit;
176  ecc.addConfigurationMode(Trk::ExtrapolationMode::StopWithPathLimit);
177  }
178  // screen output
179  ATH_MSG_DEBUG("===> forward extrapolation - collecting information <<===");
180  // call the extrapolation engine
181  Trk::ExtrapolationCode eCode = m_extrapolationEngine->extrapolate(ecc);
182  // save the success code
183  m_endSuccessful = int(eCode.code);
184  // end the parameters if there
185  if (eCode.isSuccess() && ecc.endParameters && m_writeTTree) {
186  m_endPositionX = ecc.endParameters->position().x();
187  m_endPositionY = ecc.endParameters->position().y();
188  m_endPositionZ = ecc.endParameters->position().z();
189  m_endPositionR = ecc.endParameters->position().perp();
190  m_endPhi = ecc.endParameters->momentum().phi();
191  m_endEta = ecc.endParameters->momentum().eta();
192  m_endTheta = ecc.endParameters->momentum().theta();
193  m_endP = ecc.endParameters->momentum().mag();
194  m_endPt = ecc.endParameters->momentum().perp();
195  m_endPathLength = ecc.pathLength;
196  // fill in the step parameters
197  if (!fillStepInformationT < T, P > (ecc, 1, stepSurfaces).isSuccess()) ATH_MSG_VERBOSE(
198  "Somthing went wrong with recording the step information.");
199  // for the memory cleanup
200  const T* feParameters = ecc.endParameters;
201  // do the back extrapolation
202  if (m_backExtrapolation) {
203  // let's restart at the destination
204  ecc.restartAtDestination();
205  // screen output
206  ATH_MSG_DEBUG("===> backward extrapolation - collecting information <<===");
207  // call the extrapolation engine
208  Trk::ExtrapolationCode eCodeBwd = m_extrapolationEngine->extrapolate(ecc, &startParameters.associatedSurface());
209  // safe the success code
210  m_backSuccessful = int(eCodeBwd.code);
211  // end the parameters if there
212  if (eCodeBwd.isSuccess() && ecc.endParameters) {
213  m_backPositionX = ecc.endParameters->position().x();
214  m_backPositionY = ecc.endParameters->position().y();
215  m_backPositionZ = ecc.endParameters->position().z();
216  m_backPositionR = ecc.endParameters->position().perp();
217  m_backPhi = ecc.endParameters->momentum().phi();
218  m_backEta = ecc.endParameters->momentum().eta();
219  m_backTheta = ecc.endParameters->momentum().theta();
220  m_backP = ecc.endParameters->momentum().mag();
221  m_backPt = ecc.endParameters->momentum().perp();
222  // fill in the step parameters
223  if (!fillStepInformationT < T, P > (ecc, -1, stepSurfaces).isSuccess()) ATH_MSG_VERBOSE(
224  "Somthing went wrong with recording the step information.");
225  // memory cleanup
226  delete ecc.endParameters;
227  }
228  }
229  // next test - run stepwise extrapolation from one surface to another - this is for Kalman filtering
230  if (m_stepwiseExtrapolation && stepSurfaces.size()) {
231  // the ExtrapolationCell for stepwise
232  Trk::ExtrapolationCell < T > eccSteps(startParameters);
233  eccSteps.setParticleHypothesis(static_cast<Trk::ParticleHypothesis>(m_particleHypothesis.value()));
234  eccSteps.addConfigurationMode(Trk::ExtrapolationMode::StopAtBoundary);
235  if (m_collectMaterial) eccSteps.addConfigurationMode(Trk::ExtrapolationMode::CollectMaterial);
236  std::vector < const T * > stepParameters;
237  // prepare the extrapolation Code
238  Trk::ExtrapolationCode eCodeSteps(Trk::ExtrapolationCode::Unset);
239  // path
240  double pathLength = 0.;
241  double materialInX0 = 0.;
242  // we step through the surfaces and go from one to the next
243  ATH_MSG_DEBUG("===> stepwise extrapolation - collecting information <<===");
244  for (auto& sf : stepSurfaces) {
245  ATH_MSG_VERBOSE("===> step <===");
246  // do the extrapolation - will create a memory leak, but
247  eCodeSteps = m_extrapolationEngine->extrapolate(eccSteps, sf);
248  if (!eCodeSteps.isSuccess()) {
249  m_endStepSuccessful = int(eCodeSteps.code);
250  break;
251  }
252  // accumulate the path length
253  pathLength += eccSteps.pathLength;
254  // loop over the collected information for material collection
255  for (auto& es : eccSteps.extrapolationSteps) {
256  if (es.material) materialInX0 += es.material->thicknessInX0() * es.materialScaling;
257  }
258  //
259  eccSteps.restartAtDestination();
260  // for the memory cleanup
261  stepParameters.push_back(eccSteps.endParameters);
262  }
263  // the final extrapolation
264  if (eCodeSteps.isSuccess()) {
265  // the final one is still missing
266  ATH_MSG_VERBOSE("===> final step <===");
267  eCodeSteps = m_extrapolationEngine->extrapolate(eccSteps, &feParameters->associatedSurface());
268  // what is the final extrapolation code
269  m_endStepSuccessful = int(eCodeSteps.code);
270  // declare success ?
271  if (eCodeSteps.isSuccess()) {
272  // and now the final step
273  pathLength += eccSteps.pathLength;
274  // and the material
275  // loop over the collected information for material collection
276  for (auto& es : eccSteps.extrapolationSteps) {
277  if (es.material) materialInX0 += es.material->thicknessInX0() * es.materialScaling;
278  }
279  // stepwise end parameters
280  m_endStepPositionX = eccSteps.endParameters->position().x();
281  m_endStepPositionY = eccSteps.endParameters->position().y();
282  m_endStepPositionZ = eccSteps.endParameters->position().z();
283  m_endStepPositionR = eccSteps.endParameters->position().perp();
284  m_endStepPhi = eccSteps.endParameters->momentum().phi();
285  m_endStepEta = eccSteps.endParameters->momentum().eta();
286  m_endStepTheta = eccSteps.endParameters->momentum().theta();
287  m_endStepP = eccSteps.endParameters->momentum().mag();
288  m_endStepPt = eccSteps.endParameters->momentum().perp();
289  m_endStepPathLength = pathLength;
290  m_endStepThicknessInX0 = materialInX0;
291  }
292  }
293  // memory cleanup
294  for (auto& sp : stepParameters)
295  delete sp;
296  }
297 
298  // now clean up the first end parameters
299  delete feParameters;
300  } else {
301  ATH_MSG_WARNING(
302  "Extrapolation was not Successful - code : " << eCode.toString() << " -> Printing start parameters.");
303  ATH_MSG_WARNING(" -> start x = " << m_startPositionX);
304  ATH_MSG_WARNING(" -> start y = " << m_startPositionY);
305  ATH_MSG_WARNING(" -> start z = " << m_startPositionZ);
306  ATH_MSG_WARNING(" -> start r = " << m_startPositionR);
307  ATH_MSG_WARNING(" -> start phi = " << m_startPhi);
308  ATH_MSG_WARNING(" -> start eta = " << m_startEta);
309  ATH_MSG_WARNING(" -> start theta = " << m_startTheta);
310  ATH_MSG_WARNING(" -> start pt = " << m_startPt);
311  ATH_MSG_WARNING(" -> start p = " << m_startP);
312  }
313 
314  // use the dedicated position / momentum writeter
315  if (!m_posmomWriter.empty()) m_posmomWriter->finalizeTrack();
316  m_tree->Fill();
317  }
318  return StatusCode::SUCCESS;
319 }
320 
321 template < class T, class P > StatusCode Trk::ExtrapolationEngineTest::fillStepInformationT(
322  Trk::ExtrapolationCell < T >& ecc, int fwdBwd, std::vector < const Trk::Surface* >& stepSurfaces) {
323 
324  // loop over the collected information
325  for (auto& es : ecc.extrapolationSteps) {
326  // continue if we have parameters
327  const T* parameters = es.parameters;
328  if (parameters) {
329  // collect for stepwise extrapolation
330  if (fwdBwd && m_stepwiseExtrapolation && parameters->associatedSurface().associatedDetectorElement()) {
331  stepSurfaces.push_back(&parameters->associatedSurface());
332  }
333  // there are parameters assigned, so they need to be either
334  // sensitive -> ip = 0
335  // passive -> ip = 1
336  // boundary -> ip = 2
337  unsigned int ip = 0;
338  if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectPassive)) {
339  ip = 1;
340  } else if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectBoundary)) {
341  ip = 2;
342  }
343 
344  // use the dedicated position / momentum writeter
345  if (!m_posmomWriter.empty()) {
346  m_posmomWriter->recordTrackState(parameters->position(), parameters->momentum());
347  }
348  // fill the parameters
349  m_pPositionX[ip]->push_back(parameters->position().x());
350  m_pPositionY[ip]->push_back(parameters->position().y());
351  m_pPositionZ[ip]->push_back(parameters->position().z());
352  m_pPositionR[ip]->push_back(parameters->position().perp());
353  m_pPhi[ip]->push_back(parameters->momentum().phi());
354  m_pTheta[ip]->push_back(parameters->momentum().eta());
355  m_pEta[ip]->push_back(parameters->momentum().theta());
356  m_pP[ip]->push_back(parameters->momentum().mag());
357  m_pPt[ip]->push_back(parameters->momentum().perp());
358  // check this
359  if (ip == 0 && m_collectSensitive) {
360  // fill the layer index and the local variables
361  int idx = parameters->associatedSurface().associatedLayer() ?
362  parameters->associatedSurface().associatedLayer()->layerIndex().value() : -1;
363  m_sensitiveLayerIndex->push_back(idx);
364  int type = int(parameters->associatedSurface().type());
365  m_sensitiveSurfaceType->push_back(type);
366  m_sensitiveLocalPosX->push_back(parameters->localPosition()[Trk::loc1]);
367  m_sensitiveLocalPosY->push_back(parameters->localPosition()[Trk::loc2]);
368 
369  if (parameters->associatedSurface().associatedDetectorElement()) {
370  Identifier id = parameters->associatedSurface().associatedDetectorElement()->identify();
371 
372  m_sensitiveCenterPosX->push_back(parameters->associatedSurface().associatedDetectorElement()->center().x());
373  m_sensitiveCenterPosY->push_back(parameters->associatedSurface().associatedDetectorElement()->center().y());
374  m_sensitiveCenterPosZ->push_back(parameters->associatedSurface().associatedDetectorElement()->center().z());
375  m_sensitiveCenterPosR->push_back(parameters->associatedSurface().associatedDetectorElement()->center().perp());
376  m_sensitiveCenterPosPhi->push_back(parameters->associatedSurface().associatedDetectorElement()->center().phi());
377 
378  bool isHGTD = false;
379  if (m_idHelper->is_pixel(id)) {
380  // save pixel info
381  m_sensitiveDetector->push_back(0);
382  m_sensitiveIsInnermost->push_back(m_pixel_ID->layer_disk(id)==0);
383  m_sensitiveIsNextToInnermost->push_back(m_pixel_ID->layer_disk(id)==1);
384  m_sensitiveBarrelEndcap->push_back(m_pixel_ID->barrel_ec(id));
385  m_sensitiveLayerDisc->push_back(m_pixel_ID->layer_disk(id));
386  m_sensitiveEtaModule->push_back(m_pixel_ID->eta_module(id));
387  m_sensitivePhiModule->push_back(m_pixel_ID->phi_module(id));
388  m_sensitiveSide->push_back(0);
389  } else if (m_idHelper->is_sct(id)) {
390  // save sct info
391  m_sensitiveDetector->push_back(1);
392  m_sensitiveIsInnermost->push_back(0);
393  m_sensitiveIsNextToInnermost->push_back(0);
394  m_sensitiveBarrelEndcap->push_back(m_sct_ID->barrel_ec(id));
395  m_sensitiveLayerDisc->push_back(m_sct_ID->layer_disk(id));
396  m_sensitiveEtaModule->push_back(m_sct_ID->eta_module(id));
397  m_sensitivePhiModule->push_back(m_sct_ID->phi_module(id));
398  m_sensitiveSide->push_back(m_sct_ID->side(id));
399  } else if (m_useHGTD and m_idHelper->is_hgtd(id)) {
400  // save hgtd info
401  m_sensitiveDetector->push_back(2);
402  m_sensitiveIsInnermost->push_back(0);
403  m_sensitiveIsNextToInnermost->push_back(0);
404  m_sensitiveBarrelEndcap->push_back(m_hgtd_ID->endcap(id));
405  m_sensitiveLayerDisc->push_back(m_hgtd_ID->layer(id));
406  m_sensitiveEtaModule->push_back(m_hgtd_ID->eta_module(id));
407  m_sensitivePhiModule->push_back(m_hgtd_ID->phi_module(id));
408  m_sensitiveSide->push_back(0);
409  isHGTD = true;
410  }
411 
412  if (not isHGTD) {
413  const InDetDD::SiDetectorElement* siElement = dynamic_cast<const InDetDD::SiDetectorElement*> (parameters->associatedSurface().associatedDetectorElement());
414  double phitol = 2.5; double etatol = 5.0;
415  bool nearGap = siElement->nearBondGap(parameters->localPosition(), etatol);
416  InDetDD::SiIntersect siIn = siElement->inDetector(parameters->localPosition(), phitol, etatol);
417  bool isInside = siIn.in();
418  m_sensitiveNearBondGap->push_back(nearGap ? 1 : 0);
419  m_sensitiveisInside->push_back(isInside ? 1 : 0);
420  auto& pBounds = parameters->associatedSurface().bounds();
421  bool isInsideBounds = pBounds.inside(parameters->localPosition(),phitol,etatol);
422  m_sensitiveisInsideBound->push_back(isInsideBounds ? 1 : 0);
423  } else {
424  m_sensitiveNearBondGap->push_back(-1);
425  m_sensitiveisInside->push_back(-1);
426  m_sensitiveisInsideBound->push_back(-1);
427  }
428  }
429  m_materialtInX0AccumulatedUpTo->push_back(m_materialThicknessInX0Accumulated->empty() ? 0. : m_materialThicknessInX0Accumulated->back());
430  }
431 
432  // collect the material if configured to do so
433  if (m_collectMaterial && es.material) {
434  // VERBOSE output
435  ATH_MSG_VERBOSE("Crossing material with scaling " << es.materialScaling);
436  ATH_MSG_VERBOSE(" material information " << *es.material);
437  // thickness in X0 and L0
438  float tInX0 = es.material->thicknessInX0() * es.materialScaling;
439  float tInL0 = es.material->thicknessInL0() * es.materialScaling;
440  float tzarho = es.material->zOverAtimesRhoTimesD() * es.materialScaling;
441  // that's the path
442  double path = es.materialScaling * es.material->thickness();
443  double sigma = 0.;
444  double kazl = 0.;
445  float mean = Trk::MaterialInteraction::dE_MPV_ionization(parameters->momentum().mag(),
446  (es.material->material()),
447  Trk::muon,
448  sigma,
449  kazl,
450  path);
451  // sample the energy loss - m_tRandom.Landau(fabs(mpv),sigma);
452  // fabs(energyLoss)+energyLossSigma*CLHEP::RandLandau::shoot(m_randomEngine);
453  double eloss = mean - sigma * Trk::TrkExUnitTestBase::m_landauDist->shoot();
454 
455  // the accummulated material
456  if (fwdBwd > 0) {
457  m_materialThicknessInX0 += tInX0;
458  m_materialThicknessInL0 += tInL0;
459  m_materialThicknessZARho += tzarho;
460  m_materialEmulatedIonizationLoss += eloss;
461  } else {
462  m_materialThicknessInX0Bwd += tInX0;
463  m_materialThicknessInL0Bwd += tInL0;
464  }
465  // the stepwise material
466  m_materialThicknessInX0Accumulated->push_back(m_materialThicknessInX0);
467  m_materialThicknessInX0Steps->push_back(tInX0);
468  m_materialThicknessInX0Steps->push_back(tInL0);
469  m_materialPositionX->push_back(es.materialPosition.x());
470  m_materialPositionY->push_back(es.materialPosition.y());
471  m_materialPositionZ->push_back(es.materialPosition.z());
472  m_materialPositionR->push_back(es.materialPosition.perp());
473  m_materialPositionP->push_back(parameters->momentum().mag());
474  m_materialPositionPt->push_back(parameters->momentum().perp());
475  m_materialScaling->push_back(es.materialScaling);
476  m_stepDirection->push_back(fwdBwd);
477  // check what type of material you have
478  if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectSensitive)) {
479  m_materialThicknessInX0Sensitive += tInX0;
480  } else if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectPassive)) {
481  m_materialThicknessInX0Passive += tInX0;
482  } else if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectBoundary)) {
483  m_materialThicknessInX0Boundary += tInX0;
484  }
485  // check what type of surface you have
486  if (es.surface && es.surface->type() == Trk::SurfaceType::Cylinder) {
487  m_materialThicknessInX0Cylinder += tInX0;
488  } else if (es.surface && es.surface->type() == Trk::SurfaceType::Disc) {
489  m_materialThicknessInX0Disc += tInX0;
490  } else if (es.surface && es.surface->type() == Trk::SurfaceType::Plane) {
491  m_materialThicknessInX0Plane += tInX0;
492  }
493  }
494  // delete the parameters if there are any there
495  delete parameters;
496  }
497  }
498 
499  return StatusCode::SUCCESS;
500 }