TestEm7 Bragg Peak Incorrect for High Energy Particles

I am having difficulties using TestEm7 to successfully generate correct LET/Bragg Peak curves. Specifically, I am using the ionC12.mac file which is included with the test script. What I am finding however is that as I increase the energy of the particle, the Bragg Peak shape becomes distorted and then disappears. What I generally expect the Bragg Peak to look like as I increase the energy of the particle is the following:

Below are the plots that Geant generated. In the title of each plot is the energy of the particle (all C12) along with the “setSizeX” value given in the .mac script:

C12_2_4Gev_20cm_Water.pdf (17.1 KB)

C12_4_8Gev_40cm_Water.pdf (16.9 KB)

C12_7_2Gev_60cm_Water.pdf (17.0 KB)

C12_20_4Gev_200cm_Water.pdf (17.0 KB)

As you can see, these plots do not follow the general shape one would except an LET curve to follow. Is there something simple I am missing here? I have not edited the code in any way, and the only values changing with the ionC12.mac script are “/gun/energy 20.4 GeV” and “/testem/det/setSizeX 200 cm”. Any insight/help is greatly appreciated.


Geant4 Version: 11.1.2
Operating System: Mac OS 14.2.1

I was able to solve the problem. For some reason, including the line “/testem/phys/addPhysics binary_ion” caused the issues. Once that line was removed from the .mac file the results looked correct. For instance, I was able to recreate the following plot which is in good agreement (at least to first order) with the verified results from the TAMU cyclotron:

LET_Compared.pdf (15.9 KB)

(2nd page) https://cyclotron.tamu.edu/ref/images/LET%20vs%20Range%20for%20K150%20Heavy%20Ions.pdf

If anyone can provide a bit of understanding for why the inclusion of that line causes issues it would be greatly appreciated.


It is better to work on HadronTherapy example which can be found in advanced examples. It is more advance but helpful. They simulate 62 MeV of protons in a water phantom. You can change the energy and the ion. The code gives accurate results.

Dear djulson,

At very high energies, the probability per unit path length of having inelastic nuclear reactions (i.e. fragmentation) becomes significant. This decreases the fluence of primaries as the depth increases, which in consequence “degrades” the shape of the Bragg Peak and produces a fragmentation tail downstream the Bragg peak.

If you remove “binary_ion”, then these nuclear reactions are not taking place anymore. Then, you end up having just an EM simulation, in a similar way as (e.g.) the code SRIM does.

Hope this helps,

Hi Miguel,

Thank you for the reply, and also the explanation of the “binary_ion” portion. I tried searching for a description of it but found almost nothing, so thank you for the clarification.

I want to mention however that in the “LET_compared” plot I posted previously, which matched well with the TAMU cyclotron data, I am using “binary” but not “binary_ion”. In the code (PhysicsList.cc), “binary” corresponds to the G4HadronInelasticQBBC physics list, while “binary_ion” corresponds to “G4IonPhysics”. I am confused what the difference is between these two physics lists? are the nuclear interactions still included in “binary”? I am also confused why the simulations I made without “binary_ion” match the TAMU cyclotron so well if it does not include all interactions?

Thank you for your help!

You may find useful to take a look at their .cc files to see their differences:

The short answer is: Yes, but it depends on the particle. In this example, “binary” registers hadronic inelastic interactions for a bunch of elementary particles (p, n, pions…), but not for alphas and heavier ions. For the latter, you need “binary_ion”.


Hi Miguel,

Thank you again for the reply. I am a little concerned as to why adding in G4IonPhysics causes the simulation to disagree with the LET curves provided by the cyclotron, if it is instead more accurately modeling the underlying physics? I am looking to model the energy deposited by heavy ions, and am unsure which PhysicsLists to use in my simulations. I do not know how the cyclotron data is actually recorded, and if there are any caveats to consider there that may explain the discrepancy.

Thanks again!

It seems that the problems comes from the fact that in the simulation you calculate Edep per unit depth, which is not exactly LET. The former is affected by the fluence variation as the beam penetrates, and fluctuations.

If you switched off EM fluctuations, scattering, and hadronic interactions, then you could approximately assume that your Edep per unit length is roughly equal to LET. But please note that the LET is mean energy loss by the charged particle per unith path length, a quantity which is not affected by fluence. However, the Edep you score is affected by the fluence.