gdalwarp erzeugt leere Raster, wenn sie aus dem Python-Skript aufgerufen werden, jedoch nicht über die Befehlszeile

8

Ich versuche, einige Rasterdaten (alle .tif) mit gdalwarp in Python bis zu einer gröberen Auflösung (von 0,05 bis 0,25 Grad) zusammenzufassen, aber der Befehl funktioniert nicht. Anstatt eine Ausgabe mit einem großen Wertebereich zu erhalten, sind alle Werte der Ausgabe 0. Die Auflösung und die Pixeltiefe / der Pixeltyp sind korrekt, die Werte jedoch nicht.

Hier ist die Dokumentation zum Befehl gdalwarp: http://www.gdal.org/gdalwarp.html

Ich habe zwei Eingabedateien, die ich mit einer Auflösung von bis zu 0,25 Grad zusammenfassen möchte, um mehrere Ausgaben zu erzeugen:

  • 'NDVI_raster': Die erste Eingabe ist ein 16-Bit-Raster mit Vorzeichen, das NDVI darstellt, mit Werten zwischen -10.000 und 10.000 und Knotenwerten von -15.000.

  • 'nodata_mask': Die zweite ist eine NoData-Maske, 32-Bit-Float, wobei 1 = Werte für "gute" Daten und 0 = NoData.

Mit 'NDVI_raster' als Eingabe möchte ich 7 verschiedene Ausgaben erzeugen, die jeweils eine andere Statistik darstellen. Dazu rufe ich gdalwarp siebenmal auf und setze die Resampling-Methode (-r) jedes Mal auf einen der folgenden Werte: Durchschnitt, Modus, Max, Min, Median, Q1, Q2. Ich werde die Ausgänge NDVI_ave.tif, NDVI_mode.tif usw. aufrufen. Im Moment verwende ich GDAL 1.10.1, das nur Durchschnitt und Modus zulässt, also teste ich jetzt nur mit diesen beiden Statistiken.

Mit 'nodata_mask' als Eingabe möchte ich letztendlich eine QAL (Qualitätssicherungsschicht) erstellen. Zu diesem Zweck verwende ich gdalwarp, wobei der Resampling-Modus auf "Durchschnitt" eingestellt ist, um bis zu 0,25 Grad zu aggregieren. Dies führt dazu, dass jedes Pixel das Verhältnis von guten Pixeln zu Gesamtpixeln von der Eingabe darstellt. Nennen wir den Ausgang QAL.

Folgendes ist in meinem Code enthalten (am Beispiel des Modus für die erste Eingabe):

os.system('gdalwarp -tr .25 .25 -r mode -srcnodata -15000 %s %s' % (NDVI_raster, NDVI_mode))

Und für die QS-Schicht:

os.system('gdalwarp -tr .25 .25 -r average -srcnodata -15000 %s %s' % (nodata_mask, QAL))

Das Ergebnis sind Raster mit der richtigen Auflösung, Projektion und Pixeltiefe, aber die Pixelwerte sind alle 0.

Weiß jemand, der mit Python / GDAL vertraut ist, was los ist?


Wenn ich den Befehl gdalwarp über die Befehlszeile (Linux) aufrufe, erhalte ich das gewünschte Ergebnis. Wenn ich mit os.system gdalwarp von Python aus aufrufe, erhalte ich leere Raster. Vielleicht stimmt etwas mit meinen GDAL / Python-Bindungen nicht?


Anstatt den Befehl über os.system aufzurufen, habe ich einen Unterprozess verwendet. Das Tool über diese Methode scheint ebenfalls reibungslos zu funktionieren, aber das Ergebnis ist dasselbe: ein Raster voller Nullen.


Ich habe versucht, den gdalwarp-Aufruf in ein Bash-Shell-Skript einzufügen und dieses Shell-Skript von Python aus aufzurufen, aber das Ergebnis ist eine Reihe von -1s anstelle von 0s. Seltsamerweise hatte ich es schon einmal getestet und bin mir ziemlich sicher, dass es funktioniert hat, aber der Test wurde von meinem Server gelöscht und jetzt kann ich ihn aus irgendeinem Grund nicht neu erstellen.


