Presentación al 5º Congreso Forestal Español

Algunos de los resultados obtenidos hasta el momento en relación al uso de imágenes MODIS (250 x 250 m) para detectar el decaimiento en repoblaciones de pino silvestre y pino salgareño en la Sierra de los Filabres (Almería) van a ser presentados al 5º Congreso Forestal Español, que tendrá lugar en Ávila, del 21 al 25 de septiembre de 2009. A continuación se muestra el resumen y el texto completo en pdf.

Detección del decaimiento en masas de Pinus sylvestris L. y P. nigra Arnold. procedentes de repoblación en la Sierra de los Filabres (Almería) mediante el uso del sensor MODIS

Durante las últimas décadas se han detectado en Europa y Norte América numerosos casos de decaimiento en masas forestales asociados, de forma general, con la contaminación, la presencia de plagas forestales y el cambio climático. En España, se han descrito procesos de decaimiento en abetales, encinares y pinares. Se prevé que el cambio climático resulte en una intensificación de estos procesos, por lo que la detección temprana de este fenómeno es un paso esencial para la gestión forestal sostenible. El diseño de modelos con base biológica para llevar a cabo dicha detección en Andalucía es uno de los objetivos prioritarios del proyecto GESBOME.

En este estudio se muestran resultados preliminares sobre la aplicación de imágenes MODIS (resolución espacial ~ 250 m) como herramienta para el seguimiento del estado vegetativo de masas forestales en la Sierra de los Filabres (Almería). En esta localidad, el problema del decaimiento de pinares de Pinus nigra Arnold. y Pinus sylvestris L. se viene observando desde el año 2001.

Para ello, se seleccionaron 76 puntos de control en campo y se asignaron a tres niveles de daño: sin afectar, moderado e intenso. En cada parcela se extrajeron los valores del índice de vegetación NDVI (Normalized Difference Vegetation Index) de las imágenes MODIS cada 16 días desde comienzos del año 2000 hasta agosto de 2008. Las series temporales de los tres niveles de daños se compararon entre sí utilizando un índice de referencia basado en la mediana del total de las series temporales para todas las parcelas. Las señales espectrales permitieron discriminar correctamente el nivel de daño extremos de los daños intermedios y moderado, especialmente en pinares de Pinus sylvestris, por lo que pueden utilizarse como indicadores para la detección del decaimiento. El uso de series temporales basadas en índices espectrales puede contribuir además a explorar procesos causales del decaimiento forestal mediante la comparación con variables edafoclimáticas y/o fisiográficas.

Análisis de los factores de incitación del decaimiento en masas forestales

En la estructura inicial planteada para analizar el proceso de decaimiento en la Sierra de los Filabres, Almería, se propuso la exploración de los factores de incitación mediante el análisis correlacional de distintas variables (topográficas, estructurales, espaciales) con el porcentaje de decaimiento estimado de visu en campo.

Uno de los problemas técnicos que se planteó previo a la realización de dichos análisis fue el cálculo de las superficies de insolación, problema que se ha resuelto mediante la aplicación de modelos específicos de insolacióimplementados en GRASS. En una entrada reciente en este blog se proveen los detalles técnicos para la preparación de la bases de datos points. Esta base contiene la siguiente información:

1. Porcentaje de decaimiento observado in situ en campo.
2. Especie (Pinus nigra y P. sylvestris).
3. Área basimétrica del rodal (var. estructural).
4. Densidad de pies del rodal (var. estructural).
5. Coordenadas geográficas x e y (var. espaciales).
6. Elevación (var. topográfica).
7. Pendiente (var. topográfica).
8. Orientación (var. topográfica).
9. Insolación en el solsticio de invierno (var. topográfica).
10. Insolación en el solsticio de verano (var. topográfica).
11. Índice de humedad del suelo (var. topográfica).
    

La información referente a variables topográficas (6-10) ha sido recalculada a tres resoluciones distintas (10x10 m, 30x30 m, 50x50 m). El modelo digital de elevaciones, del que se deriva el resto de información, está a una resolución original de 10x10 m. Sin embargo, dado el grano de las imágenes MODIS utilizadas (250x250 m) y el posible error de localización de los puntos puede ser conveniente re-escalar la información a escalas más groseras como 30x30 m o 50x50 m.

