Spurious peaks in double differential cross section plot for proton beam on Li7 target

Dear experts,

In this simulation, I have considered a 9 MeV proton beam on thin Li7 target. The interaction will yield 7Be. Depending on the states of 7Be, I got different peaks in the double differential cross section spectrum. Their positions are matching nearly when compared with LISE tool. However, I have following queries regarding this:

  1. The literature suggests the peak to be Gaussian in nature, however the simulation shows a flat top distribution. Can you please guide regarding this observation ?

  2. I got few other peaks which I did not expect. I have marked them by yellow circle. (Please refer to the attached image)

  3. I have marked few regions as blue i.e. the shoulder like data after the second excited peak, which should not exist.

p_onLi7_9MeV.inp (3.9 KB)



Hello Riya,

Sorry for the slight delay.

Not exclusively. There are other channels producing neutrons out of 7Li under p irradiation. Keep reading, it gets better.

Jupp, these are due to direct reaction mechanisms that leave the residual 7Be in the ground state, 1st excited state, and 2nd excited state, as you nicely label them.

I guess what is implied is the uncertainty due to an experimental resolution - if you need it, apply a Gaussian convolution offline / during post-processing.

What you see (especially for the gs and 1st excited state peaks) are not quite flat tops. You are looking at neutrons collected within in a certain “acceptance cone” of plus/minus a few degrees around 25 and 150 deg. This probes a small domain of the angular distribution of the emitted neutron, which is not isotropic (it appears as a flat neutron energy spectrum - but if you look closely is not - simply because you’re effectively looking at a small angular domain). Actually quite a bit of work went into getting a reasonable intensity and angular dependency of these first two peaks due to direct nuclear reactions.

These are due to a process whereby the p is absorbed by the 7Li to form a compound 8Be with an excitation energy of about 26.2 MeV, which then splits straight into n+3He+4He, i.e.

p + 7Li -> 8Be* -> n + 3He + 4He

This expectation is unsubstantiated. In fact, these spectral features you point out (actually a continuum below the 2nd excited level peak) is due to a similar process, where the excited compound 8Be first splits into two alpha particles, one of them in its 1st excited level (20 MeV), which subsequently breaks into n+3He. That is:

p + 7Li -> 8Be* -> 4He* + 4He -> n + 3He + 4He

Hope this helps,



Thank you @cesc for the explanation and it was very helpful. Since I am very new to this field, pardon my ignorance. I got few more queries regarding this topic.

At first I added one RESNUCLEI card with Semi_analogue mode in RADDECAY. In the tab.lis, I found 8Be is present. Next I did another simulation described below:

  1. Like you mentioned that

Is there any way to visualize/print this reaction using FLUKA? I mean such exercise will help to understand the process and learn how to analyze the output in future studies. In order to achieve this, I have added mgdraw.f and usrrnc.f to print some information on residual nuclei.


Does the above block (above image) refer to the reaction that you mentioned?

  • If yes, then how can I print 8Be* in the same block?

  • Apart from this, if we want to print energy corresponding to that Particle ID as well, what variable should I use?
    Likewise, what variable will denote the kinetic energy of the residual nuclei? (so that I can plug that variable in this line:
    WRITE(96, ‘(A, 2(a, I3))’) 'Generated a residual ',
    & 'with A = ', IBRES, ’ , Z = ', ICRES

I am attaching the inputs here for your reference.

mgdraw.f (12.3 KB)
p_onLi7_9MeV.inp (2.0 KB)
usrrnc.f (1.5 KB)

p_onLi7_9MeV001_ISOTOPE.txt (9.2 KB)

  1. Since the incident proton energy is 9 MeV, is it possible to to acquire this excitation energy of 26.2 MeV?
  1. As you mentioned the following statement, is there any option in FLUKA to add Gaussian profile on top of the output data? I know it can be done via offline using different plotting tools, I am just checking if FLUKA has any option for that.
  1. Is there any way to activate only one channel? I tried the following method, but it did not work.



The proton absorption into 7Li frees its binding energy of 17.25 MeV, to be added to the large fraction of incident energy left available after subtracting the 8Be kinetic energy imposed by momentum conservation.

Here RADDECAY is out of context, since we are actually talking about prompt products (the 8Be decay is extremely fast).

Yes. No need for usrrnc, the mgdraw output is enough.

You cannot: it’s a transient state preceding the reaction final state.

TKHEAV (IP) for the fragments and TKI (IP) for the neutron.
Technically speaking, there is no single residual nucleus, as you can read above from mgdraw, rather two nuclear fragments.

There is a dedicated card applying to DETECT scoring.


1 Like

Dear @ceruttif ,

Thank you for the explanation. I modified the mgdraw.f to print the energies.

  1. For particles, it is printing the values. For e.g. in the attached image, neutron energy = 1.7 MeV. However, for the fragments, it showed 0. Can you please guide if there is any mistake?

mgdraw.f (12.4 KB)


  1. In this output, there are 3 types of generated particle/element. One is particle, one is fragment and the third one is residual. Considering this in the output we can see:
  • (a) For primary history 30, Generated particle = 0, but, there is 1 Fragment with A =2, Z =1 i.e. Deuteron and one residual with A=6, Z=3 i.e. 6Li.

  • (b) For primary history 34, Generated particle ID = -3 i.e. Deuteron, but fragment is 0 and residue is A=6, Z=3 i.e. 6Li.

So, we observed that in these two cases, sometimes Deuteron is registered as generated particle and sometimes as fragments. What is the physical meaning of this?

  1. From this mgdraw output, if we want to identify which one refers direct reaction and which one is not ?

for e.g.

if in the output, there is particle ID printed with 0 fragments, but with a residual, then it will denote a direct reaction. Please correct me if my understanding is wrong. I have added one such example below:


  1. There is also one case when particle Id is printed, fragments have been generated, but there is no residual. This will infer complete break up of the target nucleus. Is it so ? for e.g.



  1. The fragment kinetic energy is TKHEAV(I) ( = kinetic energy of the ith fragment), and not TKHEAV(KHEAVY(I)) (since KHEAVY(I) is the species of the ith fragment and the kinetic energy is a property of the fragment and not of its species, which instead the mass number IBHEAV and the atomic number ICHEAV refer to).

  2. The deuteron in the final state appears as a particle (b) if it is emitted in the early stages of the reaction, before reaching the statistical equilibrium of the compound nucleus. It appears as a fragment (a) if it is produced by the Fermi break-up of the latter (see also 4, where the Fermi break-up of 8Be generated instead 3He, 4He as well as one neutron, which never appears as a fragment).

  3. Strictly speaking, the direct reaction mechanism implemented in FLUKA for p+7Li implies the neutron emission and the 7Be production in its ground state or first excited level, with the corresponding photon emission in the second case, as in your history #68. The neutron spectrum peak corresponding to the second excited level of 7Be comes from a different reaction mechanism via the 8Be transient state.

  4. Yes, this is the result of the Fermi break-up of 8Be discussed in 2 above ( 8Be* -> n + 3He + 4He). Note that your history #30 is conceptually analogous (8Be* -> 2H + 6Li), with 6Li appearing as residue.

1 Like

Thank you @ceruttif for the in-depth explanation.