Статическая OSM карта на python

В настоящие время, благодаря таким библиотекам как leaflet, развернуть на своем сайте динамическую карту почти не уступающую своими возможностями сервисам, предоставляемым такими гигантами, как Yandex и Google, дело пары десятков строчек.

Но когда дело доходит до создания статического изображения в формате png или jpeg, да еще и с возможностью наненсения различных маршрутов или меток, вы столкнетесь с отстутствием удобных пакетов для реализации данной задачи.

В данной статье я покажу как можно создавать статическое изображение карты, нанести на него географический маршрут и сохранить его в файл.

Использовать мы будем данные OpenStreetMap.

Нам потребуется

  • mercantile - для того, чтобы понять какие тайлы извлекать, зная координаты области;
  • pycairo - библиотека для рисования, с лучшим функционалом, чем компонента Draw от популярного Pillow.

Собираем данные

Итак, поехали. Устанавливаем необходимые библиотеки:

pip install mercantile pycairo

Сначала нам надо узнать координаты области, которую мы хотим извлечь в виде карты. Сделать это можно следующим образом:

  • Заходим на https://www.openstreetmap.org/
  • Выбираем "Экспорт" изображения, и в выплывшем слева сайдбаре, кликаем на "Выделить другую область":

  • Выбираем желаемый регион, слева увидим все необходимые координаты:

Таким образом, прямоугольник, который ограничивает зону региона, имеет следующие координаты:

  • Верхний левый угол: 100.393 восточной долготы и 56.559 северной широты;
  • Нижний правый угол: 111.182 восточной долготы и 50.986 северной широты.

Теперь мы можем использовать эти координаты и узнать какие тайлы нам потребуются, чтобы скачать эту карту. Создаем файл geo2image.py и правим его код:

import mercantile

west = 100.393
south = 50.986
east = 111.182
north = 56.559
zoom = 5 # для начала нам хватит зума поменьше, чтобы не сильно нагружать osm

tiles = list(mercantile.tiles(west, south, east, north, zoom))
print(tiles)

В ответ получим:

[Tile(x=24, y=9, z=5), Tile(x=24, y=10, z=5), Tile(x=25, y=9, z=5), Tile(x=25, y=10, z=5)]

У каждого тайла, есть его x-, y-индексы и значения зума. Благодаря простой схеме адресации в веб-меркаторе мы теперь знаем какие тайлы запрашивать у OpenStreetMap и можем собрать карту.

Собираем карту

import io
import urllib.request
import random
from cairo import ImageSurface, FORMAT_ARGB32, Context
import mercantile

west = 100.393
south = 50.986
east = 111.182
north = 56.559
zoom = 5

tiles = list(mercantile.tiles(west, south, east, north, zoom))

tile_size = (256, 256)
# создаем пустое изображение в которое как мозайку будем вставлять тайлы
# для начала просто попробуем отобразить все четыре тайла в строчку
map_image = ImageSurface(FORMAT_ARGB32, tile_size[0] * len(tiles), tile_size[1])

# создаем контекст для рисования
ctx = Context(map_image)

for idx, t in enumerate(tiles):
    server = random.choice(['a', 'b', 'c'])  # у OSM три сервера, распределяем нагрузку
    url = 'http://{server}.tile.openstreetmap.org/{zoom}/{x}/{y}.png'.format(
        server=server,
        zoom=t.z,
        x=t.x,
        y=t.y
    )
    # запрашиваем изображение
    response = urllib.request.urlopen(url)

    # создаем cairo изображние 
    img = ImageSurface.create_from_png(io.BytesIO(response.read()))

    # рисуем изображение, с правильным сдвигом по оси x
    ctx.set_source_surface(img, idx * tile_size[0], 0)
    ctx.paint()

# сохраняем собраное изображение в файл
with open("map.png", "wb") as f:
    map_image.write_to_png(f)

Получаем такую "красоту":

Естественно, выводить в строчку не очень интересно, поэтому, прежде чем рисовать, надо узнать реальное положение тайлов. Для этого нам пригодятся параметры x и y, которые есть у тайлов. Чтобы корректно распложить изображение, а также сформировать корректные размеры болванки под запись, надо собрать минимальные и максимальные x и y координаты тайлов.

import io

# ...

tiles = list(mercantile.tiles(west, south, east, north, zoom))

min_x = min([t.x for t in tiles])
min_y = min([t.y for t in tiles])
max_x = max([t.x for t in tiles])
max_y = max([t.y for t in tiles])

tile_size = (256, 256)
# теперь указываем уже корректные размеры болванки
map_image = ImageSurface(
    FORMAT_ARGB32, 
    tile_size[0] * (max_x - min_x + 1), 
    tile_size[1] * (max_y - min_y + 1)
)

