LoginSignup
52
43

More than 3 years have passed since last update.

フューチャーAdvent Calendar 2020Day 17future
@uhey22einフューチャー株式会社

地理情報・PostGISことはじめ

Last updated at Posted at 2020-12-16

フューチャー Advent Calendar 2020 17日目です。
昨日は @alfe_below さんによる「package.json dependencies メンテの仕方 最短ルート」でした。

はじめに

地理情報を扱うアプリケーションを開発するにあたって、基礎知識として頭に入れておいたほうが良さそうなことと、PostGISの簡単な使い方をまとめてみました。

Google Maps JavaScript APIやLeaflet.jsなど、フロントエンドのライブラリの使い方については触れていません。また、測量に関する細かい話も書いていません。

基礎知識

ここでは、地理情報のデータとしての表現方法や、システムの設計・開発にあたって覚えておいたほうが良さそうなことを簡単に説明します。

そもそも地理情報とは？ということについては、以下のリンク先が参考になるかと思います。

CRSとEPSGコード

(139.767, 35.681) という点があったとして、それが地球のどこを指し示すのかは原点や単位系といったパラメータの取り方によって異なるものになります。

このパラメータ群をまとめてCRS（空間参照系）と呼び、その組み合わせについてEPSG（European Petroleum Survey Group）という団体がユニークなIDを採番したものがEPSGコードです。よく使われているものの例としてEPSG:4326があり、GPSで使われています。

CRSは計算によって変換（射影）が可能で、各言語でライブラリが開発されています。

参考リンク

Well-Known Text (WKT)

駅の位置、道路の始点と終点、土地の形状など地理情報といった情報はベクトルによって表現されます。たとえば駅の位置であれば、2次元のベクトルが1個として表現できます。

Well-Known Textは、これらベクトルの表記方法です。例えば駅の位置（点）であれば、POINT(139.767 35.681) というように書くことが出来ます。

テキストではなくバイナリデータで表す場合、Well-Known Binary (WKB)という表記方法があります。

参考リンク

各種データフォーマット

地理情報オブジェクトを保存する際のフォーマットにはいくつか種類があります。

Shapefile

地理情報システムの老舗であるESRIが開発した形式です。仕様も公開されており、デファクトスタンダードと言って良いかと思います。

また政府発行の地理情報データも、だいたいShapefile形式での提供があります。

一方シェープ”ファイル”と言いながらも、拡張子違いで最低3つのファイルを同じディレクトリに配置して読み込み、という独特の扱いにくさがあります。

参考リンク

KML

Google Earth、Google Mapsで主に使われる形式です。XMLです。

参考リンク

GeoJSON

その名の通りJSON形式で地理情報データを表現します。その仕様はRFC7946で定められています。

Web Friendlyなフォーマットで、たとえばLeaflet.jsではGeoJSONオブジェクトを簡単にレンダリング出来るようになっています。
https://leafletjs.com/examples/geojson/

他にも多くのサービス、ツール類がGeoJSONに対応しています。

参考リンク

位置情報のハッシュ化

GeohashやQuadkeyといったアルゴリズムを使用すると、位置情報をハッシュ化することが可能です。

位置情報を扱う上で近傍検索というのはよくあるユースケースだと思いますが、これを緯度・経度から1件ずつ距離を計算していては計算コストが非常に高くなってしまいます。位置情報をハッシュ化しておくと、これが文字列の前方一致検索で可能になります。

参考リンク

標準地域メッシュ

国が地域を緯度・経度に基づいて分割したものです。日本ではその範囲の広さに応じて1次メッシュ、2次メッシュ、3次メッシュ…がJISで定められています。

緯度・経度から計算でメッシュコードが求められ、また近傍のメッシュは前方一致するので、ハッシュコードのような扱い方も可能です。

参考リンク

PostGIS

PostGISは、PostgreSQLで地理情報を扱えるようにするための拡張機能です。PostGISを有効化すると、いくつかのデータ型と関数が使えるようになります。

公式サイトで紹介されている次のクエリでは、ゴッサム・シティに所在するスーパーヒーローの名前を取得することが出来ます。

SELECT superhero.name
FROM city, superhero
WHERE ST_Contains(city.geom, superhero.geom)
AND city.name = 'Gotham';

ST_Contains() は、第1引数のオブジェクトが第2引数のオブジェクトを内包する場合にTrueを返します。はみ出していても被っていればTrueとしたい場合はST_WithIn() が使用可能です。他にもPostGISは多数の関数を提供しています。

GeometryとGeography

PostGISで扱えるデータ型に、GeometryとGeographyというふたつの型があります。これらの違いは、その座標をデカルト座標として扱うか否かです。

