Наша команда работает над системой геовизуализации, позволяющей отображать статистику на интерактивной карте. Одной из задач была отрисовка полигонов (в нашем случае — прямоугольников, рис.1) по центральной точке из массива данных, предоставленного заказчиком. Достаточно быстро мы столкнулись с ситуацией искажения размеров этих полигонов в зависимости от широты: чем ближе к северу находилась точка, тем меньше становились полигоны, а на юге они и вовсе начали перекрываться друг другом. Для тех, кто спешит, сразу перейду к объяснению решения задачи. Для остальных, кто настроен на более неспешное чтение, во второй части статьи приведена справочная информация о проекции Меркатора и её особенностях.

Рис.1. Полигоны
Рис.1. Полигоны

Решение задачи

В качестве языка программирования мы использовали Python. Из одной центральной точки нам нужно получить четыре новые точки – это будут вершины прямоугольника, который мы построим вокруг центра. Введём несколько понятий:

  • row и col – ряд и колонка на глобальной сетке полигонов, покрывающей всю
    Землю. Эти индексы определяют, в каком полигоне находится данная точка с
    координатами (lat, lon).

  • resolution – параметр, задающий размерность сетки. Например, если resolution
    = 2, мир делится на 4 × 4 = 16 полигонов (2² по каждой оси). Если resolution =
    15, то получаем сетку 32 768 × 32 768. Чем выше resolution, тем мельче и точнее
    разбита сетка. В нашем случае мы используем проекцию WGS 84 / Pseudo-
    Mercator
     (она же Web Mercator, EPSG:3857), где Земля моделируется как сфера.

Соответственно, формулы для преобразования широты/долготы в координаты сетки (и обратно) берём для сферической модели Земли – это упрощённый подход, принятый в Web Mercator. Такая сферическая проекция работает быстрее и достаточно точно для наших целей (разница с более точной эллипсоидной моделью здесь невелика).

Функция row_to_lat (конверсия Y-координаты в широту)

Первая функция row_to_lat преобразует индекс ряда (ось Y) в географическую широту. Параметр shift определяет, какую именно точку полигона мы получим: верхнюю границу, центр или нижнюю границу:

def row_to_lat(row: int, resolution: int, shift:
  # shift=0 -- top edge
  # shift=0.5 -- center
  # shift=1 -- bottom edge
float = 0.5) -> float:
  return degrees(atan(sinh(float(row + shift)
/ (2 ** resolution) * 2.0 * pi)))

Здесь (row + shift) / (2 ** resolution) – это доля от общей высоты мирового покрытия, а умножение на 2π подготавливает её для обратной формулы проекции Меркатора. Выражение atan(sinh(...)) и перевод из радиан в градусы (degrees) – как раз реализация обратной функции Меркатора для сферической модели. Проще говоря, мы берём долю по вертикали и получаем соответствующую географическую широту. Эта формула эквивалентна более классической записи:

φ=atan⁡ ⁣(sinh⁡(π⋅(1−2y/n)))⋅180π,φ=atan(sinh(π⋅(1−2y/n)))⋅π180​,

где n = 2^{\text{resolution}} и y = \text{row} + \text{shift} – положение полигона. В нашем коде формула реализована с некоторыми условностями индексации, но по сути она делает то же самое – находит широту по Y-координате полигона (сферический вариант проекции Меркатора).

Проверка диапазона: если row = 0 с shift=1 (верх самого верхнего полигона), функция даст широту примерно +85.05113°, а для самого низа (последний ряд с shift=0) – примерно -85.05113°. Эти ~85° северной и южной широты соответствуют предельным значениям, которые отображает проекция Псевдо- Меркатора (она не доходит до ±90°, обрезая полюса по этим широтам).

Функция col_to_lon (конверсия X-координаты в долготу)

Вторая функция col_to_lon преобразует индекс колонки (ось X) в географическую долготу. Параметр shift аналогично определяет положение внутри полигона по горизонтали: левый край, центр или правый край полигона.

def col_to_lon(col: int, resolution: int, shift:
  # shift=0 -- left
  # shift=0.5 -- center
  # shift=1 -- right
float = 0.5) -> float:
  x = float(col + shift) / (2 ** resolution) * 360.0
  if x > 180:
    x = -(360 - x)
  return x

