PostGISで環境省植生図をマージするの巻 - 日々是戯言・Geo出張版・はてなブログ分室

ということで，なんとかできたのでまとめてみます。
対象とするのは関東地方の植生図，1次メッシュコードでいうと，5239，5240，5339，5340，5439，5440で，ファイル数は現状で119。動作環境はWindowsXP，PostgreSQLは8.4.0-1，PostGISは1.4.0-2です。

まず，ogr2ogrを使って処理対象となるファイルを一つのshpファイルにまとめてしまう*1。

つぎに，空間データベースを作る。PostGISをインストール時にも作ってくれるが，コマンドで入力するとすると以下の通り。

createdb -U postgres -T template_postgis -E UTF-8 kanto_veg_db

ここで，-Uは作成するデータベースの所有者を指定するオプション（その後ろがユーザー名），-Tがテンプレート使用に関する指定（ここでは，template_postgisを使用），-Eはデータベースのエンコード指定（UTF-8を指定），そして最後が作成するデータベース名である。

次に，shpファイルをデータベースに投入するためのsqlファイルを作る。まず，

shp2pgsql -p -i -s 32654 -W CP932 -D -I veg_wgs84_utm54.shp kanto_veg > kanto_veg_p.sql

といれる。pはプリペアモードのことで，データ投入は行わず，テーブルのヘッダーだけ作る。-I は空間インデックスを作るというやつ。ちなみに，ここに出てくるkanto_vegは，データベース名ではなく，テーブル名である。
続けて，データ投入用のsqlをつくる。

shp2pqsql -a -i -s 32654 -W CP932 -D veg_wgs84_utm54.shp kanto_veg > kanto_veg_a.sql

ここで，-pの代わりに-aを使う。これは，データ部分だけを投入するということ。さらに，前もって空間インデックスは作られているので，-Iは不要みたいである。つか，付けたらエラーが出た。

で，このsql文を実行する。まず，テーブルを作る

psql -U postgres -h localhost -f kanto_veg_p.sql kanto_veg_db

ここで，最後のkanto_veg_dbはデータベース名である。次にデータ投入用sqlを実行

psql -U postgres -h localhost -f kanto_veg_a.sql kanto_veg_db

つぎ。これに対してST_Bufferをかけるわけだが，その前に投入先のテーブルを準備する必要がある。
ちなみに，ここまでの作業はコマンドラインでやってたけど，以下の二つはpgAdminIIIのSQLエディタでやってます。

まず，前出のkanto_veg_p.sqlを修正した以下のSQLを実行してテーブルを作成。

SET CLIENT_ENCODING TO UTF8;
SET STANDARD_CONFORMING_STRINGS TO ON;
BEGIN;
CREATE TABLE "kanto_veg_union" (gid serial PRIMARY KEY,
"hanrei_c" varchar(6),
"hanrei_n" varchar(254));
SELECT AddGeometryColumn('','kanto_veg_union','the_geom','32654','MULTIPOLYGON',2);
CREATE INDEX "kanto_veg_union_the_geom_gist" ON "kanto_veg_union" using gist ("the_geom" gist_geometry_ops);
END;

続けて，以下のSQLを実行

INSERT INTO kanto_veg_union (the_geom,hanrei_c,hanrei_n)
SELECT ST_Multi(ST_Buffer(ST_Collect(the_geom), 0) ) AS the_geom,hanrei_c,hanrei_n FROM kanto_veg
~~SELECT setsrid(multi(ST_Buffer(ST_Collect(the_geom), 0) ), 32654) AS the_geom,hanrei_c,hanrei_n FROM kanto_veg~~*2
GROUP BY hanrei_c,hanrei_n

えっとすみません，細かい説明はﾐﾘ(汗
ということでやってみると，１８３０秒，つまり30分ほどで完了。
とにもかくにも，こうしてマルチポリゴンで，ディゾルブがかかった植生図のできあがり！！*3
んで最後にこいつをシングルポリゴン化する。ます，投入用テーブルを作るためのsqlの作成。ここは，またコマンドラインからの操作。

shp2pgsql -p -i -S -s 32654 -D -I veg_wgs84_utm54.shp kanto_veg_union2 > kanto_veg_union2.sql

ここで-Sがあるが，シングルポリゴンでデータベースを作成しろということ。それから，kanto_veg_union2.sqlをテキストエディタで開いて，不要な列を削除する（使うのは，hanrei_cとhanrei_nのみ）。次に，データベース内に投入用テーブルを作成。コマンドは以下の通り。

psql -U postgres -h localhost -f kanto_veg_union2.sql kanto_veg_db

で，こいつにkanto_veg_unionをシングルポリゴン化して投入してやる。これは，pgAdminIIIから実行です。

INSERT INTO kanto_veg_union (the_geom,hanrei_c,hanrei_n)
SELECT (ST_Dump(the_geom)).geom AS the_geom,hanrei_c,hanrei_n FROM kanto_veg_union

このsqlの実行にかかった時間は約100秒でした。

以上で，シングルポリゴンのファイルができあがりだ！！
しかし，実は微妙に境界線が残ったりしてるけど*4，その辺は後で考えるよう(汗

はー。
つか，本当の作業はここから始まるわけですが(汗

*1:これに関しては，shp2pgsqlの-pと-sオプションを上手く使ってやるとやらなくてもすむかも知れないが，今回は先にまとめてあったので，この手順で行きます。

*2:20090926修正

*3:ちなみに，ここでST_Dumpをかけてしまうことも可能なんだけど，そうしたらポリゴンが落ちてしまったのは，id:wata909:20090919#p1の通りです。

*4:ちなみに，最後のST_Dumpの前にST_Unionをかけてみても残ってしまった。謎である。