緯度と経度は角度ですので、Geographyとして扱うのが適切です。しかしGeometry型に緯度経度を入れることも可能で、用途次第では特に問題は生じず計算コストが低いというメリットがあります。

公式ドキュメントでは、「正直どっちがいいかわからないならGeographyがいいよ」とされています。

参考リンク

とりあえず使ってみる

何はともあれPostGISを動かして見るときはDockerを使うと便利です。
postgis/postgisというイメージがDockerHubに登録されているのでこれを使います。

$ docker run -ti --rm -p 15432:5432 -e POSTGRES_PASSWORD='postgres' postgis/postgis:11-2.5

このイメージの使用方法はpostgresとほぼ同じです。
コンテナが立ち上がれば、あとはお好みのクライアントから接続してください。

DBeaverというクライアントを使用すると、Geography型のデータをOpenStreetMap上でプレビューしてくれてとても便利です。

dbeaver.PNG

InsertとSelectの例

空港を保存するテーブルを例として考えます。ここではGeography型を使います。

create table airports ( id serial primary key , name text , iata_code text , location geography(POINT, 4326) )

このテーブルに羽田空港を登録するときは、次のようなINSERT文になります。

insert
    into
    airports (name, iata_code, location)
values ( '東京国際空港' , 'HND' , ST_GeogFromText('SRID=4326;POINT(139.7776446 35.5493975)') )

Geography(Geometry)型のカラムにデータを投入する際は、ST_GeogFromText(ST_GeomFromText)のような関数を使用します。この関数によって、EWKT形式のテキストからデータをInsertすることが出来ます。

羽田空港（HND）から最も近い空港を取得するクエリは次のようになります。

select
    airports.iata_code
    , airports.location
    , ST_Distance(
        airports.location, ref_airport.location, false
    ) as distance
from
    airports
    , (
        select
            location
        from
            airports
        where
            iata_code = 'HND'
        limit 1
    ) ref_airport
order by
    distance asc
limit 1

結果は次のようになります。Geography型を使っている場合、distanceの単位はメートルとなります。地理院地図で羽田空港と成田空港の直線距離を計測したところ約60kmだったので、正しく計算が出来ていると言えます。

iata_code|location                      |distance      |
---------|------------------------------|--------------|
NRT      |POINT (140.3906614 35.7719867)|60661.39893749|

空間インデックス

PostGISではGeometry / Geography型のカラムにインデックスを作成することが出来ます。先の羽田から近い空港を探すようなクエリであれば、locationカラムに空間インデックスを貼ることでパフォーマンスの向上が見込めます。

create index gix_airports_location on airports using gist(location)

Pythonとの連携

PostgreSQLを扱うpsycopg2、ジオメトリオブジェクトを扱うShapely、Pandasのようにジオメトリオブジェクトを含むデータを扱えるGeoPandasという組み合わせを紹介します。

Shapelyを使って点を定義

from shapely.geometry import Point

p = Point(139.7776446, 35.5493975)

psycopg2でデータベースに保存

import psycopg2

conn = psycopg2.connect("postgresql://...")

with conn:
    sql = """
    insert into airports (name, iata_code, location)
    values (%s, %s, ST_GeogFromText(%s))
    """
    with conn.cursor() as cur:
        cur.execute(sql, ("東京国際空港", "HND", "SRID=4326;" + p.wkt))

conn.close()

データベースからGeoPandasのDataFrameにデータをロード

import geopandas as gpd

sql = """
    select id, iata_code, location from airports
"""
conn = psycopg2.connect(_dsn)
gdf = gpd.read_postgis(sql, conn, geom_col="location")
_logger.info(gdf)

#    id iata_code                    location
# 0   1       HND  POINT (139.77764 35.54940)
# 1   2       NRT  POINT (140.39066 35.77199)

GeoPandas.DataFrameからGeoJSONへ変換

gdf.to_json(indent=2)

# {
#   "type": "FeatureCollection",
#   "features": [
#     {
#       "id": "0",
#       "type": "Feature",
#       "properties": {
#         "iata_code": "HND",
#         "id": 1
#       },
#       "geometry": {
#         "type": "Point",
#         "coordinates": [
#           139.7776446,
#           35.5493975
#         ]
#       }
#     },
#     ...
#   ]
# }

参考リンク

おわりに

地理情報データおよび測量はとっても奥が深い世界でした。アリストテレスの時代から知見が積み重ねられてきた分野なので当然といえば当然です。

この記事がこれから位置情報を扱うアプリケーションを開発する誰かの役に立てば幸いです。

52
43
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
Sign upLogin
52
43