Здесь (col + shift) / (2 ** resolution) * 360 вычисляет долготу в диапазоне [0°, 360°). Далее, если значение выходит за 180°, оно преобразуется в отрицательное: так мы переводим, например, 190° в -170°. Это нужно, чтобы всегда получать долготы в стандартном диапазоне [-180°, 180°]. Иными словами, функция реализует формулу:

(col + shift) / (2 ** resolution) * 360 

затем приводит её к отрезку [-180°, 180°]. Эта формула эквивалентна широко
используемой в веб-картографии: \lambda = \frac{\text{xtile}}{n} \cdot 360^\circ - 180^\circ, где n = 2^{\text{resolution}}.

Функция get_polygon (формирование вершин полигона)

Теперь, используя две вспомогательные функции, составим полигон вокруг центральной точки. Функция get_polygon принимает row и col полигона, а также параметр is_space_between. Если is_space_between=true, то мы сделаем небольшой отступ внутрь полигона, чтобы между соседними полигонами оставался зазор (это косметический эффект, чтобы границы не сливались).

def get_polygon(row: int, col: int, resolution:
int, is_space_between: bool) ->;
List[List[List[int]]]:
  SPACE_BETWEEN_POLYGONS = 0.0001 if
is_space_between else 0
  y1 = row_to_lat(row, shift=1,
resolution=resolution) - SPACE_BETWEEN_POLYGONS
# left top lat
  x1 = col_to_lon(col, shift=0,
resolution=resolution) + SPACE_BETWEEN_POLYGONS
# left top lon
  y2 = row_to_lat(row, shift=0,
resolution=resolution) + SPACE_BETWEEN_POLYGONS
# right bottom lat
  x2 = col_to_lon(col, shift=1,
resolution=resolution) - SPACE_BETWEEN_POLYGONS
# right bottom lon
return [
  [x1, y2], # left bottom vertex
  [x1, y1], # left top vertex
  [x2, y1], # right top vertex
  [x2, y2], # right bottom vertex
]

Как работает эта функция? Она последовательно вычисляет координаты четырёх углов полигона:

  • Левый верхний угол (x1, y1): y1 – широта верхней границы полигона
    (shift=1 даёт верхний край), x1 – долгота левой границы полигона (shift=0 для
    левого края).

  • Правый нижний угол (x2, y2): y2 – широта нижней границы полигона (shift=0 –
    низ), x2 – долгота правой границы (shift=1 – правый край).

Если параметр is_space_between=true, мы вычитаем/добавляем небольшое значение 0.0001 градусов, слегка смещая координаты внутрь полигона, чтобы между смежными полигонами был визуальный зазор (полигон тогда не вплотную прилегает к соседнему).

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

Почему возникало искажение и как мы его устранили

Теперь, разобрав функции, становится понятнее природа исходной проблемы. Изначально полигоны искажались по размеру, потому что единицы широты и долготы «не эквивалентны» на разных широтах. Проекция Меркатора искажает масштаб: чем ближе к полюсам, тем сильнее масштаб растягивается по широте. Один и тот же промежуток в градусах охватывает на земле разное расстояние: например, 1° долготы у экватора — это ~111 км, а на широте 60° — уже около 55 км. Если бы мы строили полигоны, откладывая фиксированные Δlat и Δlon от центральной точки, то на высоких широтах такой полигон покрывал бы куда меньшую область на карте, чем аналогичный полигон на экваторе. В нашем случае при движении на север полигоны визуально уменьшались, а ближе к экватору (югу) — наоборот, накладывались друг на друга, поскольку выбранный шаг в градусах давал слишком большой охват территории.

Используя же формулы проекции Меркатора, мы адаптируем размеры полигонов к широте. Функции row_to_lat и col_to_lon фактически вычисляют, какой диапазон широт и долгот должен быть у полигона на заданной широте, чтобы на карте (в проекции) он имел нужные размеры. Например, при одной и той же resolution широта Δlat между y1 и y2 для полигона около экватора будет больше (в градусах), чем Δlat для полигона ближе к полюсу – и это компенсирует то, что на карте северные области растянуты по вертикали. Таким образом, каждый наш полигон в координатах проекции имеет одинаковые размеры на карте, чего мы и добивались. Мы как бы разделили мир на равные прямоугольники в пространстве проекции, поэтому на экране все полигоны одинаковы, хотя их реальная площадь на поверхности Земли разнится (это неизбежно для проекции Меркатора).

Надеюсь, эти пояснения помогут разобраться в решении. Далее – немного истории и интересных фактов о проекции Меркатора, которая сделала возможным такой подход к отображению данных.

Немного о проекции Меркатора

