Tipps und Tricks rund um ArcGIS - Heute: Automatisierung von Arbeitsprozessen mittels Python
Herzlich Willkommen bei der Blog-Serie „Tipps und Tricks rund um ArcGIS“
Wer kennt es nicht? - Wiederkehrende Aufgaben die viel Zeit in Anspruch nehmen und aus vielen „Klicks“ bestehen, gehören zum Alltagsgeschäft eines GIS-Experten. In unserem Szenario müssen 100 Rasterdatensätze, die in dem Koordinatensystem Gauß-Krüger vorliegen, in das Koordinatensystem ETRS89 transformiert werden. Unter normalen Umständen ist man nun einen Großteil des Tages beschäftigt, das Werkzeug Project Raster aufzurufen und jeweils 100 mal die entsprechenden Eingabedatensätze sowie die entsprechenden Parameter einzustellen, die für alle 100 Datensätze auch noch immer identisch sind.
Um sich die Bearbeitung von Massendatensätzen zu vereinfachen, kann man in ArcGIS for Desktop auf 2 Lösungswege setzen. Lösungsweg Nummer 1 ist der ModelBuilder, Lösungsweg Nummer 2 ist ein Python-Skript. In diesem Eintrag liegt der Fokus auf der Entwicklung eines Python-Skriptes. Ziel ist es, ein Skript zu erstellen, welches alle 100 Rasterdatensätze in einem Ordner automatisiert ausliest und als Ergebnis in das Koordinatensystem ETRS89 transformiert. Benötigt werden nur eine Entwicklungsumgebung – hier empfiehlt sich die kostenfreie Software PyScripter (http://sourceforge.net/projects/pyscripter/) – sowie ein PC mit Maus und Tastatur und installiertem und lizenziertem ArcGIS for Desktop… 😉
Zu Beginn eines Skriptes müssen zunächst alle benötigten Module importiert werden. Module ergänzen die Basisfunktionalität von Python um weitere, teils fachspezifische Funktionen. Mit dem Modul math erhält man z.B. Zugriff auf komplexe mathematische Formeln, während man mit dem Modul arcpy einen Zugriff auf alle Werkzeuge aus der ArcGIS Toolbox erhält.
Die erste Zeile in dem Skript lautet nun wie folgt:
import arcpy
Mit dieser Zeile gewährleisten wir, dass das Skript einen Zugriff auf alle Funktionen von ArcGIS hat.
Als nächstes geben wir die Arbeitsumgebung an, d.h. den Ordner, in welchem unsere Rasterdatensätze liegen. Dies geschieht über folgende Zeile:
arcpy.env.workspace = r"C:\temp\Ordner_mit_100_Rasterdaten"
In einem nächsten Schritt müssen alle Raster, die in dem Ordner abgelegt sind, ausgelesen und in eine Liste geschrieben werden.
rasterListe = arcpy.ListRasters()
Die Liste beinhaltet nun die Dateinamen aller Rasterdatensätze, die vom Skript im angegebenen Ordner gefunden worden sind. Als optionale Parameter könnte man in der obigen Zeile u.a. noch eine Wildcard vergeben, um nur nach bestimmten Rasterdatensätzen zu suchen.
Da wir nun alle Rasterdatensätze identifiziert und in einer Liste abgelegt haben, können wir diese nun über eine Schleife abarbeiten. Hierzu ergänzen wir im Skript folgende Zeilen Code:
for raster in rasterList:
arcpy.ProjectRaster_management(raster, r"ETRS89_" + raster, arcpy.SpatialReference(102329), "#", "#", r"DHDN_To_ETRS_1989_8_NTv2")
Der Code-Block besteht aus 2 wesentlichen Elementen: die for-Schleife sowie das Tool „Project Raster".
Die for-Schleife for raster in rasterList: arbeitet Element für Element aus der erzeugten Liste mit Rasterdaten ab. Jedes Raster wird hierbei in der Funktion
arcpy.ProjectRaster_management(raster, r"ETRS89_" + raster, arcpy.SpatialReference(102329), "#", "#", r"DHDN_To_ETRS_1989_8_NTv2")
bearbeitet und in das Zielkoordinatensystem umprojiziert. Schaut man sich die Funktion genauer an, so haben die einzelnen Parameter jeweils die folgende Bedeutung:
arcpy.ProjectRaster_management = Aufruf der Funktion
raster = Das jeweils aktuelle Raster aus der Liste.
r"ETRS89_" + raster = Ergebnis: Ursprungsrasters plus Präfix.
arcpy.SpatialReference(102329) = Zielkoordinatensystems / EPSG-Code
"#" = Platzhalter für einen optionalen Parameter.
r"DHDN_To_ETRS_1989_8_NTv2" = Transformationsmethode
Setzt man nun alle Zeilen zusammen, ergibt sich folgendes Skript, welches hier noch durch Kommentare ergänzt wurde:
Möchte man nun stattdessen Featureklassen oder Shapefiles umprojizieren, ersetzt man arcpy.ListRasters() durch arcpy.ListFeatures() und arcpy.ProjectRaster_management() durch arcpy.Project_management(). Die Transformation in ein anderes Koordinatensystem wird durch Änderung des EPSG-Codes sowie Umstellung der Transformation erreicht. Eine Liste sämtlicher verfügbarer Transformationen finden Sie unter http://arcg.is/1OjcyNO.
Dieser Blogartikel hat gezeigt, wie Sie sich mittels Python in ArcGIS die Bearbeitung von Massendatensätzen und immer wiederkehrenden Aufgaben erleichtern können.
Ich wünsche viel Spaß beim Nachbauen und Experimentieren mit Python!