Contoh-contoh Skrip Python untuk geoprocessing shapefile tanpa menggunakan arcpy

33

Saya ingin menggunakan skrip Python yang tidak didasarkan pada arcpy untuk melakukan hal-hal seperti kueri sebuah shapefile dengan atribut, membuat layer baru dari seleksi, dan menghitung area poligon dan mengkonversi poligon menjadi poin.

Adakah yang punya contoh kode menggunakan modul atau pustaka Python lainnya? Saya dapat melakukan ini dengan mudah menggunakan arcpy tetapi saya ingin menjelajahi opsi lain.

sherpa
sumber
geopanda adalah teman Anda untuk file vektor. Rasterio untuk raster.
RutgerH

Jawaban:

54

Itu aneh, seolah-olah orang tiba-tiba menemukan kekuatan Python (tanpa ArcPy yang hanya satu modul Python antara lain), lihat misalnya pertanyaan Visualisasikan shapefile dengan Python :

  • pemrosesan geospasial dengan Python memiliki sejarah yang sangat panjang, jauh lebih tua dari Arcpy (atau arcgisscripting) -> tidak ada "meniru" kemampuan ArcPy di ​​sini, seperti yang dikatakan Paul, sebagian besar sudah ada sebelum ArcPy.
  • referensi untuk modul Python adalah Python Package Index ( Pypi ) dan ada bagian khusus: Topik :: Ilmiah / Teknik :: GIS
  • Anda dapat melakukan apa saja dengan modul-modul ini dan seringkali lebih mudah dan lebih cepat daripada ArcPy karena itu adalah Python murni (tidak ada kursor ...).
  • Shapely adalah salah satu modul ini untuk memproses geometri geospasial -> menghitung luas poligon dan mengonversi poligon menjadi titik ..
  • jika Anda ingin memproses layer vektor, ada osgeo / ogr , Fiona atau Pyshp (dan yang lainnya, yang kurang digunakan) -> query sebuah shapefile dengan atribut, buat layer baru dari seleksi, hitung luas poligon, konversikan poligon ke titik
  • untuk memproses raster, standarnya adalah osgeo / gdal
  • untuk analisis spasial, ada Pysal
  • untuk 3D, Anda dapat menggunakan modul Ilmiah lainnya seperti numpy atau scipy (algoritma 3D, grid, tetapi juga statistik, geostatistik, 2D atau 3D)
  • Dan saya tidak berbicara tentang mapnik , matplotlib / basemap , Geodjango dan ...

Anda dapat menggabungkan semua (Pysal dengan rupawan, ...) dan mencampurnya dengan modul Ilmiah lainnya.

Jadi untuk contoh Skrip Python, cari Pyshp Fiona, ogr, gdal atau bentuknya di gis.stackexchange atau internet (banyak contoh, tidak hanya dalam bahasa Inggris).)
Salah satunya dalam bahasa Perancis (skrip dan gambarnya universal!):
- Python: Menggunakan vektor dan lapisan raster dalam perspektif geologis, tanpa perangkat lunak GIS
yang lain dalam bahasa Inggris:
- GIS dengan Python, Shapely, dan Fiona
dan dalam bahasa Spanyol
- Penentuan area poligon tidak beraturan menggunakan koordinat titik
di gis.stackexchange
- Profil ketinggian 10 km setiap sisi garis
- Memperbarui Atribut menggunakan Pyshp
- Bagaimana cara membuat shapefile 3D dari raster?
- Skrip Python untuk mendapatkan perbedaan elevasi antara dua titik
- dll

Script yang disajikan oleh Aaron dapat ditulis lebih sederhana dengan Fiona yang hanya menggunakan kamus Python:

import fiona
with fiona.open('sites.shp', 'r') as input:
    with open('hw1a.txt', 'w') as output:
       for pt in input:
           id = pt['properties']['id']
           cover = pt['properties']['cover']
           x = str(point['geometry']['coordinates'][0])
           y = str(point['geometry']['coordinates'][21])
           output.write(id + ' ' + x + ' ' + y+ ' ' + cover + '\n')

dan jika Anda menggunakan Shapely sebagai tambahan:

from shapely.geometry import shape
with fiona.open('sites.shp', 'r') as input:
    with open('hw1a.txt', 'w') as output:
       for pt in input:
           id = pt['properties']['id']
           cover = pt['properties']['cover']
           x = str(shape(pt['geometry']).x)
           y = str(shape(pt['geometry']).y)
           output.write(id + ' ' + x + ' ' + y+ ' ' + cover + '\n')

