CS_FluorLine_Kissel fails on some energies

[uvainio]

Hi!

I've been using xraylib a lot and I really appreciate it very much! Thank you so much for continuing the work on it, Tom.

Only recently it happened that I came across this problem. So out of millions of calls only one failed so far! So it definitely does not happen often, but nevertheless it causes of a crash, so I thought better to report it in case there is a nice fix for it. So if you call the CS_FluorLine_Kissel with certain energy near the edge Ag L2, the function call fails. Just a little bit different energy will work. I did not test systematically all elements and lines and energies down to fifth decimal, but that might be a fun exercise to run (and would take some calculation time!).

xrl.CS_FluorLine_Kissel(47, -63, 3.5282)
Traceback (most recent call last):
File ..., in xrl.CS_FluorLine_Kissel(47, -63, 3.5282)
File "C:\ProgramData\Anaconda3\lib\site-packages\xraylib.py", line 4424, in CS_FluorLine_Kissel
return _xraylib.CS_FluorLine_Kissel(Z, line, E)
ValueError: Spline extrapolation is not allowed

Kind regards,

Ulla

[tschoonj]

Hi Ulla, this looks like a bug indeed. I will try and come up with a fix, but it may take several weeks before I have the time to do so.

Best, Tom

[uvainio]

Hi Tom! Thanks a lot! Oleg actually did the test since he is able to build the branch. The fix worked and he got a reasonable value. We did not test yet to see if there are other cases or if this was the only case.

[tschoonj]

Ok great, I will merge it in for now, let me know if you run into problems.

[shirokobrod]

Hi Tom! Ulla asked me to test your fix and it was OK. I had a look at your fix and I suggest more accurate fix. If we replace in initial release version line 135 in kissel_pe.c if (EdgeEnergy_Kissel[Z][shell] > EdgeEnergy_arr[Z][shell] && E < EdgeEnergy_Kissel[Z][shell]) with if(ln_E < E_Photo_Partial_Kissel[Z][shell][0]) function will extrapolate for missing energies. It will extrapolate in case of original code and when energy is between edge energy and first table energy. I made a test. Z=47, line =-63, shell = L2. E(L2) = 3.5237 keV; E(L2)_Kissel = 3.5282 keV.; E[0] = 3.528770120566. Four cases E1 = 3.526 keV;E2=3.5282;E3=3.528268 keV;E4=3.529. E1 is between edges (original code case), E2 = E(L2)_Kissel, E(L2)_Kissel < E3 < E[0], E4 > E[0].Original code gives CS (CS_FluorLine_Kissel) CS1=7.9084769116; CS2 & CS3 failed, CS4=8.8606749103. Your fix CST1= CS1, CST2=CST3=7.9135685724, CST4=CS4. My fix CSO1=CS1, CSO2=7.9134112988, CSO3=7.9135638162, CSO4=CS4. As you can see my fix needs to change only one line of the code, more accurate and consistent with extrapolation logic. kissel_pe.exe filed because it uses CSb_Photo_Partial(47, L2_SHELL, 3.5282, &error); calculated with your fix. Regards, Oleg

[tschoonj]

Hi @shirokobrod,

I had a look at your suggestion and it looks indeed like a better fix for this issue. I was wondering if you would be willing to open a pull request with your patch? I will happily guide you through the process if necessary.

Best,

Tom