Skip to content

Commit 1f11a00

Browse files
[mini] fix bug if next state was outside (apt-sim#323)
If a step was recorded at the edge of the world such that the NextState would be "outside" the reconstruction algorithm would try to access a outside volume that doesn't exist in the lookup table. This would lead to a segfault. This PR adds a safeguard against that case. The scenario is unrealistic in a real setup, but should be avoided for stability (and the error was observed in some small tests)
1 parent f797c1c commit 1f11a00

File tree

1 file changed

+21
-20
lines changed

1 file changed

+21
-20
lines changed

src/AdePTGeant4Integration.cpp

+21-20
Original file line numberDiff line numberDiff line change
@@ -81,10 +81,7 @@ void Deleter::operator()(ScoringObjects *ptr)
8181

8282
} // namespace AdePTGeant4Integration_detail
8383

84-
AdePTGeant4Integration::~AdePTGeant4Integration()
85-
{
86-
87-
}
84+
AdePTGeant4Integration::~AdePTGeant4Integration() {}
8885

8986
void AdePTGeant4Integration::CreateVecGeomWorld(std::string filename)
9087
{
@@ -268,7 +265,6 @@ void AdePTGeant4Integration::InitScoringData(adeptint::VolAuxData *volAuxData)
268265
// Though we only record and reconstruct hits for sensitive volumes, this map needs to store every
269266
// volume in the geometry, as a step may begin in a sensitive volume and end in a non-sensitive one
270267
fglobal_vecgeom_to_g4_map.insert(std::pair<int, const G4VPhysicalVolume *>(vg_pvol->id(), g4_pvol));
271-
272268
// Now do the daughters
273269
for (int id = 0; id < g4_lvol->GetNoDaughters(); ++id) {
274270
auto g4pvol_d = g4_lvol->GetDaughter(id);
@@ -292,18 +288,20 @@ void AdePTGeant4Integration::ProcessGPUHit(GPUHit const &hit)
292288
}
293289

294290
// Reconstruct G4NavigationHistory and G4Step, and call the SD code for each hit
295-
vecgeom::NavigationState const & preNavState = hit.fPreStepPoint.fNavigationState;
291+
vecgeom::NavigationState const &preNavState = hit.fPreStepPoint.fNavigationState;
296292
// Reconstruct Pre-Step point G4NavigationHistory
297293
FillG4NavigationHistory(preNavState, &fScoringObjects->fPreG4NavigationHistory);
298294
(*fScoringObjects->fPreG4TouchableHistoryHandle)
299295
->UpdateYourself(fScoringObjects->fPreG4NavigationHistory.GetTopVolume(),
300296
&fScoringObjects->fPreG4NavigationHistory);
301297
// Reconstruct Post-Step point G4NavigationHistory
302-
vecgeom::NavigationState const & postNavState = hit.fPostStepPoint.fNavigationState;
303-
FillG4NavigationHistory(postNavState, &fScoringObjects->fPostG4NavigationHistory);
304-
(*fScoringObjects->fPostG4TouchableHistoryHandle)
305-
->UpdateYourself(fScoringObjects->fPostG4NavigationHistory.GetTopVolume(),
306-
&fScoringObjects->fPostG4NavigationHistory);
298+
vecgeom::NavigationState const &postNavState = hit.fPostStepPoint.fNavigationState;
299+
if (!postNavState.IsOutside()) {
300+
FillG4NavigationHistory(postNavState, &fScoringObjects->fPostG4NavigationHistory);
301+
(*fScoringObjects->fPostG4TouchableHistoryHandle)
302+
->UpdateYourself(fScoringObjects->fPostG4NavigationHistory.GetTopVolume(),
303+
&fScoringObjects->fPostG4NavigationHistory);
304+
}
307305

308306
// Reconstruct G4Step
309307
switch (hit.fParticleType) {
@@ -385,12 +383,15 @@ void AdePTGeant4Integration::FillG4Step(GPUHit const *aGPUHit, G4Step *aG4Step,
385383
G4TouchableHandle &aPreG4TouchableHandle,
386384
G4TouchableHandle &aPostG4TouchableHandle) const
387385
{
388-
const G4ThreeVector aPostStepPointMomentumDirection(aGPUHit->fPostStepPoint.fMomentumDirection.x(), aGPUHit->fPostStepPoint.fMomentumDirection.y(),
389-
aGPUHit->fPostStepPoint.fMomentumDirection.z());
390-
const G4ThreeVector aPostStepPointPolarization(aGPUHit->fPostStepPoint.fPolarization.x(), aGPUHit->fPostStepPoint.fPolarization.y(),
391-
aGPUHit->fPostStepPoint.fPolarization.z());
392-
const G4ThreeVector aPostStepPointPosition(aGPUHit->fPostStepPoint.fPosition.x(), aGPUHit->fPostStepPoint.fPosition.y(),
393-
aGPUHit->fPostStepPoint.fPosition.z());
386+
const G4ThreeVector aPostStepPointMomentumDirection(aGPUHit->fPostStepPoint.fMomentumDirection.x(),
387+
aGPUHit->fPostStepPoint.fMomentumDirection.y(),
388+
aGPUHit->fPostStepPoint.fMomentumDirection.z());
389+
const G4ThreeVector aPostStepPointPolarization(aGPUHit->fPostStepPoint.fPolarization.x(),
390+
aGPUHit->fPostStepPoint.fPolarization.y(),
391+
aGPUHit->fPostStepPoint.fPolarization.z());
392+
const G4ThreeVector aPostStepPointPosition(aGPUHit->fPostStepPoint.fPosition.x(),
393+
aGPUHit->fPostStepPoint.fPosition.y(),
394+
aGPUHit->fPostStepPoint.fPosition.z());
394395

395396
// G4Step
396397
aG4Step->SetStepLength(aGPUHit->fStepLength); // Real data
@@ -404,8 +405,8 @@ void AdePTGeant4Integration::FillG4Step(GPUHit const *aGPUHit, G4Step *aG4Step,
404405

405406
// G4Track
406407
G4Track *aTrack = aG4Step->GetTrack();
407-
aTrack->SetTrackID(aGPUHit->fParentID); // Missing data
408-
aTrack->SetParentID(aGPUHit->fParentID); // ID of the initial particle that entered AdePT
408+
aTrack->SetTrackID(aGPUHit->fParentID); // Missing data
409+
aTrack->SetParentID(aGPUHit->fParentID); // ID of the initial particle that entered AdePT
409410
aTrack->SetPosition(aPostStepPointPosition); // Real data
410411
// aTrack->SetGlobalTime(0); // Missing data
411412
// aTrack->SetLocalTime(0); // Missing data
@@ -468,7 +469,7 @@ void AdePTGeant4Integration::FillG4Step(GPUHit const *aGPUHit, G4Step *aG4Step,
468469
// aPostStepPoint->SetGlobalTime(0); // Missing data
469470
// aPostStepPoint->SetProperTime(0); // Missing data
470471
aPostStepPoint->SetMomentumDirection(aPostStepPointMomentumDirection); // Real data
471-
aPostStepPoint->SetKineticEnergy(aGPUHit->fPostStepPoint.fEKin); // Real data
472+
aPostStepPoint->SetKineticEnergy(aGPUHit->fPostStepPoint.fEKin); // Real data
472473
// aPostStepPoint->SetVelocity(0); // Missing data
473474
if (const auto postVolume = (*fScoringObjects->fPostG4TouchableHistoryHandle)->GetVolume();
474475
postVolume != nullptr) { // protect against nullptr if postNavState is outside

0 commit comments

Comments
 (0)