Итак, вторая попытка ответить на ваш вопрос с помощью чистого решения GDAL.
Во-первых, GDAL (Библиотека абстракции геопространственных данных) изначально была просто библиотекой для работы с растровыми геопространственными данными, тогда как отдельная библиотека OGR предназначалась для работы с векторными данными. Однако две библиотеки теперь частично объединены и обычно загружаются и устанавливаются вместе под общим именем GDAL. Таким образом, решение действительно подпадает под OGR. У вас есть это в вашем исходном коде, поэтому я думаю, что вы знали это, но это важное различие, которое нужно помнить при поиске советов и подсказок.
Чтобы прочитать данные из векторного слоя, вам нужен исходный код:
from osgeo import ogr
shapefile = ogr.Open(shapefile)
layer = shapefile.GetLayer(0)
for i in range(layer.GetFeatureCount()):
feature = layer.GetFeature(i)
name = feature.GetField("NAME")
geometry = feature.GetGeometryRef()
print i, name, geometry.GetGeometryName()
Нам нужно создать новую функцию, прежде чем мы сможем записать ее в шейп-файл (или любой другой набор векторных данных). Для создания нового объекта нам сначала необходимо: - Геометрия - Определение объекта, которое, вероятно, будет включать определения полей. Используйте конструктор Geometry ogr.Geometry (), чтобы создать пустой объект Geometry. Определите геометрию для каждого типа по-разному (точка, линия, многоугольник и т. Д.). Так, например:
point = ogr.Geometry(ogr.wkbPoint)
point.AddPoint(10,20)
или
line = ogr.Geometry(ogr.wkbLineString)
line.AddPoint(10,10)
line.AddPoint(20,20)
line.SetPoint(0,30,30) #(10,10) -> (30,30)
Для определения поля
fieldDefn = ogr.FieldDefn('id', ogr.OFTInteger)
Теперь вы можете создать свой векторный слой. В данном случае квадратный многоугольник:
#create simple square polygon shapefile:
from osgeo import ogr
driver = ogr.GetDriverByName('ESRI Shapefile')
datasource = driver.CreateDataSource('YOUR_PATH')
layer = datasource.CreateLayer('layerName',geom_type=ogr.wkbPolygon)
#create polygon object:
myRing = ogr.Geometry(type=ogr.wkbLinearRing)
myRing.AddPoint(0.0, 0.0) #LowerLeft
myRing.AddPoint(0.0, 10.0) #UpperLeft
myRing.AddPoint(10.0, 10.0) #UpperRight
myRing.AddPoint(10.0, 0.0) #Lower Right
myRing.AddPoint(0.0, 0.0) #close ring
myPoly = ogr.Geometry(type=ogr.wkbPolygon)
myPoly.AddGeometry(myRing)
print ('Polygon area =',myPoly.GetArea()) #returns correct area of 100.0
#create feature object with point geometry type from layer object:
feature = ogr.Feature( layer.GetLayerDefn())
feature.SetGeometry(myPoly)
layer.CreateFeature(feature)
#flush memory - very important
feature.Destroy()
datasource.Destroy()
Мне посчастливилось читать и писать в слои. В частности, у меня есть код, который будет считывать слой шейп-файла, содержащий полилинии, и выводить геометрию каждого объекта в текстовые файлы (используемые в качестве входных данных для старой модели).
Кажется, это может быть полезно, чтобы получить все функции из ваших слоев.
Запись в другой слой не должна быть слишком сложной. Примерно так должно работать в теории:
Отсюда вы сможете получать данные о каждой функции и записывать новые функции на новый слой.
Дэн
источник