¿Está esta información proyectada a distintas resoluciones muy correlacionada entre sí? Si es así (p.e. r = 0.95), no importará mucho a qué resolución midamos las variables. De lo contrario, habrá que investigar si el cambio de resolución tiene un efecto sobre los resultados de los análisis de correlación con la variable respuesta (i.e. decaimiento).

Las tres bases de datos (points10, points30, points50) con la información referente a variables topográficas calculada a tres resoluciones distintas (10x10 m, 30x30 m y 50x50 m respectivamente) se puede descargar en R con el siguiente código (sólo la base de datos de 30x30 contiene además la información referente a la variable índice de humedad del suelo).

load(url("http://archivos-para-subir.googlegroups.com/web/points_correlaciones.Rob?hl=es&gda=mttqMzwAAACkOfwqHAVd4YqgfIB09GDRafzhFXYXrkFb8NzylP1tVw7qabgqw0xDKbwB-h3MnSf9Wm-ajmzVoAFUlE7c_fAt"))
ls()

[1] "points10" "points30" "points50"

Representamos ahora la matriz de correlaciones.


Cómo podemos observar, el cambio de resolución en la elevación no tiene ningún efecto perceptible sobre los valores de esta variable. Sin embargo, en la pendiente y orientación sí que se observa un cambio importante en los valores de estas variables según se cambia la resolución, sobretodo cuando se cambia la resolución de 10x10 m o 30x30 m a 50x50 m.

Vamos a seleccionar una resolución de 30x30 m para las variables topográficas, ya que es más representativa de los procesos a escala de rodal que la de 10x10 m y, al mismo tiempo, no se pierde mucha información con respecto a una escala de más detalle, cosa que sí pasa con la resolución de 50x50 m. Representamos ahora gráficamente los datos.
Parece que las variables que influyen más como factores de incitación en el proceso de decaimiento son la elevación, sobretodo para el caso de Pinus sylvestris, y la insolación invernal. Relaciones con las variables topográficas reescaladas a otra resolución (10x10 m, 50x50 m) no cambian apenas los resultados obtenidos. Es necesario, no obstante, analizar en mayor profundidad estos datos.

Como utilizar el potencial de GRASS desde R

R tiene varios paquetes que pueden utilizarse para el análisis y visualización de datos espaciales (sp, maptools, rgdal, etc.). Sin embargo, R no es un SIG y, por tanto, tiene ciertas limitaciones a la hora de trabajar con información espacial. Por ejemplo, para algo tan simple como extraer la pendiente y la orientación a partir de un modelo digital de elevaciones, R necesita del paquete RSAGA, que a su vez necesita del software gratuito SAGA diseñado por Alexander Brenning, que funciona únicamente bajo Windows. Esto supone ciertas limitaciones al uso de R como SIG. Estas limitaciones no existen cuando se utiliza GRASS en combinación con R. GRASS ofrece un complemento perfecto a R para la implementación de operaciones típicamente realizadas en SIG y R permite potenciar las capacidades de GRASS para la implementación de análisis geoespaciales.

Para poder hacer uso de ambos programas en necesario tener instalados GRASS y una versión de R posterior a 2.1.0. Desde el 'shell' se puede abrir una sesión de R. Lo primero será instalar los paquetes necesarios para poder utilizar R en combinación con GRASS (sp, rgdal, maptools, spgrass6, spGDAL, spmaptools).

GRASS 6.3.0 (Filabres):~ > R

install.packages(c("sp", "rgdal", "maptools"), dependencies= TRUE)

rS <- "http://r-spatial.sourceforge.net/R" install.packages(c("spgrass6", "spGDAL", "spmaptools"), repos=rS, dependencies=TRUE)

Sin salir de la sesión de R, debemos de cargar el paquete spgrass6 y la 'location' de GRASS donde se encuentran todas las capas de información raster o vectorial.

library(spgrass6)
G <- gmeta6()
str(G)

Si ahora queremos acceder a alguna de las capas raster del 'location' en el que estamos trabajando dentro de GRASS tendremos que utilizar la función de R readCELL6sp() ò readRAST6(). Esto genera un objeto del tipo SpatialGridDataFrame.

mde <- readCELL6sp("mde")

