Posteado por: rengifo | diciembre 11, 2015

Comparando los efectos de dos diferentes interpolaciones sobre la calidad del modelo de elevación digital (MED).


Estimado lectores,

Los he tenido algo abandonados. Como siempre por razones de tiempo. En esta oportunidad  me gustaria discutir  sobre los efectos que  el uso de una técnica equivocada de interpolación tiene sobre  la calidad del MED. Para ello uso un MED de 20 metros de resolución en una zona en Noruega aleatoriamente selecionada (ver imágen 1 ).

2015-12-09 16_46_26-Microsoft PowerPoint - [DTManalyseESRIKon2 .pptm].png

Imágen 1: Zona de estudio

Para ello uso dos MED de la zona en cuestión: Un MED reproyectado de la zona  utm 32 a  la zona utm 33 usando el método de las celdas más cercanas ( nearest neighbor resampling ; NNR) y un MED  reproyectado  de igual manera , pero  usando el método bilinear (bilinear resampling; BIR). A simple vise ambos MEDs  parecen ser iguales. Ambos parecen presentar la misma distribución de altura (ver imagen 2).

DTM20(NNRBIR)

Imágen 2: MED  con 20 metros de resolución. A la izquierda MED  reproyectado usando  el  método de neareast neighbor (NNR). A la derecha MED reproyectado a la zona utm33 usando el método bilinear (BIR)

Sin  embargo el cálculo de la pendiente y de la orientación de la pendiente,  comienzan a revelar que los dos MED difieren uno del otro en la forma como las celdas fueron reorganizadas al ser  reproyectadas de una zona utm 32 a la zona utm 33 (ver imágen 3).

helding_samlet

Imágen 3: Cálculo de la pendiente de ambos MED en grados

Igualmente el cáculo de de la orientacíon de la pendiente  muestra más claramente lo que se conoce como artefactos (errores sistemáticos introducidos por un  algorítmo)  (ver imágen 4).

heldingorientering_samlet

Imágen 4: Cálculo de la  orientación de la pendiente en grados

Asi que surge la primera pregunta: Que genera esos artefactos en el MED  reproyectado con  el método NNR?  Para responder esta pregunta, importé ambos MEDs  a la base de datos  PostgreSQL /PostGIS usando el  siguiente comando:

raster2pgsql -s 32633 -d -I -C -M -R -l 4 c:\MEDNNR.tif -F -t 1×1 medNNR |psql -U postgres -d postgres

raster2pgsql -s 32633 -d -I -C -M -R -l 4 c:\MEDBIR.tif -F -t 1×1 medBIR |psql -U postgres -d postgres

Posteriormente, después de haber  importado ambos MEDs , usé siguiente SQL comando en QGIS:

SELECT rid,rast::geometry
FROM medNNR

y

SELECT rid,rast::geometry
FROM medNNR

para generar un indice espacial  de cada celdas de ambos MED para asi poder  extraer el punto central (centroid) de cada celda mediante  el uso de QGIS.

Comparando los centroids de cada celda en ambos MED, pude determinar que los artefactos en el MED reproyectado usando el método NNR, son generado por el desplazamiento de  las celdas hacia el Este. Dicho desplazamiento es generado por el  algorítmo implementado en el método NNR para efecto de reproyección y es de 3, 385 metros en este caso (imágen 5).

 

2015-12-10 14_35_16-Microsoft PowerPoint - [DTManalyseESRIKon2 .pptm]

Imágen 5: muestra el desplazamiento relativo del MED reproyectado con NNR  (rojo) con referencia al MED reproyectado con BIR (azul).

 Los artefactos conforman bloques artificiales de celdas de 10*10 metros y se diferencian uno del otro por cambios abruptos de la altura, pendiente  y la orientación de la pendiente, que  no reflejan la variación real del terreno (ver imágen  6).

2015-12-10 14_46_14-PowerPoint-lysbildefremvisning - [DTManalyseESRIKon2 .pptm].png

Imágen 6: Diferencia relativa de la altura (metros) a la izquierda  y   diferencia relativa de la pediente (grados ) a la derecha.

