GIS奮闘記

元GISエンジニアの技術紹介ブログ。主にPythonを使用。

スポンサーリンク

【Python と GIS】 GDAL を使って Shapeファイルを読み込んで属性検索してみよう!

さて、本日は GDAL について書いてみようと思います。GDAL については以下のように以前に何個か記事を書いているのですが、紹介していない機能がまだまだたくさん残っています。今回は属性検索について紹介します。

www.gis-py.com

www.gis-py.com

www.gis-py.com

www.gis-py.com

GDALとは

GDAL はフリーの GIS ライブラリです。GDAL — GDAL documentation
これがあればShapeの読み込みからジオメトリの作成まで基本的なことは何でもできるかと思います。

インストール

公式サイトの手順に従ってインストールをお願いします。https://gdal.org/download.html

サンプルデータ

以下エントリーで紹介した全国市区町村界データを使用します。

www.gis-py.com

サンプルデータの属性情報

フィールド名

  • JCODE・・・市区町村コード
  • KEN・・・都道府県名
  • SICHO ・・・支庁名・振興局名
  • GUN・・・郡名(町村部のみ)
  • SEIREI・・・政令指定都市の市名
  • SIKUCHOSON・・・市区町村名
  • CITY_ENG・・・市区町村名(英語)
  • P_NUM・・・人口
  • H_NUM・・・世帯数

実行環境

Windows10(64bit)
Python3.6.5
Jupyter Notebook

GDAL で属性検索をするには

以下のような手順を踏む必要があります。

  1. Shape ファイルの読み込み
  2. レイヤの取得
  3. フィルター設定
  4. 属性の取得

サンプルコード

KEN = '神奈川県'のフィーチャを属性検索するサンプルです。

from osgeo import ogr

daShapefile = r"D:\data\japan_ver81.shp"
dataSource = ogr.Open(daShapefile)
layer = dataSource.GetLayer(0)
layerDefinition = layer.GetLayerDefn()
layer.SetAttributeFilter("KEN = '神奈川県'")

for feature in layer:
    attribute_list = [feature.GetField(layerDefinition.GetFieldDefn(column).GetName()) for column in range(layerDefinition.GetFieldCount())]
    print(attribute_list)


結果を見ると、想定通り神奈川県のデータを取得できたことがわかりました。
f:id:sanvarie:20190713143258p:plain

処理のポイント

サンプルコードの処理を重要なポイントごとに分けると以下のようになります。

Shapeファイルの読み込み
daShapefile = r"D:\data\japan_ver81.shp"
dataSource = ogr.Open(daShapefile)
レイヤーの取得
layer = dataSource.GetLayer(0)
フィルター設定
layer.SetAttributeFilter("KEN = '神奈川県'")
属性の取得
attribute_list = [feature.GetField(layerDefinition.GetFieldDefn(column).GetName()) for column in range(layerDefinition.GetFieldCount())]


ここが少し複雑なので、この部分の処理を分解して解説します。

  • レイヤーのフィールド数の取得
for column in range(layerDefinition.GetFieldCount())
  • フィールド名の取得

上記の「レイヤーのフィールド数の取得」で取得した要素を引数にします。

layerDefinition.GetFieldDefn(column).GetName()
  • 属性の取得
feature.GetField("フィールド名")

でそのフィールドの属性を取得することができます。

なので、上記の「フィールド名の取得」で取得した要素、つまり、フィールド名を引数にしています。

feature.GetField(layerDefinition.GetFieldDefn(column).GetName())

少し長いのですが、上記の処理をリスト内包表記にしたものをサンプルコードに載せました。

本日は以上になりますが、まだまだ紹介しきれていない使い方がたくさんあるので、 また GDAL について書いてみようと思います。