Sin embargo, es importante tener en cuenta que si cargamos muchas capas simultáneamente, especialmente si son capas con un gran número de celdas, se puede agotar la memoria virtual. En este caso R dará el mensaje 'vector memory exhausted'. Cuando esto me ha ocurrido a mi, R se ha quedado colgado y he tenido que reiniciar la sesión. Así que lo que aconsejo es ir eliminando los objetos que ya no necesitemos una vez que extraigamos de ellos la información necesaria. Para conocer el uso que se va haciendo de la memoria se puede utilizar la función gc().

gc()

used (Mb) gc trigger (Mb) max used (Mb)
Ncells 279590 7.5 531268 14.2 361111 9.7
Vcells 3160949 24.2 15608917 119.1 30379099 231.8

Es decir que se están utilizando alrededor de 31 MB de memoria virtual (hay que mirar la primera columna). Ahora eliminamos la capa mde.

rm(mde)
gc()

used (Mb) gc trigger (Mb) max used (Mb)
Ncells 279513 7.5 531268 14.2 361111 9.7
Vcells 137524 1.1 12487133 95.3 30379099 231.8

Y comprobamos que disminuye considerablemente el uso de la memoria virtual.

Antes de acceder a las capas raster, vamos a leer una capa de puntos en formato *.txt (parcelas muestreadas en campo) para la cuál extraeremos la información de elevación, pendiente, altitud e insolación de las capas raster.

points <- read.table(url("http://archivos-para-subir.googlegroups.com/web/Analisis+correlacionales.txt?gda=0ywgME4AAACkOfwqHAVd4YqgfIB09GDRN3No94DJejqv6LWrffOcUtQTec6xGbbz9B_IWZmp37itXp5Ud9d8afGj09bAbTSQ47Cl1bPl-23V2XOW7kn5sQ"), header = T, sep = "\t", dec = ",")
str(points)

'data.frame': 76 obs. of 7 variables:
$ x : num 543726 543486 542903 535211 539086 ...
$ y : num 4124164 4124224 4122417 4118071 4121908 ...
$ Id : Factor w/ 76 levels "1_12","1_13",..: 1 2 4 6 7 8 9 10 11 12 ...
$ Especie : Factor w/ 2 levels "Pinus nigra",..: 1 1 1 1 1 1 1 2 1 1 ...
$ Defoliacion: num 55 55 30 55 30 55 55 60 45 55 ...
$ Area_basim : num 25.8 28.1 24.2 22.0 18.0 ...
$ Densidad : num 286 286 1592 923 1974 ...

El siguiente paso es convertir esta capa de puntos en un objeto espacial del tipo SpatialPointsDataFrame y asignarle una proyección (ha de ser la misma que la de las capas raster y sino habría que reproyectar esta capa).

coordinates(points) <- data.frame(points$x, points$y)
proj4string(points) <- CRS("+proj=utm +zone=30 +ellps=intl +units=m +no_defs")

Y ahora sí, importamos cada una de las capas raster de GRASS a R, extraemos los valores de las capas raster para los puntos incluidos en la capa vectorial points y, por último, borramos las capas raster para no agotar la memoria virtual.

mde <- readCELL6sp("mde")
extract.mde <- overlay(mde, points)
rm(mde)

slope <- readCELL6sp("slope")
extract.slope <- overlay(slope, points)
rm(slope)

aspect <- readCELL6sp("aspect")
extract.aspect <- overlay(aspect, points)
rm(aspect)

rad.summer <- readCELL6sp("beam_radiation_172")
extract.rad.sum <- overlay(rad.summer, points)
rm(rad.summer)

rad.winter <- readCELL6sp("beam_radiation_354")

extract.rad.win <- overlay(rad.winter, points)
rm(rad.winter)

wetness <- readCELL6sp("wetness_index")
extract.wetness <- overlay(wetness, points)
rm(wetness)

Y ahora juntamos toda la información extraida en un único data.frame.

points@data <- cbind(points@data, mde = extract.mde@data, slope = extract.slope@data, aspect = extract.aspect@data, rad.summer = extract.rad.sum@data, rad.winter = extract.rad.win@data, wetness = extract.wetness)

El resultado final, es por tanto, un data.frame con información sobre valores de decaimiento (incluidos en el archivo points) y toda una serie de variables físicas obtenidas a partir de capas raster (pendiente, altitud, etc) que podemos analizar estadísticamente en R, cosa que en principio no podríamos hacer en GRASS. En la siguiente entrada se muestra un análisis preliminar de estos datos en R.

Cálculo de superficies de insolación

