gdal_calc składnia kalkulatora rastrowego dla operatorów logicznych i innych funkcji

W dokumentacji gdal_calc podano kalkulator rastrowy wiersza poleceń ze składnią numpy . Później jest kilka przykładów, w których w jednym z nich:

gdal_calc.py -A input.tif --outfile = result.tif --calc = "A * (A> 0)" --NoDataValue = 0 - oznacza ustawione wartości od zera i poniżej na zero

Niestety nie ma przykładu operatorów logicznych, takich jak:

--calc = "A * (A> 0 i A> B)" - oznacza utrzymanie A, jeśli A większe zero i większe B, a resztę ustawić na zero

Na podstawie funkcji logicznych Numpy / Scipy spodziewałbym się pisać operatorów logicznych jako:

--calc = "A * logiczny_i (A> 0, A> B)"

Próbowałem tego i wydaje się, że działa, ale chciałbym mieć pewność, że jest to poprawne.

W podobny sposób, jeśli chcesz minimum A i B:

--calc = „A * (A <= B) + B * (A> B)”

Możesz po prostu napisać:

--calc = „minimum (A, B)”

Mój problem polega na tym, że nie mogę znaleźć żadnej książki kucharskiej, aby upewnić się, że dobrze to zrobię. Czy jest jakaś dobra książka kucharska z zaawansowanymi przykładami tego, co jest i nie jest możliwe w gdal_calc?

W źródle pliku gdal_calc.py obliczenia są wykonywane bezpośrednio przy użyciu eval:

myResult = eval(opts.calc, global_namespace, local_namespace)

To sugerowałoby, że każde poprawnie sformułowane wyrażenie, które również ocenia w linii poleceń, będzie działać. Zgodnie z dokumentacją możesz używać składni gdalnumerycznej z funkcjami +-/*i / lub numpy. Możesz przetestować swoje funkcje za pomocą małych atrapy tablic w interaktywnej powłoce, a następnie użyć tych samych wywołań w gdal_calc.

Należy pamiętać, że łączenie ze sobą wielu numpyfunkcji może spowodować powstanie tymczasowych macierzy w pamięci, które mogą znacznie zwiększyć wykorzystanie pamięci, szczególnie w przypadku dużych obrazów.

Możesz zajrzeć do dokumentacji numpy, aby uzyskać listę wszystkich funkcji: procedur . Te, których szukasz, są prawdopodobnie tutaj: matematyka lub tutaj: routines.logic .

Stąd pochodzą funkcje takie jak minimum, po prostu przestrzeń nazw jest już zaimportowana. Naprawdę to numpy.minimum itp

Kontynuując odpowiedź Benjamina, możesz użyć logical_or () lub logical_and (). Zobacz http://docs.scipy.org/doc/numpy/reference/routines.logic.html . Poniższy przykład działał dla mnie dobrze. To ustawia wszystkie wartości od 177 do 185 (włącznie) na 0, co jest następnie traktowane jako nodata.

gdal_calc.py -A input.tif --outfile=output.tif --calc="A*logical_or(A<=177,A>=185)" --NoDataValue=0

Miałem raster, w którym wartości mieściły się w zakresie od -1 do 3, gdzie zero jest poprawną liczbą. Miałem pewne problemy z wyrażeniem gdal_calc, więc zrobiłem to szybkie i wściekłe rozwiązanie.

#!/usr/bin/env python3

fileNameIn = "/tmp/geotiff/Global_taxonomic_richness_of_soil_fungi.tif"
fileNameOut = "/tmp/geotiff/Global_taxonomic_richness_of_soil_fungi.tiff"
dst_options = ['COMPRESS=DEFLATE',"PREDICTOR=3","TILED=YES"]
noDataValue = -3.4028234663852886e+38

from osgeo import gdal
import numpy

src_ds = gdal.Open(fileNameIn)
format = "GTiff"
driver = gdal.GetDriverByName(format)
dst_ds = driver.CreateCopy(fileNameOut, src_ds, False ,dst_options)

# Set location
dst_ds.SetGeoTransform(src_ds.GetGeoTransform())
# Set projection
dst_ds.SetProjection(src_ds.GetProjection())
srcband = src_ds.GetRasterBand(1)

dataraster = srcband.ReadAsArray().astype(numpy.float)
#Rplace the nan value with the predefiend noDataValue
dataraster[numpy.isnan(dataraster)]=noDataValue

dst_ds.GetRasterBand(1).WriteArray(dataraster)
dst_ds.GetRasterBand(1).SetNoDataValue(noDataValue)

dst_ds = None