Construyendo mapas de sombras hillshade con R

A raíz de la entrada 'Efecto hillshade o sombras de laderas' del blog de Dominic Royé, he probado una forma menos ambiciosa de obtener con R un sombreado hillshade de un mapa de elevaciones, pero con la premisa de no usar ninguna librería específica haciendo el método fácilmente trasladable a cualquier lenguaje con notación matricial como Python o Matlab. El código además está vectorizado para no usar bucles lo que hace su ejecución instantánea.

Pese a su fuerte apariencia volumétrica un hillshade no es una infografía 3D, sino un cálculo en 2D de iluminación por zonas basado en el ángulo de incidencia de los rayos solares sobre el terreno.

Con un truco vectorial que usaba en la universidad estimamos la dirección normal al terreno en cada posición del mapa. Consiste en realizar el producto vectorial de dos vectores definidos por las elevaciones de las dos parejas de celdas adyacentes en los ejes X e Y. Elegante y más sencillo que calcular el gradiente o hacer regresiones.

Contrastando la dirección obtenida con la dirección de incidencia definida para la iluminación mediante un producto escalar, tenemos de forma directa el nivel de luminosidad reflejada en cada celda.

Se muestran las ecuaciones de los cálculos descritos incluyendo la forma de vectorizar la obtención del vector normal en todo el mapa (hacer clic para ampliar):

Aplicamos el cálculo a los datos DEM de la Sierra de Guadarrama en Madrid descargados desde el Centro Nacional de Información Geográfica del Ministerio de Transportes, que representados en escala de grises tienen el siguiente aspecto (hacer clic para visualizar en alta resolución):

El resultado del hillshade para una dirección de iluminación que correspondería al atardecer resulta del siguiente modo (hacer clic para visualizar en alta resolución):

Finalmente con Photoshop fusionamos la luminosidad del hillshade con el color del mapa de elevaciones coloreado con las curvas de la paleta 'topo.colors' que calculamos en 'Mapas de color 'viridis' en Photoshop con R' (hacer clic para ver en todo su esplendor):

Aplicando el mismo procedimiento a la impresionante orografía de Tenerife tenemos (hacer clic para alta resolución):

Y haciendo lo propio con la sierra de Alicante (clic para alta resolución):

Acabamos con un hillshade de Marte obtenido con datos de la NASA, marcando el punto de mayor y menor altitud de su superficie (clic para alta resolución, NOTA: son 32Mpx):

Hacemos zoom en el Olympus Mons, la mayor montaña conocida del Sistema Solar. Se trata de un volcán extinguido con una extensión similar a Francia y una altitud respecto a su área circundante de más de 20km. Lo coloreamos con tonos marcianos y aplicamos un gradiente para simular el amanecer (procesado en Photoshop aquí, clic para mayor resolución):

~~~

Repositorio con el código R: GitHub.