いくつかのシステムでポリゴンを左回りにするか右回りにするかを見て回る

はじめに

するどいコメント

PostGISユーザがSQL Server 2017の空間機能を少し触ってみた では、TSQLの空間機能についてざらっと触った結果を報告したわけですが、そこで、ozo360さんから次のようなコメントを頂きました。

STUnionの結果が「うん、これは同じ結果になりましたね」とありますが順番が異なっているように思います。

これは、PostGISでは右回り、TSQLでは左回り、と違っているので、本当にこれでいいの?というご指摘。

これについては、しっかりと教え諭してやらないといけないと思い、ありのままに書きました。

TSQLは左回りでPostGISが右回りだ、という話ですね。
正直言います、見落としてました。
もうこれは穴が合ったら入って啓蟄になってから出てきたいレベル。

事実です、恥ずかしい。

システム内でそろってればOKなので訂正はしません

でも、結果的に訂正をする必要は無い、と判断しました。

これまで、右手系と左手系でさんざんそろえろよとかなんとか言ってたので、ポリゴンが右回りと左回りとで違っても訂正しないのが分かりにくいかも知れません。が、右手系左手系の話でも、右回り左回りの話でもそうですが、システム内で統一されていれば、どちらでもいいです。システム内で統一されてないなら、まずい(あとあと泣きますよ)、ということです。

以上です

以上で終わってもいいかなと思うのですが、余談的に、いくつかのシステムでどちらに回しているかを見ていきましょう。

試験内容

ポリゴンが右回りか左回りかに敏感でなければならない関数(ここではST_Union)の引数に左回りポリゴンと右回りポリゴンを与え、返り値が左回りか右回りかを見ます。
関数内部では(たぶん)右回りか左回りかのどちらかにそろえてくる(はず)なので、返り値で判定できる(はず)。

さらに、左回りポリゴンを右回りに、右回りポリゴンを左回りに、それぞれ変更して引数として与えた結果も見てみます。2種の引数で返り値の方向が同じなら、それはそろえている、と言ってもいいのではないかと思います。

余談:ポリゴンが右回り左回りがどうでもいい関数もある

たとえば、ST_Area()。これは、外積(のZ成分の値)の合計を2で割ると求まりますが、右回りか左回りかで正負符号が変わります(これは、右回りか左回りかを判定する方法そのものです)。ただ、面積はその絶対値を取ればいいので、正負がどうであろうとどうでもいいので、右回りであろうが左回りであろうがどうでもいいです。

SQLServerの場合

左回り+右回り

SELECT (  
  geometry::STGeomFromText('POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))',0)  
).STUnion(  
  geometry::STGeomFromText('POLYGON((1 0, 1 1, 2 1, 2 0, 1 0))',0)  
).STAsText();  
--------  
POLYGON ((1 0, 2 0, 2 1, 2 2, 1 2, 1 1, 1 0))   

右回り+左回り

SELECT (  
  geometry::STGeomFromText('POLYGON((1 1, 1 2, 2 2, 2 1, 1 1))',0)  
).STUnion(  
  geometry::STGeomFromText('POLYGON((1 0, 2 0, 2 1, 1 1, 1 0))',0)  
).STAsText();  
--------  
POLYGON ((1 0, 2 0, 2 1, 2 2, 1 2, 1 1, 1 0))  

SQLServerのまとめ

左回りにしています。

PostGISの場合

左回り+右回り

SELECT ST_AsText(  
  ST_Union(  
    'POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))'::GEOMETRY,  
    'POLYGON((1 0, 1 1, 2 1, 2 0, 1 0))'::GEOMETRY  
  )  
);  
--------  
POLYGON((1 1,1 2,2 2,2 1,2 0,1 0,1 1))  

右回り+左回り

SELECT ST_AsText(  
  ST_Union(  
    'POLYGON((1 1, 1 2, 2 2, 2 1, 1 1))'::GEOMETRY,  
    'POLYGON((1 0, 2 0, 2 1, 1 1, 1 0))'::GEOMETRY  
  )  
);  
--------  
POLYGON((1 1,1 2,2 2,2 1,2 0,1 0,1 1))  

PostGISのまとめ

右回りにしています。

MySQLの場合

左回り+右回り

SELECT ST_AsText(  
  ST_Union(  
    ST_GeomFromText('POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))'),  
    ST_GeomFromText('POLYGON((1 0, 1 1, 2 1, 2 0, 1 0))')  
  )  
);  
--------  
POLYGON((2 0,2 2,1 2,1 1,1 0,2 0))  

右回り+左回り