A de notarse que la diferencia de altura relativa muestra  que los bloques artificiales de celdas (10*1o metros)  presentan una variación de la altura sistemática  con valores negativos  en la esquina superior izquierda de cada bloque y valores positivos en la esquina inferior de cada bloque.

Esto es cierto solo para  pendiente orientadas al suroeste. Para  las pendiente orientadas al sureste, esta relación se invierte. En las zonas planas la diferencia es cero.

Al observar la diferencia relativa de la pendiente se hace aún más claro el efecto del  NNR sobre el MED. Los artefactos verticales presenta valores negativos de  pendiente, mientras que los artefactos horizontales presenta valores positivos. Esta relación al igual que en la diferencia relativa de altura se  invierte en función a la orientación de la pendiente y llega a ser cero en las zonas planas.

Mediante el uso de SAGA GIS hice un zoom de  ambos MEDs y de la  diferencia relativa de la pendiente .   Alli se puede  ver numéricamente el efecto que se  describo anteriormente (ver imágen 7).

2015-12-10 15_17_20-PowerPoint-lysbildefremvisning - [DTManalyseESRIKon2 .pptm].png

Imágen 7: Zoom en la zona los bloques  se diferencia por cambios abruptos de pendiente. Las dos imágenes  superiores muestras ambos MED . La imágen  inferior muestra una zona entre dos bloques de celdas.

Si el NNR  tiene efecto sobre la altura, la pediente y la orientación de la pediente es de esperarse que tenga  efectos sobre otros derivativos como por ejemplo  los derivativos usado en cálculos hidrológicos, especialmente cáculo de dirección del agua (aDir) y la acumulación del agua (aAcc). Asi surge la  segunda pregunta: Cual es el efecto numérico de  NNR sobre aAcc?

Primero calculé  el aDir  para derivar  aAcc usando  ArcGIS (por razones de trabajo). aAcc fue calculada usando  el algorítmo de fluido sencillo (single flow algorithm, SFA). Para ello, no se estableció ninguna restricción  para determinar la dirección del fluido del agua (ver imágen 8).

flowAcc_samlet

Imágen 8: Comparación de aAcc  calculada con SFA en ArcGIS y restrigida a valores entre  1000 y 7000 celdas para efectos de visualización.

Hay que notar que a simple vista ambos aAcc parecen diferir poco uno del otro. Uno estaría tentado a asegurar que el  aAcc calculado a partir del MED (NNR) es el correcto, ya que parece captar mejor la naturaleza del los rios.

Pero un vistazo cercano ala estadística de cada una de las direcciones  de drenajes asignadas  a los cada uno de los  aAccs por parte de ArcGIS,  revelan otra realidad (ver imágen 9).

2015-12-10 15_54_13-Microsoft PowerPoint - [DTManalyseESRIKon2 .pptm]

Imágen 9: La parte superior muestra una representación esquemática de como ArcGIS asigna los valores de la accumulación de agua y como posteriormente los codifica en direciones cardinales, comenzando con 1 en el este y terminando con 128 en noreste. La tabla en la parte inferior, muestra la substracción de los valores de la dirección del agua entre ambos MED.

Notese que en la dirección noroeste (código 8)  hay menos celdas que drenan  en esa dirección. Mientras que hacia el sur (código 4)  y  oeste (código 16) la contribución es positiva. Para evaluar el efecto de los artefactos  hice uso análisis direccional de superficies y estadística circular. Para ello usé GIS GRASS para generar gráficos polares (ver imágen 10).

LogFlowaccNNR

imágen 10: La imágen superior muestra el logaritmo de la accumulación del agua (aAcc), basada  en el MED (NNR). La imágen inferior muestra un  gráfico circular de la orientación de la pendiente de celdas que drenan hacia el noroeste (código 8), la misma fue generada con el comando d.polar en GRASS GIS .

Como se puede observar en la imágen 10 (superior) los artefactos generan lineas transversal  en la pendiente que dan la impresón que los rios se desplazan en forma de zigzag a los largo del pendiente. Si recuerdan el perfíl de  altura en la imágen 1, se darán cuenta que el terreno es de pendientes moderada. El algorítmo de fluido sencillo (single flow algorithm, SFA) que se utiliza para calcular aAcc, debería mostrar en pendiente moderadas un fluido del agua laminar y no en forma de zigzag.