Wenn Sie den Befehl gdalwarp in ein Bash-Shell-Skript einfügen und dieses Shell-Skript dann über die Befehlszeile aufrufen, erhalte ich das gewünschte Ergebnis. Das Aufrufen des gleichen Shell-Skripts aus Python ist jedoch nicht möglich. Es sieht so aus, als ob etwas mit Python nicht stimmt, aber was und wie behebe ich es?

user20408
quelle
2
Es kann kein verbindliches Problem sein, Sie haben sie nicht verwendet. Sie führen den Befehl gdalwarp aus. Sie sollten den Fehlerstatus des Befehls gdalwarp überprüfen, siehe: stackoverflow.com/questions/3791465/…
Zoltan
2
Um zu erweitern, was Zoltan sagt, rufen Sie im Wesentlichen dasselbe Programm (gdalwarp) auf. os.systemStartet tatsächlich eine Shell und führt den angegebenen Befehl aus. Zusätzlich zum Rückgabewert sollten Sie wahrscheinlich die Variablen überprüfen, die Sie übergeben möchten (z. B. Probleme beim Zitieren? Usw.)
Evil Genius
Sie sollten auch prüfen , eine Bewegung os.systemzu subprocessdem das aktuelle Modul und enthält eine Reihe von (Sicherheits-) Verbesserungen.
Kersten
@Zoltan welche Fehlermeldung? Ich habe keine Fehlermeldung erhalten, da der Befehl technisch funktioniert hat, nur die Ausgabe war nicht richtig (Werte von 0) .plz siehe mein Update. Vielen Dank!
user20408
2
Ein Problem ist, dass Sie meiner Meinung nach -ot Float32für die Qualitätsmaske etwas Ähnliches angeben müssen , da Sie mit einem ganzzahligen Raster keine Brüche darstellen können.
Nat Wilson

Antworten:

2

Sie sagen nicht, wie Sie Ihr Python / GDALWARP-Skript starten. Ich habe festgestellt, dass ein Cronjob nicht immer dieselbe Umgebung wie meine Befehlszeilenumgebung hat. Ich musste anfangen, Laufzeitumgebungen für diese Art von Skripten zu erstellen. Wenn Sie Ihr Skript beispielsweise von einem Symbol auf Ihrem Desktop aus starten, hat es möglicherweise nicht dieselbe Laufzeitumgebung wie Ihre Befehlszeilenumgebung. "PYTHONPATH" kann eine der Umgebungsvariablen sein, die Sie festlegen müssen. Außerdem müssen Sie möglicherweise Variablen für gdalwarp festlegen. Schließlich befinden sich Ihre Datendateien möglicherweise nicht am richtigen Speicherort. Möglicherweise müssen Sie einen absoluten Pfad wie / xxx / xxxx / NDVI_raster festlegen oder tilda ~ / NDVI_raster verwenden. Wie bei PYTHONPATH müssen Sie möglicherweise auch PATH und andere Umgebungsvariablen einrichten ."exportieren oder quellen " Sie die gleichen Einstellungen am Anfang Ihres Skripts.

Greg
quelle
1

Ich hatte auch dieses Problem. In meinem Fall stellte ich schließlich fest, dass ich das Image auf der Festplatte gerade mit Python- gdalBindungen erstellt hatte, das gdal.DatasetIn-Memory-Objekt jedoch nicht geschlossen hatte , sodass das Schreiben auf die Festplatte nur teilweise abgeschlossen war. Seltsamerweise ist der einzige Weg, einen gdal.Datasetin Python zu schließen, der folgende: del variable_name_of_dataset- so hässlich!

Es Cgibt eine GDALClose()Methode, die derzeit nicht von der Python GDAL-API implementiert wird, aber Folgendes Rasteriotut: https://github.com/sgillies/rasterio/blob/876b9a1e2bf04e349b485e05ebc4a8674ace3cf0/rasterio/_io.pyx#L1463

Siehe auch: Warum ein Dataset in GDAL Python schließen?

Hazzles
quelle