En GRASS es relativamente sencillo implementar un modelo para estimar la insolación a partir de los mapas de altitud, pendiente y orientación por medio del comando r.sun.

Se puede acceder a este comando a través del 'shell' o del gestor GIS. Si optamos por la segunda opción, hay que ir a Raster - Solar radiance and shadows - Solar irradiance and irradiation. Se abrirá una ventana con tres pestañas: (1) Options; (2) Input_options; (3) Output_options. Algunos de los campos que hay que rellenar son obligatorios pero muchos otros son optativos y, o bien no son completamente necesarios, o están definidos por defecto.


Antes de utilizar este comando es conveniente leerse la ayuda de la función. r.sun puede utilizarse para calcular tanto la radiación directa, difusa y reflejada [W.m-2] en un momento puntual del día (modo 1), como la radiación directa, difusa y reflejada acumulada [Wh.m-2.day-1] a lo largo de todo un día (modo 2). El segundo modo suele ser más informativo, sobretodo en estudios relacionados con la producción de cultivos o la presencia de formaciones vegetales. Esta radiación acumulada será la que calculemos en este ejemplo.

En la primera pestaña (Options) se pide información sobre:

• el día del año en el que se quiere calcular la radiación (día juliano).
• el intervalo del tiempo en el que se calcula la radiación cuando se estima la radiación acumulada durante todo el día. Esto es sólo aplicable al modo 2. 0.5 (valor por defecto) significaría que se estima la radiación cada media hora a lo largo del día y luego se calcula la integral para todas las estimas puntuales.
• el valor de declinación solar (radianes) se estima por defecto de manera automática utilizando el día del año. No obstante, si se tiene una medida más precisa de este valor, puede introducirse en la casilla correspondiente.
• la hora solar para la que se quiere calcular la radiación, sólo aplicable en el caso del modo 1.
  

Por tanto, en el caso de estar interesados en el cálculo de la radiación acumulada (modo 2), sólo se tiene que especificar el día del año y el intervalo de tiempo para el cuál se hacen las estimas parciales.

En la siguiente pestaña (Input_options) hay que introducir el nombre de los archivos que contienen las capas raster (en formato Grass) de elevación, (elevation) pendiente (slope) y orientación (aspect).

Por defecto, hay establecido un coeficiente de turbidez atmosférica (Linke atmospheric turbidity coefficient) de 3.0. Este coeficiente refleja la presencia de aerosoles y nubosidad que puede influir en la radiación difusa y reflejada. En la página de ayuda muestran una tabla en donde se sugieren coeficientes de Linke orientativos según el tipo de área de estudio y mes del año. Si se dispusiera de una capa raster con el valor de este coeficiente para toda el área de estudio, se podría introducir en la casilla que dice 'Name of the Linke atmospheric turbidity coefficient input raster map'. Como generalmente no se dispone de una información tan detallada, se deja este espacio en blanco.

En la tercera pestaña (Output_options) es necesario especificar el nombre de los archivos de salida para la radiación directa (beam_rad), difusa (diff_rad) y reflejada (refl_rad). Además, en el caso del modo 1 habrá que especificar el nombre de un archivo de salida con el ángulo de incidencia del sol (incidout). En el caso del modo 2, habrá que especificar el nombre de un archivo de salida con el tiempo de insolación (insol_time), es decir, el número de horas a lo largo del día que cada pixel ha estado expuesto a la radiación directa.


Una vez definidos todos estos parámetros, le damos a RUN. El proceso, especialmente en el caso del modo 2, puede ser muy lento ya que exije una gran cantidad de cálculos computacionales. Para el caso de la Sierra de los Filabres (Almería), el cálculo de la insolación diaria total (radiación integrada) para un área de 3839 x 7086 píxeles de 10 x 10 m duró aproximadamente 3-4 horas.

Por otro lado hay que tener en cuenta que para la definición de las características de insolación de un lugar deben estimarse los modelos de insolación representativos de un conjunto suficiente de periodos anuales. En efecto, al menos en latitudes medias o altas los contrastes estacionales son fuertes y contienen una información ambiental que no puede despreciarse. Por este motivo no es suficiente caracterizar la insolación para un valor medio anual, al menos en zonas no ecuatoriales. Por esta razón, he calculado la insolación en el solsticio de verano (día juliano 172) y en el solsticio de invierno (día juliano 354). Luego estas dos variables pueden considerarse individual o conjuntamente (mediante la media) en análisis posteriores.

