−目次

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

Basemap

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

Basemap tutorial — Basemap tutorial 0.1 documentation

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

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

python + OpenStreetMapで地図にデータをプロットする - Qiita
- 以下のコードはほとんど上記サイトのコピペ。感謝。

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

描画範囲を決定
描画範囲のBasemapオブジェクトを作成
描画範囲の地図画像を取得、Basemapオブジェクトに貼り付ける
経緯度データを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のサーバに付加をかけないよう、キャッシュしておこう。

        
              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

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

        
              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座標群を与えることもできる。

        
              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()