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.
Main cotrans function.
The main function in both IDL SPEDAS and python pyspedas is called "cotrans". In python it requires either a tplot name or a list of times and a list of data. Here is a comparison of IDL to python code for a simple cotrans operation:
IDL code: SPEDAS
trange = ['2010-02-25/00:00:00', '2010-02-25/23:59:59'] probe = 'a' name_in = "tha_state_pos" name_out = "tha_state_pos_new_geo" thm_load_state, probe=probe, trange=trange cotrans, name_in, name_out, /gei2geo get_data, name_in,data=din, limit=l_in, dl=dl_in get_data, name_out,data=dout, limit=d_l_in, dl=d_dl_in print, "Input data length: ", n_elements(din.x) print, "Input data [100-105]:" transpose(din.y[100:105, *]) print, "Output data length: ", n_elements(dout.x) print, "Output data [100-105]:" transpose(dout.y[100:105, *])
IDL results:
Input data length: 1440
Input data [100-105]:
-45141.852 13490.809 -6501.9253
-45287.578 13420.789 -6511.8516
-45432.707 13350.592 -6521.6914
-45577.234 13280.219 -6531.4453
-45721.172 13209.673 -6541.1143
-45864.516 13138.956 -6550.6982
Output data length: 1440
Output data [100-105]:
45184.977 -13345.651 -6501.9253
45271.961 -13473.378 -6511.8516
45358.098 -13601.920 -6521.6914
45443.371 -13731.304 -6531.4453
45527.785 -13861.525 -6541.1143
45611.340 -13992.556 -6550.6982
Python code: pySPEDAS
import pyspedas
import pytplot
from pyspedas.cotrans.cotrans import cotrans
trange = ['2010-02-25/00:00:00', '2010-02-25/23:59:59']
probe = 'a'
name_in = "tha_pos"
name_out = "tha_pos_new_geo"
pyspedas.themis.state(probe=probe, trange=trange, time_clip=True, varnames=[name_in])
cotrans(name_in=name_in, name_out=name_out, coord_in="gei", coord_out="geo")
din = pytplot.get_data(name_in)
dout = pytplot.get_data(name_out)
print("Input data length: " + str(len(din[0])))
y = din[1]
print("Input data [100-105]:")
print(y[100:106, :])
print("Output data length: " + str(len(dout[0])))
yy = dout[1]
print("Output data [100-105]:")
print(yy[100:106, :])
Python results:
Input data length: 1440 Input data [100-105]: [[-45141.85 13490.809 -6501.9253] [-45287.58 13420.789 -6511.8516] [-45432.707 13350.592 -6521.6914] [-45577.234 13280.219 -6531.4453] [-45721.17 13209.673 -6541.1143] [-45864.516 13138.956 -6550.698 ]] Output data length: 1440 Output data [100-105]: [[ 45184.97759641 -13345.65392237 -6501.92529297] [ 45271.96223983 -13473.37175482 -6511.8515625 ] [ 45358.09480658 -13601.92651913 -6521.69140625] [ 45443.36651621 -13731.31250082 -6531.4453125 ] [ 45527.78422625 -13861.52509347 -6541.11425781] [ 45611.33919795 -13992.55860422 -6550.69824219]]
Basic cotrans functions.
Use the basic cotrans functions to compute the direction of Earth's magnetic axis in GEO, using the IGRF model, and test the results for some other transformation vectors with random data.
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 pyspedas.cotrans.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)
Python 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]))
Transform data from MAG to GEO.
This transformation uses many different functions internally: SM -> GSM -> GSE -> GEI -> GEO.
IDL code: SPEDAS
cotrans_lib 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]] data_in = transpose(d) t = [1577112800, 1577308800, 1577598800, 1577608800, 1577998800] t1 = time_struct(t) subMAG2GEO, t1, data_in, dout print, transpose(dout)
IDL results:
IDL> print, transpose(dout)
-13.1952 -300.294 207.556
236.685 -725.359 -136.598
555.777 48.4522 -20.7720
-64.8366 -481.680 794.762
-417.047 -548.897 -427.348
Python code: pySPEDAS
from pyspedas.cotrans.cotrans_lib import submag2geo 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] mag = submag2geo(t, d) print(mag)
Python results:
[[ -13.19521536 -300.29393966 207.55586259] [ 236.68488882 -725.35935255 -136.59821768] [ 555.7768468 48.45227165 -20.77195057] [ -64.83660486 -481.68020999 794.76242518] [-417.04680107 -548.89736495 -427.34807238]]
Transform data from GEI to J2000.
IDL code: SPEDAS
cotrans_lib 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], [264.14, 61.55, -56.32]] data_in = transpose(d) t = [1577112800, 1577308800, 1577598800, 1577608800, 1577798800, 1577998800] t1 = time_struct(t) subGEI2J2000, t1, data_in, dout print, transpose(dout)
IDL results:
IDL> print, transpose(dout)
245.02820 -103.07841 250.53148
775.01595 6.5942014 -11.479825
123.39361 544.46259 -1.2287885
305.37957 -206.64889 855.51530
461.18337 -563.58286 -357.10918
264.30012 60.387870 -56.824668
Python code: pySPEDAS
from pyspedas.cotrans.cotrans_lib import subgei2j2000 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], [264.14, 61.55, -56.32]] t = [1577112800, 1577308800, 1577598800, 1577608800, 1577798800, 1577998800] mag = subgei2j2000(t, d) print(mag)
Python results:
[[ 245.02820801 -103.07809266 250.53160298] [ 775.01595209 6.59521058 -11.47942543] [ 123.39317313 544.46268498 -1.22861055] [ 305.37948491 -206.64887151 855.51536316] [ 461.1836826 -563.58260063 -357.10921253] [ 264.30009126 60.38797267 -56.82463388]]
Daisy chain transformations.
Both IDL and python contain the following functions for transformations between coordinate systems:
subGEI2GSE, subGSE2GEI
subGSE2GSM, subGSM2GSE
subGSM2SM, subSM2GSM
subGEI2GEO, subGEO2GEI
subGEO2MAG, subMAG2GEO
subGEI2J2000, subJ20002GEI
These functions can be daisy chained. In pyspedas, this can be done automatically using the function subcotrans(time , data, from, to).
The following example is using: GEI -> GSE -> GSM
IDL code: SPEDAS
cotrans_lib 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], [264.14, 61.55, -56.32]] data_in = transpose(d) t = [1577112800, 1577308800, 1577598800, 1577608800, 1577798800, 1577998800] t1 = time_struct(t) subGEI2GSE, t1, data_in, dout data_in = dout subGSE2GSM, t1, data_in, dout print, transpose(dout)
IDL results:
IDL> print, transpose(dout)
0.061859131 263.75420 252.71276
45.984592 762.30346 132.67920
-480.51590 183.48641 -217.03501
-112.06107 407.08298 830.41709
724.92896 341.54693 -125.16594
21.239977 275.98252 -10.621965
Python code: pySPEDAS
from pyspedas.cotrans.cotrans_lib import subcotrans 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], [264.14, 61.55, -56.32]] t = [1577112800, 1577308800, 1577598800, 1577608800, 1577798800, 1577998800] mag = subcotrans(t, d, 'gei', 'gsm') print(mag)
Python results:
[[ 6.17328386e-02 2.63754240e+02 2.52712677e+02] [ 4.59847513e+01 7.62303493e+02 1.32679267e+02] [-4.80515919e+02 1.83486338e+02 -2.17035056e+02] [-1.12061112e+02 4.07083036e+02 8.30417143e+02] [ 7.24928845e+02 3.41547042e+02 -1.25165948e+02] [ 2.12400389e+01 2.75982462e+02 -1.06219558e+01]]