Mencari pustaka Python (selain ArcPy) untuk geoprocessing seperti buffer? [Tutup]

16

Tidak termasuk ArcPy, apakah ada pustaka python yang dapat melakukan geoprocessing, seperti buffer / intersect, dengan shapefile?

Mingshu
sumber

Jawaban:

17

The Python GDAL / OGR Cookbook memiliki beberapa contoh kode untuk Buffer sebuah Geometri .

from osgeo import ogr

wkt = "POINT (1198054.34 648493.09)"
pt = ogr.CreateGeometryFromWkt(wkt)
bufferDistance = 500
poly = pt.Buffer(bufferDistance)
print "%s buffered by %d is %s" % (pt.ExportToWkt(), bufferDistance, poly.ExportToWkt())

dan Menghitung persimpangan antara dua Geometri

from osgeo import ogr

wkt1 = "POLYGON ((1208064.271243039 624154.6783778917, 1208064.271243039 601260.9785661874, 1231345.9998651114 601260.9785661874, 1231345.9998651114 624154.6783778917, 1208064.271243039 624154.6783778917))"
wkt2 = "POLYGON ((1199915.6662253144 633079.3410163528, 1199915.6662253144 614453.958118695, 1219317.1067437078 614453.958118695, 1219317.1067437078 633079.3410163528, 1199915.6662253144 633079.3410163528)))"

poly1 = ogr.CreateGeometryFromWkt(wkt1)
poly2 = ogr.CreateGeometryFromWkt(wkt2)

intersection = poly1.Intersection(poly2)

print intersection.ExportToWkt()

Geometri dapat dibaca dan ditulis ke shapefile dan berbagai format lainnya.

PolyGeo
sumber
14

Untuk mempermudah, Shapely: manual memungkinkan semua pemrosesan geometri PostGIS dengan Python.

Premis pertama dari Shapely adalah bahwa programmer Python harus dapat melakukan operasi geometri tipe PostGIS di luar RDBMS ...

Contoh pertama PolyGeo

