2015年6月25日木曜日

astropy.coordinatesを用いた座標の取り扱い、その7 ― 座標変換

astropy.coordinatesは、座標変換を行うための様々な機能を装備しています。天文学で用いられる主な座標系は予めastropyの中に組み込まれていますが、ユーザーが定義した特殊な座標系についても、独自の座標をastropy内で定義することによって他の座標系との間の座標変換を取り扱うことができます。独自の座標系の定義については別の機会に取り扱うとして、ここではastropy.coordinatesの中でどのように座標変換を実行できるかを見てゆくことにします。astropy.coordinatesに予め組み込まれた座標系と座標変換は、ここで確認することができます。

まず、一番単純な座標変換を見てみましょう。

>>> import astropy.units as u
>>> from astropy.coordinates import SkyCoord
>>> gc = SkyCoord(l=0*u.degree, b=45*u.degree, frame='galactic')
>>> gc.fk5  
<SkyCoord (FK5: equinox=J2000.000): (ra, dec) in deg
    (229.272514629, -1.12844288043)>

より一般的には、transform_to() で変換を表します。transform_to() 用いることで入力方法が柔軟になり、なおかつ、より細かいオプションの設定が可能となります。

>>> from astropy.coordinates import FK5
>>> gc.transform_to('fk5')  
<SkyCoord (FK5: equinox=J2000.000): (ra, dec) in deg
    (229.272514629, -1.12844288043)>
>>> gc.transform_to(FK5)  
<SkyCoord (FK5: equinox=J2000.000): (ra, dec) in deg
    (229.272514629, -1.12844288043)>
>>> gc.transform_to(FK5(equinox='J1980.0'))  
<SkyCoord (FK5: equinox=J1980.000): (ra, dec) in deg
    (229.014693505, -1.05560349378)>

transform_to() には、SkyCoordオブジェクトをそのまま指定することもできます。SkyCoordオブジェクトを指定すると、そのSkyCoordオブジェクトで指定されている座標系に変換されます。

>>> sc = SkyCoord(ra=1.0, dec=2.0, unit='deg', frame=FK5, equinox='J1980.0')
>>> gc.transform_to(sc)  
<SkyCoord (FK5: equinox=J1980.000): (ra, dec) in deg
    (229.014693505, -1.05560349378)>

また、FK5, FK4, FK4NoETermsなどの座標系においては、同じ座標系の中での変換(例えば、分点を変更する)を行うことができます。

>>> fk5c = FK5('02h31m49.09s', '+89d15m50.8s')
>>> fk5c.equinox
<Time object: scale='utc' format='jyear_str' value=J2000.000>
>>> fk5c  
<SkyCoord (FK5: equinox=J2000.000): (ra, dec) in deg
    (37.9545416667, 89.2641111111)>
>>> fk5_2005 = FK5(equinox='J2005')  # 分点を変更します
>>> fk5c.transform_to(fk5_2005)  
<SkyCoord (FK5: equinox=J2005.000): (ra, dec) in deg
    (39.3931763878, 89.2858442155)>

分点の変更は、以下のようにTimeオブジェクトを用いて行うこともできます。

>>> from astropy.time import Time
>>> fk5c = FK5('02h31m49.09s', '+89d15m50.8s',
...            equinox=Time('J1970', scale='utc'))
>>> fk5_2000 = FK5(equinox=Time(2000, format='jyear', scale='utc'))
>>> fk5c.transform_to(fk5_2000)  
<SkyCoord (FK5: equinox=2000.0): (ra, dec) in deg
    (48.0231710002, 89.386724854)>