¿Qué herramientas de alisado/generalización de tramas existen?

Question

Tengo un MDE que me gustaría suavizar o generalizar para eliminar los extremos topográficos (cortar los picos y rellenar los valles). Idealmente, también me gustaría tener control sobre el radio o el nivel de "desenfoque". Al final, necesitaré un conjunto de rásteres que vayan de ligeramente borrosos a realmente borrosos. (En teoría, el más borroso sería un raster constante de la media aritmética de todos los valores).

¿Existen herramientas o métodos que pueda utilizar (basados en Esri, GDAL, GRASS)?

¿Necesito hornear en casa mi propio Desenfoque gaussiano ¿Rutina?

¿Podría utilizar un filtro de paso bajo (por ejemplo, el filtro de ArcGIS filtrar ), y si es así, ¿tendría que ejecutarlo un montón de veces para conseguir el efecto de un radio grande?

Answer 1

He estado explorando la herramienta de SciPy señal.convolve enfoque (basado en este libro de cocina ), y estoy teniendo un buen éxito con el siguiente fragmento:

import numpy as np
from scipy.signal import fftconvolve

def gaussian_blur(in_array, size):
    # expand in_array to fit edge of kernel
    padded_array = np.pad(in_array, size, 'symmetric')
    # build kernel
    x, y = np.mgrid[-size:size + 1, -size:size + 1]
    g = np.exp(-(x**2 / float(size) + y**2 / float(size)))
    g = (g / g.sum()).astype(in_array.dtype)
    # do the Gaussian blur
    return fftconvolve(padded_array, g, mode='valid')

Utilizo esto en otra función que lee/escribe GeoTIFFs float32 a través de GDAL (sin necesidad de reescalar a 0-255 bytes para el procesamiento de la imagen), y he estado utilizando tamaños de píxel de tentativa (por ejemplo, 2, 5, 20) y tiene una salida realmente agradable (visualizado en ArcGIS con 1:1 píxel y rango constante min/max):

Nota: esta respuesta se ha actualizado para utilizar un método mucho más rápido basado en la FFT señal.fftconvolve función de procesamiento.

Answer 2

El desenfoque gaussiano es sólo una media focal ponderada. Se puede recrear con gran precisión con una secuencia de medias circulares de corta distancia (no ponderadas): es una aplicación de la Teorema del límite central .

Tienes muchas opciones. "Filtro" es demasiado limitado -sólo sirve para vecindarios de 3 x 3-, así que no se moleste con él. La mejor opción para grandes MDE es llevar el cálculo fuera de ArcGIS en un entorno que utiliza Transformadas Rápidas de Fourier: hacen los mismos cálculos focales pero (en comparación) lo hacen increíblemente rápido. (GRASS tiene un Módulo FFT . Está pensada para el procesamiento de imágenes, pero quizá puedas ponerla al servicio de tu DEM si puedes reescalarla con una precisión razonable en el rango 0..255). Si no es así, dos soluciones al menos vale la pena considerarlo:

1. Crear un conjunto de pesos de vecindad para aproximar un desenfoque gaussiano para una vecindad considerable. Utiliza pases sucesivos de este desenfoque para crear tu secuencia de MDE cada vez más suaves.

   (Los pesos se calculan como exp(-d^2/(2r)) donde d es la distancia (en celdas si se quiere) y r es el radio efectivo (también en celdas). Deben calcularse dentro de un círculo que se extienda al menos hasta 3r . Después de hacerlo, divide cada peso por la suma de todos ellos para que al final sumen 1).

2. Alternativamente, olvídate de la ponderación; simplemente ejecuta una media focal circular repetidamente. He hecho exactamente esto para estudiar cómo cambian las mallas derivadas (como la pendiente y el aspecto) con la resolución de un MDE.

Ambos métodos funcionarán bien, y después de las primeras pasadas habrá poco que elegir entre los dos, pero hay rendimientos decrecientes: el radio efectivo de n Las medias focales sucesivas (todas utilizando el mismo tamaño de vecindad) es sólo (aproximadamente) la raíz cuadrada de n veces el radio de la media focal. Por lo tanto, para grandes cantidades de desenfoque, querrá empezar de nuevo con una vecindad de radio grande. Si utiliza una media focal no ponderada, ejecute 5-6 pases sobre el MDE. Si utiliza pesos aproximadamente gaussianos, sólo necesitará una pasada: pero tendrá que crear la matriz de pesos.

En efecto, este enfoque tiene como valor límite la media aritmética del MDE.

Answer 3

Esto podría ser un comentario para La excelente respuesta de MikeT si no fuera demasiado largo y complejo. He jugado mucho con él y he hecho un plugin de QGIS llamado Filtros de convolución FFT (en fase "experimental" todavía) basado en su función. Además de suavizar, el plugin también puede afinar los bordes restando la trama suavizada de la original.

He mejorado un poco la función de Mike en el proceso:

def __gaussian_blur1d(self, in_array, size):
        #check validity
        try:
            if 0 in in_array.shape:
                raise Exception("Null array can't be processed!")
        except TypeError:
            raise Exception("Null array can't be processed!")
        # expand in_array to fit edge of kernel
        padded_array = np.pad(in_array, size, 'symmetric').astype(float)
        # build kernel
        x, y = np.mgrid[-size:size + 1, -size:size + 1]
        g = np.exp(-(x**2 / float(size) + y**2 / float(size)))
        g = (g / g.sum()).astype(float)
        # do the Gaussian blur
        out_array = fftconvolve(padded_array, g, mode='valid')
        return out_array.astype(in_array.dtype)

Las comprobaciones de validez son bastante evidentes, pero lo importante es lanzar a flote y volver. Antes de esto, la función hacía que las matrices de enteros fueran negras (sólo ceros), debido a la división por la suma de los valores ( g / g.sum() ).

Answer 4

En QGIS, obtuve buenos resultados fácilmente utilizando Caja de herramientas de Orfeo Filtrado de imágenes. Es razonablemente rápido y el modo por lotes funciona bien. Están disponibles las difusiones gaussiana, media o anisotrópica.

Tenga en cuenta que Radius se refiere al número de celdas, no a la distancia.

Este es un ejemplo que utiliza Suavizado(gaussiano) :

• Crudo:

  

• Filtrado:

  

Answer 5

Buena solución para el desenfoque gaussiano y una animación genial. En cuanto a la herramienta de filtro de Esri mencionada anteriormente, es básicamente la herramienta de "Estadísticas Focales" de Esri codificada en un tamaño de 3x3. La herramienta "Focal Statistics" te da muchas más opciones sobre la forma de tu filtro móvil, el tamaño y la estadística que quieres ejecutar. http://desktop.arcgis.com/en/arcmap/latest/tools/spatial-analyst-toolbox/focal-statistics.htm

También puede hacer un filtro "irregular" en el que usted pasa su propio archivo de texto con los pesos a utilizar para cada celda. El archivo de texto tiene tantas filas como quieras en tu área de filtrado, con valores delimitados por espacios en blanco para las columnas. Supongo que siempre debes usar un número impar de filas y columnas, para que tu celda objetivo esté en el medio.

He creado una hoja de cálculo Excel para jugar con diferentes pesos que simplemente copio/pego en este archivo. Debería conseguir los mismos resultados que arriba si ajustas las fórmulas.