207 const G4VProcess *
process = aStep->GetPostStepPoint()->GetProcessDefinedStep();
210 G4StepPoint *preStep=aStep->GetPreStepPoint();
211 G4StepPoint *postStep=aStep->GetPostStepPoint();
213 G4ThreeVector
mom = preStep->GetMomentum();
214 const G4ThreeVector&
pos = preStep->GetPosition();
217 G4Track *
track = aStep->GetTrack();
219 int trackID=
track->GetTrackID();
220 const bool trackIsAlive = (
track->GetTrackStatus() == fAlive ||
track->GetTrackStatus() == fStopButAlive);
224 m_pdg =
track->GetDefinition()->GetPDGEncoding();
225 const G4VProcess* creation =
track->GetCreatorProcess();
226 m_scIn = creation? creation->GetProcessSubType() : -1;
230 HepMC::GenVertexPtr vtx = currentGenParticle ? currentGenParticle->production_vertex() :
nullptr;
231 m_gen = currentGenParticle? 0 : -1;
233 if (currentGenParticle) {
234 while (currentGenParticle && vtx ) {
235 int pdgID=currentGenParticle->pdg_id();
237 const HepMC::GenParticlePtr genmom = vtx->particles_in().size()>0 ? vtx->particles_in().front() :
nullptr;
238 if ( genmom && pdgID!=genmom->pdg_id() )
m_gen++;
239 else if (vtx->particles_out().size()>0 && currentGenParticle!=vtx->particles_out().front())
m_gen++;
242 HepMC::GenParticlePtr genmom = vtx->particles_in_size()>0 ? *(vtx->particles_in_const_begin()) :
nullptr;
243 if ( genmom && pdgID!=genmom->pdg_id() )
m_gen++;
244 else if (vtx->particles_out_size()>0 && currentGenParticle!=*(vtx->particles_out_const_begin()))
m_gen++;
246 vtx = genmom ? genmom->production_vertex() :
nullptr;
247 currentGenParticle = genmom;
251 int parentID=
track->GetParentID();
253 if ( genIt !=
m_trackGenMap.end())
m_gen = (genIt->second >= 0) ? genIt->second+1 : genIt->second-1;
272 double stepLength = aStep->GetStepLength();
273 double radLengthInX0 = preStep->GetMaterial()->GetRadlen();
274 double l0 = preStep->GetMaterial()->GetNuclearInterLength();
275 float stepInX0 = stepLength/radLengthInX0;
277 if (stepInX0>1.
e-06) {
286 double totNA = preStep->GetMaterial()->GetTotNbOfAtomsPerVolume();
287 const G4ElementVector* eVec = preStep->GetMaterial()->GetElementVector();
288 const G4double* atVector = preStep->GetMaterial() ->GetVecNbOfAtomsPerVolume();
289 float mFactor =stepInX0* preStep->GetMaterial()->GetDensity();
291 float zOverA = 0.;
float frSum = 0.;
292 for (
unsigned int i=0;
i<eVec->size();
i++) {
293 float fEl = atVector ? atVector[
i]/totNA : 0.;
294 m_wZ += stepInX0*fEl*((*eVec)[
i]->GetZ());
297 zOverA += fEl*((*eVec)[
i]->GetZ())/((*eVec)[
i]->GetA());
300 if (fabs(frSum-1.)>0.01)
ATH_MSG_DEBUG(
"G4 material description inconsistent, sum of element fractions:"<< frSum);
310 float eloss = postStep->GetMomentum().mag()-preStep->GetMomentum().mag();
318 unsigned int nSec = truth.numberOfChildren();
319 if (nSec>0 || !trackIsAlive ) {
324 G4ThreeVector
mom = preStep->GetMomentum();
329 m_vtx_phi = postStep->GetPosition().phi();
331 int iPrimSurv = !trackIsAlive ? 0 : 1;
334 G4ThreeVector pbal(
mom);
340 pbal -= postStep->GetMomentum();
344 for (
unsigned int isec=0; isec< nSec; isec++) {
345 G4ThreeVector secMom = truth.childP(isec);
347 m_pdg_child[isec+iPrimSurv] = truth.childPdgCode(isec);
366 G4VPhysicalVolume *preVol=preStep->GetPhysicalVolume();
367 G4VPhysicalVolume *postVol=postStep->GetPhysicalVolume();
375 m_pth =
track->GetVertexMomentumDirection().theta();
376 m_pph =
track->GetVertexMomentumDirection().phi();
377 m_p =
track->GetVertexKineticEnergy();
384 m_thEnd=postStep->GetPosition().theta();
385 m_phEnd=postStep->GetPosition().phi();
386 m_dEnd=postStep->GetPosition().mag();
401 if ( !trackIsAlive ) {
406 m_pth =
track->GetVertexMomentumDirection().theta();
407 m_pph =
track->GetVertexMomentumDirection().phi();
408 m_p =
track->GetVertexKineticEnergy();
415 m_thEnd=postStep->GetPosition().theta();
416 m_phEnd=postStep->GetPosition().phi();
417 m_dEnd=postStep->GetPosition().mag();
429 if ( preVol==postVol )
return;
432 const G4ThreeVector &postPos = postStep->GetPosition();
433 const G4ThreeVector &postMom = postStep->GetMomentum();
442 if ( nextGeoID == geoID) {
445 ATH_MSG_DEBUG(
"track moves from "<<geoID<<
" to "<<nextGeoID);
449 if (aStep->GetStepLength() == 0.) {
458 const G4ThreeVector& g4pos =
track->GetPosition();
460 const HepGeom::Point3D<double> position(g4pos.x(),g4pos.y(),g4pos.z());
462 G4ThreeVector g4mom =
track->GetMomentum();
463 const HepGeom::Vector3D<double>
momentum(g4mom.x(),g4mom.y(),g4mom.z());
466 if (
track->GetTrackStatus()==fStopAndKill) {
477 m_pth =
track->GetVertexMomentumDirection().theta();
478 m_pph =
track->GetVertexMomentumDirection().phi();
479 m_p =
track->GetVertexKineticEnergy();
486 m_thEnd=postStep->GetPosition().theta();
487 m_phEnd=postStep->GetPosition().phi();
488 m_dEnd=postStep->GetPosition().mag();