−目次

Basemap

経緯度をmatplotlib.pyplotによしなに落としこんでくれるライブラリ…という認識でいいのかな。

経緯度データを地図にプロット

matplotlib, Basemapを用いる。地図画像はOpenStreetMapにお世話になる。

一応Basemap単独でも国境線と海岸線は描けるみたいだけど、より詳細な地図となると真っ白。よってBasemapをキャンバスとして、地図画像を貼りあわせて、その上にプロットすることになる。

  1. 描画範囲を決定
  2. 描画範囲のBasemapオブジェクトを作成
  3. 描画範囲の地図画像を取得、Basemapオブジェクトに貼り付ける
  4. 経緯度データをBasemapにプロットする

OpenStreetMap

地図画像を取得する。

export APIを使えば、最小最大緯度経度と縮尺を指定すると、その範囲の画像を返してくれる。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
import requests
from PIL import Image
from io import BytesIO
 
def get_osm_img(minlat, minlon, maxlat, maxlon, scale=60000, img_format='png'):
    url = 'http://www.openstreetmap.org/export/finish'
    payload = {
        'mapnik_format': img_format,
        'mapnik_scale': scale,
        'minlon': minlon,
        'minlat': minlat,
        'maxlon': maxlon,
        'maxlat': maxlat,
        'format': 'mapnik'
    }
 
    response = requests.post(url, payload)
    # print(response.headers)
    return Image.open(BytesIO(response.content))

ただし、範囲と縮尺によって地図画像があまりに大きくなりすぎると、「Map too large」エラーが返ってくる。その場合、最終行のImage.open()でエラーが出るので、縮尺をもう少し大きな数字にしてやる。

取得した画像は、再度同じクエリを投げてOpenStreetMapのサーバに付加をかけないよう、キャッシュしておこう。

1
2
3
4
5
6
7
file_name = '{}-{}-{}-{}-{}.png'.format(minlat, maxlat, minlng, maxlng, scale)
bg_img = None
try:
    bg_img = Image.open(file_name)
except FileNotFoundError as fnfe:
    bg_img = get_osm_img(minlat=minlat, minlng=minlng, maxlat=maxlat, maxlng=maxlng, scale=scale)
    bg_img.save(file_name)

Basemap

キャンバスを用意し、地図画像を貼り合わせる。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
from PIL import Image
from mpl_toolkits.basemap import Basemap
 
minlat, maxlat, minlng, maxlng = 34.65, 34.70, 135.49, 135.52
scale = 25000
 
# メルカトル図法で指定の最大最小緯度経度の地図を用意
# lat_ts, resolutionの意味は不明
bmap = Basemap(projection='merc',
               llcrnrlat=minlat, urcrnrlat=maxlat, llcrnrlon=minlng, urcrnrlon=maxlng,
               lat_ts=0, resolution=None)
 
bg_img = # OSMから取得した地図画像のPIL.Imageオブジェクト
 
bmap.imshow(bg_img, origin='upper')

プロット

xy座標を与える。pandasを使って、x座標群、y座標群を与えることもできる。

1
2
3
4
5
6
7
8
9
10
11
import pandas as pd
 
data = pd.read_csv('input.csv')
 
# Basemap用の緯度経度群データに変換
# (カラム名は読み込みcsvによって変化)
x, y = bmap(data['lng'].values, data['lat'].values)
 
bmap.scatter(x, y, c=色, alpha=透明度, s=サイズ など)
 
plt.show()

まとめ

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
import pandas as pd
import matplotlib.pyplot as plt
import requests
 
from PIL import Image
from io import BytesIO
 
from mpl_toolkits.basemap import Basemap
 
 
# Retrive static OpenStreetMap
def get_osm_img(minlat, minlng, maxlat, maxlng, scale=60000, img_format='png'):
    url = 'http://www.openstreetmap.org/export/finish'
    payload = {
        'mapnik_format': img_format,
        'mapnik_scale': scale,
        'minlon': minlng,
        'minlat': minlat,
        'maxlon': maxlng,
        'maxlat': maxlat,
        'format': 'mapnik'
    }
 
    response = requests.post(url, payload)
    print(response.headers)
    return Image.open(BytesIO(response.content))
 
 
 
 
fig = plt.figure(figsize=(15, 15))
 
minlat, maxlat, minlng, maxlng = 34.65, 34.70, 135.49, 135.52
scale = 25000
bmap = Basemap(projection='merc',
               llcrnrlat=minlat, urcrnrlat=maxlat, llcrnrlon=minlng, urcrnrlon=maxlng,
               lat_ts=0, resolution=None)
 
file_name = '{}-{}-{}-{}-{}.png'.format(minlat, maxlat, minlng, maxlng, scale)
bg_img = None
try:
    bg_img = Image.open(file_name)
except FileNotFoundError as fnfe:
    bg_img = get_osm_img(minlat=minlat, minlng=minlng, maxlat=maxlat, maxlng=maxlng, scale=scale)
    bg_img.save(file_name)
bmap.imshow(bg_img, origin='upper')
 
data = pd.read_csv('input.csv')
x, y = bmap(data['lng'].values, data['lat'].values)
bmap.scatter(x, y, c='r', alpha=0.5, s=20)
 
plt.show()