Pythonでshapefileから逆ジオコーディング

はじめに

緯度経度から住所(市区町村・番地)を求める逆ジオコーディングですが、Rを使ったバージョンは以前書きました。

y-mattu.hatenablog.com

また、上記は市までわかっているバージョンだたのですが、そもそもの県や市を求めたいときは以下も参考になります。

qiita.com

今回は市までわかってるverをPythonでやってみます。

方針

やりかたは、Rのときと同じで「知りたい緯度経度がshapefileのポリゴンに入っているかを判定し、入っているところの番地を返す」です。shapefileは geopandas モジュールで読み込みます。pip で入ります。

コード

reverse_geo という関数を書きました。引数はshapeflieのデータフレーム、経度、緯度です。

import geopandas as gpd
from geopandas.geoseries import *

def reverse_geo(shdf, lon, lat):
    """
    逆ジオコーディング
    shdfは、geopandasでshapefileを読み込んだデータフレーム
    lon, latは経度、緯度
    """
    polygon = shdf['geometry']
    point = Point(lon, lat)
    
    whichtrue = polygon.contains(point)
    whichrow = whichtrue[whichtrue == True].index

    geos = shdf.ix[whichrow]
    
    if geos.empty == True:
        return(None)

    reslist = [geos['KEN_NAME'].to_string(index=False),
               geos['GST_NAME'].to_string(index=False),
               geos['CSS_NAME'].to_string(index=False),
               geos['MOJI'].to_string(index=False)]

    res = ','.join(reslist)
    
    return(res)

やってみる

千葉県浦安市のshapefileをこちらからダウンロードします。今回はデスクトップに解凍するものとします。

shdf2 = gpd.read_file('~/Desktop/A002005212015DDSWC12227/h27ka12227.shp')

reverse_geo(shdf2, 139.884678, 35.626065) # 東京ディズニーシーのある地点
#=> '千葉県,浦安市,,舞浜'

カンマで分けているのは区切りを分けられるようにするためです。これであれば、データフレームに適用したときに後から「県、市、区、番地の列に分ける」みたいなことも簡単にできます。