В подмосковном Подольске, в микрорайоне Силикатная-2, есть один лайфхак — когда на дворе уже 9 вечера, и пиво в магазинах уже не продают — достаточно просто перейти дорогу, чтобы его купить. Через дорогу Москва — в ней желаемое до 11 найти можно.

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

Исторически на этот вопрос отвечает «обратный геокодер»(reverse geocoder). Он важная часть практически всех картографических АПИ — Google, Яндекс, и даже OSM. Но в большинстве случаев его ответ предназначается человеку, и содержит исключительно текстовое описание локации.

Это-не-технологично! И уж точно непрактично. Esosedi, кушали этот кактус пару лет, а потом просто сделали свой обратный геокодер. Главное как и зачем.

Совсем недавно на хабре искали Смерть Кащееву (nested set и вложеность административных рубрик), ходили по районам(отображение данных регионов на карте), и (не)попадали на счетчик Яндекса (прямой геокодер). А теперь разберем, что такое обратный геокодер, и зачем он нужен. А потом разберем механики его работы.

Часть 1.

Так уж исторический сложилось, что esosedi это «места на карте мира». Олдскульные квадратики и привязкой к геокоординатам. Жаль только, что геокоординаты для нормальных людей просто набор чисел. Нормальным людям нужен адрес.

Для перевода координат в адрес и используются «обратные геокодеры». Ты ему координаты, а он тебе некий адрес. Этого текста достаточно для понимания человеком расположения некого объекта. И, по человечески, никаких претензий к стандартным реализациям у меня нет. Но для «системы» это не подходит — тут и текстовый «матчинг» сам по себе не стабилен (и тормознут), так еще «нормальные» геокодеры ищут сразу по нескольким источникам, или по крайней мере по данным различных «провайдеров». Два запроса географически близких адресов могут быть найдены в разных базах, и ответы «текстами» друг с другом не состыкуются.

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

Так например встречается иногда Город Москва, который Регион Город Москва, иногда еще и Городской Округ Город Москва, в котором сам Город Москва. Хотя чисто технически — это правильно. Нормализовано по образу и подобию других городов. Проблема в том, что иногда Москва — просто Москва, да и вообще по флагам в ответах она не город. Спасибо, что «Большую Москву» еще не придумали.

Еще одна стандартная ошибки геокодера гугла — «притяжка» координат в улицам, не аккуратная работа с городскими округами, сокращения (Алтай Республиц), использование различных наименований местности. Россия может быть и Россией и РФ и нормальной Российской Федерацией.

Я решил не ждать, когда за мной прийдут Аутономусы из Республики Чечения(всего 5 вариантов написания), и два года назад сделал свой геокодер на основе модуля регионов.

Был наш — стал ваш. После случайного комментария в топике dadata геокодер стал публичным.

Работает все просто — вы ему координаты, а в ответ получаете все osmRelationId, в которые эта точка входит. С разбивкой по административным уровням. Везде циферки, везде связь с другими источниками, везде «проверяемость».

Для визуализации работы есть пример — jsfiddle.net/azvjukh8, работа которого мне так понравилась, что теперь тоже самое красуется на главной esosedi.
Алгоритм работы примера такой:

• Берем текущий центр карты и отправляем запрос в геокодер.
• Находим в ответе объект, размер которого нас устраивает.
• Показываем его на карте. (На самом деле его родителя).
• При движении карты — повторить.

С момента выпуска osme многие люди меня просили облегчить поиск нужного региона — вот наконец и решается одна из проблем модуля подсветки административного деления — сложность найти нужный тебе «номер» региона. Ранее предлагалось идти на морду data.esosedi.org и искать там свое место рекурсивно погружаясь в пучины административного деления.

Теперь все проще — jsonp ручка data.esosedi.org/geocode/v1?[lng=(ru|en)]&point=x,y[&seq=?][&callback=?] — работает 24 часа в сутки, 7 дней в неделю, с ограничением на 6* запросов в секунду средствами nginx…

И позволяет, в среднем, делать в 3 раза больше запросов чем геокодеры АПИ Карт. А для тех кто сможет распределить нагрузку по фронтам (средствами DNS) выйдет примерно 20 RPS (+ burst x4).

И самое главное — никаких домов и улиц, к которым может быть неправомерно притянута точка. А то иногда нормальная улица может находиться и в паре сотен километров от искомой точки, совершенно в другой стране. Только страны, кварталы, районы…

Часть 2

Геокодер работает на основе геометрии, которую отдает ручка регионов. Ему надо просто найти самый маленький полигон, куда входит точка, и все полигоны «выше». Но искать «выше» не надо — эта информация уже есть. Но собрал ее именно геокодер. Как уже говорилось — все это дело основано на геометрии OSM. Из OSM достаточно легко вытащить геометрию, но намного сложнее её переварить. Тем более склеить в четкую иерхическую структуру административного деления.

В принципе и у меня и первый этап не особо то получился — очень много геометрии в текущей базе просто нет, в частности многих городов, для которых кто-то почему-то (до сих пор) не создал relation.

