В настоящее время я работаю над проектом, в котором мне нужно построить топологическую сеть из геометрических элементов, которые я нахожу в шейп-файлах. До сих пор, используя проект с открытым исходным кодом Бена Рейли, мне удалось преобразовать линейные строки в ребра сети x, а также обнаружить близкие объекты (как говорят другие линейные строки) и добавить их в ближайшую точку, чтобы я мог запускать алгоритмы кратчайшего пути.
Но это хорошо для одного шейп-файла. Однако теперь мне нужно соединить объекты из разных шейп-файлов в большой граф networkx. Так, например, если точка находится внутри многоугольника, я бы подключил ее (под подключением я имею в виду добавление ребра networkx - add_edge (g.GetPoint (1), g.GetPoint (2)) с точкой в следующем шейп-файле, который также находится внутри многоугольника, который имеет похожий атрибут (скажем, ID). Обратите внимание, что многоугольники в разных shps имеют одни и те же идентификаторы, а не координаты. Точки, попадающие в многоугольники, также не имеют одинаковые координаты.
Мое решение этой проблемы состояло в том, чтобы определить точку, которая находится в многоугольнике, сохранить ее, найти точку в следующем шейп-файле, который находится в многоугольнике с тем же идентификатором, а затем добавить ребро networkx между ними.
Как узнать, находится ли точка внутри многоугольника? Ну, есть хорошо известный алгоритм: алгоритм RayCasting , который делает это. Это то, где я действительно застрял, потому что для реализации алгоритма мне нужны координаты многоугольника, и я не знаю, как получить к ним доступ прямо сейчас, даже после просмотра документации по геометрии OGR. Итак, вопрос, который я задаю, - как получить доступ к точкам многоугольника или координатам ИЛИ есть ли более простой способ определить, попадает ли точка в многоугольник? Используя python с библиотекой osgeo.ogr, я написал следующее:
if g.GetGeometryType() == 3: #polygon
c = g.GetDimension()
x = g.GetPointCount()
y = g.GetY()
z = g.GetZ()
см. изображение для лучшего понимания моей проблемы.
[РЕДАКТИРОВАТЬ] До сих пор я пытался сохранить все объекты многоугольника в списке, с которым я бы затем сравнил первую и последнюю точку строки . Но пример Паоло связан с использованием ссылки на объект Point и ссылки на объект Polygon, которая не будет работать со ссылкой на объект линии, поскольку не вся линия находится внутри многоугольника, а, скорее, первая или последняя точка его строки.
[EDIT3] Создание нового точечного объекта Geometry из координат первой и последней точки линейной линии и последующее использование этой точки для сравнения с объектами геометрии многоугольника, сохраненными в списке, работает нормально:
for findex in xrange(lyr.GetFeatureCount()):
f = lyr.GetFeature(findex)
flddata = getfieldinfo(lyr,f,fields)
g = f.geometry()
if g.GetGeometryType() == 2:
for j in xrange(g.GetPointCount()):
if j == 0 or j == g.GetPointCount():
point = ogr.Geometry(ogr.wkbPoint)
point.AddPoint(g.Getx(j),g.GetY(j))
if point.Within(Network.polygons[x][0].GetGeometryRef()):
print g.GetPoint(j)
Я не знаком с networkx, но если я правильно понял ваш вопрос, вы можете использовать shapely и OGR lib для нахождения точки в полигоне из shapefile. Вот один пример того, как работает поиск одной точки (2000,1200) с любым полигоном из одного шейп-файла. Для результата он печатает координаты этого многоугольника.
Я надеюсь, что это помогает.
источник