# ...

Исправляем цикл, чтобы он правильно складывал мозайку:

for t in tiles:  # тут теперь просто по тайлам идем
    server = random.choice(['a', 'b', 'c'])
    url = 'http://{server}.tile.openstreetmap.org/{zoom}/{x}/{y}.png'.format(
        server=server,
        zoom=t.z,
        x=t.x,
        y=t.y
    )
    response = urllib.request.urlopen(url)
    img = ImageSurface.create_from_png(io.BytesIO(response.read()))

    # а тут теперь указываем корректный сдвиг
    ctx.set_source_surface(
        img, 
        (t.x - min_x) * tile_size[0], 
        (t.y - min_y) * tile_size[1]
    )
    ctx.paint()

Запускам и вуаля:

А если масштаб 7-ку поставить:

# ...
zoom = 7
# ... 

Получаем красоту:

Как вы, наверное, уже заметили, область, занимаемая изображением, отличается от той зоны, которая была непосредственно указана координатами в самом начале. Произошло это из-за того, что данные хранятся в виде квадратных тайлов, чьи размеры строго заданы в 256x256 пикселей, и которые резать специально для вас никто, к сожалению, не будет.

Загоняем код в функцию

Давайте загоним весь код в функцию, и пусть функция возвращает нам избражение в виде cairo.ImageSurface, а также настоящую область, которую занимает это изображение.

Получим такую штуку:

def get_map(west, south, east, north, zoom):
    tiles = list(mercantile.tiles(west, south, east, north, zoom))

    min_x = min([t.x for t in tiles])
    min_y = min([t.y for t in tiles])
    max_x = max([t.x for t in tiles])
    max_y = max([t.y for t in tiles])

    tile_size = (256, 256)
    map_image = ImageSurface(
        FORMAT_ARGB32,
        tile_size[0] * (max_x - min_x + 1),
        tile_size[1] * (max_y - min_y + 1)
    )

    ctx = Context(map_image)

    for t in tiles:
        server = random.choice(['a', 'b', 'c'])
        url = 'http://{server}.tile.openstreetmap.org/{zoom}/{x}/{y}.png'.format(
            server=server,
            zoom=t.z,
            x=t.x,
            y=t.y
        )
        response = urllib.request.urlopen(url)
        img = ImageSurface.create_from_png(io.BytesIO(response.read()))

        ctx.set_source_surface(
            img,
            (t.x - min_x) * tile_size[0],
            (t.y - min_y) * tile_size[0]
        )
        ctx.paint()

    return {
        'image': map_image,
        # определим реально занимаемую тайлами область,
        # для этого будем конвертить тайл в gps координаты,
        # а потом уже искать минимум и максимум
        'bounds': {
            "west": min([mercantile.bounds(t).west for t in tiles]),
            "east": max([mercantile.bounds(t).east for t in tiles]),
            "south": min([mercantile.bounds(t).south for t in tiles]),
            "north": max([mercantile.bounds(t).north for t in tiles]),
        }
    }

Собираем данные для маршрута

Итак, давайте теперь попробуем нарисовать какой-нибудь маршрут. Для этого нам надо его откуда-то взять, опять идем на https://www.openstreetmap.org/, открываем консоль, открываем вкладку Network и мониторим запросы, кликаем правой кнопкой на карту и выбираем "Маршрут отсюда" и "Маршрут сюда" (не выбирайте сильно длинный маршрут, иначе получите в ответ какой-нибудь плохой статус)

В консоли увидим запрос, копируем его, он будет выглядеть примерно так:

https://router.project-osrm.org/route/v1/driving/104.2565836,52.333758;101.601701597635,56.16934105?overview=false&geometries=polyline&steps=true

Исправляем параметры, чтобы ответ приходил в формате geojson:

https://router.project-osrm.org/route/v1/driving/104.2565836,52.333758;101.601701597635,56.16934105?overview=full&geometries=geojson

т.е.

?overview=false&geometries=polyline&steps=true

меняем на

?overview=full&geometries=geojson

Исправленный запрос загоняем в браузер, и ждем результата в виде json, ну либо если вы "крутой хакер" дергаем curl-ом:

curl 'https://router.project-osrm.org/route/v1/driving/104.2565836,52.333758;101.601701597635,56.16934105?overview=full&geometries=geojson' > route.json

Сохраняем в файл.

Формат json данных пути имеет вид:

{
  "routes": [
    {
      "geometry": {
        "coordinates": [
          [
            104.25675,
            52.333886
          ], // ...
        ]
      }
      "legs": [
        {
          "summary": "",
          "weight": 83136.3,
          "duration": 41389.6,
          "steps": [],
          "distance": 614894.9
        }
      ],
      "weight_name": "routability",
      "weight": 83136.3,
      "duration": 41389.6,
      "distance": 614894.9
    }
  ],
  "waypoints": [],
  "code": "Ok"
}

