MySQL8.0で北極点と南極点2点間に線を引いてみた

本エントリの概要
前提
極点の表し方
2点間で線を引く
線はどこにひかれたのか
ちょっといたずら
まとめ

本エントリの概要

この記事は RDBMS-GIS(MySQL,PostgreSQLなど) Advent Calendar 2022 の6日目として作成しています。
ふと思いついてやってみたけどよくわからないことになったメモです。

前提

手元にあったMySQL8.0.17で検証しています。
SRIDは4326を使います。
SRID=4326の場合、AXIS["Lat",NORTH],AXIS["Lon",EAST]の順番なので、緯度→経度の順番で書きます。
今回はざっくりと算出したいので細かい精度はこだわりません。

極点の表し方

まずは2点の場所について考えます。

北極点WKT：POINT(90 0)
南極点WKT：POINT(-90 0)

極点に経度はないので仮に0とします。

SELECT ST_Distance_Sphere(
    ST_GeomFromText('POINT(90 0)',4326),
    ST_GeomFromText('POINT(-90 0)',4326)
);

結果：20015114.352233686 (単位はm(メートル)です)
地球一周が4万km弱なのでだいたいあってますね。
ちなみに経度は-180～180の間で設定できますので、こんなことをいろいろ試してみましたが、距離を測るといずれも同じ結果でした。

SELECT ST_Distance_Sphere(
    ST_GeomFromText('POINT(90 180)',4326),
    ST_GeomFromText('POINT(-90 0)',4326)
);
SELECT ST_Distance_Sphere(
    ST_GeomFromText('POINT(90 0)',4326),
    ST_GeomFromText('POINT(-90 180)',4326)
);

ここまでは特に問題ありません。

2点間で線を引く

では、北極点と南極点を結ぶ線を描いてみます。
WNTはLineString(90 0,-90 0)とします。
この線の長さを測ってみます。

SELECT ST_Length(ST_GeomFromText('LineString(90 0,-90 0)',4326));

結果：20003917.356955886
ヨシ、約2万kmです。
2点間距離の結果と少し違いますが、まぁそういうものだろうということで差は一旦忘れます。

線はどこにひかれたのか

「さて、この線が地球上のどこにひかれているのでしょうか。」（これが本エントリの主題です）
北極点と南極点の最短ルートであること、いずれかの緯線に沿っていることは確信するところですが、経度はどこを通ってもここでの計算上は同じになるはずです。
「とりあえずまぁ本初子午線を通ってるんじゃまいか」という仮説を立ててST_Contains()で緯度経度0 0を通っているか確かめます。

SELECT ST_Contains(
    ST_GeomFromText('LineString(90 0,-90 0)',4326),
    ST_GeomFromText('Point(0 0)',4326)
);

結果：1 （引数1に引数2が含まれている場合1、含まれない場合0）
ヨシ。やっぱり。
念のためほかの経度を通っていないかも確認します。

SELECT
ST_Contains(
    ST_GeomFromText('LineString(90 0,-90 0)',4326),
    ST_GeomFromText('Point(0 0)',4326)
) AS Lon0,
ST_Contains(
    ST_GeomFromText('LineString(90 0,-90 0)',4326),
    ST_GeomFromText('Point(0 90)',4326)
) AS Lon90,
ST_Contains(
    ST_GeomFromText('LineString(90 0,-90 0)',4326),
    ST_GeomFromText('Point(0 180)',4326)
) AS Lon180,
ST_Contains(
    ST_GeomFromText('LineString(90 0,-90 0)',4326),
    ST_GeomFromText('Point(0 -90)',4326)
) AS Lonm90;

結果：Lon0のみ1、それ以外は0
ヨシ。やっぱり経度0のみです。

さて、さっき「経度は-180～180の間で設定できます」って試したのをこちらでもやってみます。

SELECT ST_Length(ST_GeomFromText('LineString(90 180,-90 180)',4326));

結果：20003917.356955886
同じ長さです。
この場合も経度0を通っているのか確認してみます。

SELECT ST_Contains(
    ST_GeomFromText('LineString(90 180,-90 180)',4326),
    ST_GeomFromText('Point(0 0)',4326)
);

結果：0
あれ、通っていません。
ということは経度180を指定しているのが効いているのでしょうか。

SELECT ST_Contains(
    ST_GeomFromText('LineString(90 180,-90 180)',4326),
    ST_GeomFromText('Point(0 180)',4326)
);

結果：1

そのようです。
念のためほかの経度も確認します。

SELECT
ST_Contains(
    ST_GeomFromText('LineString(90 180,-90 180)',4326),
    ST_GeomFromText('Point(0 0)',4326)
) AS Lon0,
ST_Contains(
    ST_GeomFromText('LineString(90 180,-90 180)',4326),
    ST_GeomFromText('Point(0 90)',4326)
) AS Lon90,
ST_Contains(
    ST_GeomFromText('LineString(90 180,-90 180)',4326),
    ST_GeomFromText('Point(0 180)',4326)
) AS Lon180,
ST_Contains(
    ST_GeomFromText('LineString(90 180,-90 180)',4326),
    ST_GeomFromText('Point(0 -90)',4326)
) AS Lonm90;

結果：Lon180のみ1、それ以外は0
問題なさそうです。

ちょっといたずら

北極点と南極点で違う経度を指定して、まずは距離を測ります。

SELECT ST_Length(ST_GeomFromText('LineString(90 0,-90 180)',4326));

結果：20003917.35695591
今度は距離が目視できないレベルで変化しました。
が、約2万kmであることは変わらないです。

さてこれはどこを通っているのでしょうか。

SELECT
ST_Contains(
    ST_GeomFromText('LineString(90 0,-90 180)',4326),
    ST_GeomFromText('Point(0 0)',4326)
) AS Lon0,
ST_Contains(
    ST_GeomFromText('LineString(90 0,-90 180)',4326),
    ST_GeomFromText('Point(0 90)',4326)
) AS Lon90,
ST_Contains(
    ST_GeomFromText('LineString(90 0,-90 180)',4326),
    ST_GeomFromText('Point(0 180)',4326)
) AS Lon180,
ST_Contains(
    ST_GeomFromText('LineString(90 0,-90 180)',4326),
    ST_GeomFromText('Point(0 -90)',4326)
) AS Lonm90;

結果：Lon0とLon180が1、それ以外は0
え？あれ？いったいなにが起こっているの。
経度0と180を緯度を10度ずつずらしながら検査してゆきます。

SET @line1 = ST_GeomFromText('LineString(90 0,-90 180)',4326);
SELECT
ST_Contains(@line1, ST_GeomFromText('Point(80 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(70 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(60 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(50 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(40 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(30 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(20 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(10 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(0 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-10 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-20 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-30 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-40 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-50 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-60 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-70 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-80 0)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(80 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(70 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(60 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(50 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(40 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(30 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(20 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(10 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(0 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-10 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-20 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-30 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-40 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-50 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-60 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-70 180)',4326)) AND
ST_Contains(@line1, ST_GeomFromText('Point(-80 180)',4326));

結果：1
一つでも0が混ざれば結果は0になるはずですので、すべて通っていることになります。
これはつまり「経度0と経度180の両方に線があるので地球1周分のはずが、長さは約2万kmである」ということになりますね。
なぞです。

まとめ

極点を扱う場合は経度に指定する数値を何かに固定して運用するのがよいのかもしれません。
これを解析するだけの能力は私にはないようですので、「変なことをするのが悪い」ということにしたいと思います。
手元にMySQL以外のGIS環境がないので他のシステムでどうなるのかはわかりません。

untitled .engineer

技術系のブログ（仮）