PyOsmiumを用いたOpenStreetMapデータ処理

# はじめに
本記事はOpenStreetMapより取得したデータを処理するライブラリの1つであるPyOsmiumについて、その基本的な使用法及び概念についてまとめたものです。
ご意見等ありましたらどうぞよろしくお願いいたします。

# OpenStreetMapとは
OpenStreetMap(OSM)は、世界中の人々による共同作業で作られた、自由に利用・編集可能な世界地図プロジェクトです。オープンソースの地図データが公開されており地図のスタイルは柔軟に変更可能です（[様々なスタイルの例](https://leaflet-extras.github.io/leaflet-providers/preview/)）。現在、多くのサイト（例：[Yahooマップ](http://maps.loco.yahoo.co.jp/maps)）で活用されています。

### ライセンス・クレジット表記
OSMはライセンスとしてOpen Database License (ODbL)を使用しています（[参考1](https://en.wikipedia.org/wiki/Open_Database_License), [参考2](https://opendatacommons.org/licenses/odbl/)）。そして、OSMを公に使用する場合にはガイドラインに従って以下の2条件を守る必要があります([引用元](https://www.openstreetmap.org/copyright))。
> ・Provide credit to OpenStreetMap by displaying our copyright notice.
> ・Make clear that the data is available under the Open Database License.
>> 日本語訳
>> ・OpenStreetMapの著作権表示を表示し、OpenStreetMapへのクレジットを提供する。
>> ・データがOpenDatabaseLicenseの下で利用可能であることを明確にする。

**クレジット表示**は[詳細な表示ガイドライン](https://osmfoundation.org/wiki/Licence/Attribution_Guidelines)を参考に媒体によって適切な方法を選ぶことができます。**ODbLの下で利用可能であることの明確化**には[著作権ページ](https://www.openstreetmap.org/copyright)をリンクすることで対応することができます。ライセンス関連のF&Qは[こちら](https://osmfoundation.org/wiki/Licence/Licence_and_Legal_FAQ)に沢山記載されています。


# PyOsmiumとは
PyOsmiumはOSMデータを処理するためのPythonライブラリです。これは C++ ライブラリであるOsmiumのラッパーであり、OSMデータを高速かつ効率的に順次処理できます。その他のPythonを用いたOSMデータ処理のライブラリとして[Pyrosm](https://pyrosm.readthedocs.io/en/latest/)が挙げられます。リファレンスを読んだところ、PyrosmはPyOsmiumに比べ、ライブラリ単体でできることが多い（データのダウンロード、図示）一方で、[メモリ効率はPyOsmiumのほうが良い](https://stackoverflow.com/questions/68845934/is-there-a-way-to-analyze-oms-open-street-maps-pbf-files-for-large-areas-germ)ようです。

### Osmium対象ファイル形式
Osmiumの対象ファイル形式は以下です（[参考](https://docs.osmcode.org/osmium/latest/osmium-file-formats.html)）。
- XML 形式(.osm.xml; .osh.xml; .osc.xml)
- PBF バイナリ形式 (.osm.pbf; .pbf)
- OPL 形式 (.osm.opl; .opl)
- O5M/O5C 形式 (.o5m; .o5c) (読み取り専用)
- デバッグ形式 (.osm.debug) (書き込み専用)

### OSMデータの取得
OSMデータの取得には以下に示す例を含め多くの方法があります。
- [QGISを用いたデータ取得](https://learnosm.org/ja/osm-data/osm-in-qgis/)
- [OpenStreetMapのサイトから取得](https://www.openstreetmap.org/export)
- [配布サイト](https://download.geofabrik.de/asia/japan.html)からのダウンロード
- Pyrosmを用いたダウンロード

### OSMデータモデル
OSMファイルの解析はそれぞれのオブジェクトに対して行います。解析対象のオブジェクトは以下の5種です。この中で、基礎的なモデルとなっているのは**Node**、**Way**、**Relation**であり、AreaはWayまたはRelationに含まれています（[参考1](https://wiki.openstreetmap.org/wiki/JA:%E8%A6%81%E7%B4%A0), [参考2](https://osmcode.org/osmium-concepts/)）。

```plantuml
@startuml
<style>
root {
  Padding 0
  Margin 0
}
header {
  FontColor black
  FontSize 30
}
</style>
left header OSMデータモデル概観(Osmium Concepts Manualより引用)
title <img:https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/1974522/91deaeba-9563-3052-4b46-d4c7536e74ff.png>
@enduml
```

#### Node
地上にある特定の地点。
```xml
  <node id="63686086" version="11" timestamp="2023-10-18T17:19:41Z">
    <tag k="contact:website" v="http://www.kyohaku.go.jp/"/>
    <tag k="name" v="京都国立博物館"/>
    <tag k="name:en" v="Kyoto National Museum"/>
    <tag k="name:fr" v="Musée national de Kyoto"/>
    <tag k="name:ja" v="京都国立博物館"/>
    <tag k="name:ja_rm" v="Kyōto Kokuritsu Hakubutsukan"/>
    <tag k="tourism" v="museum"/>
    <tag k="wikidata" v="Q147286"/>
    <tag k="wikipedia" v="ja:京都国立博物館"/>
  </node>
```
#### Way
Nodeの順序のある集合。連続した線分やAreaの境界を表し、Areaの境界を表す場合そのWayの先頭と末尾のNodeは同じ座標である。
```xml
  <way id="785043062" version="1" timestamp="2020-03-27T09:51:34Z">
    <nd ref="7334478940"/>
    <nd ref="7334478937"/>
    <tag k="highway" v="service"/>
    <tag k="source" v="Bing"/>
  </way>
```
#### Relation
2つ以上のデータ要素（NodeやWayや他のRelation）の関係を記述する多目的のデータ構造（詳細な記事は[こちら](https://qiita.com/nyampire/items/13f6dd0daad5cc072ea8)）。
```xml
 <relation id="56688" user="kmvar" uid="56190" visible="true" version="28" changeset="6947637" timestamp="2011-01-12T14:23:49Z">
  <member type="node" ref="294942404" role=""/>
  ...
  <member type="node" ref="364933006" role=""/>
  <member type="way" ref="4579143" role=""/>
  ...
  <member type="node" ref="249673494" role=""/>
  <tag k="name" v="Küstenbus Linie 123"/>
  <tag k="network" v="VVW"/>
  <tag k="operator" v="Regionalverkehr Küste"/>
  <tag k="ref" v="123"/>
  <tag k="route" v="bus"/>
  <tag k="type" v="route"/>
 </relation>
```
#### Area
1つ以上の直線で囲まれ、適切なタグでマークされた2次元オブジェクト。最も単純な領域は、同じ始点と終点を持つ閉じたWay。
#### Changeset
OSMデータベース内の関連する変更のセットが記述されている

# PyOsmiumを用いたデータ分析
PyOsmiumの基本的な使用法について、以下に示す複数の例を元に解説します。
PyOsmiumのリファレンスは[こちら](https://docs.osmcode.org/pyosmium/latest/reference.html)、さらに詳しい使用法は[こちら](https://lonvia.github.io/geopython17-pyosmium/)。

## PyOsmiumのインストール
```
pip install osmium
```

## 基本的な使い方
OsmiumによるOSMファイルの解析は、ハンドラに基づいて構築されています（現状、常にosmium.SimpleHandler）。各関数はファイルから読み込まれたオブジェクトに対応しており、処理したいオブジェクトの種類ごとに関数を実装する必要があります。ハンドラーへのOSMファイルの適用は`apply_file()`を用います。以下のコードではOSMファイルに含まれるそれぞれのオブジェクトの個数を数え上げる処理を記載しています。
```python
import osmium

class ObjectCounterHandler(osmium.SimpleHandler):
    def __init__(self):
        osmium.SimpleHandler.__init__(self)
        self.object_counter = {'node':0, 'way':0, 'relation':0, 'area':0, 'changeset':0}

    def count_object(self, key):
        self.object_counter[key] += 1

    def node(self, n):
        self.count_object("node")

    def way(self, w):
        self.count_object("way")

    def relation(self, r):
        self.count_object("relation")
    
    def area(self, a):
        self.count_object("area")

    def changeset(self, c):
        self.count_object("changeset")


if __name__ == '__main__':
    h = ObjectCounterHandler()
    h.apply_file("./kansai-latest.osm.pbf")
    print(h.object_counter) 
```

## データの取得
OSMファイルからデータを得たい場合、ハンドラ内でオブジェクトを変数に格納し後で値を取得するのではなく、**ハンドラ内で値を取得し変数に格納する**必要があります。クラス内では、あくまでOsmiumはファイルからオブジェクトを読み込み、ハンドラ関数に渡し、メモリを解放するため、後でオブジェクトを参照しても値を得ることはできません。
以下のコードではNodeオブジェクトのタグがホテルであるものの名前を取得する方法について記載しています。
オブジェクトから取得できる内容は[こちら](https://docs.osmcode.org/pyosmium/latest/ref_osm.html)、オブジェクトのタグ一覧は[こちら](https://wiki.openstreetmap.org/wiki/JA:Map_Features)。

### データ取得が不可能な例
```python
class HotelHandler(osmium.SimpleHandler):
    def __init__(self):
        osmium.SimpleHandler.__init__(self)
        self.hotels = []

    def node(self, o):
        if o.tags.get('tourism') == 'hotel':
            self.hotels.append(o)       # オブジェクトの保持は意味が無い
            

h = HotelHandler()
h.apply_file("./kansai-latest.osm.pbf")

hotel_names = []
for o in h.hotels:
    if 'name' in o.tags: # Illegal access to removed OSM object
        hotel_names.append(o.tags['name'])

print(sorted(hotel_names))
```

### データ取得が可能な例
```python
import osmium

class HotelHandler(osmium.SimpleHandler):
    def __init__(self):
        osmium.SimpleHandler.__init__(self)
        self.hotels = []

    def node(self, o):
        if o.tags.get('tourism') == 'hotel' and 'name' in o.tags:
            self.hotels.append(o.tags['name']) # 値を保持する


h = HotelHandler()
h.apply_file("./kansai-latest.osm.pbf")

print(sorted(h.hotels))
```

## 地理的計算
Osmiumより得たデータについて、外部ライブラリと連携して地理的な計算をすることが可能です。また、ハンドラがWayやAreaのジオメトリ処理をする場合はNodeについて内部的にキャッシュする必要があり、`locations=True`と明示的に示す必要があります。また`idx`では使用するキャッシュの種類を指定することができ、デフォルトの`flex_mem` はほとんどの用途に適しています。
以下のコードでは、Osmium形式をWell-Known Binary形式（WKB）に変更し、`shapely`を用いてWayの全長の計算を記載しています。

```python
import osmium
import shapely.wkb as wkblib

wkbfab = osmium.geom.WKBFactory()

class WayLenHandler(osmium.SimpleHandler):
    def __init__(self):
        osmium.SimpleHandler.__init__(self)
        self.total = 0

    def way(self, w):
        wkb = wkbfab.create_linestring(w) # wkb形式とする
        line = wkblib.loads(wkb, hex=True)
        self.total += line.length

if __name__ == '__main__':
    h = WayLenHandler()
    h.apply_file("./kansai-latest.osm.pbff", locations=True, idx='flex_mem')
    print("Total length: %f" % h.total)
```

## Geojsonによるデータ出力
得たデータをGeojson形式にて出力したい場合は`osmium.geom.GeoJSONFactory()`を用いることができます。
以下のコードでは、タグ名が`highway`のWayについてGeojson形式の緯度経度情報、タグ名、名称について出力する処理を記載しています。

```python
import osmium
import json
import os

gjfac = osmium.geom.GeoJSONFactory()

class WayHandler(osmium.SimpleHandler):
    def __init__(self):
        osmium.SimpleHandler.__init__(self)
        self.total = 0
        self.something = []

    def way(self, w):
        if 'highway' in w.tags:
            wkb = gjfac.create_linestring(w)
            self.something.append({"geojson":wkb, "name":w.tags.get('name', '<unknown>'), "type":w.tags.get('highway', '<unknown>')})

if __name__ == '__main__':
    h = WayHandler()
    h.apply_file('./kansai-latest.osm.pbf', locations=True, idx='flex_mem')
    count = len(str(len(h.something)))
    for i, json_string in enumerate(h.something):
        json_data = json.loads(json_string['geojson'])
        json_data["name"] = json_string['name']
        os.makedirs(json_string['type'], exist_ok=True)
        with open(f"{json_string['type']}/{str(i).zfill(count)}-way.geojson", "wt") as f:
            json.dump(json_data, f, indent=2)
```


# まとめ
本記事ではPyOsmiumの使用法について記載しました。
概念やデータ構造が難しく、少し使うだけでもかなり苦戦しました。
使ってないのでわかりませんが、Pyrosmの方がとっつきやすいライブラリな気がします。