Problem 2561 - Version 11.1.2 G4PhantomParameterisation producing nans in sensitive detector when running setregularstructureID(1).
Summary: Version 11.1.2 G4PhantomParameterisation producing nans in sensitive detector...
Status: RESOLVED FIXED
Alias: None
Product: Geant4
Classification: Unclassified
Component: geometry/navigation (show other problems)
Version: 11.1
Hardware: All All
: P4 major
Assignee: Pedro.Arce
URL:
Depends on:
Blocks:
 
Reported: 2023-09-05 08:58 CEST by James
Modified: 2023-11-04 07:25 CET (History)
1 user (show)

See Also:


Attachments
Example code reproducing the nans. (48.21 KB, application/x-zip-compressed)
2023-09-05 08:58 CEST, James
Details
G4THitsMap.hh (25.20 KB, text/plain)
2023-09-22 04:52 CEST, James
Details
G4THitsMap.hh (25.20 KB, text/plain)
2023-09-22 04:59 CEST, James
Details
Geant4 Installation excluding data (27.27 MB, application/x-zip-compressed)
2023-10-07 15:15 CEST, James
Details

Note You need to log in before you can comment on or make changes to this problem.
Description James 2023-09-05 08:58:41 CEST
Created attachment 823 [details]
Example code reproducing the nans.

I have had this question active on the forums for a while and additionally have reached out to group affiliations who work on geant4 so I feel I have done due diligence before submitting this report. 

I have uploaded my code for this issue in the attachments. It is modified from example RE02 to use G4phantomparameterisation instead of nested parameterisation. I did this to try and isolate the cause of my issue which is that nans (not a number)s are produced in the sensitive detector when setregularstructureID(1) is used. When it is not used the simulations run correctly but are horrifically slow for any large voxelised object (order of one event per second).

The intent of my code is to measure dosimetry in low energy x-ray imaging, around 40keV. The nans seem to be produced primarily in low density materials as a result of secondary particle (ionised electron) interactions. Given the issue only occurs when using the navigation with setregularstructureiD(1) it makes me think the problem is there. 
The test case in the is just a plane of air voxels. Cut values do no seem to affect the presence of nans, unless by completely stopping all electron secondaries being produced.
Additionally when I tried rolling back to version 10.6.3 I found that nans were no longer produced but instead I got navigation errors occurring from stuck tracks. Given all these points I feel my last option is some kind of bug in the code, or in the very least some kind of obscurity not known to me or the many I have tried contacting.
I really appreciate this being looked into.
Kind regards,
James
Comment 1 James 2023-09-18 09:49:07 CEST
Not sure how many tickets there are in the queue but I'd be very interested to hear if this problem has been looked at / validated.
Comment 2 Pedro.Arce 2023-09-21 14:31:09 CEST
I am running your code. Where do you produce the nans? I do not see it. Maybe is in your output files? I cannot reproduce that because I do not have your G4THitsMap.hh or H5C.hh files
Comment 3 James 2023-09-22 04:52:12 CEST
Created attachment 825 [details]
G4THitsMap.hh
Comment 4 James 2023-09-22 04:59:45 CEST
Created attachment 826 [details]
G4THitsMap.hh
Comment 5 James 2023-09-22 05:00:07 CEST
Hi Pedro,
Thank you for getting back to me. I just downloaded the zip files and confirmed to myself they reproduce the error immediately on running. The nans are produced during the event action, I have a print statement during EndOfEventAction that will print most of the time when nans occur. They also save into a hdf5 file during the EndOfRunAction into the same directory the executable is in. 

I use the same G4THitsMap.hh file that make in the include directory of the Geant4-11.1 installation. There are no modifications made to it from there but I have attached it as well. 
The hdf5 files I have used are from the hdf5cpp package. They are most easily opened using imagej free software, if for whatever reason you cannot install that package I will see about removing that dependency. 

Note that these nans have occurred for me on both linux and windows installations.
Comment 6 James 2023-10-07 15:15:44 CEST
Created attachment 827 [details]
Geant4 Installation excluding data

I have attached a copy of my geant4 installation excluding the data files. Please update with any progress made in reproducing the problem.
Comment 7 James 2023-10-11 08:59:25 CEST
Hope all is well, let me know if I should resubmit my report.
Comment 8 Pedro.Arce 2023-11-03 17:32:46 CET
Sorry to answer so late. I found your problem. It is just adding a safeguard:



G4int G4EnergySplitter::SplitEnergyInVolumes(const G4Step* aStep )
{
  theEnergies.clear();

  G4double edep = aStep->GetTotalEnergyDeposit();

  if( edep == 0. ) {
    return 0;
  }
Comment 9 James 2023-11-04 01:16:58 CET
Hey Pedro, 
Apologies but I'm not quite across how to implement this, what should this safeguard go and what is it doing?
Comment 10 James 2023-11-04 07:25:11 CET
Disregard the previous comment sorry, I realised the if-return statement should be added to the G4EnergySplitter.cc source file and my geant4 installation recompiled. Would this change also be in future editions of geant4 or would I need to changed it each time I redownload?

I'd still be curious to hear about what was happening here if you have time. I understand that the new if statement just exits out of the SplitEnergyInVolumes function if no energy was deposited in that step but I'm not sure how the rest of that function could create nans otherwise. 

Regardless though thank you very much for looking into it, I really appreciate it.