はじめに

通常運転の記事をFOSS4G Advent Calendar 2017のネタとしたものです。

本題

PostGISユーザがMySQL 8の空間拡張機能を少し触ってみたで、MySQLにおいては、WKT表現では座標系のaxis orderはその座標系の定義によることと、内部表現ではWKTと違い横-縦に固定していることとを述べました。

次の結果を見てみると、後ろ32桁が、16桁の前半後半で入れ替わっていることが分かります。

mysql> SELECT HEX(ST_GeomFromText('POINT(1 2)', 4326));
+----------------------------------------------------+
| HEX(ST_GeomFromText('POINT(1 2)', 4326))           |
+----------------------------------------------------+
| E610000001010000000000000000000040000000000000F03F |
+----------------------------------------------------+
1 row in set (0.00 sec)

mysql> SELECT HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 4326)));
+----------------------------------------------------+
| HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 4326))) |
+----------------------------------------------------+
| 0101000000000000000000F03F0000000000000040         |
+----------------------------------------------------+
1 row in set (0.00 sec)

他の座標参照系ではどうなるかを確認してみたいと思います。また、それが本当に座標系の定義に合致しているかを見てみたいと思います (疑り深い)。

前提

EPSG:4326においては、WKTと内部表現のaxis orderが違います。

POINT(1 2)のWKBは0101000000000000000000F03F0000000000000040です。うち1を示す部分は000000000000F03Fで、2を示す部分は0000000000000040です。

1, 2の順で渡した場合 (かつリトルエンディアンである場合)、末尾が40になっているときはWKTと内部表現のaxis orderが違い、3Fのときは同じ、となります。

ざーっと見ていきます

mysql> SELECT HEX(ST_GeomFromText('POINT(1 2)', 0));
+----------------------------------------------------+
| HEX(ST_GeomFromText('POINT(1 2)', 0))              |
+----------------------------------------------------+
| 000000000101000000000000000000F03F0000000000000040 |
+----------------------------------------------------+
1 row in set (0.00 sec)

mysql> SELECT HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 0)));
+-------------------------------------------------+
| HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 0))) |
+-------------------------------------------------+
| 0101000000000000000000F03F0000000000000040      |
+-------------------------------------------------+
1 row in set (0.00 sec)
mysql> SELECT HEX(ST_GeomFromText('POINT(1 2)', 4326));
+----------------------------------------------------+
| HEX(ST_GeomFromText('POINT(1 2)', 4326))           |
+----------------------------------------------------+
| E610000001010000000000000000000040000000000000F03F |
+----------------------------------------------------+
1 row in set (0.00 sec)

mysql> SELECT HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 4326)));
+----------------------------------------------------+
| HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 4326))) |
+----------------------------------------------------+
| 0101000000000000000000F03F0000000000000040         |
+----------------------------------------------------+
1 row in set (0.00 sec)
mysql> SELECT HEX(ST_GeomFromText('POINT(1 2)', 3857));
+----------------------------------------------------+
| HEX(ST_GeomFromText('POINT(1 2)', 3857))           |
+----------------------------------------------------+
| 110F00000101000000000000000000F03F0000000000000040 |
+----------------------------------------------------+
1 row in set (0.00 sec)

mysql> SELECT HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 3857)));
+----------------------------------------------------+
| HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 3857))) |
+----------------------------------------------------+
| 0101000000000000000000F03F0000000000000040         |
+----------------------------------------------------+
1 row in set (0.00 sec)

SRID0 (未定義)のとき、3857 (Webメルカトル)のときは、WKTと内部表現のaxis orderが同じです。
4326はWKTと内部表現のaxis orderが違います。

本当かどうか確認する

http://www.epsg-registry.org/ に、EPSGの空間参照系検索サイトがあります。

左上の"retrieve by code"をクリックして、コード検索を行います。

トップページの左上部分

EPSG:4326は北-東のオーダー

Codeに4326を入れてみます。
検索結果は、次のように2件出ますが、typeがareaのものは今回は不要なので、WGS 84のみ見てください。

2件がヒットした様子

座標系の定義がざーっと表示されます。

今回は、axis orderのみを見るので、"Ellipsoidal CS"を開いてみます。

EPSG:4326の座標軸の定義

"Axes"の表を見ると、1番目が"Lat"で北方向、2番目が"Lon"で東方向、となっているのが分かります。

以上から、EPSG:4326のaxis orderは縦(南北)-横(東西)の順になることが、これで分かります。

EPSG:3857は東-北のオーダー

同じようにして、EPSG:3857も見てみます。
同じように検索して表示させるのですが、注目するカテゴリ名は、"Cartesian CS"となります。EPSG:4326の時は"Ellipsoidal CS"でしたので、名称が異なりますが、とりあえず気にしないで下さい。

EPSG:3857の座標軸の定義

以上から、EPSG:3857のaxis orderは横(東西)-縦(南北)の順になることが、これで分かります。

もうちょっとやってみる

今度は、JGD2000平面直角座標系1系で同じことをやってみましょう。

mysql> SELECT HEX(ST_GeomFromText('POINT(1 2)', 2443));
+----------------------------------------------------+
| HEX(ST_GeomFromText('POINT(1 2)', 2443))           |
+----------------------------------------------------+
| 8B0900000101000000000000000000F03F0000000000000040 |
+----------------------------------------------------+
1 row in set (0.00 sec)

mysql> SELECT HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 2443)));
+----------------------------------------------------+
| HEX(ST_AsWKB(ST_GeomFromText('POINT(1 2)', 2443))) |
+----------------------------------------------------+
| 0101000000000000000000F03F0000000000000040         |
+----------------------------------------------------+
1 row in set (0.00 sec)

ということで、平面直角座標系はWKTと内部表現のaxis orderが同じ、となってしまいます。

間違ってますやん

EPSGレジストリを見てみます。

平面直角座標系1系のaxis order

ご覧の通り、North-Eastの順となっているので、EPSG:4326と同様、axis orderは縦(南北)-横(東西)の順になります。

ということは、WKTと内部表現のaxis orderが違わなければいけません。

とても困ります

axis orderが想定と違うととても困ります。

MapServerのWMS 1.3.0で、axis orderを、左(東)-上(北)固定、としていたのを、座標系の定義による、と変更したことがありました。このアップデートの直後の頃は、QGIS等のクライアントは従来の解釈である「左(東)-上(北)固定」でした。

たとえばEXTENTminx,miny,maxx,maxyで、x, yが出てきます。これをどう解釈するかが、サーバとクライアントで異なったのです。

サーバ側が何も考えずにMapServerをアップグレードしたため、クライアント側から「おかしくなった」という報告が出て、でも原因特定がすぐにはできなかった、という話もあったようです。

おわりに

MySQLのaxis orderについてもうちょっと見ました。また、座標系の定義を提供してくれるサイトで座標系のaxis orderを確認してみました。

そして、axis orderが座標系の定義と異なるものがあることを示しました。

どこに言っていけばいいのか、それとも開発チームは実はこの問題を把握しているのか、承知しておりません。もし把握されていないなら、これはとても困る問題です。