Проекция Меркатора — одна из самых известных и значимых картографических проекций в истории. Разработанная в XVI веке фламандским картографом Герардом Меркатором, она произвела революцию в морской навигации, позволив мореплавателям прокладывать прямые курсы на картах без постоянной корректировки компаса. Несмотря на искажения, эта проекция остаётся важной и сегодня – особенно в электронных картографических сервисах.

До XVI века мореходы пользовались портоланами — лоциями с сеткой румбов (направлений ветров), удобными для плавания в пределах Средиземного моря, но малопригодными для океанских путешествий. С началом эпохи Великих географических открытий возникла острая необходимость в картах, учитывающих кривизну Земли. Герард Меркатор (1512–1594), изучив труды античных учёных (например, Птолемея) и современные ему математические достижения, стремился создать карту, которая сохраняла бы углы между направлениями – это было критически важно для навигации.

В 1569 году Меркатор опубликовал свою знаменитую карту мира под названием «Nova et Aucta Orbis Terrae Descriptio ad Usum Navigantium Emendate Accommodata» («Новое и дополненное изображение Земли, исправленное для навигационного употребления»). В проекции этой карты заложены основные свойства, сделавшие её незаменимой для мореходов:

  • Меридианы изображены параллельными вертикальными линиями (хотя в
    реальности они сходятся к полюсам).

  • Параллели (линии широт) растягиваются по мере удаления от экватора, чтобы
    сохранить углы и форму небольших объектов (проекция Меркатора является
    конформной, т.е. сохраняет угловые соотношения).

  • Прямая линия на карте соответствует курсу с постоянным азимутом, то есть
    локсодромии – линии, пересекающей все меридианы под одним и тем же углом.
    Такая линия на глобусе обычно представляет спираль к полюсу, но на карте
    Меркатора она изображается прямой. Именно это свойство упрощало навигацию:
    мореплавателю достаточно было держать постоянный курс по компасу, и на карте
    его путь отображался как прямая линия.

Плюсы проекции Меркатора:

  • Сохранение направлений. Проекция Меркатора удобна для морской навигации,
    так как локально сохраняет углы. Курс, проложенный по компасу, отображается на
    карте прямой линией, что сильно упрощает ориентирование.

  • Простота использования карты. До появления проекции Меркатора штурману
    приходилось постоянно корректировать курс (поскольку кратчайший путь — дуга
    большого круга — постоянно меняет азимут). На карте Меркатора же достаточно
    начертить прямую от пункта А до пункта Б и следовать по постоянному азимуту.
    Это была огромная практическая ценность в эпоху парусного флота.

Минусы проекции Меркатора:

  • Сильное искажение площадей. Поскольку масштаб по широте растёт по мере
    приближения к полюсам, разные регионы отображаются непропорционально.
    Например, Гренландия выглядит почти размером с Африку, хотя на самом деле
    она примерно в 14 раз меньше!

  • Локсодромия не является кратчайшим путём. Прямая линия маршрута на карте Меркатора соблазнительно проста, но на практике это не самый короткий путь между точками на Земле. Большой круг (ортодромия) — истинно кратчайший маршрут — на карте Меркатора выглядит изогнутой линией. На длинных дистанциях разница существенна: например, на широте 60° путь по локсодромии длиннее великокругового примерно на 8.5%. Тем не менее, в эпоху парусников моряки предпочитали чуть более долгий, но простой курс по локсодромии, чтобы не возиться с постоянными поправками направления.

  • Непригодность для полярных областей. Проекция Меркатора по своей природе не может отобразить полюса: она «растягивает» карту к 90° широты до бесконечности. Практически все онлайн-карты обрезают изображение около ±85°. Вы никогда не увидите точку Северного полюса в Google Maps или Яндекс.Картах – она всегда остаётся за краем проекции. Это серьёзный недостаток для отображения арктических и антарктических регионов.

Проекция Меркатора стала стандартом для морских карт и использовалась веками. В XX веке появились альтернативные проекции (например, равнопропорциональная проекция Галла–Петерса), пытающиеся исправить искажения площадей, но меркаторовская по-прежнему остаётся популярной в навигации и веб-картах (Google Maps, Яндекс.Карты и др.). Изобретение Меркатора стало важнейшим шагом в развитии картографии и навигации. И хотя у этой проекции есть недостатки, её практическая ценность для мореплавания, авиации и геоинформационных сервисов остаётся неоспоримой.

Комментарии (0)