Я пытаюсь преобразовать географические и геомагнитные координаты. Я нашел следующий скрипт на Прологе, но я недостаточно понимаю его, чтобы преобразовать его самостоятельно. Целевой язык — Java, но подойдет любой понятный язык (C, Python, VB и т. д.).
http://idlastro.gsfc.nasa.gov/ftp/pro/astro/geo2mag.pro
Если бы кто-нибудь мог помочь с преобразованием этого скрипта или объяснить, что именно он делает (эти операции с массивами меня сбивают с толку), я был бы очень признателен.
Спасибо
В зависимости от приложения высота может быть важной переменной в этом преобразовании координат, поскольку геомагнитные координаты представляют собой отображение дипольного магнитного поля Земли.
В Python вы можете легко преобразовать географические координаты в геомагнитные (и наоборот) с помощью SpacePy (http://sourceforge.net/projects/spacepy/).
Поскольку вы ищете исходный код для преобразования в Java, SpacePy реализует библиотеку Fortran International Radiation Belt Environment Modeling (IRBEM), исходный код которой доступен (http://irbem.svn.sourceforge.net/viewvc/irbem/web/index.html)
В Python, если другие ищут быстрое решение:
import spacepy.coordinates as coord
from spacepy.time import Ticktock
import numpy as np
def geotomag(alt,lat,lon):
#call with altitude in kilometers and lat/lon in degrees
Re=6371.0 #mean Earth radius in kilometers
#setup the geographic coordinate object with altitude in earth radii
cvals = coord.Coords([np.float(alt+Re)/Re, np.float(lat), np.float(lon)], 'GEO', 'sph',['Re','deg','deg'])
#set time epoch for coordinates:
cvals.ticks=Ticktock(['2012-01-01T12:00:00'], 'ISO')
#return the magnetic coords in the same units as the geographic:
return cvals.convert('MAG','sph')
Я сделал это для кода Python и попытался проверить на этом веб-сайте http://wdc.kugi.kyoto-u.ac.jp/igrf/gggm/index.html . я нашел это
Я использовал значение для Киото, Япония (35N, 135,45 Вт). Веб-страница рассчитана (25,18, -155,80). У меня получилось (25.33580652, -155.82724011). Так что я не совсем уверен, что это может быть полезно...
import numpy as np
from numpy import pi, cos, sin, arctan2, sqrt, dot
def geo2mag(incoord):
"""geographic coordinate to magnetic coordinate:
incoord is numpy array of shape (2,*)
array([[glat0,glat1,glat2,...],
[glon0,glon1,glon2,...])
where glat, glon are geographic latitude and longitude
(or if you have only one point it is [[glat,glon]])
returns
array([mlat0,mlat1,...],
[mlon0,mlon1,...]])
"""
# SOME 'constants'...
lon = 288.59 # or 71.41W
lat = 79.3
r = 1.0
# convert first to radians
lon, lat = [x*pi/180 for x in lon,lat]
glat = incoord[0] * pi / 180.0
glon = incoord[1] * pi / 180.0
galt = glat * 0. + r
coord = np.vstack([glat,glon,galt])
# convert to rectangular coordinates
x = coord[2]*cos(coord[0])*cos(coord[1])
y = coord[2]*cos(coord[0])*sin(coord[1])
z = coord[2]*sin(coord[0])
xyz = np.vstack((x,y,z))
# computer 1st rotation matrix:
geo2maglon = np.zeros((3,3), dtype='float64')
geo2maglon[0,0] = cos(lon)
geo2maglon[0,1] = sin(lon)
geo2maglon[1,0] = -sin(lon)
geo2maglon[1,1] = cos(lon)
geo2maglon[2,2] = 1.
out = dot(geo2maglon , xyz)
tomaglat = np.zeros((3,3), dtype='float64')
tomaglat[0,0] = cos(.5*pi-lat)
tomaglat[0,2] = -sin(.5*pi-lat)
tomaglat[2,0] = sin(.5*pi-lat)
tomaglat[2,2] = cos(.5*pi-lat)
tomaglat[1,1] = 1.
out = dot(tomaglat , out)
mlat = arctan2(out[2],
sqrt(out[0]*out[0] + out[1]*out[1]))
mlat = mlat * 180 / pi
mlon = arctan2(out[1], out[0])
mlon = mlon * 180 / pi
outcoord = np.vstack((mlat, mlon))
return outcoord
if __name__ == '__main__':
mag = geo2mag(np.array([[79.3,288.59]]).T).T
print mag # should be [90,0]
mag = geo2mag(np.array([[90,0]]).T).T
print mag # should be [79.3,*]
mag = geo2mag(np.array([
[79.3,288.59],
[90,0]
]).T).T
print mag # should be [ [90,0]. [79.3,*] ]
# kyoto, japan
mag = geo2mag(np.array([[35.,135.45]]).T).T
print mag # should be [25.18, -155.80], according to
# this site using value for 1995
# http://wdc.kugi.kyoto-u.ac.jp/igrf/gggm/index.html
