Przykłady skryptu Python dla geoprzetwarzania plików kształtów bez użycia arcpy

33

Chciałbym użyć skryptu w języku Python, który nie jest oparty na arcpy, do wykonywania zadań takich jak wyszukiwanie pliku kształtu według atrybutów, tworzenie nowej warstwy z zaznaczenia oraz obliczanie obszarów wielokąta i konwertowanie wielokątów na punkty.

Czy ktoś ma jakieś przykłady użycia innych modułów lub bibliotek Python? Jestem w stanie to zrobić z łatwością za pomocą arcpy, ale chciałem zbadać inne opcje.

Sherpas
źródło
geopandas jest twoim przyjacielem dla plików wektorowych. Rasterio dla rastra.
RutgerH

Odpowiedzi:

54

To dziwne, jakby ludzie nagle odkryli moc Pythona (bez ArcPy, który jest tylko jednym z modułów Pythona), zobacz na przykład pytanie Wizualizuj plik kształtów w Pythonie :

  • przetwarzanie geoprzestrzenne w Pythonie ma bardzo długą historię, znacznie starszą niż Arcpy (lub arcgisscripting) -> nie „naśladuje” możliwości ArcPy tutaj, jak mówi Paul, większość była już tam przed ArcPy.
  • odnośnikiem do modułów Python jest Indeks pakietów Pythona ( Pypi ) i jest dedykowana sekcja: Temat :: Nauka / Inżynieria :: GIS
  • możesz zrobić wszystko z tymi modułami i często jest to łatwiejsze i szybsze niż ArcPy, ponieważ jest to czysty Python (bez kursorów ...).
  • Shapely jest jednym z tych modułów do przetwarzania geometrii geoprzestrzennych -> oblicz obszary wielokąta i konwertuj wielokąty na punkty.
  • jeśli chcesz przetwarzać warstwy wektorowe, istnieje osgeo / ogr , Fiona lub Pyshp (i inne, rzadziej używane) -> przeszukuj plik kształtów według atrybutów, twórz nową warstwę z zaznaczenia, obliczaj obszary wielokąta, konwertuj wielokąty na punkty
  • w przypadku przetwarzania rastrów standardem jest osgeo / gdal
  • do analizy przestrzennej istnieje Pysal
  • do 3D można użyć innych modułów naukowych, takich jak numpy lub scipy (algorytmy 3D, siatki, ale także statystyki, geostatystyka, 2D lub 3D)
  • I nie mówię o mapniku , matplotlib / basemap , Geodjango i ...

Możesz połączyć wszystkie (Pysal z foremnym, ...) i mieszać je z innymi modułami naukowymi.

Tak więc dla przykładów skryptu Python wyszukaj Pyshp Fiona, ogr, gdal lub shapely w gis.stackexchange lub w Internecie (wiele przykładów, nie tylko w języku angielskim).)
Jeden z nich w języku francuskim (skrypty i liczby są uniwersalne!):
- Python: Korzystanie z warstw wektorowych i rastrowych w perspektywie geologicznej, bez oprogramowania GIS
i inne w języku angielskim:
- GIS z Python, Shapely i Fiona
oraz w języku hiszpańskim
- Określanie obszarów nieregularnych wielokątów za pomocą współrzędnych wierzchołków
w gis.stackexchange
- Profil wysokości 10 km po każdej stronie linii
- Aktualizowanie atrybutów za pomocą Pyshp
- Jak utworzyć plik kształtu 3D z rastra?
- Skrypt Python do uzyskiwania różnicy wysokości między dwoma punktami
- itp

Skrypt zaprezentowany przez Aarona można napisać prościej za pomocą Fiony, która używa tylko słowników 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')

a jeśli dodatkowo używasz zgrabnie:

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')

Istnieją również dwie książki:

Python Geospatial Development of Eric Westra.

wprowadź opis zdjęcia tutaj

Nauka analizy geoprzestrzennej z Pythonem Joela Lawhead'a

wprowadź opis zdjęcia tutaj

Python jest również używany jako język skryptowy w innych aplikacjach GIS, takich jak QGIS (Quantum GIS), GRASS GIS, gvSIG lub OpenJump lub modelerach 3D, takich jak Paraview (i Blender również!). I możesz używać większości modułów geoprzestrzennych we wszystkich tych aplikacjach (zobacz Wizualizacja danych QGIS za pomocą Blendera )

17 obrotów
źródło
O czym mówisz w Pythonie;)
Nathan W
Wygląda na to, że Fiona zgłasza błąd DLL w systemie Windows.
multigoodverse
Jak zainstalowałeś Fionę? nie ma dla mnie problemu
gen
19

Gorąco polecam stronę USU Geoprocessing with Python using Open Source GIS, aby zacząć. Używają przede wszystkim biblioteki GDAL / OGR podczas ćwiczeń. Instalacja GDAL / OGR może być trudnym wyzwaniem, więc ten wpis na blogu może być dla Ciebie pomocny: Instalowanie GDAL (i OGR) dla Pythona w systemie Windows . Sprawdź także alternatywy dla korzystania z Arcpy na GIS.SE.

Poniższy przykład skryptu geoprzetwarzania typu open source (z witryny USU) służy do wyodrębnienia danych atrybutów i zapisania ich w pliku tekstowym:

# 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()
Aaron
źródło
4
.Destroyto niesamowita nazwa metody: p
Jason
5

Możesz być zainteresowany GDAL / OGR .

GDAL służy do przetwarzania rastrów, podczas gdy OGR jest używany do wektorów. Obie są bibliotekami typu open source.

Jeśli chcesz usunąć pewną zależność od ArcPy, możesz naśladować niektóre możliwości, odczytując informacje do tablicy i uruchamiając własne obliczenia w czystym języku Python.

Ostatnio zrobiłem to, wybierając punkty w wielokącie, jak widać tutaj . Wykorzystuje algorytm odlewania promieniowego, aby ustalić, czy punkt leży w wielokącie, biorąc pod uwagę współrzędne wierzchołków wielokąta.

Paweł
źródło
1
prosimy o uwzględnienie wystarczającej części esencji rozwiązania, którą można uchwycić i zrozumieć przed odwiedzeniem i przeczytaniem strony. Z czasem ta strona prawdopodobnie nie będzie pod tym adresem, co czyni tę odpowiedź niezbyt przydatną. :)
matt wilkie
1

Nigdy nie korzystałem z tego osobiście, ale inni w biurze lubią używać foremnie: https://pypi.python.org/pypi/Shapely

Jason
źródło
Czy jest szansa, że ​​możesz zamieścić przykładowe kody przy użyciu foremnego?
sherpas
5
Tylko linki nie są pomocne na dłuższą metę, ponieważ nieuchronnie się psują. Podaj wystarczającą ilość informacji o miejscu docelowym, aby: a) jego nowy dom mógł zostać ponownie odkryty, oraz b) esencję rozwiązania można uchwycić i zrozumieć przed odwiedzeniem i przeczytaniem strony.
matt wilkie