folium（位置情報の可視化）の使い方メモ

## 概要
地理情報をPythonで扱うのに便利なfoliumというパッケージについてまとめる。

### foliumについて
foliumは、[Leaflet](https://leafletjs.com/)というJavaScriptライブラリを上手いことPythonでラップしてくれているライブラリ。

1. Mapオブジェクトを作成
1. Mapオブジェクトに様々な形を追加
1. Mapオブジェクトをhtmlとして保存

という手順で、簡単にインタラクティブな地図を含んだHTMLファイルを作成することができる。JupyterLab上では、単純にMapオブジェクトを表示させるだけで画面内で地図が表示される。

このエントリは、JupyterLab上で書いている。当エントリ内のコードを順番にJupyterLabで実行すれば、同じ結果が得られる（はず）。

### バージョン情報
```{python}
import sys
import folium
print(f"""Python
{sys.version}

folium
{folium.__version__}""")
```

```
Python
3.7.3 (default, Mar 27 2019, 16:54:48) 
[Clang 4.0.1 (tags/RELEASE_401/final)]

folium
0.10.1
```

## 準備
以下で使う定数を定義しておく。

```python
# 虎ノ門ヒルズの緯度経度
office_lat = 35.66687568
office_lng = 139.74947495
```

## 基本
### とりあえず地図を表示
Jupyter上で地図を表示するには、単純にfolium.Mapを評価すればよい。もしくは、Map.saveでHTMLファイルに吐き出せる。

```{python}
fmap1 = folium.Map(
    location=[office_lat, office_lng],
    tiles = "OpenStreetMap",
    zoom_start = 20, # 描画時の倍率 1〜20
    width = 800, height = 800 # 地図のサイズ
) 
fmap1 # もしくは fmap1.save("1.html")
```
![スクリーンショット 2020-02-03 16.08.01.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/7973/b53e902a-f897-573f-76ea-13e10ef3257b.png)

`tiles` は、

- OpenStreetMap
- Mapbox Bright (Limited levels of zoom for free tiles)
- Mapbox Control Room (Limited levels of zoom for free tiles)
- Stamen (Terrain, Toner, and Watercolor)
- Cloudmade (Must pass API key)
- Mapbox (Must pass API key)
- CartoDB (positron and dark_matter)

から選べる。


### マーカーを置く
マーカーを置く。マーカーにはポップアップを登録できる。

```{python}
fmap2 = folium.Map(
    location=[office_lat, office_lng], 
    zoom_start=20
)
folium.Marker([office_lat, office_lng], popup="Datawiseはここにあります").add_to(fmap2)
fmap2
```
![スクリーンショット 2020-02-03 16.09.19.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/7973/ff61af24-28ea-1e0c-7f64-b91660a3ba50.png)


### 線を書いてみる
PolyLineに(lat, lng)のリストを渡すと、折れ線を描く

```{python}
import itertools as it

# 虎ノ門ヒルズを中心にした四角形の頂点
sq = [
    (office_lat + dy * pow(10, -3), office_lng + dx * pow(10, -3))
    for dx, dy in it.product([-1, 1], [-1, 1])
]
fmap3 = folium.Map(location=[office_lat, office_lng], zoom_start=20)
folium.PolyLine(locations=sq).add_to(fmap3)
fmap3
```

![スクリーンショット 2020-02-03 16.10.25.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/7973/c1bd94b7-50ec-ee5a-5a76-8cd083c58890.png)

### 面を塗ってみる
Polygonで多角形が描ける。

```{python}
sw, nw, se, ne = sq
fmap4 = folium.Map(location=[office_lat, office_lng], zoom_start=20)
folium.Polygon(
    locations=[sw, se, ne, nw], # 多角形の頂点
    color="red", # 線の色
    weight=10, # 線の太さ
    fill=True, # 塗りつぶす
    fill_opacity=0.5 # 透明度（1=不透明）
).add_to(fmap4)

fmap4
```

![スクリーンショット 2020-02-03 16.11.31.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/7973/a3c6d9d4-53b9-ece0-d4b9-3197f214033c.png)

原理的には、これくらいの道具があれば、地図上に適当なデータの情報を可視化できる・・が、実務的には多角形を頂点のリストで管理するなどはダルい（かもしれない）。

地図上の形状を表現するための規格として[GeoJSON](http://geojson.io/)が定義されており、foliumでも扱うことができる。

## GeoJSONオブジェクトを使う
GeoJSONを扱うためにはgeojsonライブラリを使う。

### GeoJSONの注意点
GeoJSONは、仕様をよく読まないと死ぬ。自分は以下の点に気がつくまでに数時間を費やした・・。これらを同僚に伝えたいが為にこのエントリを書いたと言っても過言ではない。

- GeoJSONでは、座標は(latitude, longitude)ではなく(longitude, latitude)で表される（foliumと逆）
- geojson.Polygonに渡すリストは、最初と最後が同じ値でなければならない
- geojson.Polygonのcoordinatesは「(longitude, latitude)のリスト」ではなく、「(longitude, latitude)のリスト"のリスト"」を受け取る[^1]

詳しくは、以下のコードで・・

```python
import geojson as gj

# (lat、lng)のリスト
# ポイント1. 最初と最後の要素が同じ値
lat_lng = [sw, se, ne, nw, sw] 

# ポイント2. (lng, lat)に変換する
def swap(p):
    return p[1], p[0]
lng_lat = list(map(swap, lat_lng))

# ポイント3. （lng、lat）のリストのリストにする
lng_lat2 = [lng_lat]

poly5 = gj.Polygon(lng_lat2)
fmap5 = folium.Map(location=[office_lat, office_lng], zoom_start=20)
folium.GeoJson(poly5).add_to(fmap5)
fmap5
```

![スクリーンショット 2020-02-03 16.14.29.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/7973/74cb721b-bc2f-408e-0c4b-f758b70f3d64.png)

詳しくは調べてないが、上手くいかない時にはHTMLで保存→ブラウザで開いてJSのエラーなどを見ると何か情報があるかもしれない。

### 複数のポリゴンを描画する
GeoJSONでは、複数のオブジェクトの集まりをFeatureCollectionで表現する。FeatureCollectionを使えば、一気に何個も描画できる。（実は、folium.GeoJsonは、渡されるGeoJSONとしてFeatureCollectionを期待しているのだが、それ以外が渡された時には内部でFeatureCollectionに変換している）

```python
def slide(poly, i): 
    """
    ポリゴンを、ちょっとズラす関数
    """
    vtx = poly["coordinates"][0] # gj.Polygonのcoodinateは、頂点のリスト"のリスト"
    vtx2 = [
        (lng + i * pow(10, -3), lat + i * pow(10, -3))
        for lng, lat in vtx
    ]
    return gj.Polygon([vtx2]) # gj.Polygonのcoodinateは（略）





fmap6 = folium.Map(location=[office_lat, office_lng], zoom_start=16)
polys6 = [slide(poly5, i) for i in range(-2, 3)]
fc6 = gj.FeatureCollection(polys6)
folium.GeoJson(fc6).add_to(fmap6)
fmap6
```

![スクリーンショット 2020-02-03 16.24.42.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/7973/703db39c-24f7-eea6-f5ab-50a81b76a137.png)


上で上手くいっているように見えるが、実はGeoJSONの仕様としては間違ったことをしている。[FeatureCollectionの仕様](https://s.kitazaki.name/docs/geojson-spec-ja.html#id7)上は、FeatureCollectionのfeaturesは、type="feature"を持つオブジェクトのリストでなければならない。一方でfc6の"features"を見てみると

```python
fc6["features"]
```
```
[{"coordinates": [[[139.746475, 35.663876], [139.748475, 35.663876], [139.748475, 35.665876], [139.746475, 35.665876], [139.746475, 35.663876]]], "type": "Polygon"},
 {"coordinates": [[[139.747475, 35.664876], [139.749475, 35.664876], [139.749475, 35.666876], [139.747475, 35.666876], [139.747475, 35.664876]]], "type": "Polygon"},
 {"coordinates": [[[139.748475, 35.665876], [139.750475, 35.665876], [139.750475, 35.667876], [139.748475, 35.667876], [139.748475, 35.665876]]], "type": "Polygon"},
 {"coordinates": [[[139.749475, 35.666876], [139.751475, 35.666876], [139.751475, 35.668876], [139.749475, 35.668876], [139.749475, 35.666876]]], "type": "Polygon"},
 {"coordinates": [[[139.750475, 35.667876], [139.752475, 35.667876], [139.752475, 35.669876], [139.750475, 35.669876], [139.750475, 35.667876]]], "type": "Polygon"}]
```

見づらいが、fc6の"features"は、type="Polygon"を持つオブジェクトの配列となっている。（上の例が上手く行っているのは、foliumがよしなにやってくれているから、だったと思う、確か）

以下のように書くと、正しいGeoJSONオブジェクトが得られる。

```python
fmap7 = folium.Map(location=[office_lat, office_lng], zoom_start=16)
fc7 = gj.FeatureCollection(
    features=[
        gj.Feature(
            geometry=p,
            id=i
        ) for i, p in enumerate(polys6)
    ]
)
folium.GeoJson(fc7).add_to(fmap7)
fmap7
```

![スクリーンショット 2020-02-03 16.16.20.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/7973/c4960676-a465-4e7f-81ff-6ac499c60285.png)

### Polygonの書式を変更する
Polygonの書式を変更するには、folium.GeoJsonにstyle_functionを渡す。各Featureにstyle_functionを適用した結果が、そのFeatureを描画する際の書式として利用される。

```python
fmap8 = folium.Map(location=[office_lat, office_lng], zoom_start=16)
folium.GeoJson(
    fc7,
    style_function=lambda feature: {
        "fillColor": "red",
        "color": "black",
        "weight": 10 / (feature["id"] + 1),
        "fillOpacity": feature["id"] * 0.2
    }).add_to(fmap8)
fmap8
```

![スクリーンショット 2020-02-03 16.17.16.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/7973/87deef42-56bb-62af-d7d5-a49bdefba527.png)

[^1]: この日本語は分かりずらいので補足します。GeoJSONの中では、地図上の点は長さ2のリスト[longitude, latitude]で表現されます。多角形は複数の点（=各頂点）で表現されるので、多角形の外周は、[longitude, latitude]を複数含んだリストで表現されます。これで多角形を表せる、と思いきや、GeoJSONのPolygonは、穴の空いた多角形も表現できるようになっています。穴の空いた多角形は、1番最初の要素が外周を表すリストで、2番目以降が内部の穴の周を表すリストであるようなリストで表現されます。GeoJSONの仕様上、穴の無い場合にも、外周を表すリスト（これは、[longitude, latitude]という長さ2のリストを複数含むリスト）を1つだけ含むリストを用意する必要があります。例えば、[0,0], [0,1], [1,1], [1,0]の4点を頂点とする四角形を定義するためには、`[[0,0], [0,1], [1,1], [1,0]]` という頂点のリストではなくて、`[[[0,0], [0,1], [1,1], [1,0]]]`という頂点のリストのリストを用意する必要があります。