SELECT ST_AsText(  
  ST_Union(  
    ST_GeomFromText('POLYGON((1 1, 1 2, 2 2, 2 1, 1 1))'),  
    ST_GeomFromText('POLYGON((1 0, 2 0, 2 1, 1 1, 1 0))')  
  )  
);  
--------  
POLYGON((2 0,2 2,1 2,1 1,1 0,2 0))  

MySQLのまとめ

左回りにしています。

SpatiaLiteの場合

左回り+右回り

SELECT ST_AsText(  
  ST_Union(  
    ST_GeomFromText('POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))'),  
    ST_GeomFromText('POLYGON((1 0, 1 1, 2 1, 2 0, 1 0))')  
  )  
);  
--------  
POLYGON((1 1, 1 2, 2 2, 2 1, 2 0, 1 0, 1 1))  

右回り+左回り

SELECT ST_AsText(  
  ST_Union(  
    ST_GeomFromText('POLYGON((1 1, 1 2, 2 2, 2 1, 1 1))'),  
    ST_GeomFromText('POLYGON((1 0, 2 0, 2 1, 1 1, 1 0))')  
  )  
);  
--------  
POLYGON((1 1, 1 2, 2 2, 2 1, 2 0, 1 0, 1 1))  

SpatiaLiteのまとめ

右回りにしています。

内部表現生成時に方向を強制するようにはなっていない

2019年1月8日になって、内部表現は左回りまたは右回りを強制しているかもしれない、と思いました。

WKTからパースした時点でどちらか一方に強制されても、上記のような結果が出るはずです。

WKTからパースした時点なのかST_Union()関数にかけた時点なのか、ちょっと気になったので確かめました。

右回り左回り両方のWKTから内部表現を生成し、この時点で方向を強制されているなら、ST_AsText()の出力に影響を与えるはずですので、観察してみます。

SQLServerの場合

左回り

SELECT (  
  geometry::STGeomFromText('POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))',0)  
).STAsText();  
--------  
POLYGON ((1 1, 2 1, 2 2, 1 2, 1 1))  

右回り

SELECT (  
  geometry::STGeomFromText('POLYGON((1 1, 1 2, 2 2, 2 1, 1 1))',0)  
).STAsText();  
--------  
POLYGON ((1 1, 1 2, 2 2, 2 1, 1 1))  

PostgreGIS

左回り

SELECT ST_AsText(  
  'POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))'::GEOMETRY  
);  
--------  
POLYGON((1 1,2 1,2 2,1 2,1 1))  

右回り

SELECT ST_AsText(  
  'POLYGON((1 1, 1 2, 2 2, 2 1, 1 1))'::GEOMETRY  
);  
--------  
POLYGON((1 1,1 2,2 2,2 1,1 1))  

## MySQLの場合  

### 左回り  