Por último, es posible obtener una representación del mapa de insolación utilizando el siguiente código en la consola de comandos de GRASS.

g.region rast=mde
d.mon x1
r.colors beam_radiation_172 col=grey
d.his i_map=beam_radiation_172 h_map=mde

Análisis de modelos digitales de elevación en GRASS

En GRASS se pueden calcular los mapas de pendiente y orientación a partir de un Modelo Digital de Elevaciones (Digital Elevation Model (DEM) en inglés) desde el administrador GIS pulsando:

Raster - Terrain analysis - Slope and aspect

O bien desde la consola de comandos (shell) con el comando r.slope.aspect.

r.slope.aspect el=dem as=aspect sl=slope

Los dos mapas son procesados a la vez. Los ángulos del mapa de orientación son medidos en el sentido opuesto a las manecillas del reloj a partir del Este. Las pendientes son calculadas por defecto en grados.

Si se quieren visualizar estos dos archivos basta con usar los comandos d.mon y d.rast para abrir un monitor y visualizar (display) en él la capa raster respectivamente. Alternativamente el comando d.rast.leg agrega una leyenda sencilla al monitor.

d.mon x1
d.rast aspect
d.mon x2
d.rast.leg slope

Existen módulos adicionales para trabajar con los DEMs: las áreas con depresiones pueden ser llenadas con r.fill.dir, las líneas de flujo pueden ser calculadas con r.flow. El análisis de cuencas puede ser realizado con r.watershed y con r.terraflow en archivos muy grandes. Todos estos comandos pueden ser utilizados también desde el administrador GIS (en la pestaña raster).

¡Bienvenido a GRASS!

