Cotrans
Both the IDL spedas and the python pyspedas contain routines for coordinate transformations in the following systems: GSE, GSM, SM, GEI, GEO, MAG, J2000.
Below, we compare IDL code to python code using some of these cotrans functions, including the IGRF-13 model.
IDL code: SPEDAS
; Compile contrans library cotrans_lib ; Define some data d = [[245.0, -102.0, 251.0], [775.0, 10.0, -10], [121.0, 545.0, -1.0], [304.65, -205.3, 856.1], [464.34, -561.55, -356.22]] ; Define times t = [1577112800, 1577308800, 1577598800, 1577608800, 1577998800] t0 = time_string(t) t1 = time_struct(t0) print, t0 ; Compute direction of Earth's magnetic axis in GEO, using the IGRF model. cdipdir_vect,transpose(t1.year[*]),transpose(t1.doy[*]),gd1,gd2,gd3 print, gd1,gd2,gd3 ; Compute sun direction in GEI system csundir_vect,transpose(t1.year[*]),transpose(t1.doy[*]),transpose(t1.hour[*]),transpose(t1.min[*]),transpose(t1.sec[*]),gst,slong,sra,sdec,obliq print,gst,slong,sra,sdec,obliq ; Compute GEI to GSE transformation. tgeigse_vect,transpose(t1.year[*]),transpose(t1.doy[*]),transpose(t1.hour[*]),transpose(t1.min[*]),transpose(t1.sec[*]),transpose(d[0, *]),transpose(d[1, *]),transpose(d[2, *]),xgse,ygse,zgse print,xgse,ygse,zgse ; Compute GSE to GSM transformation. tgsegsm_vect,transpose(t1.year[*]),transpose(t1.doy[*]),transpose(t1.hour[*]),transpose(t1.min[*]),transpose(t1.sec[*]),transpose(d[0, *]),transpose(d[1, *]),transpose(d[2, *]),xgsm,ygsm,zgsm print,xgsm,ygsm,zgsm
IDL results:
IDL> print, t0
2019-12-23/14:53:20
2019-12-25/21:20:00
2019-12-29/05:53:20
2019-12-29/08:40:00
2020-01-02/21:00:00
IDL> print, gd1,gd2,gd3
0.0486864 0.0486846 0.0486811 0.0486811 0.0486776
-0.156116 -0.156111 -0.156101 -0.156101 -0.156091
0.986538 0.986539 0.986541 0.986541 0.986542
IDL> print,gst,slong,sra,sdec,obliq
5.50119 0.944179 3.24176 3.97097 0.994296
4.73823 4.77856 4.83825 4.84031 4.92061
4.74045 4.78439 4.84933 4.85157 4.93861
-0.408904 -0.408101 -0.405626 -0.405513 -0.399712
0.409047 0.409047 0.409047 0.409047 0.409047
IDL> print,xgse,ygse,zgse
0.0618591 45.9846 -480.516 -112.061 738.670
245.080 773.652 182.717 321.559 318.591
270.862 -13.1524 -217.683 867.127 -103.484
IDL> print,xgsm,ygsm,zgsm
245.000 775.000 121.000 304.650 464.340
-83.8585 11.7033 545.000 -118.216 -460.356
257.629 -7.93938 0.929342 872.399 -479.899
Python code: pySPEDAS
from cotrans_lib import * d = [[245.0, -102.0, 251.0], [775.0, 10.0, -10], [121.0, 545.0, -1.0], [304.65, -205.3, 856.1], [464.34, -561.55, -356.22]] t = [1577112800, 1577308800, 1577598800, 1577608800, 1577998800] a = cdipdir_vect(t) b = csundir_vect(t) gse = tgeigse_vect(t, d) gsm = tgsegsm_vect(t, d)
Results:
a Out[2]: (array([0.04868638, 0.04868462, 0.0486811 , 0.0486811 , 0.04867761]), array([-0.15611589, -0.15611097, -0.15610113, -0.15610113, -0.15609089]), array([0.98653812, 0.98653899, 0.98654072, 0.98654072, 0.98654251])) b Out[3]: (array([5.50118781, 0.94417896, 3.24175902, 3.9709706 , 0.9942959 ]), array([4.73823063, 4.77856341, 4.83825495, 4.84031354, 4.92060693]), array([4.74044499, 4.78438591, 4.84932932, 4.85156629, 4.93861416]), array([-0.40890358, -0.40810141, -0.40562582, -0.40551314, -0.39971157]), array([0.4090472 , 0.40904719, 0.40904716, 0.40904716, 0.40904714])) gse Out[4]: (array([ 6.17328386e-02, 4.59847513e+01, -4.80515919e+02, -1.12061112e+02, 7.38670168e+02]), array([245.07961052, 773.65200071, 182.71689404, 321.55872916, 318.59101371]), array([ 270.86155264, -13.15235508, -217.68322895, 867.12698798, -103.48369801])) gsm Out[5]: (array([245. , 775. , 121. , 304.65, 464.34]), array([ -83.85842685, 11.70325822, 545.00012517, -118.21614302, -460.35592829]), array([ 257.6291215 , -7.93937951, 0.92928139, 872.39913086, -479.89947926]))