Ada juga dua buku:

Pengembangan Geospasial Python dari Eric Westra.

masukkan deskripsi gambar di sini

Belajar Analisis Geospasial dengan Python of Joel Lawhead

masukkan deskripsi gambar di sini

Python juga digunakan sebagai bahasa scripting dalam aplikasi SIG lainnya seperti QGIS (Quantum GIS), GRASS GIS, gvSIG atau OpenJump atau pemodel 3D seperti Paraview (dan Blender juga!). Dan Anda dapat menggunakan sebagian besar modul geospasial di semua aplikasi ini (lihat Memvisualisasikan data QGIS dengan Blender )

gen
sumber
Apa hal Python ini yang Anda ucapkan;)
Nathan W
Fiona tampaknya melempar kesalahan DLL pada Windows.
multigoodverse
Bagaimana Anda menginstal Fiona? tidak masalah bagi saya
gen
19

Saya sangat merekomendasikan situs USU Geoprocessing dengan Python menggunakan Open Source GIS untuk membantu Anda memulai. Mereka terutama menggunakan perpustakaan GDAL / OGR selama latihan. Menginstal GDAL / OGR bisa menjadi sedikit tantangan, jadi entri blog ini mungkin dapat membantu Anda: Menginstal GDAL (dan OGR) untuk Python di Windows . Lihat juga Alternatif untuk menggunakan Arcpy di GIS.SE.

Contoh skrip geoprocessing opensource berikut (dari situs USU) digunakan untuk mengekstrak data atribut dan menulisnya ke file teks:

# import modules
import ogr, os, sys

# set the working directory
os.chdir('f:/data/classes/python/data')

# open the output text file for writing
file = open('hw1a.txt', 'w')

# get the shapefile driver
driver = ogr.GetDriverByName('ESRI Shapefile')

# open the data source
datasource = driver.Open('sites.shp', 0)
if datasource is None:
  print 'Could not open file'
  sys.exit(1)

# get the data layer
layer = datasource.GetLayer()

# loop through the features in the layer
feature = layer.GetNextFeature()
while feature:

  # get the attributes
  id = feature.GetFieldAsString('id')
  cover = feature.GetFieldAsString('cover')

  # get the x,y coordinates for the point
  geom = feature.GetGeometryRef()
  x = str(geom.GetX())
  y = str(geom.GetY())

  # write info out to the text file
  file.write(id + ' ' + x + ' ' + y + ' ' + cover + '\n')

  # destroy the feature and get a new one
  feature.Destroy()
  feature = layer.GetNextFeature()

# close the data source and text file
datasource.Destroy()
file.close()
Harun
sumber
4
.Destroyadalah nama metode yang luar biasa: p
Jason
5

Anda mungkin tertarik dengan GDAL / OGR .

GDAL digunakan untuk memproses raster sedangkan OGR digunakan untuk vektor. Keduanya adalah perpustakaan open source.

Jika Anda ingin menghapus ketergantungan pada ArcPy, Anda dapat meniru beberapa kemampuan dengan membaca informasi ke sebuah array dan menjalankan perhitungan Anda sendiri dengan Python murni.

Baru-baru ini saya melakukan ini dengan memilih titik dalam poligon, seperti yang terlihat di sini . Ini menggunakan algoritma pengecoran sinar untuk menentukan apakah suatu titik terletak di dalam poligon, mengingat koordinat titik poligon.

Paul
sumber
1
harap sertakan cukup esensi dari solusi yang dapat dipahami dan dipahami sebelum mengunjungi dan membaca halaman. Pada saatnya halaman itu mungkin tidak ada di alamat yang memberikan jawaban ini tidak terlalu berguna. :)
matt wilkie
1

Saya tidak pernah menggunakan ini secara pribadi, tetapi orang lain di kantor suka menggunakan rupawan: https://pypi.python.org/pypi/Shapely

Jason
sumber
Apakah Anda bisa memposting beberapa kode sampel menggunakan rupawan?
sherpa
5
Tautan hanya jawaban tidak membantu dalam jangka panjang, karena mereka pasti akan rusak. Harap sertakan informasi yang cukup tentang tujuan bahwa a) rumah baru dapat ditemukan kembali, dan b) inti dari solusi dapat dipahami dan dipahami sebelum pergi mengunjungi dan membaca halaman.
matt wilkie