Второй «узкий» момент заключается в том, что многие relation имеют неверный admin_level, или не имеют его вообще. Ну а флаг is_in(кто у кого в подчинении) вообще редкий зверь. В итоге самый надежный вариант — произвести геометрическую «пробивку» геометрии.

На вход поступает полигон. Далее требуется найти некую точку внутри него и найти все другие полигоны, в которые эта точка входит. Отсортировать их по размеру и пересортировать в конечную административную вложенность — более мелкий обьект всегда входит в более крупный. С этим все четко.



Итого задача разбивается на две:

1. Быстро находить объекты, в которые входит точка.
2. Быстро производить более точную проверку «PointInPolygon».

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

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

Лично я решил проблему просто — я беру любую точку которая может входить в полигон, например геометрический центр, и проверяю ее «стандартным» алгоритмом.

Стандартный алгоритм (Ray casting, O(n)) работает просто — из бесконечности слева в точку пускается луч, после чего считаются его пересечения с гранями полигона.

Если пересекли нечетное количество раз — мы внутри. И точка нам подходит.

Если четное — снаружи.

И под конечную точку используем середину отрезка от первой точки пересечения до второй — она точно в полигоне.

Всего на свете есть много алгоритмов для определения принадлежности точки полигону, с совершенно разной сложностью. Часть из них требует препроцесинг, часть — нет. Не тривиальные варианты или переводят полигон в более «математические» формы, снижая алгоритмическую сложность, или просто снижают N(количество вершин). Например любой, самый сложный, полигон можно просто разрезать на части. Для этого хорошо подойдет geojson-vt.

Второй интересный момент — поиск минимального обьекта, который может пересекать точку. Чем меньше обьект — тем быстрее проверка. И точнее. Для этой цели хорошо работают Quad или R-деревья. Но это не особо лампово, посему давайте сделаем это через Z коды.

Z код (N, Morton, Peano) — это число, полученное из 2Д координат (X,Y), такое что обьекты выше или левее всегда будут уметь его меньше, а правее и ниже — больше. Nested set.

Фактически мы можем просто найти все обьекты, у которых один из уголов умеет Z координату больше искомой, и отсортировать по Z (у mapbox есть хорошая библитека, которая использует этот эффект — earcut).

Хотя, если бы это работало — название у этих кодов было бы другое — Z код идет зигзагами. Когда мы попытаемся искать «ниже и правее» — его уведет налево. Спасает тут то, что под маской Z(оrro) скрывается bit interleaving. Из Z-кода через наложение маски всегда можно получить X или Y координату, и дополнительно ограничить поиск по ним… Ну и совершенно академическии можно найти точку «выхода» (LITMAX) из региона поиска и точку «входа» (BIGMIN). Этого, кстати, в earcut не хватает.



В принципе использование «spatial filling curves» — это один из самых эффективных способов для 2D поиска, если у вас обычная база с циферками. Например MySQL. Хотя, меня до сих пор не покидает чувство, что PostGIS…

Остается только загрузить из базы реальную геометрию, и произвести более точные вычисления. Тут возможны различные оптимизации, главная из которых лежит в основе данных OSM — не надо проверять полигон как полигон. Полигон это relation — набор отдельных сегментов. OSM это почти нормальная «топологическая» база данных. Можно проверять только маленькие сегменты, часто общие для двух вероятных ответов.

Насколько это все быстро работает? — На рабочем ноуте при concurrency>10 получается порядка 600RPS. На боевых серверах с одной стороны быстрее, с другой стороны пинги до Германии — 60мсек.

Итого на свете стало на еще один обратный геокодер больше. Он не такой как все — он «проверяемый»(термин википедии), выдает не только «текстовое» описание локации, но и ID, в том числе из разных баз, на которые можно опираться. Он не ищет дома, он не ищет улицы — только административную вложенность. (И только ту, что попала в текущую базу. Екатеринбурга, например, так нет — только городской округ.)
Если вам нужен реально «нормальный» гекодер — возьмите Nominatim, который работает поверх OSM, или нормальные геокодеры нормальных АПИ Карт.

… но тогда пропадет «дух охоты»....

Если жив он в тебе, %username%, то жертва твоя — пример на фидле — jsfiddle.net/azvjukh8.

Кому-то(dadata) такой геокодер может быть полезен для определения района улицы. Кому-то(калькулятор доставки) для фильтрации замкадышей.

Кому-то(мне) он просто позволяет приятно проводить вечера.

Документация по формату ответа доступна на https://github.com/esosedi/regions.

Еще можно почитай немного теории, да практики:

• en.wikipedia.org/wiki/Point_in_polygon
• en.wikipedia.org/wiki/Winged_edge
• www.gdmc.nl/oosterom/gap.pdf
• www.geos.ed.ac.uk/~gisteac/gis_book_abridged/files/ch27.pdf

Или собери гекодер сам:

• github.com/twain47/Nominatim
• http://wiki.openstreetmap.org/wiki/RU:%D0%AD%D0%BA%D1%81%D0%BF%D0%BE%D1%80%D1%82

PS: Ограничение на продажу алкогольной продукции сняли 18 декабря 2014 года (5 января 2015). Теперь по Области действуют федеральные нормативы. Но до магазинов в Подольске, как и до меня, эта новость не дошла.