En la imágen 10 (inferior) se puede observar que el promedio de la orientación del las celdas de drenaje es en dirección suoreste (azul). pero tambien se puede ver (rojo) que el efecto zigzag genera en una propagación de la acumulación de agua hacia el sur y el oeste que no es real. Esto se hace evidente si comparamos la imágen 10 con el aAcc calculado a partir del MED (BIR) (ver imágen 11)

LogFlowaccBIR

imágen 11: La imágen superior muestra el logaritmo de la accumulación del agua (aAcc), basada en el MED (BIR). La imágen inferior muestra un gráfico circular de la orientación de la pendiente de celdas que drenan hacia el noroeste (código 8), la misma fue  generada con el comando d.polar de GRASS GIS.

Observando la imágen 11 se puede observar que la dirección promedio coincide (azul) con la imágen 10, pero  no presenta mucha propagación de la accumulación del agua hacia el oeste o al sur (rojo).  La imágen 11 muestra un fluido laminar en la pendiente, que es la manera como el SFA debería comportarse en en este tipo de  terreno de pendientes moderadas.

De esta manera se puede concluir que:

  • El método NNR (stándar en ArcGIS)  para la reproyección  de   modélos  de elevación  no es apto para reproyectar datos continuos, ay que introduce errrores de altura mediante el desplazamiento del as celdas hacia el este.
  • El método NNR  reorganiza  la extructura interna del MED y bloque artificiales de celdas de 10*10 metros, llamados aretefactos. Es bloques artificiales se diferencian uno del otro por marcadas variaciones de altura y pendiente y da origen errores en la calculación de la dirección y acumulación de agua a partir de un MED.
  • La cantidad de celdas, hacia las cuales se drenan varia  dependiendo de la altura,  la pendiente y la orientación de la pendiente.
  • Desplazamiento lateral de la accumulación del agua en terreno moderado indica que los artefactos infuencia la forma  como el algorítmo de fluido sencillo (single flow algorithm, SFA).
  • Con excepción de los casos donde ya existe una red hidrológica calculada con anterioridad, que se pueda pueda usar para forzar el SFA a seguir una dirección específica, no se  recomienda usar MED reproyectado con el método NNR para ningún  análisis de terreno. En su lugar, deberia usarse el MED reproyectado con el método bilinear (BIR).

 

Espero haberle despertado la curiosidad. Abajo lagunas referencias de la literatura que use para  escribir el papel original de donde saqué las imágenes e ideas.

Les deseo una feliz navidad y un próspero año nuevo !!

 

Referencias

 

ESRI ArcGIS Resources: how flow accumulation works, http://resources.arcgis.com/EN/HELP/MAIN/10.1/index.html#//009z00000062000000 last accessed 03.11.2015.

Conrad, O., Bechtel, B., Bock, M., Dietrich, H., Fischer, E., Gerlitz, L., Wehberg, J., Wichmann, V., and Böhner, J. (2015). System for Automated Geoscientific Analyses (SAGA) v. 2.1.4, Geosci. Model Dev., 8, 1991-2007, doi:10.5194/gmd-8-1991-2015.

Neteler, M. and Mitasova,H., 2008, Open Source GIS: a GRASS GIS Approach. Springer, New York.

Parker, J.A., Kenyon R.V. and Troxel, D.E. Comparison of interpolating methods for image resampling, IEEE transactions of medical imaging, vol, MI-2, No.1 ,1983, pp.31-39.

R Development Core Team (2008). R: A language and environment for statistical computing. R Foundation for Statistical Computing, Vienna, Austria. ISBN 3-900051-07-0, URL http://www.R-project.org.

Anuncios

Responder

Introduce tus datos o haz clic en un icono para iniciar sesión:

Logo de WordPress.com

Estás comentando usando tu cuenta de WordPress.com. Cerrar sesión / Cambiar )

Imagen de Twitter

Estás comentando usando tu cuenta de Twitter. Cerrar sesión / Cambiar )

Foto de Facebook

Estás comentando usando tu cuenta de Facebook. Cerrar sesión / Cambiar )

Google+ photo

Estás comentando usando tu cuenta de Google+. Cerrar sesión / Cambiar )

Conectando a %s

Categorías

A %d blogueros les gusta esto: