国土数値情報から行政区界を作成する


この記事とだいたい一緒です。

* [国土交通省のデータを使って、緯度経度から市区町村までを取り出す](http://qiita.com/masuidrive/items/21a282c7bf54fd6a4985)

ここでは、住所コードからその市区町村の形状を描いてみます。流れとしては以下のようになります。

1. [国土数値情報ダウンロードサービス](http://nlftp.mlit.go.jp/ksj/gml/datalist/KsjTmplt-N03.html)から都道府県ごとのZIPファイルをダウンロード
2. ZIPを展開し、Shape形式のデータをPostGISに登録
3. PostGISを使って形状データを市区町村ごとに結合
4. GeoJSONで出力して描画

まずは平成25年度のデータをダウンロードします。平成26年度のデータが公開されたらそちらを使った方が良いでしょう。クリックを繰り返すだけの単調作業なので、時間があるときに実施しておきます。47都道府県のZIPファイルを、どこか適当なディレクトリにまとめて保管します。

## EC2インスタンスを起動

AWS の EC2 で作業します。インスタンスのイメージIDは [Amazon Linux AMI](http://aws.amazon.com/jp/amazon-linux-ami/) のページに記載されています。その時の一番新しいタイプを使うと便利なことは多いですが、ライブラリのバージョン依存に悩む場合は、ちょっと古いタイプを試してみましょう。`--user-data` を使って自動アップデートを止める方法は下記の記事が参考になります。

* [Amazon Linux AMI 2014.03 Released! 早速ハマったこと (Perl・glibc編)](http://dev.classmethod.jp/cloud/aws/amazon-linux-ami-2014-03-trap/)

PostGIS を `yum` でインストールする場合、EPEL の *gdal-libs* に依存します。
しかし、2014.03 のレポジトリだと EPEL の *gdal-libs* が依存する *poppler* のバージョンが合いません。 *gdal-libs* は *poppler-0.12* に依存しますが、最新のインスタンスには *poppler-0.22.5* がインストールされています。ひとつ古い 2013.09 は *poppler-0.12.4* がインストールされていますので、ソースコードからビルドするのが煩雑な場合はこちらを使うのが良いでしょう。

ということで、以下の構成で [awscli](https://pypi.python.org/pypi/awscli) のコマンドラインツールを使って EC2 を起動します。パッケージレポジトリのバージョンを固定しておきたいので、`--user-data` で起動時に自動アップグレードしないよう設定します。


| 種別 | 内容 |
|----:|:----|
| AMI | 2013.09 の 64bit |
| リージョン | アジアパシフィック東京 |
| イメージID | ami-3561fe34 |
| インスタンスタイプ | t1.micro |
| SSH鍵 | aws-tokyo-sandbox |
| セキュリティグループ | SSHEnabledGroup |


```bash:EC2インスタンスの起動
aws ec2 run-instances \
    --image-id ami-3561fe34 \
    --instance-type t1.micro \
    --count 1 \
    --key-name ec2-tokyo-sandbox \
    --security-groups SSHEnabledGroup \
    --user-data `cat << EOF | base64
#cloud-config
repo_upgrade: none
EOF`
```

SSH でログインしたらレポジトリの設定を変更します。

```bash:レポジトリ設定
curl -O http://yum.postgresql.org/9.3/redhat/rhel-6-x86_64/pgdg-redhat93-9.3-1.noarch.rpm
sudo rpm -ivh pgdg-redhat93-9.3-1.noarch.rpm
sudo vim /etc/yum.conf  # releasever=latest をコメントアウト
sudo vim /etc/yum.repos.d/pgdg-93-redhat.repo  # $releasever ⇒ latest
```

次に、PostGIS をインストールします。インストールが成功すると `shp2pgsql` コマンドが使えるようになっているはずです。

```bash:PostGISのインストール
sudo yum install postgresql93-server postgresql93-contrib -y
sudo yum install --enablerepo epel postgis2_93 postgis2_93_devel postgis2_93_utils -y
```

PostGIS がインストールできたら PostgreSQL を起動します。
作業を簡単に進めるために、アクセス権限を変更しておきましょう。SSHでログインしている場合はユーザー名を指定するだけにします。


```bash:データベース設定
sudo service postgresql-9.3 initdb
sudo service postgresql-9.3 start
sudo su - postgres
vi /var/lib/pgsql/9.3/data/pg_hba.conf    # 下記参照
```

```conf:/var/lib/pgsql/9.3/data/pg_hba.conf
local        all        all        trust
```

設定ファイルを変更したら PostgreSQL を再起動します。

```bash
sudo service postgresql-9.3 restart
```

データベースを作成し、PostGIS 拡張を有効にします。

```bash:データベース作成
createdb -U postgres geo
cat << EOF | psql -U postgres -d geo
CREATE EXTENSION postgis;
CREATE EXTENSION postgis_topology;
CREATE EXTENSION fuzzystrmatch;
CREATE EXTENSION postgis_tiger_geocoder;
SELECT postgis_full_version();
EOF
```

これでデータベースの設定ができました。

## ZIPファイルの中身をPostGISに登録

国土数値情報ダウンロードサービスからダウンロードしたZIPファイルをEC2のインスタンスに転送します。転送できたら、ZIPファイルを展開して `shp2pgsql` でSQLファイルを作成し、PostGISに登録します。

```bash
for f in N03-*.zip
do
    unzip -d ${f%.zip} $f
done

for f in N03-*/*.shp
do
    t=`basename $f`
    sql=${t%.shp}.sql
    shp2pgsql -W "Shift_JIS" $f | psql -U postgres -d geo
done
```

ファイル名に応じたテーブルが生成されますので、47個のテーブルが生成されます。

データを見てみると、いくつかの都道府県のデータには注意が必要であることが分かります。 ([北海道地方の支庁・地方事務所一覧](http://ja.wikipedia.org/wiki/%E6%94%AF%E5%BA%81#.E5.8C.97.E6.B5.B7.E9.81.93.E5.9C.B0.E6.96.B9) - Wikipedia)

```sql:北海道の支庁と振興局が混在している
SELECT distinct n03_002
FROM "n03-13_01_130401"
WHERE n03_002 IS NOT NULL
ORDER BY n03_002;

    n03_002
----------------
 上川支庁
 十勝支庁
 宗谷支庁
 後志支庁
 日高支庁
 根室振興局
 根室支庁
 檜山支庁
 渡島支庁
 留萌支庁
 石狩支庁
 空知支庁
 網走支庁
 胆振支庁
 釧路支庁
 釧路総合振興局
(16 rows)
```

```sql:香川県の郡・政令市名が香川県のレコードがある
\x

SELECT n03_001, n03_002, n03_003, n03_004, n03_007
FROM "n03-13_37_130401"
WHERE n03_001 = n03_003;

-[ RECORD 1 ]---
n03_001 | 香川県
n03_002 |
n03_003 | 香川県
n03_004 | 丸亀市
n03_007 | 37202
```

```sql:所属未定地の確認
SELECT n03_001, n03_002, n03_003, n03_004, n03_007 FROM "n03-13_13_130401" WHERE n03_007 IS NULL;

-- 結果は省略
```


たとえば、住所コードで集約したときに名称が不一致になるレコードを抽出すると、以下のレコードをリストアップできます。名寄せが必要なレコードには _UPDATE_ 文を発行しておきましょう。

| n03_001 | n03_002 | n03_003 | n03_004 | n03_007 | レコード数 |
|:--------|:--------|:--------|:--------|:--------|---------:|
 北海道  | 根室振興局     |         | 根室市  | 01223   |  76
 北海道  |                |         | 根室市  | 01223   | 464
 北海道  | 釧路支庁       | 厚岸郡  | 浜中町  | 01663   | 139
 北海道  | 釧路総合振興局 | 厚岸郡  | 浜中町  | 01663   |   1
 三重県  |         | 北牟婁郡 | 紀北町       | 24543   | 416
 三重県  |         | 北牟婁郡 | 紀北町入会地 | 24543   | 110
 香川県  |         | 香川県  | 丸亀市  | 37202   |   1
 香川県  |         |         | 丸亀市  | 37202   |  14

市区町村には同一名称もありますので、集約させるキーには注意しましょう。また、テーブル定義で *UNIQUE* 制約を付ける場合は慎重に検討しましょう。

* [同一名称の市区町村一覧](http://ja.wikipedia.org/wiki/%E5%90%8C%E4%B8%80%E5%90%8D%E7%A7%B0%E3%81%AE%E5%B8%82%E5%8C%BA%E7%94%BA%E6%9D%91%E4%B8%80%E8%A6%A7) - Wikipedia

## PostGISを使って形状データを市区町村ごとに結合

都道府県ごとのテーブルにあるレコードを住所コードで集約します。PostGIS では複数の形状データを `ST_Union` で結合できます。
今回の場合は、住所コード (n03_007 フィールド) で _GROUP BY_ し、 _geom_ フィールドを結合させます。変換先のテーブルは以下のように定義します。

```sql:変換先のテーブル定義
CREATE TABLE address (
    code varchar(5) PRIMARY KEY CHECK (length(code) IN (2, 5)),
    geom geometry(MULTIPOLYGON, 4326) NOT NULL
);
```

変換クエリはSQLファイルに保存し、 `psql` コマンドの "-v" オプションを使って変数をループさせます。都道府県ごとに別々のテーブルになっているためです。"table" にはテーブル名、"code" には都道府県コードを指定します。

```sql:transform.sql
INSERT INTO address
SELECT n03_007, ST_Union(ST_SetSRID(geom, 4326))
FROM :"table"
WHERE n03_007 IS NOT NULL
GROUP BY n03_007;

INSERT INTO address
SELECT :'code', ST_Union(ST_SetSRID(geom, 4326))
FROM :"table"
WHERE n03_007 IS NULL
GROUP BY n03_001;
```

```bash:都道府県ごとにクエリを実行
for i in {01..47}
do
    psql -U postgres -d geo -f transform.sql -v table="n03-13_${i}_130401" -v code=$i
done
```

## GeoJSONで出力して描画

PostGIS の `ST_AsGeoJSON` を使うと形状データを GeoJSON で出力できます。
市区町村ごとにファイルを作成してみましょう。

```bash
echo "SELECT code, ST_AsGeoJSON(geom) FROM address" | psql -U postgres -d geo -t -A --field-separator=$'\t' > address.txt
awk '{ print $2 >> $1".geojson" '} address.txt
```

出来上がった GeoJSON ファイルを [geojson.io](http://geojson.io/) に貼付けてみると、市区町村の形状を確認できます。自前で描画したい場合は Leaflet や D3 を使って JavaScript を書けば良いでしょう。

住所コードに対応する名称は総務省の[全国地方公共団体コード](http://www.soumu.go.jp/denshijiti/code.html)で定義されています。利用する住所データに対応する年度のコード表を使うとマッピングできます。なお、いくつか注意点があります。

* 年度が異なる場合は市区町村の合併などに適応できない場合があります。
* 総務省のデータには６桁目のチェックディジットがあります。
* 総務省のデータには「郡」がありません。
* 総務省のデータは都道府県市区町村名はカタカナですが、政令指定都市名はひらがなで併記されています。

PostgreSQL などのデータベースでデータを管理しておくと、JOIN を使った結合処理を簡単に記述できる点は便利だと言えます。表示するラベルを工夫したい場合はデータを出力するクエリを頑張って記述しましょう。