ATLAS Offline Software
ExtrapolationEngineTest.icc
Go to the documentation of this file.
1 /*
2  Copyright (C) 2002-2017 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((Trk::ParticleHypothesis) m_particleHypothesis, 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((Trk::ParticleHypothesis) m_particleHypothesis);
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((Trk::ParticleHypothesis) m_particleHypothesis);
234  eccSteps.addConfigurationMode(Trk::ExtrapolationMode::StopAtBoundary);
235  // if (m_collectSensitive) eccSteps.addConfigurationMode(Trk::ExtrapolationMode::CollectSensitive);
236  // if (m_collectPassive) eccSteps.addConfigurationMode(Trk::ExtrapolationMode::CollectPassive);
237  // if (m_collectBoundary) eccSteps.addConfigurationMode(Trk::ExtrapolationMode::CollectBoundary);
238  if (m_collectMaterial) eccSteps.addConfigurationMode(Trk::ExtrapolationMode::CollectMaterial);
239  std::vector < const T * > stepParameters;
240  // prepare the extrapolation Code
241  Trk::ExtrapolationCode eCodeSteps(Trk::ExtrapolationCode::Unset);
242  // path
243  double pathLength = 0.;
244  double materialInX0 = 0.;
245  // we step through the surfaces and go from one to the next
246  ATH_MSG_DEBUG("===> stepwise extrapolation - collecting information <<===");
247  for (auto& sf : stepSurfaces) {
248  ATH_MSG_VERBOSE("===> step <===");
249  // do the extrapolation - will create a memory leak, but
250  eCodeSteps = m_extrapolationEngine->extrapolate(eccSteps, sf);
251  if (!eCodeSteps.isSuccess()) {
252  m_endStepSuccessful = int(eCodeSteps.code);
253  break;
254  }
255  // accumulate the path length
256  pathLength += eccSteps.pathLength;
257  // loop over the collected information for material collection
258  for (auto& es : eccSteps.extrapolationSteps) {
259  if (es.material) materialInX0 += es.material->thicknessInX0() * es.materialScaling;
260  }
261  //
262  eccSteps.restartAtDestination();
263  // for the memory cleanup
264  stepParameters.push_back(eccSteps.endParameters);
265  }
266  // the final extrapolation
267  if (eCodeSteps.isSuccess()) {
268  // the final one is still missing
269  ATH_MSG_VERBOSE("===> final step <===");
270  eCodeSteps = m_extrapolationEngine->extrapolate(eccSteps, &feParameters->associatedSurface());
271  // what is the final extrapolation code
272  m_endStepSuccessful = int(eCodeSteps.code);
273  // declare success ?
274  if (eCodeSteps.isSuccess()) {
275  // and now the final step
276  pathLength += eccSteps.pathLength;
277  // and the material
278  // loop over the collected information for material collection
279  for (auto& es : eccSteps.extrapolationSteps) {
280  if (es.material) materialInX0 += es.material->thicknessInX0() * es.materialScaling;
281  }
282  // stepwise end parameters
283  m_endStepPositionX = eccSteps.endParameters->position().x();
284  m_endStepPositionY = eccSteps.endParameters->position().y();
285  m_endStepPositionZ = eccSteps.endParameters->position().z();
286  m_endStepPositionR = eccSteps.endParameters->position().perp();
287  m_endStepPhi = eccSteps.endParameters->momentum().phi();
288  m_endStepEta = eccSteps.endParameters->momentum().eta();
289  m_endStepTheta = eccSteps.endParameters->momentum().theta();
290  m_endStepP = eccSteps.endParameters->momentum().mag();
291  m_endStepPt = eccSteps.endParameters->momentum().perp();
292  m_endStepPathLength = pathLength;
293  m_endStepThicknessInX0 = materialInX0;
294  }
295  }
296  // memory cleanup
297  for (auto& sp : stepParameters)
298  delete sp;
299  }
300 
301  // now clean up the first end parameters
302  delete feParameters;
303  } else {
304  ATH_MSG_WARNING(
305  "Extrapolation was not Successful - code : " << eCode.toString() << " -> Printing start parameters.");
306  ATH_MSG_WARNING(" -> start x = " << m_startPositionX);
307  ATH_MSG_WARNING(" -> start y = " << m_startPositionY);
308  ATH_MSG_WARNING(" -> start z = " << m_startPositionZ);
309  ATH_MSG_WARNING(" -> start r = " << m_startPositionR);
310  ATH_MSG_WARNING(" -> start phi = " << m_startPhi);
311  ATH_MSG_WARNING(" -> start eta = " << m_startEta);
312  ATH_MSG_WARNING(" -> start theta = " << m_startTheta);
313  ATH_MSG_WARNING(" -> start pt = " << m_startPt);
314  ATH_MSG_WARNING(" -> start p = " << m_startP);
315  }
316 
317  // use the dedicated position / momentum writeter
318  if (!m_posmomWriter.empty()) m_posmomWriter->finalizeTrack();
319  m_tree->Fill();
320  }
321  return StatusCode::SUCCESS;
322 }
323 
324 template < class T, class P > StatusCode Trk::ExtrapolationEngineTest::fillStepInformationT(
325  Trk::ExtrapolationCell < T >& ecc, int fwdBwd, std::vector < const Trk::Surface* >& stepSurfaces) {
326 
327  // loop over the collected information
328  for (auto& es : ecc.extrapolationSteps) {
329  // continue if we have parameters
330  const T* parameters = es.parameters;
331  if (parameters) {
332  // collect for stepwise extrapolation
333  if (fwdBwd && m_stepwiseExtrapolation && parameters->associatedSurface().associatedDetectorElement()) {
334  stepSurfaces.push_back(&parameters->associatedSurface());
335  }
336  // there are parameters assigned, so they need to be either
337  // sensitive -> ip = 0
338  // passive -> ip = 1
339  // boundary -> ip = 2
340  unsigned int ip = 0;
341  if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectPassive)) {
342  ip = 1;
343  } else if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectBoundary)) {
344  ip = 2;
345  }
346 
347  // use the dedicated position / momentum writeter
348  if (!m_posmomWriter.empty()) {
349  m_posmomWriter->recordTrackState(parameters->position(), parameters->momentum());
350  }
351  // fill the parameters
352  m_pPositionX[ip]->push_back(parameters->position().x());
353  m_pPositionY[ip]->push_back(parameters->position().y());
354  m_pPositionZ[ip]->push_back(parameters->position().z());
355  m_pPositionR[ip]->push_back(parameters->position().perp());
356  m_pPhi[ip]->push_back(parameters->momentum().phi());
357  m_pTheta[ip]->push_back(parameters->momentum().eta());
358  m_pEta[ip]->push_back(parameters->momentum().theta());
359  m_pP[ip]->push_back(parameters->momentum().mag());
360  m_pPt[ip]->push_back(parameters->momentum().perp());
361  // check this
362  if (ip == 0 && m_collectSensitive) {
363  // fill the layer index and the local variables
364  int idx = parameters->associatedSurface().associatedLayer() ?
365  parameters->associatedSurface().associatedLayer()->layerIndex().value() : -1;
366  m_sensitiveLayerIndex->push_back(idx);
367  int type = int(parameters->associatedSurface().type());
368  m_sensitiveSurfaceType->push_back(type);
369  m_sensitiveLocalPosX->push_back(parameters->localPosition()[Trk::loc1]);
370  m_sensitiveLocalPosY->push_back(parameters->localPosition()[Trk::loc2]);
371 
372  if (parameters->associatedSurface().associatedDetectorElement()) {
373  Identifier id = parameters->associatedSurface().associatedDetectorElement()->identify();
374 
375  m_sensitiveCenterPosX->push_back(parameters->associatedSurface().associatedDetectorElement()->center().x());
376  m_sensitiveCenterPosY->push_back(parameters->associatedSurface().associatedDetectorElement()->center().y());
377  m_sensitiveCenterPosZ->push_back(parameters->associatedSurface().associatedDetectorElement()->center().z());
378  m_sensitiveCenterPosR->push_back(parameters->associatedSurface().associatedDetectorElement()->center().perp());
379  m_sensitiveCenterPosPhi->push_back(parameters->associatedSurface().associatedDetectorElement()->center().phi());
380 
381  bool isHGTD = false;
382  if (m_idHelper->is_pixel(id)) {
383  // save pixel info
384  m_sensitiveDetector->push_back(0);
385  m_sensitiveIsInnermost->push_back(m_pixel_ID->layer_disk(id)==0);
386  m_sensitiveIsNextToInnermost->push_back(m_pixel_ID->layer_disk(id)==1);
387  m_sensitiveBarrelEndcap->push_back(m_pixel_ID->barrel_ec(id));
388  m_sensitiveLayerDisc->push_back(m_pixel_ID->layer_disk(id));
389  m_sensitiveEtaModule->push_back(m_pixel_ID->eta_module(id));
390  m_sensitivePhiModule->push_back(m_pixel_ID->phi_module(id));
391  m_sensitiveSide->push_back(0);
392  } else if (m_idHelper->is_sct(id)) {
393  // save sct info
394  m_sensitiveDetector->push_back(1);
395  m_sensitiveIsInnermost->push_back(0);
396  m_sensitiveIsNextToInnermost->push_back(0);
397  m_sensitiveBarrelEndcap->push_back(m_sct_ID->barrel_ec(id));
398  m_sensitiveLayerDisc->push_back(m_sct_ID->layer_disk(id));
399  m_sensitiveEtaModule->push_back(m_sct_ID->eta_module(id));
400  m_sensitivePhiModule->push_back(m_sct_ID->phi_module(id));
401  m_sensitiveSide->push_back(m_sct_ID->side(id));
402  } else if (m_useHGTD and m_idHelper->is_hgtd(id)) {
403  // save hgtd info
404  m_sensitiveDetector->push_back(2);
405  m_sensitiveIsInnermost->push_back(0);
406  m_sensitiveIsNextToInnermost->push_back(0);
407  m_sensitiveBarrelEndcap->push_back(m_hgtd_ID->endcap(id));
408  m_sensitiveLayerDisc->push_back(m_hgtd_ID->layer(id));
409  m_sensitiveEtaModule->push_back(m_hgtd_ID->eta_module(id));
410  m_sensitivePhiModule->push_back(m_hgtd_ID->phi_module(id));
411  m_sensitiveSide->push_back(0);
412  isHGTD = true;
413  }
414 
415  if (not isHGTD) {
416  const InDetDD::SiDetectorElement* siElement = dynamic_cast<const InDetDD::SiDetectorElement*> (parameters->associatedSurface().associatedDetectorElement());
417  double phitol = 2.5; double etatol = 5.0;
418  bool nearGap = siElement->nearBondGap(parameters->localPosition(), etatol);
419  InDetDD::SiIntersect siIn = siElement->inDetector(parameters->localPosition(), phitol, etatol);
420  bool isInside = siIn.in();
421  m_sensitiveNearBondGap->push_back(nearGap ? 1 : 0);
422  m_sensitiveisInside->push_back(isInside ? 1 : 0);
423  auto& pBounds = parameters->associatedSurface().bounds();
424  bool isInsideBounds = pBounds.inside(parameters->localPosition(),phitol,etatol);
425  m_sensitiveisInsideBound->push_back(isInsideBounds ? 1 : 0);
426  } else {
427  m_sensitiveNearBondGap->push_back(-1);
428  m_sensitiveisInside->push_back(-1);
429  m_sensitiveisInsideBound->push_back(-1);
430  }
431  }
432  m_materialtInX0AccumulatedUpTo->push_back(m_materialThicknessInX0Accumulated->empty() ? 0. : m_materialThicknessInX0Accumulated->back());
433  }
434 
435  // collect the material if configured to do so
436  if (m_collectMaterial && es.material) {
437  // VERBOSE output
438  ATH_MSG_VERBOSE("Crossing material with scaling " << es.materialScaling);
439  ATH_MSG_VERBOSE(" material information " << *es.material);
440  // thickness in X0 and L0
441  float tInX0 = es.material->thicknessInX0() * es.materialScaling;
442  float tInL0 = es.material->thicknessInL0() * es.materialScaling;
443  float tzarho = es.material->zOverAtimesRhoTimesD() * es.materialScaling;
444  // that's the path
445  double path = es.materialScaling * es.material->thickness();
446  double sigma = 0.;
447  double kazl = 0.;
448  float mean = Trk::MaterialInteraction::dE_MPV_ionization(parameters->momentum().mag(),
449  (es.material->material()),
450  Trk::muon,
451  sigma,
452  kazl,
453  path);
454  // sample the energy loss - m_tRandom.Landau(fabs(mpv),sigma);
455  // fabs(energyLoss)+energyLossSigma*CLHEP::RandLandau::shoot(m_randomEngine);
456  double eloss = mean - sigma * Trk::TrkExUnitTestBase::m_landauDist->shoot();
457 
458  // the accummulated material
459  if (fwdBwd > 0) {
460  m_materialThicknessInX0 += tInX0;
461  m_materialThicknessInL0 += tInL0;
462  m_materialThicknessZARho += tzarho;
463  m_materialEmulatedIonizationLoss += eloss;
464  } else {
465  m_materialThicknessInX0Bwd += tInX0;
466  m_materialThicknessInL0Bwd += tInL0;
467  }
468  // the stepwise material
469  m_materialThicknessInX0Accumulated->push_back(m_materialThicknessInX0);
470  m_materialThicknessInX0Steps->push_back(tInX0);
471  m_materialThicknessInX0Steps->push_back(tInL0);
472  m_materialPositionX->push_back(es.materialPosition.x());
473  m_materialPositionY->push_back(es.materialPosition.y());
474  m_materialPositionZ->push_back(es.materialPosition.z());
475  m_materialPositionR->push_back(es.materialPosition.perp());
476  m_materialPositionP->push_back(parameters->momentum().mag());
477  m_materialPositionPt->push_back(parameters->momentum().perp());
478  m_materialScaling->push_back(es.materialScaling);
479  m_stepDirection->push_back(fwdBwd);
480  // check what type of material you have
481  if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectSensitive)) {
482  m_materialThicknessInX0Sensitive += tInX0;
483  } else if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectPassive)) {
484  m_materialThicknessInX0Passive += tInX0;
485  } else if (es.stepConfiguration.checkMode(Trk::ExtrapolationMode::CollectBoundary)) {
486  m_materialThicknessInX0Boundary += tInX0;
487  }
488  // check what type of surface you have
489  if (es.surface && es.surface->type() == Trk::SurfaceType::Cylinder) {
490  m_materialThicknessInX0Cylinder += tInX0;
491  } else if (es.surface && es.surface->type() == Trk::SurfaceType::Disc) {
492  m_materialThicknessInX0Disc += tInX0;
493  } else if (es.surface && es.surface->type() == Trk::SurfaceType::Plane) {
494  m_materialThicknessInX0Plane += tInX0;
495  }
496  }
497  // delete the parameters if there are any there
498  delete parameters;
499  }
500  }
501 
502  return StatusCode::SUCCESS;
503 }
Generated on Thu Nov 7 2024 21:15:08 for ATLAS Offline Software by   doxygen 1.8.18