Нас будет интересовать то, что находится в массиве coordinates.

Итак, у нас есть карта в виде png, и есть список координат. Давайте нарисуем маршрут на карте.

Рисуем маршрут на изображении карты

Чтобы корректно рисовать что-то на карте, нам надо научится мэпить координаты в пиксели. Более того, чтобы это корректно делать правильнее будет преобразовать границы, определяющие размер нашей карты в координаты веб-меркатора. Определить размеры нашей карты, а затем расчитать коэффициент сжатия по x и y.

Рисуем:

west = 100.393
south = 50.986
east = 111.182
north = 56.559
zoom = 7

out = get_map(west, south, east, north, zoom)

# рассчитываем координаты углов в веб-меркаторе
leftTop = mercantile.xy(out['bounds']['west'], out['bounds']['north'])
rightBottom = mercantile.xy(out['bounds']['east'], out['bounds']['south'])

# расчитываем коэффициенты
kx = out['image'].get_width() / (rightBottom[0] - leftTop[0])
ky = out['image'].get_height() / (rightBottom[1] - leftTop[1])

# загружаем координаты
with open("route.json") as f:
    route = json.load(f)
coordinates = route['routes'][0]['geometry']['coordinates']

# а теперь порисуем
context = Context(out['image'])
for c in coordinates:
    # конвертируем gps в web-mercator
    x, y = mercantile.xy(c[0], c[1])
    # переводим x, y в координаты изображения
    x = (x - leftTop[0]) * kx
    y = (y - leftTop[1]) * ky
    # проводим линию
    context.line_to(x, y)

# заливаем наш путь
context.set_source_rgba(1, 0, 0, 0.5)  # красный, полупрозрачный
context.set_line_width(10)  # ширина 10 пикселей
context.stroke() 

# сохраняем результат
with open("map_with_route.png", "wb") as f:
    out['image'].write_to_png(f)

Получаем такую картинку

По-моему, это замечательно - таким простым незамысловатым способом у нас получается сформировать статическую карту.

Обрезаем по исходным координатам

Теперь, когда мы знаем как переходить от gps координат в координаты изображения, модифицируем нашу функцию, чтобы она обрезала скачанное изображение по указаным границам:

def get_map(west, south, east, north, zoom):
    tiles = list(mercantile.tiles(west, south, east, north, zoom))

    min_x = min([t.x for t in tiles])
    min_y = min([t.y for t in tiles])
    max_x = max([t.x for t in tiles])
    max_y = max([t.y for t in tiles])

    tile_size = (256, 256)
    # создаем пустое изображение в которое как мозайку будем вставлять тайлы
    # для начала просто попробуем отобразить все четыре тайла в строчку
    map_image = ImageSurface(
        FORMAT_ARGB32,
        tile_size[0] * (max_x - min_x + 1),
        tile_size[1] * (max_y - min_y + 1)
    )

    ctx = Context(map_image)

    for t in tiles:
        server = random.choice(['a', 'b', 'c'])  # у OSM три сервера, распределяем нагрузку
        url = 'http://{server}.tile.openstreetmap.org/{zoom}/{x}/{y}.png'.format(
            server=server,
            zoom=t.z,
            x=t.x,
            y=t.y
        )
        response = urllib.request.urlopen(url)
        img = ImageSurface.create_from_png(io.BytesIO(response.read()))

        ctx.set_source_surface(
            img,
            (t.x - min_x) * tile_size[0],
            (t.y - min_y) * tile_size[0]
        )
        ctx.paint()

    # расчитываем коэффициенты
    bounds = {
        "left": min([mercantile.xy_bounds(t).left for t in tiles]),
        "right": max([mercantile.xy_bounds(t).right for t in tiles]),
        "bottom": min([mercantile.xy_bounds(t) .bottom for t in tiles]),
        "top": max([mercantile.xy_bounds(t).top for t in tiles]),
    }

    # коэффициенты скалирования по оси x и y
    kx = map_image.get_width() / (bounds['right'] - bounds['left'])
    ky = map_image.get_height() / (bounds['top'] - bounds['bottom'])

    # пересчитываем размеры по которым будем обрезать
    left_top = mercantile.xy(west, north)
    right_bottom = mercantile.xy(east, south)
    offset_left = (left_top[0] - bounds['left']) * kx
    offset_top = (bounds['top'] - left_top[1]) * ky
    offset_right = (bounds['right'] - right_bottom[0]) * kx
    offset_bottom = (right_bottom[1] - bounds['bottom']) * ky

    # обрезанное изображение
    map_image_clipped = ImageSurface(
        FORMAT_ARGB32,
        map_image.get_width() - int(offset_left + offset_right),
        map_image.get_height() - int(offset_top + offset_bottom),
    )

    # вставляем кусок исходного изображения
    ctx = Context(map_image_clipped)
    ctx.set_source_surface(map_image, -offset_left, -offset_top)
    ctx.paint()

    # возвращаем уже просто изображение, так как размеры уже соответствуют исходному изображения
    return map_image_clipped

