У лучшего в мире управляющего по имени Пенультимо родилась очередная гениальнейшая идея, peализовать которую вам и предстоит. Он верит, что поток туристов на Исла-де-Эдукадос повысится, если он сможет рассказать всему миру, как же много замечательных дорожных знаков с длинными надписями eсть у них на острове. Вам предлагается придумать алгоритм, позволяющий подсчитать суммарное количество букв на всех знаках Въезд в город Х на острове, а затем применить полученные знания для подсчёта аналогичной метрики для Республики Беларусь. Обратите внимание язык, используемый для обозначения населённых пунктов, а также тот факт, что въездов в город может быть несколько. Пенультимо также приветствует инициативность, так что можете исследовать этот вопрос для отдельных областей, провести сравнение с количеством людей, проживающих в области, а также провести любые другие исследования, которые покажутся Вам интересными.
Под катом покажу точное решение этой и других похожих задач, например: Сколько АЗС находится в пределах Москвы?
Общий метод решения
Если взглянуть на карту OpenStreetMap, то на ум сразу приходит следующая идея: а давайте для каждого города получим его границы и дороги, находящиеся внутри его границ, а потом найдем их пересечения, на которых будут стоять знаки! Как будем искать пересечения: берем отрезок границы, потом отрезок дороги и смотрим, пересекаются ли они (типичная геометрическая задача). И так пока не кончатся все отрезки и города.
Каждый элемент имеет свой ID, а также свои теги.
- Узел это просто точка на карте, имеющая кроме ID также широту и долготу
- Линия это отсортированный список узлов, который представляет контур здания или отдельную улицу
- Отношение это список, состоящий из узлов, линий или других отношений, представляющий логические или географические связи между объектами на карте
OverPass
OverPass Это API для получения данных из OpenStreetMap. Оно имеет свой язык составления запросов, про него подробно можно почитать В этой статье.
Для того чтобы было проще и удобнее составлять запросы, есть инструмент Overpass-turbo, где результат выполнения запроса можно посмотреть в удобном и интерактивном виде.
Использование OverPass API в Python
Для обработки данных из OSM в Питоне можно использовать пакет Overpy в качестве оболочки.
Для отправки запросов и получения данных нужно проделать следующее:
import overpyapi = 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 overpyapi = 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].lonYb=Data.ways[0].nodes[1].latNodeID=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 = 0A = [0, 0]B = [0, 0]C = [0, 0]D = [0, 0]amount = 0progressbar = 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 = 0print('Total', Total) #и вывод всего
Заметки по коду
Я достаточно долго составлял запрос, подбирая разные типы дорог чтобы было меньше считать и чтобы не пропустить знаки. В итоговом запросе просто убраны те дороги, на которых нет знаков, например residential, service, footway, track и т. п.
Список городов я спарсил с википедии и сохранил в формате.тхт
Код выполняется достаточно долго, у меня даже один раз возникло желание переписать его на С++, но решил оставить так как есть. У меня он выполнялся два дня, все из-за проблем с
18981 буква
Что хочу сказать насчет правильности цифры: все упирается в качество данных самой OSM, то есть там есть места где, например, одну дорогу пересекает две линии границы, или где-нибудь на развязке граница проведена чуть-чуть не так, и в результате имеем лишнее/недостающее пересечение. Но это особенность конкретно этой не имеющей практического смысла задачи, в остальном OSM сила.
Вторая задача
Сейчас давайте посчитаем количество заправок в пределах Москвы:
area[name="Москва"];( node["amenity"="fuel"](area); way["amenity"="fuel"](area); relation["amenity"="fuel"](area););out body;>;out skel qt;
import overpyapi = 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 заправок: