24.08.2020 14:37
Duke_nukeum 14 2900 Источник
Недавно в рамках одного собеседования мне понадобилось решить задачу, условие которой приведено ниже:
У лучшего в мире управляющего по имени Пенультимо родилась очередная гениальнейшая идея, peализовать которую вам и предстоит. Он верит, что поток туристов на Исла-де-Эдукадос повысится, если он сможет рассказать всему миру, как же много замечательных дорожных знаков с длинными надписями eсть у них на острове. Вам предлагается придумать алгоритм, позволяющий подсчитать суммарное количество букв на всех знаках «Въезд в город Х» на острове, а затем применить полученные знания для подсчёта аналогичной метрики для Республики Беларусь. Обратите внимание язык, используемый для обозначения населённых пунктов, а также тот факт, что въездов в город может быть несколько. Пенультимо также приветствует инициативность, так что можете исследовать этот вопрос для отдельных областей, провести сравнение с количеством людей, проживающих в области, а также провести любые другие исследования, которые покажутся Вам интересными.

Под катом покажу точное решение этой и других похожих задач, например: «Сколько АЗС находится в пределах Москвы?»

Общий метод решения


Если взглянуть на карту OpenStreetMap, то на ум сразу приходит следующая идея: а давайте для каждого города получим его границы и дороги, находящиеся внутри его границ, а потом найдем их пересечения, на которых будут стоять знаки! Как будем искать пересечения: берем отрезок границы, потом отрезок дороги и смотрим, пересекаются ли они (типичная геометрическая задача). И так пока не кончатся все отрезки и города.

Про архитектуру данных OSM
Вообще, есть три главных компонента: узлы, линии и отношения.
Каждый элемент имеет свой ID, а также свои теги.

  • Узел — это просто точка на карте, имеющая кроме ID также широту и долготу
  • Линия — это отсортированный список узлов, который представляет контур здания или отдельную улицу
  • Отношение — это список, состоящий из узлов, линий или других отношений, представляющий логические или географические связи между объектами на карте


OverPass


OverPass — Это API для получения данных из OpenStreetMap. Оно имеет свой язык составления запросов, про него подробно можно почитать В этой статье.

Для того чтобы было проще и удобнее составлять запросы, есть инструмент Overpass-turbo, где результат выполнения запроса можно посмотреть в удобном и интерактивном виде.

Использование OverPass API в Python


Для обработки данных из OSM в Питоне можно использовать пакет Overpy в качестве оболочки.
Для отправки запросов и получения данных нужно проделать следующее:

import overpy

api = overpy.Overpass()
Data = api.query("""
*ваш запрос*
""")

где в переменной(?) Data и содержится все, что выдал нам сервер.

Как обработать эти данные? Предположим, что мы ввели следующий запрос на получение границ Минска:

relation["type"="boundary"]["boundary"="administrative"]["name:be"="Мінск"];
//Дословно: отношение вида административная граница города Минска
>; out skel qt;

На выходе имеем файл XML (можно выбрать Json), имеющий следующую структуру:

<*шапка файла*>
<далее идет подобное перечисление всех узлов>
  <node id="277218521" lat="53.8605688" lon="27.3946601"/>
  <node id="4623647835" lat="53.8603938" lon="27.3966685"/>
  <node id="4713906615" lat="53.8605343" lon="27.3998220"/>
  <node id="4713906616" lat="53.8605398" lon="27.3966820"/>
  <node id="4713906617" lat="53.8605986" lon="27.3947987"/>
  <node id="277050633" lat="53.8463790" lon="27.4431241"/>
  <node id="277050634" lat="53.8455797" lon="27.4452681"/>
  <node id="4713906607" lat="53.8460017" lon="27.4439797"/>
<далее указаны пути и ID узлов, из которых они состоят>
<way id="572768148">
    <nd ref="5502433452"/>
    <nd ref="277218520"/>
    <nd ref="4713906620"/>
    <nd ref="277218521"/>
    <nd ref="4713906617"/>
    <nd ref="4623647835"/>
    <nd ref="4713906616"/>
</way>
<way id="29079842">
    <nd ref="277212682"/>
    <nd ref="277051005"/>
    <nd ref="4739822889"/>
    <nd ref="4739822888"/>
    <nd ref="4739845423"/>
    <nd ref="4739845422"/>
    <nd ref="4739845421"/>
</way>

Давайте получим некоторые данные:

import overpy

api = overpy.Overpass()
Data = api.query("""
relation["type"="boundary"]["boundary"="administrative"]["name:be"="Мінск"];
>; out skel qt;
""")
Xa=Data.ways[0].nodes[0].lon #получаем долготу первого узла первой линии
Ya=Data.ways[0].nodes[0].lat #получаем широту
Xb=Data.ways[0].nodes[1].lon
Yb=Data.ways[0].nodes[1].lat
NodeID=Data.ways[0]._node_ids[0] #получаем ID первого узла первой линии
print(len(Data.nodes)) #получаем общее количество узлов
print(NodeID)
print(Xa,Ya)
print(Xb,Yb)

С точки зрения работы с OpenStreetMap в питоне, это всё, что понадобится, чтобы достать данные.

Перейдем непосредственно к задаче


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

Заголовок спойлера 
import overpy


