I’d like to score the LETd of carbon ions in water.
I used GETLET routine comscw.f. to weight, and it went well when i scored all the particle.
But when I used the AUXSCORE card to focus on the LETd of carbon ions(setting Z=6), LETd did not dropped after the peak as it is supposed to.
I referred to your answer at " Strange LET shape for 800 keV protons in LiF", at the end of which you have mentioned that:
" *The latter are numerically prone to numerical troubles at the end of the range as the fluence drops to zero (given the 1/0 behavior)."
It seems that i have met the similar problem, and i wonder how could i deal with it.
Any suggestions and help are important to me, and Thank you very much!
As it was also explained in the post you have linked (Strange LET shape for 800 keV protons in LiF), LETd defined in literature is derived from fluence (see equation 2). To evaluate the numerator in the equation you should use the fluscw.f routine so as to weight the scored fluence with the requested power of the stopping power. The comscw.f is used instead to weight deposited energy, stars or residual nuclei.
Please have a look and do not hesitate to come back to us in case of further questions.
Following your advice, I tried to use use the fluscw.f routine, but met with some problems.
I used the USRBIN with option of “ALL-PART” to score fluence and tried to weight with fluscw.f routine. I mainly wrote the fluscw.f file based on the comscw file I used before. It compiled well, but when it ran, the problems came up.
The attachment is my fluscw.f file and the err. file, maybe you could have a look and give me some suggestions. And thanks again for your help!
although the routine compiles, there may be some issues with it.
Note that the input arguments of fluscw are almost enough for your needs, and you certainly do not need to use variables such as DTRACK accessible via the trackr common. Kinetic energy or momentum are provided via the PLA variable and the current region number is NREG. Maybe it was intended or maybe it is just a typo, but with MEDFLK(MREG,2) you get the index of the material assigned to a region during the transport of decay radiation: it should be MEDFLK(MREG,1) for prompt.
Please find attached the routine with some adjustments: I left some comments in case it is not fully clear.
I suggest to have a second careful look at the post liked above. The numerical problems occur in the post processing when you have to divide the estimate of the numerator of LET_d and LET_t by the denominator which also goes to zero when the fluence goes to zero. This issues have nothing to do with the way you estimate the numerator of LET_d and LET_t with the FLUKA routine, they are somehow intrinsic to the definition.
As a side note, the function GETLET returns the (un)restricted LET in keV/(\mu m g/cm3) so the post-processing will be easier if you multiply by the density and convert to GeV/cm directly in the routine (see the attachment).
Last comment (I promise ;)): if you want to separate the contribution from the primary carbon ions you can select BEAMPART in USRBIN
It is very thoughtful of you to post the fluscw routine file, and your comments also did great help to me . It went well, and i got the result i wanted.