単純なデータ処理をやってみる

農業・食品産業技術総合研究機構近畿中国四国農業研究センター
寺元郁博

はじめに

PostGIS では、幾何関数やジオメトリ型を、PostgreSQL のクエリの中に混ぜ込むことができます。 WHERE 句に幾何関数を使ったり、INSERT, UPDATE 等で使用することができます。ビューで扱うこともできます。

今回はごく簡単なデータ処理の例を示してみます。広島県内の市区町村ごとの駅の数を数え上げて、多いところ少ないところを色分けして表示します。数え上げまでを PostGIS で行います。

本例での前提は次の通りです。

データベース名は pgis
格納するデータは2種類
- 国土数値情報 (行政区域) のうち広島県
- 国土数値情報 (鉄道)
テーブルは3つ
- admarea
- station
- rail (使用しません)
全テーブルのジオメトリカラム名は the_geom (shp2pgsql のデフォルト)

鉄道データの取得と格納

国土数値情報ダウンロードサービス (JPGIS準拠) に行き、「鉄道（線）」をダウンロードします。このデータは、全国でひとつのファイルになっています。年度を選択してダウンロード、unzip して下さい。

行政区域と同じように、鉄道データもシェープに変換します。

N02-08_EB02.shp (路線) と N02-08_EB03.shp (駅) の2つのシェープファイルが生成されます。前者は今回は使いませんが、せっかくなので両方とも格納してしまいましょう。 (注:掲載の都合上、1行の処理を改行して記述しています)

% shp2pgsql -s 4612 -D -i -I -W cp932 N02-08_EB02.shp rail \
   > N02-08_EB02.sql
% psql pgis -f N02-08_EB02.sql

% shp2pgsql -s 4612 -D -i -I -W cp932 N02-08_EB03.shp station \
   > N02-08_EB03.sql
% psql pgis -f N02-08_EB03.sql

これで、rail, station という2つのテーブルが新たに作られました。

さっそくビューを作ってみる

CREATE VIEW v_stncnt_1 AS
  SELECT
    gid,
    (SELECT COUNT(*) FROM station WHERE
     ST_WithIn(station.the_geom,admarea.the_geom))::int4 AS stncnt,
    the_geom FROM admarea;

3行目 (gid カラム) は shp2pgsql で自動でついてきたものです。 QGIS で表示する場合にはユニークなIDカラムが必要とされていますので、これを用いています。

4行目と5行目を抜き出すと、次のようになります。

SELECT COUNT(*) FROM station 
 WHERE ST_WithIn(station.the_geom,admarea.the_geom))::int4 AS stncnt,

ST_WithIn(A,B)は、「AがBに含まれるかどうか」を判定する関数です。この場合、station.the_geom (駅) が、admarea.the_geom (市区町村ポリゴン) に含まれるか、が WHERE 条件となっています。つまり、ここでは、station テーブルから市区町村ポリゴンに含まれるものを抜き出してカウントする、ということを行っています。

4行目と5行目を考慮に入れない場合は、次のようになります。

SELECT gid, the_geom FROM admarea;

これは、単にadmarea (国土数値情報 (行政区域) テーブル) の内容を表示するに過ぎません。

このビューの全体を見直すと、admarea テーブルの全てのレコードについて、そのポリゴンに含まれる駅の数を集計した結果とポリゴン自体との組み合わせを返すビューであることが分かります。

駅の中心点を追加する

今回の扱うサンプルデータでは、このままではちょっと困った問題が発生します。それは、station の中に市境を跨いでいるものがあることです。証拠写真を示します。赤線が駅を、黒線が市境を、それぞれ示します。

ST_WithIn(A,B) は、A が B に完全に含まれている必要があります。市境をまたぐ場合には、どちらの市でも、ST_WithIn() が TRUE になりません。この駅はどの市にも属さないことになります。

そこで、今回のお題で使用するデータは、むりやりどちらかの市に所属させるため、線データから点データに加工して行うことにします。

新しいジオメトリカラムを追加する

SELECT AddGeometryColumn('', 'station', 'c_geom', 4612, 'POINT', 2);

第1引数から第3引数は、順に対象スキーマ名 (空文字列の場合デフォルト名)、テーブル名、追加するカラム名を意味します。
第4引数の 4612 はおまじないと考えて下さい。
第5引数はポイント型データとすることを意味します。線データなら 'LINESTRING', ポリゴンデータなら 'POLYGON' とします。
第6引数は次元数で2次元 (X,Y) とすることを意味します。

新しいカラムができたので、これに空間インデックスを貼ります。

CREATE INDEX c_geom USING GiST;

以上で、新しいカラムの追加ができます。 CREATE TABLE でジオメトリカラム以外のカラムを作成した後に、このやり方でジオメトリカラムを追加します。

新しいジオメトリカラムに幾何計算結果を放り込む

元データでは駅を線で表現しているのですが、駅ごとに、その駅を表現する線から重心を計算して、これをその駅の点データとするようにします。重心を計算するには、関数 ST_Centroid() を用います。

UPDATE station c_geom = ST_Centroid(the_geom);

なお、お気づきのことと思いますが、UPDATE などでも、ジオメトリデータは普通に使えます。 INSERT でも普通のデータと同じように扱えます。

新しく作ったジオメトリカラム c_geom を QGIS で見ると、次のようになります。

改めて表示してみる

今度は、点データを使ったビューを作成します。

CREATE VIEW v_stncnt AS
    SELECT
        gid,
        (
         (SELECT count(*) AS count 
          FROM station 
          WHERE ST_WithIn(station.c_geom, admarea.the_geom)
         )
        )::integer AS stncnt,
        the_geom
    FROM admarea;

このビューを表示してみましょう。

QGISの「レイヤー」メニューから「PostGISレイヤの追加」を選択します。
「pgis」に接続して、「テーブル」が「v_stncnt」になっているアイテムを押さえて「追加」をクリックします。
とりあえず出ました。

次に、stncnt(第2カラム)の値を色の濃さで表示してみましょう。

QGISの「レイヤ」ペイン(左端にあります)に「v_stncnt」レイヤの上でマウス右ボタンでサブメニューを出して「プロパティ」を選択します。
「シンボル」カテゴリを選択して、「凡例タイプ」を「連続色」にします。
「分類フィールド」を「stncnt」として、最小値色、最大値色を任意に設定します。ここでは、白と赤にしました。
できあがりです。

あと、せっかくなので、駅データと重ね合わせてみましょう。

QGISの「レイヤー」メニューから「PostGISレイヤの追加」を選択します。
「pgis」に接続して、「テーブル」が「station」、「ジオメトリカラム」が「c_geom」になっているアイテムを押さえて「追加」をクリックします。
駅もあわせて表示した結果です。赤の濃い市区町村ほど、確かに駅が多いのが分かります。

本当ならビューより別テーブルに叩き込む

PostGIS を使用する際に期待されているのは、大量なデータのハンドリングです。ここまで書いておいて言うのも何ですが、ビューで見るよりも、新しいテーブルに叩き込んだ方が良い場合が多いです。

上に示した例では広島県内のみでしたが、日本全国に拡大させるとなると、数十倍の計算コストがかかります。移動・拡大縮小等で表示範囲を変更するたびにクエリを実行しようとするので、時間ばかりかかることになります。日本全国を扱うなら、INSERT INTO などで別テーブルに叩き込んでしまう方が、本当は良いでしょう。

最後に

PostGIS に関する、ごくごく簡単な解説をやってみました。「おまじない」とだけ述べた部分をはじめ、もっと空間データを扱うために必要なところが多くあると思います。でも、まずはなにより、フリーでここまでできるということや、空間データは決して遠くにあるものではないということを、少しでも感じていただけたなら幸いです。

(2010年7月9日公開)

PostGIS 単純なデータ処理をやってみる