I’m running simulations in FLUKA to compare the photopeak signal between two cases:
The decay of Cs-137 (which includes the emission of the 662 keV gamma line).
A monochromatic gamma source at 662 keV.
For the comparison, in the monochromatic case, I scaled the result to account for the intensity of the Cs-137 gamma line (approximately 85%). However, I’m seeing a 5% difference in the photopeak, with the signal from the Cs-137 decay being higher than the one from the monochromatic source. This difference does not appear to be covered by statistical uncertainty.
In the Cs-137 decay simulation, I applied a 100 GeV cut on electrons in the DETECT card to ensure that energy deposition from electrons is not included. The spectra produced are almost identical, except for this 5% factor in the photopeak.
My question is: Is there any scaling factor I might be missing, or am I doing something wrong in one of the simulations? I’ve checked that the other simulation parameters are consistent.
I checked and can confirm that at production each disintegration of 137Cs results in 0.85 photons of 0.662 MeV.
Thus the difference you obtain must be attributed to the transport of particles. I noticed your input file defines an extended source, which is currently commented out. Is it possible you turned on this option when simulating the decay source? The points in this extended source are closer to the detector than the point source, and could result in higher detection count rates.
To organize your input and making sure your simulation are exactly the same between the cases, i suggest using pre-processor variables and if/else statement rather than relying on commenting/uncommenting. For example, you could define a preprocessor “decay” and then use if/else to change the type of beam, raddecay, and scoring.
Please let me know if the problem was the extended source, or if the problem is deeper than that.
Dear Stefano,
I am sure that in both simulations the source was point-like. The cylindrical source is for another test that has given an higher discrepancy, in accordance with your hypothesis. But that 5% is from the same source the only difference is disintegration vs monochromatic.
Thank you very much for your answer
Unfortunately, I am not able to reproduce this difference between the two deposition spectra. I have tried several combinations. Could you send me a simplified input where you see the difference, either using preprocessors to switch between decay and monochromatic or two different inputs?
Once I have it, I think I can debug that quite quickly.
Dear Stefano,
I apologize again for the disorganized inputs. I’m attaching the two inputs (without the commented parts) that I used, which gave me the discrepancy I mentioned. As you can see, the number of primary events in the monochromatic emission case is already scaled by 0.85, so I would expect to see the same photopeak efficiency.
no need to apologize, my inputs look even worse when I try to debug and I can’t figure it out.
The good news is that in this form I can clearly see what is happening. Scaling the number of primaries in Fluka will not have an impact on the score itself, since the value given in the output is always expressed per primary. Instead, it will affect the statistical convergence of the score, i.e. the error bars.
What you want to do is run your simulation with a number of primaries that allows you to get good statistics (if anything the Cs137 should be run a little bit longer such that the number of produced photons is approximately the same) and the factor of 0.85 should instead be applied during post processing.
For example, if you use Flair (which looks like you do from your files), follow this procedure:
Go to the Run/Process tab and process your data so that you have combined files with uncertainties from both files
set up USR-1D plots importing the two .tab files you created with the Process
There is an entry labeled “Y” which allows you to specify the scaling factor in Y. This is where you would enter 0.85 for the 662 case, and either leave blank or put 1 for the Cs source. See picture
I would advise to use multiple runs with a decreased count rate, rather than a single one. This will give you the error bars when you process it
there was an issue with running your inputs because of the non-empty what(3) of EMF-CUT. This appears to have been possibly an issue with duplicating an existing EMF-CUT with a different SDUM. Make sure this card is appropriately set on your end.
I hope this solves your issue, otherwise let me know (I am still not entirely sure why Cs137 would be higher in this case, as the 662 source should still be higher, albeit with more statistics. This is also what I observe when I run your inputs)
Hi Stefano,
this version for processing the data is very interesting, thank you! However, I need to extract the value of the photopeak efficiency, meaning the signal only from the bin at 662 keV in both cases.
To do this, I’ve created a Python script that reads the resulting spectrum and multiplies each Y value by the number of primary events. This way, I get the number of photons “seen by the detector,” right? At this point, if this is correct I should be able to compare the two numbers.
Let me explain why I run only 1 cycle:
it allows me to have an unaveraged spectrum, meaning the spectrum I get from that single initial beam.
regarding the uncertainties, I can calculate them using the square root of the counts that is the same way I use for experimental uncertainties
absolutely, in the post processing you multiply each bin by a constant C and as long as C for for the monochromatic case is 0.85* the C for the Cs case, the two spectra will agree insofar as the 662 keV peak (there is of course also some x-rays at very low energies in the case of the isotopic source, but that will not affect this measurement).
Yes, you could multiply by the number of primaries, and that would give you the number of scored particles in your FLUKA run (although be careful if you do use multiple runs, as you will need this multiplier). More commonly, the number of primaries in FLUKA does not have a physical value, and C is chosen as to give a physical meaning to the scaled FLUKA scoring (this might not be true in all cases, especially in event-by-event investigations). Nevertheless, your approach is reasonable, and your 1st reason is correct.
Regarding the “what the detector sees” part, this would need to be informed by the detector design. In particular, the result you just simulated tell you how much energy is deposited in the detector, the process of converting that into a signal introduces noise, so it is customary to average over more than a single bin. Nevertheless, in an ideal scenario, this is indeed correct.
Assuming Poisson statistics uncertainty, i.e. sqrt(n), is not a bad approach when you can consider each event to occur independently with a probability p (e.g. the probability of a photon undergoing photoelectric interaction), which in your case, considering you are using isotropic sources, is perfectly reasonable. However, if in the future you perform a simulation where the randomness of the events is put in question, utilizing the FLUKA method for error estimation will give you better results, provided sufficient cycles.
In summary, yes what you are doing appears to be correct. Is it the case that you are still getting a discrepancy in the number of gamma rays?
Dear Stefano,
Thank you very much for your patience and the detailed explanation. At least I feel more confident about what I’m doing now. Unfortunately, when I compare the counts, I still find this 5% discrepancy, with a statistical uncertainty of 1%, which makes the measurements incompatible. Could this be a problem with FLUKA at this point? Perhaps there’s something in the simulation of the decay that we’re missing, but it seems strange that this would be a statistically driven issue.
This is very strange. I have executed your inputs with no modification save for fixing the EMF-CUT what(3) for the transport option (only specified for Cs137, which could cause indeed a bias, but as you see below it still works).
For the bin of interest around 662kev I get:
2.603E-4 with 1.502% uncertainty for the monochromatic source
And
2.188E-4 with 1.512% uncertainty for the Cs137 case
If I then scale the former by 0.85 I find that the relative difference between these two numbers is ~1% with a 2% relative error. I also ran for longer, and can get a much better agreement, although the 0.85 is itself an approximation (it’s slightly higher the real intensity).
I am assuming you get different numbers, can you send them to me? If you get the same, then it is the Python script you write that is not correct.
Hi Stefano,
thanks again and sorry, but I’ve been busy with other tasks these past few days.
Could you send me your inputs at this point? I can’t seem to understand the difference.
here are the inputs I ran to get the results I sent above. The 662kev file is exactly the one you sent me, just press run. The 137Cs file had an error in the EMF-CUT for SDUM transport at what(3), this can happen if you duplicate another EMF-CUT and change the SDUM , which leaves some hidden uninitialized parameters. All you have to do is create a new EMF-CUT and put the parameters you want.
In Flair, the files you sent me and the one I am sending back will appear exactly identical, but do be mindful about that error in EMF-CUT (although you must have it fixed, otherwise the simulation will simply not run).
Let me know if these are of help. It would also be helpful if you send me your fluka and flair versions, and if you want I can take a look at the python processing script as well.
Dear Stefano,
I’m attaching the very simple script I’m using to extract the information. After downloading your inputs, I got the same values as well, so I’m not sure what the issue was, but everything seems to be working fine now.
Thank you very much for your help; you’ve been incredibly helpful!