```sql  
SELECT ST_AsText(  
  ST_GeomFromText('POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))')  
);  
--------  
POLYGON((1 1,2 1,2 2,1 2,1 1))  

右回り

SELECT ST_AsText(  
  ST_GeomFromText('POLYGON((1 1, 1 2, 2 2, 2 1, 1 1))')  
);  
--------  
POLYGON((1 1,1 2,2 2,2 1,1 1))  

SpatiaLiteの場合

左回り

SELECT ST_AsText(  
  ST_GeomFromText('POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))')  
);  
--------  
POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))  

右回り

SELECT ST_AsText(  
  ST_GeomFromText('POLYGON((1 1, 1 2, 2 2, 2 1, 1 1))')  
);  
--------  
POLYGON((1 1, 1 2, 2 2, 2 1, 1 1))  

このセクションでのまとめ

内部表現生成時に方向を強制するようにはなっていません。

全体のまとめ

右回り+左回りであっても、左回り+右回りであっても、MySQLとSQLServerは左回りになり、PostGISとSpatiaLiteは右回りになりました。

右回りの結果が出るか左回りの結果が出るかは、システム毎に違っていて、システム内では統一されている、と言えそうです。

また、左回り右回りの強制は、方向の強制を必要とする関数に渡った時点で強制されるものであり、内部表現生成時においては、まだ方向の強制はされていないようです。

ポリゴンを右回りにするか左回りにするかは、(たぶん)システム毎に決まってい(ると思いま)す(自信は無い)。

はじめに

RDBMS-GIS(MySQL,PostgreSQLなど) Advent Calendar 2018は、RDBMSの地理空間機能ならOKということで、でもたぶんMySQLとかPostGISとかが例に挙がると思いますが、まさかのT-SQLです。

ちょっと触ってみたので、ここにメモを残します。行くべきか引き返すべきかは各自判断して下さい。

思ったことを実行していく

ジオメトリの生成

ジオメトリの生成はgeometry::STGeomFromText()でやるのだそうです。

SELECT geometry::STGeomFromText('POINT(135 35)', 4326);  
go  

--  
0x  

(1 行処理されました)  

実行できたけど、ちゃんとできたかよく分からん。

STAsTextの実行

ST_AsText()でWKTが出るはずなのでやってみるけどその前にSTAsText()ですからね、注意して下さい。

SELECT STAsText(geometry::STGeomFromText('POINT(135 35)', 4326));  
go  
メッセージ 195、レベル 15、状態 10、サーバー ...\SQLEXPRESS、行 1  
'STAsText' は 組み込み関数名 として認識されません。  
1>  

無いんだそうだ…いや、そんなわけない。

SELECT geometry::STGeomFromText('POINT(135 35)', 4326).STAsText();  
go  

----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------  
POINT (135 35)                                                                                                                                                                                                                                          

(1 行処理されました)  

うーん、geometry型のメソッドという扱いですね。そういう設計なんですね(たじたじ)。

横着して全部小文字にしたら失敗した

前後するのですが、次のようなクエリを投げてエラーが返ってきました。

SELECT geometry::STGeomFromText('POINT(135 35)', 4326).stastext();  
go  
メッセージ 6506、レベル 16、状態 10、サーバー ...\SQLEXPRESS、行 1  
型 'Microsoft.SqlServer.Types.SqlGeometry' のメソッド 'stastext' がアセンブリ 'Microsoft.SqlServer.Types' に見つかりません でした  

問題はstastext()です。STAsText()としていないためのエラーです。大文字小文字の区別があるので気を付けてください。

UNION (結合)をやってみる

STUNion()は、geometryのメソッドです。集約関数にはUnionAggregate()があるようですが、試していません。

SELECT (  
  geometry::STGeomFromText('POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))',0)  
).STUnion(  
  geometry::STGeomFromText('POLYGON((1 0, 2 0, 2 1, 1 1, 1 0))',0)  
).STAsText();  
go  

----------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------  
POLYGON ((1 0, 2 0, 2 1, 2 2, 1 2, 1 1, 1 0))                                                                                                                                                                                                     

(1 行処理されました)  

後ろにメソッドを付けていくのが、PostGISユーザから見ると許しがたいけど、慣れるとjQueryとかやってるみたいで意外と気持ちいいです(個人の感想です)。

同じようなクエリをPostGISでも実行してみました。

db=# SELECT ST_AsText(  
  ST_Union(  
    'POLYGON((1 1, 2 1, 2 2, 1 2, 1 1))'::geometry,  
    'POLYGON((1 0, 2 0, 2 1, 1 1, 1 0))'::geometry  
  )  
);  

               st_astext                  
----------------------------------------  
 POLYGON((1 1,1 2,2 2,2 1,2 0,1 0,1 1))  
(1 行)  

うん、これは同じ結果になりましたね。

ジオグラフィで距離を計測する。

まず、ジオグラフィはgeographyクラスであることに注意して下さい。そしてSTGeomFromText()で生成することには、さらなる注意が必要です。ST_GeomFromText()でない、と。

距離計測は、geographyのメソッドにSTDistanceがあります。2回目なのでもう慣れましたね?

SELECT (  
  geography::STGeomFromText('POINT(135 35)', 4326)  
).STDistance(  
  geography::STGeomFromText('POINT(136 36)', 4326)  
);  
go  

------------------------  
      143321.57798743498  

(1 行処理されました)  

実行できました。

PostGISでは次のようになります。

db=# SELECT ST_Distance(  
  'SRID=4326;POINT(135 35)'::GEOGRAPHY,  
  'SRID=4326;POINT(136 36)'::GEOGRAPHY  
);  
   st_distance     
-----------------  
 143321.57818178  
(1 行)  

ま、大丈夫ですね。

あと、ここで、よいこのみんな、もとい、あやしいおともだちは、Axis Orderが気になりますよね?さっきのクエリでは経度緯度で通ってるけど、逆にしたらどうなるだろう、って。

SELECT (  
  geography::STGeomFromText('POINT(35 135)', 4326)  
).STDistance(  
  geography::STGeomFromText('POINT(36 136)', 4326)  
);  
go  

メッセージ 6522、レベル 16、状態 1、サーバー ...\SQLEXPRESS、行 1  
ユーザー定義のルーチンまたは集計 "geography" を実行中に .NET Framework エラーが発生しました:  
System.FormatException: 24201: 緯度の値は -90 ~ 90 度の範囲で指定する必要があります。  
System.FormatException:  
   場所 Microsoft.SqlServer.Types.GeographyValidator.ValidatePoint(Double x, Double y, Nullable`1 z, Nullable`1 m)  
   場所 Microsoft.SqlServer.Types.Validator.BeginFigure(Double x, Double y, Nullable`1 z, Nullable`1 m)  
   場所 Microsoft.SqlServer.Types.ForwardingGeoDataSink.BeginFigure(Double x, Double y, Nullable`1 z, Nullable`1 m)  
   場所 Microsoft.SqlServer.Types.CoordinateReversingGeoDataSink.BeginFigure(Double x, Double y, Nullable`1 z, Nullable`1 m)  
   場所 Microsoft.SqlServer.Types.WellKnownTextReader.ParsePointText(Boolean parseParentheses)  
   場所 Microsoft.SqlServer.Types.WellKnownTextReader.ParseTaggedText(OpenGisType type)  
   場所 Microsoft.SqlServer.Types.WellKnownTextReader.Read(OpenGisType type, Int32 srid)  
   場所 Microsoft.SqlServer.Types.SqlGeography.ParseText(OpenGisType type, SqlChars taggedText, Int32 srid)  
   場所 Microsoft.SqlServer.Types.SqlGeography.GeographyFromText(OpenGisType type, SqlChars taggedText, Int32 srid)  
。  

ということで、経度-緯度 (東西-南北)のオーダーです。

Tokyoは?

SELECT (  
  geography::STGeomFromText('POINT(135 35)', 4301)  
).STDistance(  
  geography::STGeomFromText('POINT(136 36)', 4301)  
);  
go  

------------------------  
      143305.61766686899  

(1 行処理されました)  

PostGISでは次のようになりました。

db=# SELECT ST_Distance('SRID=4301;POINT(135 35)'::GEOGRAPHY,  
'SRID=4301;POINT(136 36)'::GEOGRAPHY);  
   st_distance     
-----------------  
 143305.61785684  
(1 行)  

うん、たぶん、回転楕円体がちゃんと変更されてるっぽですね。

空間参照系テーブルはどこ?

若干さがしましたが、なんとか発見。sys.spatial_reference_systemsです。

SELECT * FROM sys.spatial_reference_systems;  

これは、sqlcmdでなくManagement Studioで実行することをお勧めします。エントリが「それなりに」多いですから。

これが空間参照系テーブルなのか?

見て頂くと分かると思いますが、4000番台だけです。

TSQLはSTTransform()を持っていないので、使用局面は、geography#STDistance()での準拠回転楕円体面等の情報を取得するだけでないかと思われます。だとすると、4000番台はWGS84地理座標系は4326、JGD2000地理座標系は4612といったように、地理座標系が集まるところなので、必要ないといえば必要ないです。

いや、そういうわけにはいかない。

JGD2011は6667ですね。4000番台じゃないので、これはまずいのではないか?

ドキドキしながら、次のクエリを実行します。

SELECT (  
  geography::STGeomFromText('POINT(135 35)', 6667)  
).STDistance(  
  geography::STGeomFromText('POINT(136 36)', 6667)  
);  
go  

メッセージ 6522、レベル 16、状態 1、サーバー ...\SQLEXPRESS、行 1  
ユーザー定義のルーチンまたは集計 "geography" を実行中に .NET Framework エラーが発生しました:  
System.ArgumentException: 24204: SRID (spatial reference identifier) が無効です。指定する SRID は、sys.spatial_reference_systems カタログ ビューに表示されたサポートされている SRID の 1 つと一致する必要があります。  
System.ArgumentException:  
   場所 Microsoft.SqlServer.Types.SqlGeography.set_Srid(Int32 value)  
   場所 Microsoft.SqlServer.Types.SqlGeography..ctor(GeoData g, Int32 srid)  
   場所 Microsoft.SqlServer.Types.SqlGeography.GeographyFromText(OpenGisType type, SqlChars taggedText, Int32 srid)  
。  

だめでした。

気付いた「クセ」

ここまでのをPostGISユーザから見た差異を上げてみます。

  • 関数名は大文字小文字の区別あり
  • geometry, geographyの静的メソッドでインスタンス生成
  • 関数名等のプリフィックスはST (ST_でない)
  • 第1引数にgeometry, geographyを取る関数は、インスタンスのメソッドとして実行
  • STTransform()が無い等機能が非常に少ない
  • 集約関数は静的メソッドとして別途用意 (PostGISのST_Union(geometry set)はTSQLではUnionAggregate())
  • Axis Orderは東西-南北の順で固定(MySQLとは違う)
  • 空間参照系テーブルはsys.spatial_reference_systems
  • SRIDは実質的にはgeographyのみ使用
  • spatial_refernce_systemsに登録されているのはSRID 4000番台のみ
  • JGD 2011地理座標系は非対応

おわりに

まあそういうことです。

本記事のライセンス


この記事は クリエイティブ・コモンズ 表示 4.0 国際 ライセンスの下に提供されています。