ユーザーガイド#
ここでは jgdtrans の使い方を紹介します。各 API の詳細は API Reference をみてください。
インストール#
PyPI よりインストールできます。
pip install jgdtrans
依存は typing-extensions のみです。また Python >=3.9 であることが必要です。
par ファイルの読み込み#
jgdtrans は、以下の国土地理院が2023年時点で公開している全ての
Gridded Correction Parameter ファイル(拡張子にちなんで
par ファイルと呼びます)の読み込みに対応しています。
POS2JGD (geonetF3 and ITRF2014)
ここでは、 par ファイルの読み込み方法を紹介します。
jgdtrans は par ファイルを提供しません。
par ファイルを利用する場合は、国土地理院よりダウンロードしてください [1]。
load() と loads()
によって par ファイルを読み込みます。この
API は Transformer を返します。
format 引数を用いて、par ファイルのフォーマットを指定してください。
>>> import jgdtrans
>>> with open('SemiDyna2023.par') as fp:
... tf = jgdtrans.load(fp, format='SemiDynaEXE')
>>> tf
Transformer(unit=5, parameter=<dict (21134 length) at 0x123456789>, description='for [...]')
ヘッダ、 フォーマット、パラメータには、それぞれ、
Transformer.description 、
Transformer.format 、
Transformer.parameter にてアクセスできます。
>>> tf.description
'for SemiDynaEXE Ver.1.0.0\n[...]\nMeshCode dB(sec) dL(sec) dH(m)'
>>> tf.format
'SemiDynaEXE'
>>> tf.parameter
{
36230600: Parameter(
latitude=-0.05475,
longitude=0.04,
altitude=0.07721,
),
# and go on
}
Transformer.parameter のエントリは、
par ファイルのパラメータ部の一行に対応します。キーはメッシュコード(meshcode)、値は Parameter オブジェクトです。
TKY2JGD の高度、 PatchJGD(H) の経緯度は 0.0 になります。
>>> with open('TKY2JGD.par') as fp:
... tf = jgdtrans.load(fp, format='TKY2JGD')
>>> tf.parameter[46303582]
Parameter(latitude=12.79799, longitude=-8.13354, altitude=0.0)
Transformer.to_dict() は、
Transformer オブジェクトを dict に変換して返します。
>>> tf.to_dict()
{
'description': 'for SemiDynaEXE Ver.1.0.0\n[...]\nMeshCode dB(sec) dL(sec) dH(m)',
'format': 'SemiDynaEXE',
'parameter': {
36230600: {
'latitude': -0.05475,
'longitude': 0.04,
'altitude': 0.07721,
},
# and go on
},
}
from_dict() は、前述の dict オブジェクトを受け取る
Transformer のコンストラクタです。
>>> data = {
... 'description': 'my SemiDynaEXE',
... 'format': 'SemiDynaEXE',
... 'parameter': {
... 36230600: {
... 'latitude': -0.05475,
... 'longitude': 0.04,
... 'altitude': 0.07721,
... },
... # and go on
... },
... }
>>> tf = jgdtrans.from_dict(data)
>>> tf
Transformer(unit=5, parameter=<dict (21134 length) at 0x987654321>, description='my [...]')
これらを用いることで、 par 形式データを、 python dict や他シリアライズ形式(例えば JSON)に相互に変換できます。
座標変換#
ここでは、 Transformer オブジェクトを使った座標変換方法を紹介します。
順方向の変換は Transformer.forward() によって、
逆方向の変換は Transformer.backward() によって行います。いずれも国土地理院の web API よりも精度が高いです。
>>> tf.forward(36.10377479, 140.087855041, 2.34)
Point(latitude=36.103773017086695, longitude=140.08785924333452, altitude=2.4363138578103)
>>> tf.backward(36.103773017086695, 140.08785924333452, 2.4363138578103)
Point(latitude=36.10377479, longitude=140.087855041, altitude=2.34)
戻り値は Point オブジェクトです。経緯度、高度には各属性でアクセスできます。
>>> point = tf.forward(36.10377479, 140.087855041, 2.34)
>>> point.latitude
36.103773017086695
>>> point.longitude
140.08785924333452
>>> point.altitude
2.4363138578103
Point オブジェクトはアンパックをサポートしています(Point
は長さ 3 の Sequence[flaot] です)。
>>> origin = Point(36.10377479, 140.087855041, 2.34)
>>> result = tf.forward(*origin)
>>> tf.backward(*result)
Point(latitude=36.10377479000002, longitude=140.087855041, altitude=2.34)
引数(backward)によって順逆変換を切り替える Transformer.transform()
も提供しています。動作は前述の Transformer.forward() や Transformer.backward()
と同じです。つまり、任意の Point オブジェクト point に対して、以下の恒等式が成り立ちます。
>>> tf.transform(*point, backward=False) == tf.forward(*point)
True
>>> tf.transform(*point, backward=True) == tf.backward(*point)
True
Transformer.backward() は厳密解を提供しませんが、真の解からの誤差を Transformer.MAX_ERROR
以下であることを保証します。すなわち、経緯度の誤差は、国土地理院が公開しているパラメータの誤差
\(10^{-9}\) [deg] 以下となり、高度の誤差は、 \(10^{-5}\) [m] 以下となります(と思います)。
Transformer.backward() は国土地理院の web API と互換性はありません。互換な逆変換は
Transformer.backward_compat() で行えます。
TKY2JGD や PathJGD の高度や PatchJGD(H) や HyokoRev の経緯度のように、パタメータが 0.0 の場合、その成分上の変換は(順逆ともに)恒等変換になります。つぎは TKY2JGD の例です。
>>> tf.forward(36.10377479, 140.087855041, 2.34)
Point(latitude=36.106966279935016, longitude=140.08457686562787, altitude=2.34)
Transformer は、DMS 形式(60 進度)による I/O もサポートしています。
Point には DMS 形式をあつかうためのメソッドがあります。
Point.from_dms() は DMS 形式によるコンストラクタで、
Point.to_dms() は Point オブジェクトを DMS 形式に変換します。
>>> Point.from_dms("360613.58925", "1400516.27815", 2.34)
Point(36.10377479166667, 140.08785504166664, 2.34)
>>> Point(36.10377479166667, 140.08785504166664, 2.34).to_dms()
("360613.58925", "1400516.27815", 2.34)
Point.to_dms() の戻り値は
tuple[str, str, float] であり Point ではないことに注意してください。
これらを組み合わせることで、DMS 形式による Transformer の I/O が可能です。
>>> point = Point.from_dms("360613.58925", "1400516.27815", 2.34)
>>> tf.forward(*point).to_dms()
("360613.58286751213", "1400516.293278404", 2.4363138577557226)
メッシュに関連する実装の導入#
jgdtrans を利用するにあたり、この章を読む必要はありません。
この章では、メッシュに関連する実装の理解を助けるために、各用語の定義とその実装との関連を説明します。
各実装は mesh モジュールに公開されています。
ここでの用語は本実装に固有であり、オリジナルを含む他実装と関連はありません(meshcode を除く)。
パラメータは、日本周辺で定義された正方格子の上に定義されています。この正方格子のことをメッシュ(mesh)と呼びます。格子の座標系のことを
mesh coordinate 、格子点のことを mesh node (もしくは node)と呼びます。それぞれ MeshCoord と
MeshNode として実装しています。
mesh node には、対応するある非負整数が存在し、それをメッシュコード(meshcode)と呼びます。経緯度高度の各パラメータは、メッシュコードによって、
mesh node に紐づけられます。特に par ファイルでは、そのように紐づけられています。
メッシュは、格子間隔がおよそ 1 [km] のメッシュと、 5 [km]
のメッシュの 2 種類が存在します。格子間隔のことを mesh unit (もしくは unit)と呼び、
1 ならびに 5 リテラルで表現しています。
Mesh coordinate (MeshCoord)は unit 1 の座標系として実装しています。
unit 5 の座標系は、簡単のために、5つ飛びの unit 1 の座標系として実装しています。したがって、5つ飛びでない
node は、
unit 5 の時に扱うことができません。例えば、
node だけでは unit を決定できないので、 node から cell (すぐあとで説明します)の解決では、ユーザーが明示的に unit
を指定しなければなりません。
変換の補正を計算するには、その unit における単位胞が必要です。単位胞 のことを mesh cell(もしくは cell )と呼び、
MeshCell として実装しています。
Point にも関連するメソッドが存在します。