Sintaxis de la calculadora ráster gdal_calc para operadores lógicos y otras funciones

En la documentación de gdal_calc se indica Calculadora ráster de línea de comando con sintaxis numpy . Más adelante hay algunos ejemplos donde en uno de ellos:

gdal_calc.py -A input.tif --outfile = result.tif --calc = "A * (A> 0)" --NoDataValue = 0 - significa establecer valores de cero e inferiores a nulo

Lamentablemente, no hay ejemplos de operadores lógicos como:

--calc = "A * (A> 0 y A> B)" - significa mantener A si A mayor cero y mayor B y establecer el resto en nulo

Basado en las funciones lógicas de Numpy / Scipy , esperaría escribir operadores lógicos como:

--calc = "A * lógico_y (A> 0, A> B)"

Intenté esto y parece funcionar, pero me gustaría estar seguro de que es correcto.

De manera similar si desea un mínimo de A y B:

--calc = "A * (A <= B) + B * (A> B)"

Solo puedes escribir:

--calc = "mínimo (A, B)"

Mi problema es que no puedo encontrar ningún libro de cocina para asegurarme de hacerlo bien. ¿Hay algún buen libro de cocina con ejemplos avanzados de lo que es y no es posible con gdal_calc?

En la fuente de gdal_calc.py, el cálculo se realiza directamente usando eval:

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

Eso sugeriría que cualquier expresión bien formada que también se evalúe en la línea de comando funcionará. De acuerdo con la documentación, puede usar la sintaxis gdalnumeric con +-/*y / o numpyfunciones. Puede probar sus funciones utilizando pequeñas matrices ficticias en el shell interactivo, luego use las mismas llamadas en gdal_calc.

Tenga en cuenta que encadenar varias numpyfunciones es probable que produzca arreglos temporales en memoria que pueden aumentar sustancialmente el uso de memoria, especialmente cuando se trata de imágenes grandes.

Puede consultar la documentación numpy para obtener una lista de todas las funciones: rutinas . Es probable que los que buscas estén aquí: matemáticas o aquí: rutines.logic .

Aquí es de donde provienen funciones como el mínimo, es solo que el espacio de nombres ya está importado. Realmente, es numpy.minimum, etc.

Siguiendo con la respuesta de Benjamin, puede usar logical_or () o logical_and (). Ver http://docs.scipy.org/doc/numpy/reference/routines.logic.html . El siguiente ejemplo funcionó muy bien para mí. Esto establece todos los valores entre 177 y 185 (inclusive) en 0, que luego se trata como nodata.

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

Tenía un ráster donde los valores oscilaban entre -1 y 3, donde cero es un número válido. Tuve algunos problemas para hacer una expresión gdal_calc, así que hice esta solución rápida y furiosa.

#!/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