from shapely.geometry import Point, LineString, Polygon, mapping
from shapely.wkt import loads  
pt = Point(1198054.34,648493.09)
# or
pt = loads("POINT (1198054.34 648493.09)")
bufferDistance = 500
poly = pt.buffer(bufferDistance)
print poly.wkt
'POLYGON ((1198554.3400000001000000 648493.0899999999700000, 1198551.9323633362000000 
# GeoJSON
print mapping(poly)
{'type': 'Polygon', 'coordinates': (((1198554.34, 648493.09), (1198551.9323633362, 648444.0814298352), (1198544.7326402017, 648395.544838992), ....}

Contoh poligon dari PolyGeo:

poly1 = Polygon([(1208064.271243039,624154.6783778917), (1208064.271243039,601260.9785661874), (1231345.9998651114,601260.9785661874),(1231345.9998651114,624154.6783778917),(1208064.271243039,624154.6783778917)])    
poly2 = loads("POLYGON ((1199915.6662253144 633079.3410163528, 1199915.6662253144 614453.958118695, 1219317.1067437078 614453.958118695, 1219317.1067437078 633079.3410163528, 1199915.6662253144 633079.3410163528)))"

intersection = poly1.intersection(poly2)
print intersection.wkt
print mapping(intersection) -> GeoJSON

Premis kedua adalah bahwa kegigihan, serialisasi, dan proyeksi peta fitur adalah signifikan, tetapi masalah ortogonal. Anda mungkin tidak perlu seratus pembaca dan penulis format GIS atau banyak proyeksi State Plane, dan Shapely tidak membebani Anda dengan mereka.

Jadi Anda menggabungkannya dengan modul Python lain untuk membaca atau menulis shapefile dan memanipulasi proyeksi sebagai osgeo.ogr, Fiona atau PyShp .
Mencari di Gis StackExchange, Anda dapat menemukan banyak contoh tapi saya beri Anda satu lagi untuk menggambarkan kombinasi dari shapely dan Fiona dan penggunaan persimpangan fungsi shapely () dan buffer () (Ini bisa dilakukan dengan PyShp).

Diberi dua bentuk polyline:

masukkan deskripsi gambar di sini

Hitung persimpangan (fungsi persimpangan () dari rupawan)

from shapely.geometry import Point, Polygon, MultiPolygon, MumtiPoint, MultiLineString,shape, mapping
import fiona
# read the shapefiles and transform to MultilineString shapely geometry (shape())
layer1 = MultiLineString([shape(line['geometry']) for line in fiona.open('polyline1.shp')])  
layer2 = MultiLineString([shape(line['geometry']) for line in fiona.open('polyline2.shp')])
points_intersect = layer1.intersection(layer2)

Simpan hasilnya sebagai shapefile baru

# schema of the new shapefile
schema = {'geometry': 'MultiPoint','properties': {'test': 'int'}}
# write the new shapefile (function mapping() of shapely)
with fiona.open('intersect.shp','w','ESRI Shapefile', schema) as e:
  e.write({'geometry':mapping(points_intersect), 'properties':{'test':1}})

Hasil:

masukkan deskripsi gambar di sini

Buffer individual points (fungsi buffer () dari rupawan)

 # new schema
 schema = {'geometry': 'Polygon','properties': {'test': 'int'}}
 with fiona.open('buffer.shp','w','ESRI Shapefile', schema) as e:
     for point in points:
          e.write({'geometry':mapping(point.buffer(300)), 'properties':{'test':1}})

Hasil

masukkan deskripsi gambar di sini

Buffer geometri MultiPoint

schema = {'geometry': 'MultiPolygon','properties': {'test': 'int'}}
points.buffer(300)
with fiona.open('buffer2.shp','w','ESRI Shapefile', schema) as e:
     e.write({'geometry':mapping(points.buffer(300)), 'properties':{'test':1}})

masukkan deskripsi gambar di sini

gen
sumber
9

Shapely memberikan akses python ke GEOS yang dapat melakukan buffer / memotong / dll. GEOS adalah perpustakaan yang paling banyak digunakan program OSGeo untuk melakukan operasi-operasi itu.

HeyOverThere
sumber
9

Berikut adalah daftar perangkat lunak geoprosesing Python saya.

  • Bentuknya, python
  • OGR, python
  • QGIS, pyqgis, python
  • SagaGIS, python
  • Rumput, python
  • spatialite, pyspatialite, python
  • PostreSQL / PostGIS, Psycopg, python
  • Proyek R, rpy2, python
  • Whitebox GAT, python -GeoScript, jython
Klewis
sumber
1

Perpustakaan geo-processing 'go to' saya adalah 'Remote Sensing and GIS Library' (RSGISLib). Mudah dipasang dan digunakan dan dokumentasinya sangat bagus. Ini memiliki fungsi untuk pemrosesan vektor dan raster - Saya sangat jarang harus mendekati gui. Itu dapat ditemukan di sini: http://rsgislib.org .

Contoh dalam hal ini adalah:

rsgislib.vectorutils.buffervector(inputvector, outputvector, bufferDist, force)

Perintah untuk buffer vektor dengan jarak yang ditentukan.

Dimana:

  • inputvector adalah string yang berisi nama vektor input
  • outputvector adalah string yang berisi nama vektor output
  • bufferDist adalah float yang menentukan jarak buffer, dalam unit peta
  • force adalah bool, menentukan apakah akan memaksa penghapusan vektor output jika ada

Contoh:

from rsgislib import vectorutils
inputVector = './Vectors/injune_p142_stem_locations.shp'
outputVector = './TestOutputs/injune_p142_stem_locations_1mbuffer.shp'
bufferDist = 1
vectorutils.buffervector(inputVector, outputVector, bufferDist, True)
Nathan Thomas
sumber