###########################
def line_intersection(line1, line2): #функция поиска пересечений
    ax1 = line1[0][0]
    ay1 = line1[0][1]
    ax2 = line1[1][0]
    ay2 = line1[1][1]
    bx1 = line2[0][0]
    by1 = line2[0][1]
    bx2 = line2[1][0]
    by2 = line2[1][1]
    v1 = (bx2 - bx1) * (ay1 - by1) - (by2 - by1) * (ax1 - bx1)
    v2 = (bx2 - bx1) * (ay2 - by1) - (by2 - by1) * (ax2 - bx1)
    v3 = (ax2 - ax1) * (by1 - ay1) - (ay2 - ay1) * (bx1 - ax1)
    v4 = (ax2 - ax1) * (by2 - ay1) - (ay2 - ay1) * (bx2 - ax1)
    return (v1 * v2 < 0) & (v3 * v4 < 0)


#######################################
citytmp = []
city = []
Borderway = []
Roadway = []
Total = 0
A = [0, 0]
B = [0, 0]
C = [0, 0]
D = [0, 0]
amount = 0
progressbar = 0 
ReadyData = open('Готовые данные.txt','w')
with open("Города Беларуси.txt", "r", encoding='utf8') as file:
    for i in range(115):
        citytmp.append(file.readline())
citytmp = [line.rstrip() for line in citytmp]
for i in range(115):
    city.append('"' + citytmp[i] + '"')
city[0]='"Дзісна"'

api = overpy.Overpass()
for number in range(0,115):#главный цикл обработки, перебирает города
    borderstring = """(
relation["type"="boundary"]["boundary"="administrative"]["name:be"=""" + city[number] + """][place=town]; 
relation["type"="boundary"]["boundary"="administrative"]["name:be"=""" + city[number] + """][place=city];
);
>; out skel qt;"""
    roadstring = """(
area[place=town]["name:be"=""" + city[number] + """]; 
way["highway"][highway!=service]["highway"!="footway"]["highway"!="track"]["highway"!="path"]
    ["highway"!="cycleway"]["highway"!="pedestrian"]["highway"!="steps"]["highway"!="residential"](area);
area[place=city]["name:be"=""" + city[number] + """]; 
way["highway"][highway!=service]["highway"!="footway"]["highway"!="track"]["highway"!="path"]
    ["highway"!="cycleway"]["highway"!="pedestrian"]["highway"!="steps"]["highway"!="residential"](area);
);
out body; >; out skel qt;"""
    print('Getting data about', city[number],'...')
        road = api.query(roadstring)
        border = api.query(borderstring)
    print('got data!, city:', city[number]) #получаем данные
    for w in range(len(border.ways)): #перебирает линии границ
        for i in range(len(border.ways[w]._node_ids)):#перебирает узлы линий границ 
            progressbar = i / len(border.ways[w]._node_ids) * 100
            print(progressbar, "%;", w, "of", len(border.ways), "parts ready; city-", city[number])
            A[0] = border.ways[w].nodes[i].lon
            A[1] = border.ways[w].nodes[i].lat
            if i == len(border.ways[w]._node_ids) - 1:
                break
            B[0] = border.ways[w].nodes[i+1].lon
            B[1] = border.ways[w].nodes[i+1].lat
            for j in range(len(road.ways)):
                for k in range(len(road.ways[j]._node_ids)):
                    C[0] = road.ways[j].nodes[k].lon
                    C[1] = road.ways[j].nodes[k].lat
                    if k == len(road.ways[j]._node_ids) - 1:
                        break
                    D[0] = road.ways[j].nodes[k+1].lon
                    D[1] = road.ways[j].nodes[k+1].lat
                    if line_intersection((A, B), (C, D)) == 1:
                        amount += 1
                        print(road.ways[j]._node_ids[k])
    print(amount)
    Total += amount * len(city[number])
    ReadyData.write(str(city[number]))
    ReadyData.write(str(amount))
    ReadyData.write('\n')
    amount = 0
print('Total', Total) #и вывод всего


Заметки по коду


Я достаточно долго составлял запрос, подбирая разные типы дорог чтобы было меньше считать и чтобы не пропустить знаки. В итоговом запросе просто убраны те дороги, на которых нет знаков, например residential, service, footway, track и т. п.

Список городов я спарсил с википедии и сохранил в формате.тхт

Код выполняется достаточно долго, у меня даже один раз возникло желание переписать его на С++, но решил оставить так как есть. У меня он выполнялся два дня, все из-за проблем с диктатурой белорусским интернетом и перегруженностью сервера OverPass. Чтобы решить вторую проблему, нужно составить один запрос ко всем городам, но я еще не придумал как это нормально сделать.

Мой ответ на задачу

18981 буква



Что хочу сказать насчет правильности цифры: все упирается в качество данных самой OSM, то есть там есть места где, например, одну дорогу пересекает две линии границы, или где-нибудь на развязке граница проведена чуть-чуть не так, и в результате имеем лишнее/недостающее пересечение. Но это особенность конкретно этой не имеющей практического смысла задачи, в остальном OSM — сила.

Вторая задача


Сейчас давайте посчитаем количество заправок в пределах Москвы:
area[name="Москва"];
(
  node["amenity"="fuel"](area);
  way["amenity"="fuel"](area);
  relation["amenity"="fuel"](area);
);
out body;
>;
out skel qt;

Код 
import overpy

api = overpy.Overpass()
Data = api.query("""
area[name="Москва"];
(
  node["amenity"="fuel"](area);
  way["amenity"="fuel"](area);
  relation["amenity"="fuel"](area);
);
out body;
>;
out skel qt;
""")
print(len(Data.nodes)) #получаем общее количество узлов


Результат — 489 заправок: