Dear Thomas,
Thank you for your questions.
All issues ended up connected, and stemming from a deeper fact (mentioned as a side note in Stupre.f and stuprf.f routines / FLUKA stacks - #3 by ghugo):
STUPRE/ STUPRF versus USDRAW calls ordering:
- “pure” electromagnetic interactions: for γ projectiles, STUPRE is called before USDRAW, while for e+/e- projectiles, STUPRE is called after USDRAW.
(- photonuclear and electronuclear interactions: STUPRF is called after USDRAW.)
Even more confused? Let’s go.
MDSTACK
/STUPRE
versusMGDRAW
results discrepancies
My issue in MDSTCK and STUPRE is that the particles and energies of the secondaries are not consistent with those provided in MGDRAW. Also, I still do not understand why the energies are not the same as those for MGDRAW.
In MDSTACK
/ STUPRE
, you are printing secondaries information from GENSTK
. You should use EMFSTK
instead.
Using the following code (relying on EMFSTK
) in MDSTACK
/ STUPRE
(for electron or photon interactions indifferently) fixes your secondaries printouts:
do IT=NPSTRT, NPEMF
WRITE(99,*) "MDSTCK IT =", IT,
& ", ETEMF(IT) = ",ETEMF(IT),
& ", PMEMF(IT) = ", PMEMF(IT),
& ", UEMF(IT) = ", UEMF(IT),
& ", VEMF(IT) = ", VEMF(IT),
& ", WEMF(IT) = ", WEMF(IT),
& ", WTEMF(IT) = ", WTEMF(IT)
end do
* Off-topic: note the use of particle direction, better suited than energy
* to uniquely identify a particle (e.g. Rayleigh scattering).
Why? What is wrong with GENSTK
?
GENSTK
is always populated correctly at the time MGDRAW
is called, to be able to study any interaction.
** “pure” EM interactions with electron primary: MDSTACK
/ STUPRE
are called after MGDRAW
. Hence, when they are called, the secondaries are already properly loaded into GENSTK
.
This is why your print of secondaries information for Moller or Bremstrahlung interactions from MDSTACK
and STUPRE
is correct, and matches your print from MGDRAW
.
** “pure” EM interactions with photon primary: MDSTACK
/ STUPRE
are called before MGDRAW
. At that time, the secondaries are not already loaded into GENSTK
.
This is why, for example for Rayleigh scattering or photoelectric interactions, your print of secondaries information from MDSTACK
and STUPRE
is not correct, while your print from MGDRAW
is ok. The secondaries you are printing from MDSTACK
/ STUPRE
with KPART(IP),TKI(IP),WEI(IP)
(stored in GENSTK
) in your message, are actually the secondaries from the previous electron interaction.
- Let’s now discuss what you would actually like to do for your study. Ie, being able to know, within
USIMBS
, whether a secondary is same type as its parent.
In MGDRAW, LLOUSE variable (from TRACKR) seems to not be recognized
The use of TRACKR
variables (for instance LLOUSE
) is licit within MGDRAW.
However, you are using LLOUSE
as an array (LLOUSE(IP) = TRUE
), while LLOUSE
is a scalar.
This points to the next issue: you are trying to store information per secondary, while here we are still at the primary hierarchy level.
Which variables should I be using to flag the secondaries in MGDRAW/USRDRAW?
Instead, you can simply store the parent particle ID, for example by adding:
ISPUSR(1) = JTRACK
This information will be automatically propagated to each secondary before they are tracked.
At tracking time (e.g. in USIMBS
), you will then be able to simply compare the tracked particle type (which will be in JTRACK
) with its parent type (ISPUSR(1)
). Hence, you will be able to know whether each secondary is same particle type as its parent or not.
- Wait, are we not done yet? Why an extra point?
This relates to how the information is automatically propagated with ISPUSR(1)
from the primary hierarchy level to all of its secondaries. This propagation happens to be done in STUPRE
(for “pure” EM interactions, otherwise it is in STUPRF
):
IESPAK (1,NPNW) = ISPUSR (1)
Now, hold on, this is where we close the circle.
For photon “pure” EM interactions, STUPRE
is called before MGDRAW
.
Hence, if you only set ISPUSR(1) = JTRACK
within MGDRAW
, you will set it too late (the propagation will already have taken place in STUPRE
). You will thus not be able to access the correct parent for any photon “pure” EM interaction.
Instead, all the other interactions will be OK, because for them, STUPRE
or STUPRF
is called after MGDRAW
(hence the propagation happens after ISPUSR(1)
is set).
To solve this, you thus need to also add ISPUSR(1) = JTRACK
at the top of STUPRE
.
The good news is that eventually, these are just 2 lines of code to add
Hope this helps,
Cheers,
Gabrielle