PRBEM / IRBEM

IRBEM-LIB provides routines to compute magnetic coordinates for any location in the Earth's magnetic field, to perform coordinate conversions, to evaluate geophysics/space-physics models, and to propagate orbits in time.
https://prbem.github.io/IRBEM/
Other
33 stars 14 forks source link

Fix coordinate transform initialization on call of gsm2sm1 #17

Closed drsteve closed 3 years ago

drsteve commented 3 years ago

Calling a coordinate transformation from GSM to SM, prior to any other transformations returns bad values. This came up in spacepy/spacepy#534

I was able to reproduce the error by calling the IRBEM functions directly from the .so using Python's ctypes interface. I was also able to reproduce the error using the coords_transform method in IRBEM's own Python tools. (see the spacepy issue for the direct calls to the compiled library, but can also be done using the IRBEM module as follows):

import datetime as dt
import IRBEM
coords = IRBEM.Coords(IRBEMdir='~/github/IRBEM/', IRBEMname='onera_desp_lib_linux_x86_64.so')
coords.coords_transform(dt.datetime(2002, 2, 2, 12), [1., 2., 4.], 2, 4)

the return is then

array([[nan, nan, nan]])

gsm2sm1 was missing the calls to initize and init_DTD that are present in the other coordinate conversion routines. These used to be present, and were removed in d0af792b3643a82bf0b01a3013ad1034c7d02e72 during a raft of formatting updates.

This PR adds the initialization calls back. Command line testing as above now yields a non-NaN answer:

array([[1.92990225, 2.        , 3.64355284]])

To verify correctness, we can compare with the high-accuracy Python-based SpacePy routines (in spacepy/spacepy#296 ; note these have almost complete test coverage including comparison with LANLGeoMag and IRBEM):

import spacepy.coordinates as spc
import spacepy.time as spt
cc = spc.Coords([1,2,4], 'GSM', 'car', use_irbem=False)
cc.ticks = spt.Ticktock('2002-02-02T12:00:00')
cc.convert('SM', 'car')

and the result is

Coords( [[1.9292696248256644, 2.0, 3.643887857045691]] , 'SM', 'car')

Pull Request Checklist

drsteve commented 3 years ago

Thanks for this PR, I'll approve and merge it.

Thanks for the quick turnaround!

If you can go ahead and update the IRBEM version in spacepy, it would be great, to fix this issue but also for other updates.

Yeah, I'll make sure that we update, although it's not just a straight dump of files since we supply a custom magnetic field model and I always need to check that none of the new code comments in IRBEM prevent f2py from parsing the source. So it may not be immediate, but it's on my list since we have an open bug report that this will fix!

In particular I think the current version in spacepy lacks the IGRF coefficients after 2020.

We're at rev620 from the old SourceForge SVN. That's the commit where I merged in the IGRF13 update, so we do have the new IGRF. Somehow when I updated the IRBEM source and the revision number in the source folder I failed to update the date in the source folder name... so it says 2019 but it's actually rev620 (from 14 months ago), so has the new IGRF. Of course, we still want to bring in the other updates and fixes, but I wanted to clarify that the IGRF update was including in the latest SpacePy (see spacepy/spacepy#295).