Das deutsche Python-Forum

Hallo.
Ich habe ein Raster (AusgZust) mit Höhenwerten des Ostsee-Gebietes vor der postglazialen Hebung der Skanden (d. h. 10.000 BP). Für den Zeitraum der folgenden 5000 Jahre nehme ich eine exponentiell abnehmende Hebungsrate an, d.h. der Hebungsbetrag im Jahr a beträgt h_max * Exp(r*a). Dieser muss dann immer auf die Höhenwerte des vorangegangenen Jahres (a-1) addiert werden.
Um die Hebung darzustellen, würde ich gerne mindestens 10 Zeitschritte (zwischen 10.000 und 5.000 BP bzw. zwischen 0 und 5000) berechnen. Da ich eine rekursive Gleichung habe, mache ich das mit einer for-Schleife:

Code: Alles auswählen

import arcpy
from arcpy import env
from arcpy.sa import *

env.workspace = "..."


Zielordner = "..."
AusgZust = Raster(Zielordner+"RasterCalc_Ausgangszustand")
Hebung_max = Raster(Zielordner+"Hebungsrate_max")
r = -0.000766147
Zeit = range(0,5000)

for a in Zeit:
	outRaster = AusgZust + Hebung_max * Exp(r*a)
	outRaster.save(Zielordner+"Zustand_vor_"+str(10000-a)+"_Jahren")
	AusgZust = Raster(Zielordner+"Zustand_vor_"+str(10000-a)+"_Jahren")

Das funktioniert zwar, das Problem ist nur, dass das eeewig viel (Arbeits-)speicher braucht, weil ich so 5000 Raster berechne, von denen ich eigentlich nur 10 brauche.

Meine Idee war deshalb, vor den Speicherbefehl eine if-Schleife zu setzen, damit nur bestimmte Zeitschritte abgespeichert werden, z.B.

Code: Alles auswählen

for a in Zeit:
	outRaster = AusgZust + Hebung_max * Exp(r*a)
	if a%500 == 0: outRaster.save(Zielordner+"Zustand_vor_"+str(10000-a)+"_Jahren")
	AusgZust = outRaster

Das funktioniert aber leider nicht...
Hat jemand eine Idee, wie ich das lösen kann?

@a.ha: Wieso sollte es weniger *Arbeits*speicher benötigen wenn Du nur alle 500 Zwischenergebnisse speicherst? Damit sparst Du nur Plattenplatz.

Ich sehe da so auf Anhieb nur Einsparungspotential wenn man statt zwei Arrays(?) für die Daten vor und nach dem Addieren nur eines nimmt auf das die Werte „in place” addiert werden. Das würde jedenfalls bei Numpy-Arrays funktionieren, ich kenn mich aber nicht mit Arcpy aus und was für Objekte das da zurück gibt.

Rechenzeit könnte man sparen wenn man versucht die Rekursion in eine geschlossene Formel zu bringen, so dass man die Jahre die einen interessieren direkt berechnen kann, ohne alle Jahre dazwischen auch berechnen zu müssen.

Hm, ja, das stimmt allerdings...
Mit Numpy-Arrays kenn ich mich jetzt wiederum gar nicht aus. Eine funktionierende if-Schleife wäre aber auch wegen dem Speicherplatz ganz hilfreich, sonst müsste ich von 5000 abgespeicherten Rastern (die werden einzeln gespeichert und nicht als Array) eben 4990 wieder löschen, was auch nochmal Zeit braucht...
Gibt es eine Möglichkeit, die Bedingung für meine if-Schleife so zu formulieren, dass ich bestimmte Zeitschritte bzw. Werte meines Läufers (z.B. 500,1000,1500,...) abfragen kann?
Ich dachte an

Code: Alles auswählen

if a in range(0,5000,500): outRaster.save(...)

oder so...
Und: brauche ich in Python eine else-Bedingung oder läuft die for-Schleife danach auch so einfach weiter?

@a.ha: solche Fragen kannst Du doch ganz einfach selbst ausprobieren. Zur if-Schleife solltest Du auch hier lesen.

Code: Alles auswählen

>>> Zeit = range(0,5000)
>>> for a in Zeit:
...   if a % 500 == 0: print a
...
0
500
1000
1500
2000
2500
3000
3500
4000
4500

ich weiß nicht, wo Du Probleme hast, aber das if tut einwandfrei.

Statt Strings mit + zusammenzustückeln solltest Du Dir Stringformatierung anschauen.

@Sirius3: Nein, die if-ABFRAGE (danke, ich hab das wirklich so gelernt, aber das war auch kein Informatik-Kurs, sondern Umweltsystemmodellierung - vielleicht legen die nicht so großen Wert auf Fach auf korrekte Fachtermini...) funktioniert nicht, ich bekomme dann eine Fehlermeldung, dass das erste Raster nicht überschrieben werden könnte - was ich nicht verstehe, denn das dürfte ja theoretisch gar nicht passieren.
Ich hab mich deswegen durch diverse Python-Seiten geklickt, aber nirgendwo explizit die Info gefunden, ob man eine else-Bedingung braucht oder nicht und dachte das ist vielleicht das Problem.
Danke für die Info mit der String-Formatierung.

Das deutsche Python-Forum

Arcpy, Berechnung von Rastern in for-Schleife

Arcpy, Berechnung von Rastern in for-Schleife

Re: Arcpy, Berechnung von Rastern in for-Schleife

Re: Arcpy, Berechnung von Rastern in for-Schleife

Re: Arcpy, Berechnung von Rastern in for-Schleife

Re: Arcpy, Berechnung von Rastern in for-Schleife