В результате получаем:

А вот пример карты с масштабом 10 (можно открыть в полном размере по клику):

Итого, весь код программы имеет вид:

import io
import json
import urllib.request
import random
from cairo import ImageSurface, FORMAT_ARGB32, Context

import mercantile


def get_map(west, south, east, north, zoom):
    tiles = list(mercantile.tiles(west, south, east, north, zoom))

    min_x = min([t.x for t in tiles])
    min_y = min([t.y for t in tiles])
    max_x = max([t.x for t in tiles])
    max_y = max([t.y for t in tiles])

    tile_size = (256, 256)
    # создаем пустое изображение в которое как мозайку будем вставлять тайлы
    # для начала просто попробуем отобразить все четыре тайла в строчку
    map_image = ImageSurface(
        FORMAT_ARGB32,
        tile_size[0] * (max_x - min_x + 1),
        tile_size[1] * (max_y - min_y + 1)
    )

    ctx = Context(map_image)

    for t in tiles:
        server = random.choice(['a', 'b', 'c'])  # у OSM три сервера, распределяем нагрузку
        url = 'http://{server}.tile.openstreetmap.org/{zoom}/{x}/{y}.png'.format(
            server=server,
            zoom=t.z,
            x=t.x,
            y=t.y
        )
        response = urllib.request.urlopen(url)
        img = ImageSurface.create_from_png(io.BytesIO(response.read()))

        ctx.set_source_surface(
            img,
            (t.x - min_x) * tile_size[0],
            (t.y - min_y) * tile_size[0]
        )
        ctx.paint()

    # расчитываем коэффициенты
    bounds = {
        "left": min([mercantile.xy_bounds(t).left for t in tiles]),
        "right": max([mercantile.xy_bounds(t).right for t in tiles]),
        "bottom": min([mercantile.xy_bounds(t) .bottom for t in tiles]),
        "top": max([mercantile.xy_bounds(t).top for t in tiles]),
    }

    # коэффициенты скалирования по оси x и y
    kx = map_image.get_width() / (bounds['right'] - bounds['left'])
    ky = map_image.get_height() / (bounds['top'] - bounds['bottom'])

    # пересчитываем размеры по которым будем обрезать
    left_top = mercantile.xy(west, north)
    right_bottom = mercantile.xy(east, south)
    offset_left = (left_top[0] - bounds['left']) * kx
    offset_top = (bounds['top'] - left_top[1]) * ky
    offset_right = (bounds['right'] - right_bottom[0]) * kx
    offset_bottom = (right_bottom[1] - bounds['bottom']) * ky

    # обрезанное изображение
    map_image_clipped = ImageSurface(
        FORMAT_ARGB32,
        map_image.get_width() - int(offset_left + offset_right),
        map_image.get_height() - int(offset_top + offset_bottom),
    )

    # вставляем кусок исходного изображения
    ctx = Context(map_image_clipped)
    ctx.set_source_surface(map_image, -offset_left, -offset_top)
    ctx.paint()

    return map_image_clipped


west = 100.393
south = 50.986
east = 111.182
north = 56.559
zoom = 7

map_image = get_map(west, south, east, north, zoom)

# рассчитываем координаты углов в веб-меркаоторе
leftTop = mercantile.xy(west, north)
rightBottom = mercantile.xy(east, south)

# расчитываем коэффициенты
kx = map_image.get_width() / (rightBottom[0] - leftTop[0])
ky = map_image.get_height() / (rightBottom[1] - leftTop[1])

# загружаем координаты
with open("route.json") as f:
    route = json.load(f)
coordinates = route['routes'][0]['geometry']['coordinates']

# а теперь порисуем
context = Context(map_image)
for c in coordinates:
    # gps в web-mercator
    x, y = mercantile.xy(c[0], c[1])
    # переводим x, y в координаты изображения
    x = (x - leftTop[0]) * kx
    y = (y - leftTop[1]) * ky
    context.line_to(x, y)

# заливаем наш путь
context.set_source_rgba(1, 0, 0, 0.5)  # красный, полупрозрачный
context.set_line_width(10)  # ширина 10 пикселей
context.stroke()

# сохраняем результат
with open("map_with_route_clipped.png", "wb") as f:
    map_image.write_to_png(f)

На этом все, до новых встреч!

  openstreetmap, python

  Смотреть все посты