GRASS (por su siglas en inglés, the Geographic Resources Analysis and Support System, http://grass.itc.it) es un Sistema de Información Geografica (SIG) gratuito y de código abierto, integrado con un sistema de procesamiento de imágenes y subsistemas de visualización. GRASS permite trabajar con imágenes raster y vectorial. Tiene interfaces paraPostgreSQL, MySQL, DBF y bases de datos conectadas con ODBC. Además puede ser conectado a UMN/MapServer, R, gstat, Matlab, Octave, Povray y otros software.

Hasta hace poco siempre había trabajado con otros software comerciales como ArcGIS, Idrisi, PCI o ERDAS. GRASS siempre me había parecido muy complicado y poco amigable. Sin embargo, sus características como software libre y gratuito lo convierten en una herramienta muy versátil y potente, por lo que merece la pena hacer el esfuerzo de aprendizaje. En esta entrada explico como empezar a trabajar en GRASS y los conceptos más básicos para aquellos que, como yo hasta ayer, quieran dar el salto pero no acaban de animarse.

Instalación y manuales

GRASS está disponible desde el ITC-irst en Italia (http://grass.itc.it). GRASS está diseñado para correr en varios ambientes UNIX, tales como GNU/Linux, SUN-Solaris, Irix o MacOS X. También puede correr bajo MSWindows/NT/2000/XP con Cygwin (un simulador de Linux), aunque posiblemente las próximas versiones de GRASS se puedan correr bajo Windows directamente. Libros, tutoriales, manuales, cursos en línea y mas documentación se encuentra listada en el 'Proyecto para la Documentación de GRASS' (http://grass.itc.it/gdp). En español, hay al menos un par de tutoriales de introducción a GRASS (http://mpa.itc.it/markus/osg05/). Ambos son traducciones de un tutorial de Markus Neteler para un taller de trabajo y no difieren apenas en cuanto a contenidos entre sí. En mi opininión, estos tutoriales no son particularmente buenos, pero por lo menos sirven para entender algunos conceptos básicos del funcionamiento de GRASS, por lo que merece la pena echarles un vistazo (muchos de los contenidos de esta entrada han sido extraidos de estos tutoriales).

Estructura de los proyectos en GRASS

Cuando alguien está familiarizado con el uso de SIG y sensores remotos pero nunca antes había trabajado con GRASS, lo más importante es saber cómo se estructuran los proyectos y cómo abrir una sesión en GRASS. Esto es lo que normalmente más cuesta entender. La información de GRASS se encuentra almacenada en un directorio conocido como 'database' (también llamado 'GISDBASE') . Este directorio se debe de crear antes de empezar a trabajar con GRASS (yo lo he llamado 'home/luis/grassdata'). Dentro de este 'database', los proyectos están organizados por áreas de proyectos en subdirectorios llamados 'locations'. Un 'location' está definido por un sistema de coordenadas, una proyección cartográfica y unas fronteras geográficas. Los subdirectorios y archivos que definen un 'location' son creados automáticamente cuando GRASS es iniciado la primera vez con un nuevo 'location' (esto es lo que vamos a hacer un poco más abajo). Cada 'location' puede tener múltiples 'mapsets'. Un motivo para mantener múltiples 'mapsets' es almacenar mapas relacionados a algún aspecto específico del proyecto o subregiones. Otro motivo es permitir el trabajo simultáneo de usuarios a capas almacenadas dentro del el mismo 'location', por ejemplo, equipos trabajando en elmismo proyecto. Para trabajos en equipo, una base de datos de GRASS centralizada puede ser creada en un sistema de archivos de red (por ejemplo, NFS). Además de acceder el propio 'mapset', cada usuario puede leer capas de otros 'mapsets', aunque cada usuario solo puede modificar o borrar las capas en su propio 'mapset'.

Cuando se crea un nuevo 'location', GRASS automáticamente crea un 'mapset' especial llamado PERMANENT donde la información principal del proyecto puede ser almacenada. Datos en el 'mapset' PERMANENT pueden ser solamente agregados, modificados o borrados por el propietario del 'mapset' PERMANENT. Sin embargo, los datos pueden ser consultados, analizados y copiados por otros usuarios a sus respectivos 'mapsets'. El 'mapset' PERMANENT es útil para entregar información general (ejemplo, un modelo de elevación digital), accesible pero protegida contra escritura de todos los usuarios que se encuentran trabajando en el mismo 'location' del propietario de la base de datos. En principio si no se tiene pensado trabajar en equipos con varios usuarios y con una gran cantidad de información, se puede trabajar directamente en el 'mapset' PERMANENT.

Primeros pasos en GRASS: Cómo crear un nuevo 'location'

Una gran diferencia entre GRASS y otros SIG es que GRASS requiere los parámetros de la proyección antes de que el usuario pueda trabajar en un 'location'. La ventaja de esto es que las cosas quedan bien definidas y se evita un desorden con la combinación de proyecciones, cosa que en mi experiencia es bastante común cuando se trabaja con otros software que no siguen esta estructura de proyectos. Para crear un nuevo 'location' primero tenemos que abrir GRASS. En Linux (yo trabajo actualmente con Ubuntu Intrepid Ibex 8.10) hay que ir a la terminal (konsole) y escribir:

grass63

o simplemente:

grass

Esto abrira la ventana de bienvenida de GRASS.

Lo primero que hay que hacer es seleccionar el 'database' donde se va a guardar el nuevo 'location'. Este directorio se ha tenido que crear antes de empezar a trabajar en GRASS. Ahora selecciona el botón 'Projection values' para crear un nuevo 'location' (Define new location with...), el cual te llevará a una nueva ventana de texto.

En la pantalla ingresa un nuevo nombre para el 'location' (no puede contener espacios en blanco) y luego continua presionando la tecla ESC y la tecla RETURN de tu teclado. Te aparecerá otra ventana de texto diciendo que esa 'location' no existe, cuáles son las 'locations' disponibles y preguntándote si quieres crear una nueva 'location' con el nombre asignado. A esto último se dice que sí y se pulsa RETURN.

Esto te lleva a otra ventana de texto en dónde te explican la información que vas a necesitar para crear un nuevo 'location' (sistema de coordenadas, zona UTM y parámetros de la proyección, coordenadas del área de trabajo, etc). Te preguntará si tienes toda esa información. Se dice que sí y se pulsa RETURN.

Una nueva ventana te pedirá que definas el sistema de coordenadas para tu 'location'. En mi caso defino C (UTM) y pulso RETURN.

Lo siguiente es ingresar una pequeña descripción o título para el 'location'. En esta misma ventana te pregunta si quieres especificar un datum. Por defecto se responde que sí. Se pulsa RETURN para continuar.

Para especificar el datum se puede pulsar 'list' y te viene una lista de todos los datum disponibles. Para salir de este listado se pulsa 'q'. Se escribe el nombre del datum y se pulsa RETURN.

El último paso es definir las coordenadas del área de estudio. En mi caso por ejemplo, las hice coincidir con las del modelo digital de elevaciones con el que voy a trabajar. En otros casos pueden ser más o menos amplias y abarcar toda una región geográfica como Sierra Nevada, Andalucía o la Península Ibérica. Hay que tener en cuenta que esta es la región que define el 'location', pero que luego es posible redefinir el área de estudio en distintos 'mapset' dentro de un 'location'. Igualmente se pueden redefinir los parámetros de un 'location' una vez que ha sido creado.

Finalmente se pide que se defina la resolución de los píxeles. Una cosa buena de GRASS es que se puede cambiar el tamaño de los píxeles con relativa facilidad. Al contrario que con otros software, no se remuetrea una capa raster, sino que se reproyecta todo un proyecto ('location') a la resolución deseada. Esto agiliza mucho el procesamiento de la información sin necesidad de generar infinidad de nuevos archivos con información parcialmente redundante. Se pulsa RETURN para acabar.

Ya tenemos el proyecto creado. Ahora podemos comenzar GRASS desde la página de bienvenida (paso 1) seleccionando el 'database' donde has creado el nuevo 'location'. (/home/luis/grassdata), el 'location' que acabas de crear (Filabres) y el 'mapset' (PERMANENT). Pulsa 'Enter GRASS'.

¡Y ya podemos empezar a trabajar en GRASS!

Interfaz gráfica y consola de comandos ('shell')

En GRASS se puede trabajar de dos maneras (alternativa o simultáneamente): a través de la interfaz gráfica o a través de la consola de comandos ('shell'). Mientras que la primera es más intuitiva, la segunda da mayor versatilidad y potencia a la hora de programar funciones y procesos. Ambas opciones están disponibles todo el tiempo.

La interfaz gráfica consta de tres ventanas.

La ventana principal es el administrador GIS (GIS Manager, izquierda). Desde esta ventana se pueden realizar la mayoría de las funciones relacionadas con el procesamiento de capas raster y vectorial, incluyendo la importación y exportación de capas en otros formatos. También se pueden modificar las características de configuración del 'location'.

Una segunda ventana (Output - GIS.m, arriba derecha) nos muestra los resultados de las operaciones realizadas a través del administrador GIS así como el código implementado (este código se puede usar directamente en la consola de comandos).

Finalmente hay una tercera ventana (Map Display X, abajo derecha) que muestra los mapas seleccionados en el administrador GIS. Otra de las ventajas de GRASS es que en vez de desplegar los mapas en la pantalla, envía el contenido del monitor a un archivo *.PNG por lo que ocupa menos memoria y es mucho más rápido a la hora de visualizar datos, especialmente cuando se trata de capas raster con una gran cantidad de información.

La consola de comandos ('shell') funciona de manera similar a como funcionan la mayoría de las aplicaciones en sistemas Unix. A través del 'shell' podemos obtener toda la funcionalidad que se obtiene a través de la interfaz gráfica de GRASS, pero además podemos programar scripts y conectar con otros software como R.


Importar datos a GRASS

Ya hemos creado nuestra 'location'. Ahora vamos a importar una capa raster a nuestro 'location' para trabajar con ella desde GRASS. GRASS es compatible con la mayoría de los formatos típicos utilizados en SIG y teledetección. En este caso importaremos una capa raster ascii de ArcGIS (GRID). Para ello vamos a:

File - Import raster map . ESRI grid

Alternativamente desde la consola de comandos podríamos haber escrito el comando:

r.in.arc

Cualquiera de las dos formas nos llevarían a la siguiente ventana, en dónde debemos de definir el archivo de entrada y el nombre del archivo de salida que se guardará directamente en nuestro 'location' y 'mapset' definidos. Una vez definidos los argumentos necesarios para ejecutar la función, pulsamos RUN.


Como podemos observar en la parte inferior de la ventana aparece el código que ejecuta todo el proceso incluyendo la definición de los archivos de entrada y salida. Esta línea de texto puede ser ejecutada en la consola de comandos para realizar esta misma operación sin necesidad de tener que abrir un GUI (Graphic's User Interface).

Ahora ya podemos trabajar con esta capa raster en GRASS para realizar distintas operaciones. Para más información sobre cómo cargar capas en GRASS y visualizarlas se puede consultar el manual de Markus Neteler traducido por Guillermo Martínez.