OK, więc druga próba odpowiedzi na twoje pytanie za pomocą czystego rozwiązania GDAL.
Po pierwsze, GDAL (Geospatial Data Abstraction Library) była pierwotnie tylko biblioteką do pracy z rastrowymi danymi geoprzestrzennymi, podczas gdy oddzielna biblioteka OGR była przeznaczona do pracy z danymi wektorowymi. Jednak obie biblioteki są teraz częściowo połączone i są ogólnie pobierane i instalowane razem pod wspólną nazwą GDAL. Tak więc rozwiązanie naprawdę wchodzi w zakres OGR. Masz to w swoim początkowym kodzie, więc myślę, że o tym wiesz, ale ważne jest, aby pamiętać o tym, szukając wskazówek i wskazówek.
Aby odczytać dane z warstwy wektorowej, kod początkowy jest w porządku:
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()
Musimy utworzyć nową funkcję, zanim będziemy mogli zapisać ją do pliku kształtu (lub dowolnego innego zbioru danych wektorowych). Aby utworzyć nową operację, najpierw potrzebujemy: - Geometrii - Definicji operacji, która prawdopodobnie będzie zawierać definicje pól. Użyj konstruktora Geometrii ogr.Geometry (), aby utworzyć pusty obiekt Geometry. Zdefiniuj geometrię w inny sposób dla każdego typu (punkt, linia, wielokąt itp.). Na przykład:
point = ogr.Geometry(ogr.wkbPoint)
point.AddPoint(10,20)
lub
line = ogr.Geometry(ogr.wkbLineString)
line.AddPoint(10,10)
line.AddPoint(20,20)
line.SetPoint(0,30,30) #(10,10) -> (30,30)
Do definicji pola
fieldDefn = ogr.FieldDefn('id', ogr.OFTInteger)
Teraz możesz utworzyć warstwę wektorową. W tym przypadku kwadratowy wielokąt:
#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()
Miałem trochę szczęścia czytając i pisząc na warstwach. W szczególności mam kod, który odczyta warstwę pliku kształtu zawierającą polilinie i wyśle geometrię każdej operacji do plików tekstowych (używanych jako dane wejściowe dla starego modelu).
Wydaje się, że przydatne może być pobranie każdej funkcji z twoich warstw.
Pisanie na innej warstwie nie powinno być zbyt skomplikowane. Coś takiego powinno działać teoretycznie:
Odtąd powinieneś być w stanie pobrać dane z każdej funkcji i zapisać nowe funkcje na nowej warstwie.
Dan
źródło