Windows 環境で pgRouting を使用した経路探索を行う（日本全国対応版）

[**pgRouting**](https://pgrouting.org/) は経路探索機能を PostGIS/PostgreSQL データベースに追加する拡張機能であるが、日本全国データに対して経路探索を実行しようとすると、主にハードウェアスペックの問題でいくつかうまくいかない点があったため、最終的な解決方法について以下にまとめることにした。

## はじめに

今回使用した環境は以下の通り。

- ハードウェア
    - **CPU**： AMD Ryzen 9 3950X 16-Core Processor
    - **Memory**： Corsair DDR4 128 GB (32 GB x4)
    - **Disk**： Samsung SSD 980 PRO 2TB

- ソフトウェア
    - Microsoft Windows 10 Pro
    - PostgreSQL 13.2
        - 拡張機能： PostGIS 3.1.1
        - 拡張機能： pgRouting 3.1.3
    - osm2pgsql version 1.4.2

## PostgreSQL および拡張機能のインストール

PostgreSQL の Windows 版公式インストーラ（[**こちら**](https://www.enterprisedb.com/downloads/postgres-postgresql-downloads) から入手可能）を使用して PostgreSQL をインストールする。同時にインストールされる Application Stack Builder を使用して PostGIS および pgRouting 拡張機能もインストールすることができる。

## データベースの準備

PostgreSQL をインストールすると追加される `psql` プログラムを起動し、以下の SQL を実行してデータベースを作成する。

```sql
-- pgr_test データベースを作成
CREATE DATABASE pgr_test;
```

作成したデータベースで PostGIS および pgRouting 拡張機能を有効化する。

```sql
-- pgr_test データベースに接続
\c pgr_test;

-- pgrouting 拡張機能を有効化
--   CASCADE を指定することにより依存関係にある postgis 拡張機能も同時に有効化される
--   CASCADE に対応しない PostgreSQL 9.6 未満では postgis -> pgrouting の順にそれぞれ有効化が必要
CREATE EXTENSION pgrouting CASCADE;
```

## 地図データのインポート

今回は OpenStreetMap の地図データを使用する。OpenStreetMap の全世界地図データ（[**Planet.osm**](https://wiki.openstreetmap.org/wiki/Planet.osm)）はきわめて巨大だが、これを国や地域ごとに細かく分割したものが複数のサードパーティにより配布されている。今回は https://download.geofabrik.de/asia/japan.html から最新の **japan-latest.osm.pbf**（本稿執筆時点で約 1.5 GB のデータ）をダウンロードして使用する。

OpenStreetMap の地図データの PostgreSQL へのインポートには `osm2pgsql` と呼ばれるツールを使用する。Windows 向けバイナリは https://osm2pgsql.org/doc/install.html#installing-on-windows からダウンロードできる。 `osm2pgsql.exe` をコマンドプロンプトから以下のように実行する。

```powershell
osm2pgsql.exe -d pgr_test -U postgres -W -H 127.0.0.1 -C 96000 --latlong ..\..\japan-latest.osm.pbf
```

`-d`： PostgreSQL のデータベース名
`-U`： PostgreSQL のユーザー名
`-W`： PostgreSQL のパスワードプロンプトを表示
`-H`： PostgreSQL のホスト名または IP アドレス
`-C`： ハードウェアに搭載されているメモリの量に合わせてなるべく大きく設定（単位は MB）。
`--latlng`： 地図データを WGS84 地理座標系としてインポート（必須ではないが、指定しておけば緯度経度を使用して PostGIS クエリを実行する際に測地系の変換が不要になるので指定することを推奨）

※ その他のオプションについては `osm2pgsql.exe --help` で確認できる。

上記のコマンドを実行することで、私の環境では約 10 分で日本全国データのインポートが完了した。なお、`-C` に指定したキャッシュサイズが足りない場合はインポートに失敗する（日本全国データの場合、少なくとも 24 GB あれば足りることは確認済み）。メモリが足りない場合は、`--slim` オプションを指定すると、本来メモリ上に展開する作業データがデータベースに書き込まれるようになり、メモリをあまり消費せずにインポートが可能だが、速度は格段に遅くなる（ハードウェアスペックによっては 1 日以上かかることもある）。

> **MEMO**： Amazon EC2 (r5.xlarge) から Amazon RDS (db.t3.medium) に日本全国データをインポートする所要時間は約 30 分。一方で、Amazon EC2 (t3.medium) から `--slim` オプションを指定して Amazon RDS (db.t3.medium) に日本全国データをインポートする所要時間は約 15 時間。

## トポロジーの構築

ここまでを実行することで、PostGIS クエリの実行（ジオコーディングなど）は既に可能となっているが、経路探索を行うためには、さらにトポロジーに関する情報（どの道路とどの道路がつながっているか、というようなデータ）が必要となるため、これの構築を行う。

トポロジーを構築するには大きく分けて 2 つの方法がある。

1. OpenStreetMap 地図データと `osm2pgrouting` を使用して、新規テーブルにトポロジーを構築する。
2. `pgr_createTopology()` を実行して、既存テーブルの新規カラムにトポロジーに関する情報を追加する。

インターネットを検索すると上記 1. の方法が紹介されていることが多いが、`osm2pgrouting` を実行するのに必要なメモリがあまりにも多く、日本全国データのインポートには適さなかった（詳細は後述）ため、今回は 2. の方法でトポロジーを構築する。トポロジーに関する情報は、`osm2pgsql` を実行した際に作成された道路テーブルの新規カラムに追加されるので、まずは新規カラムを作成し、その後 [**pgr_createTopology**](https://docs.pgrouting.org/3.1/en/pgr_createTopology.html) 関数を実行することでトポロジーの構築を行う。

```sql
-- planet_osm_roads テーブルに新規カラムを追加
ALTER TABLE planet_osm_roads ADD COLUMN "source" integer;
ALTER TABLE planet_osm_roads ADD COLUMN "target" integer;

-- トポロジーの構築
--   第１引数： osm2pgsql を実行した際に作成された道路テーブルの名前
--   第２引数： 道路の端点同士が接続されているとみなす最大の距離（単位は測地系により異なる）
--   第３引数： 第１引数に指定した道路テーブルにおいてジオメトリ情報が格納されているカラムの名前
--   第４引数： 第１引数に指定した道路テーブルにおいて主キーが格納されているカラムの名前
SELECT pgr_createTopology('planet_osm_roads', 0.00001, 'way', 'osm_id');
```

私の環境では `pgr_createTopology()` の実行に約 7 分を要した。

> **MEMO**： 第 2 引数 (tolerance) の設定値とそれに対応する `pgr_analyzeGraph()` の実行結果は下表の通りであった。

> |                         Tolerance | 0.000001 | 0.00001 |  0.0001 |   0.001 |    0.01 |     0.1 |       1 |
|------------------------------------:|---------:|--------:|--------:|--------:|--------:|--------:|--------:|
|                   Isolated segments |   31,207 |  31,123 |  16,661 |   1,400 |      32 |       1 |       0 |
|                           Dead ends |  108,571 | 108,423 |  80,920 |  19,802 |   1,251 |      12 |       0 |
| Potential gaps found near dead ends |   74,605 |  74,860 |  67,275 |  16,355 |     948 |      10 |       0 |
|              Intersections detected |  510,781 | 510,504 | 499,439 | 441,886 | 281,768 | 125,735 |  24,225 |
|                     Ring geometries |   11,463 |  11,463 |  11,463 |  11,463 |  11,463 |  11,463 |  11,463 |

## 経路探索を実行

`pgr_createTopology()` を実行すると、トポロジーにおける各ノードの情報を格納した新規テーブル planet_osm_roads_vertices_pgr が作成される。経路探索はこのノード間に対して行われるため、まずは出発地と目的地のノードを決定する必要がある。今回は東京都日本橋から青森県大間崎までの経路探索を行ってみよう。

```sql
-- 東京都日本橋（東経 139.7745 北緯 35.6840）から最も近いノードの id を取得
SELECT id FROM planet_osm_roads_vertices_pgr
  ORDER BY the_geom <-> ST_GeomFromText('POINT (139.7745 35.6840)', 4326)
  LIMIT 1  -- 取得結果は id=271112（地図データにより異なる）

-- 青森県大間崎（東経 140.9130 北緯 41.5465）から最も近いノードの id を取得
SELECT id FROM planet_osm_roads_vertices_pgr
  ORDER BY the_geom <-> ST_GeomFromText('POINT (140.9130 41.5465)', 4326)
  LIMIT 1  -- 取得結果は id=616975（地図データにより異なる）
```

続けて、これらのノード間の経路を [**ダイクストラ法**](https://ja.wikipedia.org/wiki/%E3%83%80%E3%82%A4%E3%82%AF%E3%82%B9%E3%83%88%E3%83%A9%E6%B3%95) により探索する。なお、ここではコストを単純に道路の長さとしている（ `ST_Length(way) AS cost` ）ため、一般道よりも高速道路を優先的に利用するなどの調整は行われない。planet_osm_roads テーブルの `highway` カラムには道路種別に関する情報が含まれているため、これらをもとにコストを計算するようにすれば、より実際の所要時間が短い経路を探索できるであろう。ちなみに、`osm2pgrouting` を使用した場合は、道路種別に応じてあらかじめ調整されたコストがデータベースにインポートされるようになっている。

```sql
SELECT seq, node, edge, cost FROM pgr_dijkstra(
  'SELECT osm_id AS id, source, target, ST_Length(way) AS cost FROM planet_osm_roads',
  271112, -- 東京都日本橋の最寄りのノード ID
  616975, -- 青森県大間崎の最寄りのノード ID
  false
);
```

探索された経路を可視化すると下図の赤線のようになる。地方をまたいだ経路探索ができることが確認できた。
![Capture.PNG](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/176935/19a39639-8355-14d2-3433-86d60a68a5d0.png)

## （参考）osm2pgrouting の使用

先に述べたように `osm2pgrouting` を使用したトポロジーの構築は、特定の狭いエリアに限定された経路探索（特定の都道府県内など）であれば問題ないが、日本全国データを使用する場合には適さないことがわかった。以下に実際に試した内容を記録として残しておく。

#### osm2pgrouting のインストール

Windows 版の `osm2pgrouting` は PostGIS のバンドルに同梱されている。https://postgis.net/windows_downloads/ からダウンロードしたパッケージの bin ディレクトリ内にある `osm2pgrouting.exe` が実行プログラムである。

#### 地図データの変換

`osm2pgrouting` は `osm2pgsql` とは異なり、圧縮されたデータ形式である PBF ファイルをそのまま読み込むことができない。そのため、`osmconvert` と呼ばれるツールを使用してファイル形式を変換する必要がある。`osmconvert` は https://wiki.openstreetmap.org/wiki/Osmconvert から Windows 版バイナリを入手できる。これをコマンドプロンプトで以下のように実行する。

```powershell
osmconvert.exe japan-latest.osm.pbf --drop-author --drop-version --out-osm -o=japan-latest.osm
```

上記のコマンドを実行することで、2 分程度で **japan-latest.osm**（XML 形式）ファイルが得られる。なお、変換前のファイルが 1.5 GB であるのに対し、変換後のファイルは 20 GB 程度である。`--drop-author` および `--drop-version` オプションは、経路探索に不要なデータを変換中に除去するためのオプションで、これを指定しないとファイルサイズはさらに増え 33 GB 程度になる。

#### 環境変数 PATH の設定

`osm2pgrouting` を動作させるため、次のディレクトリに PATH を通しておく。

- `%PROGRAMFILES%\PostgreSQL\13\lib`
- `%PROGRAMFILES%\PostgreSQL\13\bin`

#### osm2pgrouting の実行

コマンドプロンプトで以下のように実行する。

```powershell
osm2pgrouting.exe -f ..\..\..\japan-latest.osm -c mapconfig_for_cars.xml -d pgr_test -U postgres -W password
```

このコマンドは、`--drop-author` `--drop-version` オプションを指定してデータサイズを削減した osm ファイルを使用している場合でさえも、ハードウェアに搭載されている 128 GB のメモリをすべて消費してしまい、インポート処理中に以下のようなエラーが発生した（このエラーの原因は詳しくは分からないが、おそらくメモリ不足が関連していると推測される）。そのままインポート処理は続行したものの、実際に構築されたトポロジーは不完全な状態で、経路探索はほとんどできない状態であった。

```
[***|                                               ] (7%) Total processed: 1940000      Vertices inserted: 20685       Split ways inserted 24282
[***|                                               ] (7%) Total processed: 1960000      Vertices inserted: 19912       Split ways inserted 22623
[***|                                               ] (7%) Total processed: 1980000
ERROR:  relation "__ways26876" does not exist
While processing FROM 1960000th          to: 1980000th way
count1980000 While processing FROM 1960000th     to: 1980000th way
```

> **MEMO**： 上記は 7% まで処理が進んだ時点でエラーが発生したが、実行するたびにエラーが発生するタイミングや、`relation "__ways26876"` の部分の数字が異なるため、地図データの破損などに起因しているとは考えにくい。

なお、日本全国データではなく四国４県（香川県、徳島県、愛媛県、高知県）の OpenStreetMap 地図データを使用して `osm2pgrouting` を実行したところ、こちらは成功し、四国内部であれば問題なく経路探索が実行できる状態となった。

この結果から、日本全国のトポロジーを構築するには、以下のいずれかが必要と考えられる。

- 本記事で紹介したように `pgr_createTopology()` 関数を実行する。
- メモリを増設する（ちなみに Windows 10 **Home** ではメモリは 128 GB までしか認識されない）。
- スワップ領域を使用してメモリを拡張する。
- OpenStreetMap 地図データを複数のエリア（地方や都道府県ごと）に分割し、複数回に分けて `osm2pgrouting` を実行する。
    - この方法はドキュメントに詳細が記述されていないため未確認。
    - `osm2pgrouting` の実行時に `--addnodes` オプションを指定する必要があるらしい。
    - 分割したエリア間をまたぐ経路探索ができるのかどうかが気になるポイント。
- `osm2pgrouting` の代わりに `osm2po` を使用する。
    - `osm2po` は `osm2pgrouting` よりも省メモリ動作であるらしい。
    - インストールまでは行ったが、使い方が理解できていない。