Регистрация

[Python, Data Mining, OpenStreetMap, Открытые данные] Как найти количество всех букв на всех знаках вида «въезд в город Х» в стране? Точный способ ответить на такие вопросы

Ответить на тему
Автор Сообщение
news_bot ®

Стаж: 6 лет 9 месяцев
Сообщений: 27286

Создавать темы news_bot ® написал(а)
24-Авг-2020 20:31
Цитировать

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

Под катом покажу точное решение этой и других похожих задач, например: «Сколько АЗС находится в пределах Москвы?»
Общий метод решения
Если взглянуть на карту OpenStreetMap, то на ум сразу приходит следующая идея: а давайте для каждого города получим его границы и дороги, находящиеся внутри его границ, а потом найдем их пересечения, на которых будут стоять знаки! Как будем искать пересечения: берем отрезок границы, потом отрезок дороги и смотрим, пересекаются ли они (типичная геометрическая задача). И так пока не кончатся все отрезки и города.

Про архитектуру данных OSM

SPL
Вообще, есть три главных компонента: узлы, линии и отношения
Каждый элемент имеет свой 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 в питоне, это всё, что понадобится, чтобы достать данные.
Перейдем непосредственно к задаче
Для ее решения написан код на Питоне, увидеть его можно под спойлером. Просьба сильно не ругать за качество кода, это первый такой большой проект на нем.

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

SPL
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. Чтобы решить вторую проблему, нужно составить один запрос ко всем городам, но я еще не придумал как это нормально сделать.

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

SPL
18981 буква

Что хочу сказать насчет правильности цифры: все упирается в качество данных самой OSM, то есть там есть места где, например, одну дорогу пересекает две линии границы, или где-нибудь на развязке граница проведена чуть-чуть не так, и в результате имеем лишнее/недостающее пересечение. Но это особенность конкретно этой не имеющей практического смысла задачи, в остальном OSM — сила.
Вторая задача
Сейчас давайте посчитаем количество заправок в пределах Москвы:
area[name="Москва"];
(
  node["amenity"="fuel"](area);
  way["amenity"="fuel"](area);
  relation["amenity"="fuel"](area);
);
out body;
>;
out skel qt;

Код

SPL
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 заправок:

 ===========
 Источник:
habr.com
===========
Похожие новости: Теги для поиска: #_python, #_data_mining, #_openstreetmap, #_otkrytye_dannye (Открытые данные), #_tegi_nikto_ne_chitaet (теги никто не читает), #_python, #_overpass, #_overpy, #_openstreetmap, #_osm, #_python, #_data_mining, #_openstreetmap, #_otkrytye_dannye (
Открытые данные
)
Профиль  ЛС 
Показать сообщения:     
Ответить на тему

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

Текущее время: 22-Ноя 18:45
Часовой пояс: UTC + 5