2016-09-03

GISポリゴンデータに基づく面積計算 (_AZMEA_の利用)

GISポリゴンデータに基づいて地物の面積を求めることはよくあると思います。GISソフトウェア等の機能として面積計算はありふれたものですが、データの精度は元の地図の縮尺やデジタイズの方法等によって異なるとともに、同じデータであっても投影先の座標系によって地物の形状は異なるため、使用するデータや投影方法によって、同一の地物について求められる面積が異なることもあるということは認識しておく必要があるでしょう。

既存のGISデータから求めた面積を研究論文・報文や調査業務報告書などで用いるのならば、データの出典を明記するのはもちろんのこと、元のデータの座標系が緯度経度である場合には、どの座標系に投影した形状の面積なのかも説明するべきであると考えます。

緯度経度の座標系で作成された日本国内のGISデータに基づいて地物の面積を求める場合、平面直角座標系またはUTM座標系に変換することが多いかも知れません。

しかし、平面直角座標系の複数の系、UTM座標系の複数のゾーンにまたがった広い範囲に分布する地物の面積を求める場合には、個々の地物が所在する位置に応じて変換先座標系の系番号、ゾーン番号を選択するといったプロセスが必要になります。そのようなプロセスを含むワークスペースを作成することもできますが、合理的な説明ができれば面積計算にあたっての投影法は任意で良いという場合には、_AZMEA_ (Dynamic Reprojection Equal Area) という名前の座標系を利用するのが便利です。

_AZMEA_ は、任意の座標を基準点とする「ランベルト正積方位図法」 (Lambert azimuthal equal-area projection) によって定義される座標系 (単位: メートル) であり、座標系変換を行う Reprojector または CsmapReprojector トランスフォーマーで変換先の座標系として _AZMEA_ を指定すると、ワークスペースの実行時に、個々の地物について元の座標系におけるバウンディングボックス中心を基準点とするランベルト正積方位図法の座標系が定義され、それに変換されます。これにより、系番号やゾーン番号を考慮しなくても、広い範囲に分布する地物の面積を統一的な手法で求めることができます。

ここでは、国土数値情報「行政区域」データに基づいて、_AZMEA_ を利用して全国の市区町村の面積を求めるワークスペース例を掲げます。

【ソースデータ】
国土数値情報ダウンロードサービスサイトからダウンロードした「行政区域」データ(平成28年1月1日時点, 全国, Esri Shapefile 形式ポリゴンデータ)

FME 2016.1.1.0 build 16609

FMEワークスペース例 1


















[ESRISHAPE] リーダー: 国土数値情報行政区域(全国)データを読み込む。
Tester: 市区町村コード(N03_007属性値)を持っているフィーチャーとそうでないものに分岐する。
Aggregator: 同じ市区町村コードを持っているフィーチャーを集約する。
Reprojector: フィーチャーごとに _AZMEA_ に座標系を変換する。
AreaCalculator: フィーチャーごとに面積を計算する。
Inspector: 変換結果を FME Data Inspector で確認する。

国土数値情報の行政区域データ(ポリゴン)は、飛び地や離島などがある市区町村について、それらが個別のレコードとなっている、いわゆるシングルパートポリゴンです。また、埋立地などでどの市区町村に属するか決定していない「所属未定地」のレコードには、市区町村コードが記録されていません。

そのため、このワークスペース例ではまず、Tester によって市区町村コードを持っているフィーチャーと持っていないフィーチャー(所属未定地)に分け、市区町村コードを持っているフィーチャーについて Aggregator によって市区町村ごとに集約しました。これによって、飛び地や離島がある市区町村についても、市区町村単位で全域の面積を求めることができます。

Reprojector の Destination Coordinate System (変換先の座標系)パラメーターに _AZMEA_ を設定することにより、前述のとおり、フィーチャー(市区町村または所属未定地)ごとに、それぞれバウンディングボックス中心を基準点とする「ランベルト正積方位図法」による座標系に変換されます。

Reprojector パラメーター設定画面

















AreaCalculator は、入力フィーチャーの面積を計算し、その値を格納する新たな属性を付加します。ここで、Multiplier パラメーターによって、求められた値に乗じる値を指定することができます。_AZMEA_ によって変換された座標値の単位はメートルなので、Multiplier パラメーターに 1e-6 (= 0.000001) を設定すれば、平方キロメートル単位の値に換算されます。


このワークスペースを実行して Data Inspector に変換後のポリゴンを表示させると分かりますが、_AZMEA_ は、個々のフィーチャーのバウンディングボックス中心を基準点、つまり原点 (0, 0) とする座標系なので、変換後のポリゴンは原点周辺で重なり合うことになります。地図としては用を成しません。

面積を求めることだけが目的ならばこのままでも構いません。しかし、元の座標系(この場合は JGD2011 緯度経度)におけるジオメトリとしての処理も必要な場合には、元の座標系における位置、形状を復元する必要があります。

もうひとつ Reprojector を追加して元の座標系に再変換するのは簡単ですが、座標系変換にはどうしても計算誤差が伴い、完全に復元できるとは限りません。正確に元に戻すには、あらかじめ元の座標系名とジオメトリを属性データとして抽出しておき、面積を求めた後で、それらの属性データに基づいて元のジオメトリと座標系の定義を復元するという方法があります。

上記ワークスペースにそれらの手順を追加すると、次のようになります。

FME ワークスペース例 2


















(座標系変換前)
CoordinateSystemExtractor: フィーチャーに設定されている座標系名を属性 (_coordsys) として抽出する。
GeometryExtractor: ジオメトリ(形状のデータ)を属性 (_geometry) として抽出する。

(面積計算後)
GeometryReplecer: 属性 (_geometry) に基づいてジオメトリを復元する。
CoordinateSystemSetter: 属性 (_coordsys) に基づいて元の座標系の定義を復元する。

GeometryExtractor/Replacer では、Geometry Encoding パラメーターでいくつかの属性値(文字列)としてのジオメトリの表現方法が選択できますが、そのうち FME Binary は、ジオメトリの抽出/復元に伴う誤差は全くなく、かつ、最も効率が良いので、この例のようにジオメトリを復元することを目的とする場合には最も適しています。

GeometryReplacer は、フィーチャーのジオメトリ(形状のデータ)だけを置き換えるものであり、座標系の定義は復元しないことに注意してください。そのため、このワークスペース例では、CoordinateSystemExtractor によって元の座標系名をあらかじめ抽出しておき、CoordinateSystemSetter によってその座標系の定義を復元しています。

なお、GeometryExtractor/Replacer は、FME Binary のほか、GeoJSON, GML, OGC Well Known Binary, OGC Well Known Text などのメジャーなジオメトリの文字列表現をサポートしており、用途に応じて使い分けることができます。

0 件のコメント:

コメントを投稿