T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Análisis Geo-Estadístico y de Exposición Solar
Xander Bakker
Septiembre de 2000 Santafé de Bogotá, D.C., Colombia
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
INTERPRETACION DIGITAL Y SENSORES REMOTOS
for
WINDOWS
Advertencia:
El presente documento fue elaborado como un curso de análisis geo-estadístico utilizando el sistema de información geográfica ILWIS versión 2.23 para Windows. Para el buen desarrollo de este curso cu rso los participantes participantes deben tener conceptos básicos bási cos de sistemas de información geográfica, haber recibido el curso de introducción a ILWIS o tener experiencia con el sistema de información geográfica ILWIS versión para Windows obtenido por medio de métodos autodidactas. HeRindser Ltda ni el autor de este documento documento asumen responsabilidad por resultados erróneos obtenido obteni do en el análisis análisi s geo-estadístico efectua efe ctuado do con base en este documento.
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
INTERPRETACION DIGITAL Y SENSORES REMOTOS
for
WINDOWS
Contenido : 1. INTRODUCCIÓN . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.1 Definición de Interpolación . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 Interpolación de puntos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.1 1.2.1 Interp Interpolac olació iónn de veci vecinos nos más más cercano cercanos s . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.2 1.2.2 Interpo Interpolac lación ión de prom promedi edioo móvi móvil l . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2.3 1.2.3 Interp Interpolac olació iónn de tendenci tendenciaa de superfic superficie ie .. . . . . . . . . . . . . . . . . . . . . . . . . 1.2.4 1.2.4 Interp Interpolac olació iónn de superfi superfici ciee móvil móvil .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.3 Historia del Kriging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.4 1.4 Aplicaci cacio ones nes del Krigin ging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
1 1 1 1 3 5 5 6 6
2. TEORÍA T EORÍA DE GEO-ESTADÍSTICA . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 2.1 Distribución es espacial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7 2.2 (Auto-)Correlación espacial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 11 2.2.1 2.2.1 Interp Interpret retaci ación ón de “Moran's “Moran's I” y “Geary's “Geary's c”: c”: . . . . . . . . . . . . . . . . . . . . . . 12 2.3 Semivariogramas . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13 2.4 Variograma representado como superficie . . . . . . . . . . . . . . . . . . . . . . . . . 15 2.5 Variograma de cruces . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16 2.6 Métodos de Kriging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 2.6. 2.6.11 Simp Simple le Kri Krigi ging ng .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 2.6. 2.6.22 Ordi Ordinar naryy Krig Krigin ing g . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 17 2.6. 2.6.33 Ani Anisotro sotropi picc Krigi Kriging ng .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 2.6. 2.6.44 CoKr CoKriiging ging .. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 18 2.6. 2.6.55 Univ Univer ersal sal Krig Krigin ing g . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 19 2.7 2.7 Diag Diagrrama ama de flu flujo del del Krig Krigin ing g . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20 3. EJERCICIO DE KRIGING . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1 Reconocer la información . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.1.1 Despliegue y consulta de la información . información . . . . . . . . . . . . . . . . . . . . . . . . 3.1.2 Análisis de distribución espacial . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2 3.2 Gener enerac aciión de var variogram grama as . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3.2.1 Variograma representado como superficie (Anisotropía) . . . . . . . . . . . 3.2.2 Correlación espacial bi-direccional . . . . . . . . . . . . . . . . . . . . . . . . . . . .
21 21 21 22 23 23 24
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
INTERPRETACION DIGITAL Y SENSORES REMOTOS
3.3
for
WINDOWS
Modelar el semivariograma . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 25
Adicionar modelos . modelos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Anisotropía zonal . zonal . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Correlación espacial omnidireccional . . . . . . . . . . . . . . . . . . . . . . . . . . Adicionar modelos . modelos . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Verificación Verificación y selección de los modelos (Goodness of Fitt Fitt R 2 ) . . . . . . . Interpolación por medio del Kriging . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Verificación del resultado (Confidence interval) . . . . . . . . . . . . . . . . . . . .
25 26 27 27 28 33 34
4. EXPOSICIÓN SOLAR . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1 La trayectoria solar . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1.1 Declinación solar . solar . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1.2 Latitud geográfica . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.1.3 Ángulo horario . horario . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2 Ocultamiento topográfico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2.1 Definición Definición de exposición solar . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2.2 Ángulo de incidencia . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.2.3 Irradiación . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3 Ejercicio práctico . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.1 Operaciones de conectividad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.2 Intervisibilidad . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.3 Conectividad en ILWIS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4.3.4 Iniciar el proceso . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
36 36 36 37 37 38 38 39 39 40 40 41 41 43
3.4 3.5
3.3.1 3.3.2 3.3.3 3.3.4 3.3.5
5. RESPUESTAS . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 50 Bibliografía . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53 Anexo: Declinación solar . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 54
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
1. INTRODUCCIÓN La geo-estadística se aplica a un conjunto de puntos con sus atributos con el fin de comprender su distribución y comportamiento espacial y determinar el tipo de análisi s que se puede aplicar. Tema central de este curso es el proceso de Kriging , entendido como una interpolación interpola ción de atributos puntuales pu ntuales por medio de una función funci ón estadística. estadística. Dentro Dentro de este curso se explicarán expli carán los fundamentos de la geo-estadística geo-estadística y se pondrán en práctica dichos conocimientos. Además se ha incluido una sección teórico-práctica que trata del análisis de exposición solar. 1.1 Definición de Interpolación
Es el procedimiento de estimar las propiedades en sitios no muestreados dentro de un área cubierta por valores conocidos en localidades vecinas. El estimar los valores de propiedades en lugares fuera del área cubierta por observaciones se ll ama extrapolació extrapolación. n. La calidad de la interpolación interpolaci ón depende de la confiabilidad, confiabil idad, certeza, número y distribución de los puntos conoci dos usados en el cálculo y en la precisión del modelam modelamiento iento de la función utilizada. Los valores desconocidos son calculados con esta función. El escoger el modelo apropiado, es esencial para obtener buenos resultados. (Valenzuela (Valenzuela 1989, Aronoff 1989 ) 1.2 Interpolación de puntos
La interpolación de puntos sirve para generar un superficie continua de valores a partir de valores registrados en ubicaciones puntuales. Existen en la literatura muchas diferentes formas para efectuar este tipo de interpolación y cada una de estos tipos tip os tiene su aplicación específica. ILWIS soporta varios de estas interpolaciones: vecino más cercano, promedio móvil, tendencia de superficie y superficie móvil. (ILWIS ( ILWIS 1999 ) 1.2.1 Interpolación de vecinos más cercanos La interpolación de vecinos más cercanos crea los denominados polígonos de Thiessen (o "polígonos Voronoi") los cuales definen áreas de influencia individual por cada punto de un conjunto de puntos. Es una manera de extender la información puntual asumiendo que la mejor información, para ubicaciones sin observaciones, es el valor del punto más cercano. (Aronoff (Aronoff 1989). 1989). A
A B
C
B
D
C
D E
F
E
G
F
G H
H
Figura 1: Polígonos de Thiessen
1
Introducción
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
La funcionalidad es la siguiente: para cada pixel las distancias euclidianas son calculadas hacia todos los puntos. A cada pixel se le asigna el valor del punto con la menor distancia . Este tipo de interpolación no se aplica para generar mapas de superficie continua, pero más bien para determinar áreas de influencia de, por ejemplo; escuelas con base en distancias o tiempos de viaje. (Ver gráfico abajo) (Bakker y Perez 2000 ) En la figura 2 se muestra un mapa de polígonos de Thiessen calculado sobre una superficie compleja mediante una operación de esparcimiento. (Datos: (Datos: CVC, Cuenca Tulua-Morales, Proyecto SIGPAFC, asesorado por X. Bakker, HeRindser LTDA) LTDA)
condiciones de esparcimiento
Figura 2: Polígonos de Thiessen sobre superficie compleja
Una función de esparcimiento evalúa esparcimiento evalúa fenómenos que se acumulan con la distancia. Se utiliza para evaluar tiempo o costo de transporte sobre una superficie compleja. Los movimientos pueden ser limitados por barreras parciales que reducen la velocidad de movimiento o por barreras absolutas que impiden completamente el movimiento. (Valenzuela 1989, Aronoff 1989 ). ). ILWIS permite por medio del programa “Distance” generar estos tipos de mapa, especificando especific ando un mapa de pesos (valor dificultad, dificulta d, tiempo o costo para atravesar un pixel). Dentro de esta operación se puede opcionalmente crear un mapa de polígonos de Thiessen.
2
Introducción
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
1.2.2 Interpolación de promedio móvil Una de los tipos de interpolación más aplicada para determinar un valor en ubicaciones “no-visitadas” es la del promedio móvil. Funciona por medio de una ventana circular donde el valor resultante resultante es determinado determinado por la influencia influ encia de puntos conocidos conoc idos dentro de esta ventana. La influencia i nfluencia es ponderada con base en la distancia inversa o decrecimiento lineal (En las figuras 3 y 4: i es la influencia y d es una función de distancia). (Burrough (Burrough 1985, ILWIS 1999 ). ). En ILWIS se pueden usar dos funciones para ponderar la influencia como función de distancia. La primera es la “Distancia Inversa” (ver figura 3) y la segunda es el “Decrecimiento Lineal” (ver figura 4). La figura 3 muestra gráficamente la influencia de cada punto como función de distancia (Distancia (Distancia Inversa ). ). Los puntos más cercanos tendrán la mayor influencia en el valor resultante. Al incrementar la variable n los puntos a mayor distancia casi no tienen influencia en el resultado.
Figura 3: Influencia como función de distancia (distancia inversa)
La figura 4 muestra gráficamente la influencia influenc ia de cada punto como función de distancia (Decrecimiento Lineal ). ). Los puntos más cercanos tendrán la mayor influencia en el valor resultant. Al incrementar la variable n los puntos a mayor distancia ejercen influencia en el resultado. Con el valor 1 para la variable n el decrecimiento será una linea recta. Figura 4: Influencia como función de distancia (decrecimiento lineal)
3
Introducción
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
En la figura 5 se muestra el resultado de una interpolación de puntos, usando el método de promedio móvil con una función de Distancia Inversa con n=2.5. Información de precipitación en el área de estudio de la Bioreserva Amboseli, Kenia.
Figura 5: Interpolación por medio del promedio móvil
Los puntos incluidos para la interpolación también están ubicados fuera del área de estudio, para evitar que los errores de interpolación ocurran dentro del área de interés. Posteriormente se recortó el área útil. 4
Introducción
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
1.2.3 Interpolación de tendencia de superficie Una forma relativamente sencilla para describir variaciones graduales sobre grandes distancias es por medio del modelamiento con una regresión “polynomial” . La idea es ajustar una superficie polynomial, usando “least “least squares fit ” mediante los puntos conocidos. Las funciones que se pueden usar varían desde la más sencilla “plano “plano ” para o áreas de poca variación hasta la más compleja “6 “6 orden superficie polynomial ” para superficies con un alto grado de variación. La superficie que resulta trata de ajustarse a la máxima cantidad de puntos. (Burrough (Burrough 1985, ILWIS 1999 ). ). En ILWIS la interpolación interpolac ión de tendencia de super superficie ficie únicamente puede ser desarrollada con mapa de puntos con dominio de tipo valor o con base en un atributo numérico.
Figura 7: Regresión polynomial
Figura 6: Regresión lineal
Las figuras muestran dos funciones para ajustar la superficie a los puntos conocidos. La figura 6 muestra una función lineal, y la figura 7 muestra una función polynomial. Horizontalmente está graficada grafica da la distancia distanc ia y verticalmente la la varianza. varianza. En el análisis análi sis de tendencias la linea de regresión se ajusta de tal forma que la suma de los cuadrados de diferencia sea la mínima ( b) y no la distancia mínima entre el punto y la linea ( a) ni tampoco la desviación estándar en distancia (c). (Burrough (Burrough 1985 ) 1.2.4 Interpolación de superficie móvil La interpolación interpolaci ón de superficie móvil es una combinación de la interpolación de tendencias de superficie y el promedio móvil. Usa las mismas funciones que la anterior interpolación, interpolación, pero adicionalmente adicion almente utiliza un método de distancia como c omo el empleado empleado en la interpolación de promedio móvil. Se pueden usar las dos funciones de distancias (Distancia Inversa y Decrecimiento Lineal). Tanto en la interpolación de superficie móvil como en la de tendencias de superficie se recomienda determinar la distribución espacial con el programa “Pattern “Pattern Analysis ” y la correlación espacial por medio del programa “Spatial “ Spatial Correlation ”. ”. 5
Introducción
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
1.3 Historia del Kriging
En la mayoría de los casos la interpolación de puntos (ej. promedio móvil) genera un resultado satisfactorio, sin embargo deja unas preguntas sin responder: a) b) c)
Que Que tan tan grand grandee deber debería ía ser la dist distanc ancia ia de de los los punto puntoss que que ejerc ejercen en inf influe luenci nciaa en el valor del resultado? Exis Existe tenn mej mejor ores es maner aneras as para para dete determ rmin inar ar el resu resultltad ado? o? Cual Cual es el erro errorr aso asoci ciad adoo con con el resu resultltad adoo de de la la int inter erpo pola laci ción ón??
En 1971 estas preguntas motivaron al geo-matemático Francés Georges Matheron y el ingeniero de minería de Sud-África D.G. Krige a desarrollar un método para la interpolación interpolac ión de puntos con el objeto de ser aplicado en la l a industria minera. El método de interpolación es óptimo en el sentido que la función de ponderación se selecciona perfeccionando la función de interpolación con el fin de proveer el mejor estimado lineal imparcial (Best (Best Linear Unbiased Estimate ). ). La aplicabilidad aplicabi lidad del método se determina determina con la aceptación de que la distribución espacial espaci al de cualquie cualquierr propiedad geológica, hidrológica hidrológi ca o de suelos, conocida como variable regionalizada, es demasiado irregular para ser modelada usando una función matemática que promedie y que además puede ser representada con mejores resultados con una superficie empíricamente determinada. Antes de efectuar la interpolación se tiene que explorar y modelar estos aspectos por medio de funciones estadísticas que determinen las variables que van a ser utili zadas en el proceso. Por medio de la correlación espacial se selecciona selecci ona la función estadística que mejor se adapte a los datos disponibles. (Burrough (Burrough 1985, ILWIS 1999 ) 1.4 Aplicaciones del Kriging
En documentos existentes el proceso del Kriging tiene muchas aplicaciones en campos de exploración y explotación, modelar aspectos hidrológicos hasta la generación de mapas de caracterización de d e suelos. Es importante notar que que el proceso del Kriging ddentro entro de ILWIS únicamente se puede efectuar con base en atributos numéricos.
6
Introducción
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
2. TEORÍA T EORÍA DE GEO-ESTADÍSTICA ILWIS cuenta con las herramientas para determinar la distribución espacial y la correlación de atributos. Todos estos procesos son pre-requisito pre-requisito para la operación del Kriging. En este capítulo se explica la utilidad de cada una de estas operaciones. 2.1 Distribución espacial
La distribución espacial (Pattern (Pattern Analysis ) es una técnica que ofrece al usuario más información acerca la distribución di stribución del de l fenómeno a examinar. Se pueden distinguir distinguir 4 tipos de distribución espacial; Complete Spatial Randomness CSR (aleatorio), Cluster (agrupado), Regular (regular) y Paired (en parejas). Ver figuras 8 a 15. Es posible identificar visualmente la distribución espacial, sin embargo se recomienda graficar la distancia contra la probabilidad para encontrar 1 punto (columna Prob1Pnt graficado Prob1Pnt graficado en rojo) y para encontrar todos los puntos (columna ProbAllPnt graficado ProbAllPnt graficado en verde).
Figura 8: Distribución aleatoria
Figura 10: Distribución agrupada
Figura 9: Distancia versus probabilidad (aleatoria)
Figura 11: Distancia versus probabilidad (agrupada)
7
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Figura 12: Distribución regular
Figura 13: Distancia versus Probabilidad (regular)
Figura 14: Distribución en pareja
Figura 15: Distancia versus probabilidad (en pareja)
Notese la diferencia en el comportamiento de las curvas de probab ilidad tanto tan to en el caso de encontr encontrar ar un punto como en el de encontrarlos todos. En la distribución distribuci ón aleatoria y en la distribución regular la certeza (probabilidad (probabilidad igual a 1) de encontrar 1 punto ocurre a menor distancia que en las otras distribuciones. La diferencia entre la distribución aleatoria y la distribución regular se puede notar en la inclinación de la linea roja: la distribución aleatoria tiene un incremento gradual, mientras que la distribución regular causa una linea casi vertical. La curva roja (probabilidad de encontrar 1 punto) es muy similar para las distribuciones agrupadas y en pareja. Sin embargo, se pueden distinguir fácilmente por la curva verde (probabilidad (probabili dad de encontrar encontrar todos los puntos): para la l a distribución distribuci ón en pareja esta aumenta gradualmente, mientras para la distribución agrupada se incrementa con un comportamiento escalar. Una operación de Kriging Kriging requiere de una distribución distribuci ón aleatoria de los puntos. 8
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
En ILWIS existen dos técnicas para determinar dete rminar la distribución espacial con base en un mapa de puntos: Nearest Neighbour Distance NND (Distancia de vecino más cercano) Este método método consiste en analizar la distribución espacial, espacia l, utilizando uti lizando características característica s de la distancia entre pares de puntos. La distribución se analiza calculando las distancias entre cada punto y su primer, segundo, tercer, cuarto, quinto y sexto vecino más cercano. El resultado se puede comprobar contra las distancias esperadas en una CSR (distribuci CSR (distribución ón completamente aleatoria). En el caso de que los pares de puntos sean más cercanos se evidenciará una distribución agrupada. En caso de que estén mas distanciados, se tomará tomará como una distribución regular. regula r.
Figura 16: Tabla de distribución espacial (NND)
Reflexive Nearest Neighbour (RNN) En este método dos puntos son considerados RRN RRN cuando cuando dentro de un par de puntos cada uno es el vecino mas cercano del otro. Este cálculo se hace para el primer vecino más cercano hasta el sexto vecino más cercano (ILWIS (ILWIS 1997 ). ).
Figura 17: Información adicional acerca RNN
9
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
A continuación se presenta información adicional acerca de RNN y NND que ILWIS genera para cuatro mapas de puntos con las cuatro cuatro diferentes distribuciones: distribuci ones: aleatoria (1), agrupada (2), regular (3) y en pareja (4): Reflexive Nearest Neighbours
Reflexive Nearest Neighbours
Order Obs. values Assumed CSR ================================== 1 62.00 62.15 2 40.00 32.91 3 30.00 24.31 4 18.00 20.15 5 22.00 17.60 6 12.00 15.82
Order Obs. values Assumed CSR ================================== 1 64.00 62.15 2 32.00 32.91 3 22.00 24.31 4 24.00 20.15 5 16.00 17.60 6 8.00 15.82
Distance to Nearest Neighbours
Distance to Nearest Neighbours
Order Obs. values Assumed CSR ================================== 1 5052.37 4899.90 2 7778.14 7349.84 3 10258.35 9187.30 4 12084.69 10718.03 5 13870.36 12058.64 6 15288.82 13264.02
Order Obs. values Assumed CSR ================================== 1 1092.95 3590.43 2 1742.72 5385.65 3 2272.11 6732.06 4 2743.68 7853.72 5 3127.95 8836.06 6 3468.96 9719.31
Cuadro 1: RNN y NND para distribución aleatoria
Cuadro 2: RNN y NND para distribución agrupada
Reflexive Nearest Neighbours
Reflexive Nearest Neighbours
Order Obs. values Assumed CSR ================================== 1 2.00 61.53 2 16.00 32.58 3 0.00 24.07 4 2.00 19.95 5 30.00 17.42 6 2.00 15.66
Order Obs. values Assumed CSR ================================== 1 100.00 62.15 2 34.00 32.91 3 30.00 24.31 4 14.00 20.15 5 18.00 17.60 6 14.00 15.82
Distance to Nearest Neighbours
Distance to Nearest Neighbours
Order Obs. values Assumed CSR ================================== 1 10000.00 4494.67 2 10000.00 6742.00 3 10167.36 8427.50 4 11742.91 9831.63 5 14378.82 11061.37 6 16367.65 12167.06
Order Obs. values Assumed CSR ================================== 1 759.91 4918.68 2 8138.23 7378.03 3 8632.07 9222.53 4 11342.81 10759.13 5 11812.38 12104.88 6 14382.59 13314.88
Cuadro 3: RNN y NND para distribución regular
Cuadro 4: RNN y NND para distribución en pareja
10
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
2.2 (Auto-) Correlación espacial
Muchas variables tienen diversos d iversos valores en diferentes ubicaciones, estas pueden ser consideradas como ubicaciones geográficas aleatorias y ser analizadas usando una (auto-)correlación espacial. espacial . Ejemplos Ejemplos de este tipo de fenómeno son datos pluviométricos, concentración de elementos ele mentos tóxicos, elevación en puntos de triangulación, triangulac ión, etc. Por medio medio de la autocorrelación autocorrelaci ón espacial se puede determinar la variación espacial espacia l del fenómeno a estudiar. Normalmente, parejas de puntos con menor distancia intermedia presentan una mayor correlación en la variable, mientras que puntos con mayor distancia intermedia presentan más variación. El “autocorrelogram “autocorrelogram “ cuantifica esta relación entre distancia y varianza y permite obtener una mejor comprensión del fenómeno a estudiar. ILWIS permite el cálculo de correlación espacial omnidireccional (en todas las direcciones) direc ciones) y bidireccional (en una dirección especificada con su dirección perpendicular). perpendicular). Requiere la especificación del “Lagspacing “Lagspacing ” (intervalo de distancia). En caso de seleccionar la opción bi-direccional bi-direcci onal se debe indicar la dirección, tolerancia y opcionalmente opcional mente el ancho de banda. (Ver figura 18).
En caso de que exista un punto en el origen de la figura 18, se evalúan todos aquellos puntos que se encuentren encuentren entre la dirección (linea azul) menos la tolerancia (lineas rojas) y la dirección más la tolerancia. Cuando se encuentre un punto su varianza y otros atributos serán registrados en su correspondiente clase de “Lag spacing” (columna Distance de Distance de la figura 19).
Figura 18: Correlación espacial bi-direccional
El resultado de esta operación produce una tabla con las siguientes columnas:
Figura 19: tabla que resulta de la correlación espacial bi-direccional
11
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
Distance NrPairs I c AvgLag1 NrPairs1 SemiVar1 AvgLag2 NrPairs2 SemiVar2
WINDOWS
Indica el centro de cada clase de distancia (“Lag spacing”) El número de pares de puntos encontrados en las clases de distancia. ”Moran's I ” representa la correlación espacial espa cial de los puntos dentro de cada clase de distancia. (La relación entre el producto de la diferencia de los valores entre puntos con la diferencia global). “Geary's c ” representa la varianza espacial de los puntos dentro de cada clase de distancia (compara las diferencias en cuadrado con la media de todos los valores). Para cada clase de distancia se determina la distancia promedio entre los pares de puntos (1 es para la dirección principal de la operación bidireccional). Números de pares encontrados en cada clase de distancia (dirección principal) Semivarianza experimental para cada clase de distancia (direcc (dirección ión principal) Para cada clase de distancia se determina la distancia promedio entre los pares de puntos (2 es para la dirección perpendicular perpendicular a la dirección principal de la operación bi-direccional). bi-direccional). Números de pares encontrados en cada clase de distancia (dirección perpendicular) Semivarianza experimental para cada clase de distancia (dirección perpendicular)
En caso de efectuar la operación omnidireccional, únicamente se genera una columna para el número de pares (NrPairs), distancia promedio (AvgLag) y semivarianza experimental (SemiVar). La operación bi-direccional bi-di reccional únicamente úni camente se hace cuando se supone un comportamiento de Anisotropía en los datos. El ángulo de Anisotropía puede ser determinado por medio de una operación de “Variograma representado como superficie” (ver 2.4). Es importante que el “Lag spacing” sea seleccionado de tal forma que las clases de distancia tengan mínimo 30 puntos en cada clase. 2.2.1 Interpretación de “Moran's I” y “Geary's c” Los valores de “Moran’s “Moran’s I ” y “Geary’s “Geary’s c ” pueden interpretarse así: 0 < c < 1 c > 1 c = 1
Autocorrelación positiva Autocorrelación negativa Distribución aleatorio de los valores
12
I > 0 I < 0 I = 0
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
2.3 Semivariogramas
El modelamiento modelamiento de las semivariogramas es un aspecto muy importante en la preparación de la operación de Kriging. En este proceso se determina que modelo estadístico se ajusta mejor al comportamiento (correlación espacial) de los datos y cuales son las variables a ser usadas durante el Kriging. A partir de la tabla generada en la operación de correlación espacial se despliega la columna “Distance “Distance ” o la columna “AvgLag “AvgLag ” horizontalmente contra la columna col umna “SemiVar “SemiVar ” verticalmente. Al seleccionar el rango de distancia, no debe especificarse más que la mitad de la distancia máxima, puesto que se asume que después de una determinada distancia no existe e xiste correlación espac espacial ial entre los puntos. En caso que se use la l a opción bidireccional en la correlación espacial, se debe hacer este proceso para ambas direcciones. ILWIS soporta los siguientes variogramas teóricos para modelar el comportamiento de la correlación espacial (ver figura 20):
Figura 20: Modelos disponibles en ILWIS para modelar el semivariograma
En el ejemplo anterior se utilizaron las siguientes variables: Power model : Gaussian model: Spherical model : Circular model : Exponential model: Rational Quadratic model: Wave model:
Nugget = 215, Slope = 0.000015 y Power = 1.6 Nugget = 350, Sill = 800, Range = 35km Nugget = 215, Sill = 750, Range = 39km Nugget = 150, Sill = 650, Range = 40km Nugget = 300, Sill = 700, Range = 35km Nugget = 260, Sill = 550, Range = 35km Nugget = 0, Sill = 350, Range = 4km 13
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Notése que todos los modelos utilizan uti lizan las tres variables vari ables “Nugget”, “Sill” “Sil l” y “Range”, menos menos el Power model . El siguiente gráfico contiene la explicación de cada variable:
Figura 21: semivariograma esférico
Nugget
Sill Range Slope
Es el valor inicial de la semivarianza (a la distancia cero) donde empieza el semivariogram (corte con el eje vertical). En teoría este valor debería ser cero, puesto que a distancia cero no debería existir varianza en el atributo. Sin embargo, embargo, en la práctica existen casos donde se requiere requ iere especificar especifica r un valor diferente a cero. Es el valor de semivarianza donde la curva se deja de incrementar. Corresponde a la distancia donde la curva se deja de incrementar. Para el “Power model” representa la pendiente de la curva (ver figura 20). En caso de usar un “Power “Power ” de 1 se obtiene una linea recta. Si en el ejemplo de la figura 20 se hubiera usado un Power igual a 1, la pendiente hubiera sido la siguiente:
∆S e m i v a r = ∆ Distancia Power
800 65000
≈ 0 .012
Para el “Power model” representa la curvatura de la linea (ver figura 20).
14
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
2.4 Variograma representado como superficie
La operación “variogram “va riogram surface” utilice un mapa de puntos (o mapa raster) y calcula una superficie de valores de semivarianza, donde cada pixel en la superficie representa una clase de distancia distanci a direccional. direccional . Esta superficie permite determinar visualm visualmente ente la dirección de una posible anisotropía. anisotropía. En caso de establecerse la anisotropía anisotropía en los datos, debe medirse el ángulo principal; luego realizar la correlación espacial usando este ángulo y modelar el semivariograma, estableciendo los parámetros necesarios para el proceso de Anisotropic Kriging .
Figura 22: Generación de un variograma representado como superficie
La figura 22 contiene 4 gráficos: el primero representa un mapa de puntos donde todos los pares de puntos han sido identificados mediante lineas de color (son 6). El segundo representa estas lineas (separation (separation vectors ) dibujadas desde el orígen en ambas direcciones (es decir dirección CB y también dirección BC) con su correspondiente longitud (distancia entre puntos). El tercer gráfico muestra las clases de distancias (Lag spacing) como un mapa raster con c on las lineas li neas sobrepuestas. sobrepue stas. El cuarto muestra el resultado resultado del variograma: un mapa raster con valores de semivarianza en los pixeles correspondientes a las clases de distancia direccional. Para obtener buenos resultados se recomienda tener unos 30 pares de puntos por cada pixel. Normalmente el semivariograma representado como superficie debe ser circular (baja varianza cerca del orígen {color azul} y alta varianza {color rojo} hacia afuera). En caso que el comportamiento sea más similar a una elipse se puede concluir concluir que los datos tienen anisotropía ani sotropía y se debe usar el proceso “Anisotropic “Anisotropic Kriging ”. ”. En este caso debe medirse el ángulo de poca varianza y su tolerancia. Figura 23: “variogram surface”
15
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
2.5 Variograma de cruces
El variograma de cruces calcula valores val ores de semivarianza experimental para dos variables de entrada y un variograma de cruces para la combinación de estas dos variables. Esta aplicación aplicaci ón es útil cuando se tiene un conjunto de puntos donde la variable de interés inte rés tiene pocas observaciones (debido al costo o dificultad para medir la variable), y otra variable que si tiene observaciones. observa ciones. En caso que exista exista entre entre estas dos variables una un a correlación determinante (positiva o negativa), se puede usar el comportamiento de esta co-variable para interpolar interpol ar la variable de d e interés. Este método método de interpolación es llamado ll amado CoKriging .
Figura 24: tabla que resulta de la operación “Cross variogram”
En este caso se debe modelar el semivariograma de la semivarianza de la variable de interés (columna SemiVarA), SemiVarA), la co-variable (columna SemiVarB ) y la semivarianza de cruce entre las dos. (columna CrossVarAB ). ). (ver figura 25).
Figura 25: modelamiento de los tres semivariogramas
16
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Los valores del semivariograma de cruce pueden incrementarse o reducirse con la distancia dista ncia dependiendo de la correlación entre las variables variables A y B. B. Los valores de los semivariogramas de las variables A y B son positivos. Cuando se modelan estos tres variogramas la relación definida por Cauchy-Schwartz debe cumplirse con el fin de garantizar resultados confiables de la operación CoKriging .
Para todos los valores de distancia entre cero y la distancia limitante en la operación de CoKriging . donde: es la variable del variograma de cruce, es el valor del semivariograma para la variable A, es el valor del semivariograma para la co-variable B. La distancia limitante es la máxima (de un pixel) con respecto a un punto donde este punto todavía ejerza influencia en el resultado (del pixel). 2.6 Métodos de Kriging
ILWIS soporta diferentes tipos de procesos de Kriging: Simple Kriging, Ordinary Kriging, Anisotropic Kriging, CoKriging y Universal Kriging . Cada uno de ellos tiene aplicabilidad especifica así: 2.6.1 Simple Kriging El Simple Kriging como Kriging como su nombre lo indica ind ica es la l a forma más sencilla de Kriging. Kriging. Se debe especificar especifica r el modelo del del semivariograma semivariograma y sus variables variab les (Nugget, Sill y Range o Nugget, Nugge t, Slope y Power). Todos los puntos punto s serán incluidos inclui dos en este proceso, no puede especificarse especificarse el número de puntos a ser usado para determinar los valores de salida. Opcionalmente (como en todos los procesos de Kriging) puede generarse un mapa de errores. 2.6.2 Ordinary Kriging El Ordinary Kriging es Kriging es muy similar a Simple Kriging , con la diferencia que en este caso si puede especificarse el número de puntos que serán usados para determinar los valores en el mapa raster de salida. Este se puede hacer indicando la distancia limitante, el número mínimo mínimo y máximo máximo de puntos que se requieren para determinar los valores en el mapa de salida. Este método es más verosímil, puesto que qu e normalmente después de una cierta distancia no existe ninguna correlación entre los atributos de los puntos.
17
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
2.6.3 Anisotropic Kriging Como Como fue explicado en el numeral 2.4 (variograma representado como superficie) superfici e) en caso que exista anisotropía en los datos (comportamiento distinto en diferentes ángulos), se debe usar el Anisotropic Kriging . Para este proceso se debe modelar el semivariograma para ambas direcciones (resultado de la correlación correla ción espacial bi-direccional). Además, del ángulo de anisotropía, Anisotropic Kriging tal Kriging tal como Ordinary Kriging Kriging permite permite al usuario definir la distancia limitante y el mínimo y máximo número de puntos. En caso de obtener dos diferentes tipos de semivariogramas (ej: (ej: esférico y gausiano) gausi ano) se habla de Anisotropía zonal (no soportada por ILWIS, ver figura 26).
Figura 26: Anisotropía zonal
Es importante destacar que en el modelamiento de los semivariogramas, no se trata de ajustar todos los puntos. Los puntos más importantes están a menor distancia. La distancia en donde no se puede ajustar el semivariograma corresponde a la distancia limitante. En la figura 26 la distancia limitante es de 30 a 35 km. 2.6.4 CoKriging El proceso de CoKriging, como fue explicado en el numeral 2.5 (variograma de cruces), es útil en caso de contar con un conjunto de puntos donde la variable de interés tiene pocas observaciones, observacion es, pero existe otra variable que si tiene observaciones (debido al costo o dificultad para medir la variable).
18
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
En caso de que exista entre estas dos variables una correlación fuerte (positiva o negativa), se puede usar el comportamiento de la variable observada para interpolar la variable de interés. Este método de interpolación es llamado CoKriging . 2.6.5 Universal Kriging El Universal Kriging es Kriging es una variación del Ordinary Kriging . La diferencia es que este usa una tendencia local. El Universal Kriging consiste Kriging consiste en una superficie de tendencias que cambian gradualmente gradualmente sobre la cual se superpone sup erpone la variación variac ión a ser interpolada. interpol ada. Cuando la media de una variable regionalizada no es constante para toda la superficie se denomina “no estacionaria”. Este tipo de variables es muy frecuente y tiene dos componentes: el promedio o valor esperado de la variable regionalizada un residual determinado por la diferencia entre el valor medido y el valor esperado Puede usarse el proceso de Universal Kriging en Kriging en vez de quitar este tipo de variación debido a la tendencia local.
19
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
2.7 Diagrama de flujo de Kriging El siguiente es un diagrama de flujo para determinar el proceso de Kriging que debe usarse y para verificar si el conjunto de datos es apto para su aplicación:
Los rombos representan la toma de una decisión deci sión y los cuadros representan representan operaciones dentro de ILWIS. En varios casos no debe usarse Kriging: Cuando no se tiene una buena distribución Si se tienen pocas observaciones y no se cuenta con una co-variable Cuando el resultado del “Cross variogram” no cumple la relación Cauchy-Schwartz. Si se detecta una “Anisotropía “Anisotropía zonal ” es posible efectuar Simple o Simple o Ordinary Kriging . Luego de efectuar un proceso de Kriging se recomienda hacer una verificación de los resultados por po r medio medio de un proceso llamado l lamado “intervalo “intervalo de confiabilidad ” (ver numeral 3.4).
20
Teoría de geo-estadísticas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
3. EJERCICIO DE KRIGING Para este ejercicio se tomaron como base los datos existentes en la versión 2.2 de ILWIS elaborados por A. Gieske de la Universidad de Botswana y la Universidad Libre de Amsterdam. Dichos datos fueron editados para mejorar el nivel de aprehensión del ejercicio y del mismo proceso de Kriging . Cambie al sub-directorio “C14 “C14 " del directorio “C:\ILWIS22\DATA\ “C:\ILWIS22\DATA\”” . Se asume que ILWIS esta instalado en el directorio “C:\ILWIS22\ “ C:\ILWIS22\”” y los datos para este curso fueron copiados en dicho directorio. 3.1 Reconozca la información
El primer paso es conocer los datos y visualmente. Determine su distribución espacial y atributos. Establezca si cuenta con observaciones suficientes. Estos pasos definen si el conjunto de datos es apto para ser utilizado en el proceso de Kriging . 3.1.1 Despliegue y consulte la información Despliegue los puntos para evaluar de manera visual su distribución espacial. Haga doble clic sobre el mapa de puntos “C14 “C14 ". ". En la ventana del diálogo oprima el botón “Symbol “Symbol”” y en la ventana de la definición del símbolo seleccione el “Fill color ” rojo, luego oprima dos veces el botón “OK “ OK”. ”. Los puntos puntos se despliegan despli egan dentro de una venta ventana na con círculos rojos. Los puntos están e stán bien bien distribuidos dentro dentro de la ventana. Que tipo de distribución es, según su opinión? Aleatoria Agrupada Regular En pareja Cambie el despliegue de los símbolos con el fin de poder medir una distancia entre los puntos: Seleccione el menú “Layers “Layers”” en la ventana del mapa, luego seleccione la opción Options” y “1 pnt C14". C14". En la ventana de las opciones del despliegue “Display Options” Plus”, luego oprima dos veces oprima el botón “Symbol “Symbol”” y cambie el símbolo s ímbolo a “+ “+ Plus”, el botón “OK “OK”. ”. Con el botón
mida las distancias distancias mínima mínima y máxima máxima existente entre los puntos:
Distancia mínima: ............. .............
Distancia máxima: ................ ................
21
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Que tipo de dominio tiene el el ma mapa (ver propiedades del ma mapa)?
...............
Cuantos puntos tiene el mapa (ver propiedades del mapa)?
...............
Cuantos pares de puntos se pueden crear con estos puntos?
...............
Despliegue la tabla de atributos del mapa: Haga doble clic sobre la tabla “C14 “C14 ". ". La tabla se abre y se puede observarse que tiene además de los códigos del dominio un atributo (columna “C14 “C14 "). "). Los valores que aparecen el la columna “C14 “ C14 " representan porcentajes de contenido de 14 C (p.m.c. = porcentaje de “carbono “carbono 14 ”). ”). Cuales son los valores mínimos y máximos de la columna “C14 “ C14 "? "? Valor mínimo: ............. ............. Valor máximo: ................ ................ Se debe establecer cual es la varianza de los datos de atributos. atributos . Esto se hace dentro de la tabla de atributos: Entre al menú “Column “Column”” y seleccione seleccione la l a opción “Statistics “Statistics”. ”. Seleccione la función “Variance ” y asegure que la columna “C14 “C14 " este seleccionada como columna. Oprima el botón “OK “OK”. ”. Anote el valor de varianza:................ La varianza es un indicador para la variable “Sill” cuando el semivariograma no tiene “Nugget”. Cierre la tabla de atributos y cierre el mapa de puntos. 3.1.2 Análisis de distribución espacial En este proceso se determina la distribución espacial del mapa de puntos: En el menú “sensitivo “sensitivo al contexto ” del mapa de puntos “C14 “C14 " seleccione la opción Statistics” y luego “Patter n Analysis Analysi s”. Por defecto el sistema tiene seleccionado “Statistics” “Pattern selecci onado en mapa de puntos “C14 “C14 " como mapa de entrada. Escriba “distrib “distrib ” como nombre de la tabla de salida, active la opción “Show “Show”” y oprima el botón “OK “ OK”. ”.
22
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Con el fin de interpretar los valores generados en esta tabla, se debe graficar la distanci a (abscisa) y la probabilidad de encontrar 1 punto y todos los l os puntos (ordenada): Graph”. En la caja En el menú “Options “Options”” de la tabla seleccione la opción “Show “Show Graph”. de diálogo seleccione la columna “Distance “Distance ” como eje X (izquierda) y la columna “Prob1Pnt ” como eje Y (derecha). Oprima el botón “OK “OK”. ”. Acepte los valores por Graph” y oprime el botón “OK defecto de la ventana “Edit “Edit Graph” “ OK”. ”. La curva de distancia versus la probabilidad se despliega en color rojo. En la ventana del gráfico, entre al menú “Graph “Graph”” y seleccione la opción “Add “Add Graph form Columns”. Columns”. Seleccione “Distance “Distance ” como eje X, “ProbAllPnt “ProbAllPnt ” como eje Y y Graph” oprima el botón “OK”. Acepte los valores por defecto de la ventana “ Edit Graph” y oprime el botón “OK “ OK”. ”. Se actualiza el despliegue asignando a la curva de distancia versus probabilidad de encontrar todos todos los puntos el color verde. Compare Compare este gráfico con los l os figuras figur as de las l as páginas 7 y 8 de este manual. Con que distribución corresponde el comportamiento del gráfico? Aleatoria Agrupada Regular En pareja Este corresponde con lo que usted visualmente había concluido? 3.2 Generación de variogramas
Siguiendo el diagrama de flujo de la página 20, el proceso anterior determinó que el conjunto de puntos tiene una distribución espacial adecuada con buenas observaciones (todos los puntos tiene información acerca el atributo de interés C 14 ). En este sección generamos el variograma vario grama representado como superficie superfic ie para determinar si el el conjunto de datos presenta anisotropía. ani sotropía. En este momento momento queda descartado el proceso p roceso de “CoKriging” por no tener otra (co-)variable y además se cuenta con buenas observaciones. 3.2.1 Variograma representado como superficie (Anisotropía) El proceso de la generación generaci ón de un variograma vari ograma representado representado como superficie nos permite verificar la existencia de anisotropía en los datos: En el menú “sensitivo “sensitivo al contexto ” del mapa de puntos “C14 “C14 " seleccione la opción Statistics” y luego “Variogram Surface”. Asegure que “Statistics” “Variogram Surface”. que el mapa “C14 “C14 " y la columna “C14 " están seleccionados. 23
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
spacing” (intervalo de distancia) Especifique un ‘Lag ‘Lag spacing” dista ncia) de 5000 metros y mantenga el número de “Lags “Lags”” en 10. Especifique el e l nombre “supvario “supvario ” como como mapa de salida, sali da, acepta los valores por defecto de rango y precisión, prenda la opción “Show “Show”” y oprima el botón “OK “OK”. ”. En la caja del diálogo de opciones de despliegue oprima el botón “OK “OK”” para visualizar el mapa. Se crea un mapa de 19 por 19 pixeles. Para una mejor comprensión se debe adicionar una grilla para indicar el origen: Annotation” y luego la opción En el menú “Layers “Layers”” seleccione la opción “Add “Add Annotation” Lines”. Desactive la opción “Coordinates “Grid Lines”. “ Coordinates”, ”, cambie la distancia (“Grid (“Grid distance”) distance”) a 50000 y oprima el botón “OK “ OK”. ”. La grilla muestra el origen del mapa (el centro). Esto permite al usuario medir el ángulo de posible anisotropía. La georeferenciación es interna y no está disponible para el usuario. En caso de no tener anisotr anisotropía opía el mapa tendrá un comportamiento comportamiento circular, con los colores azules (baja varianza) en el centro y colores rojo en las partes afuera (alta varianza). Es evidente que el comportamiento no es circular, pero más bien de tipo elíptico indicando indicand o la presencia de anisotropía en los datos. Para mediar el ángulo de la anisotropía se debe usar la herramienta herramienta de medición y medir la línea central de la la elipse por los pixeles de baja varianza (color azul). Anot Anotee el áng ángul uloo prin princi cipa pall de la la anis anisot otro ropí pía: a:
.... ...... .... .... .... ..° °
Anote también un án ángulo de de tolerancia: ia:
............°
Este valor del ángulo áng ulo principal princi pal de anisotropía aniso tropía se debe debe usar en la correlación espacial bi-direccional. bi-direccional. Figura 29: mapa supvario
3.2.2 Correlación espacial espa cial bi-direccional Ahora debe determinar la correlación espacial bi-direccional: En el menú “sensitivo “sensitivo al contexto ” del mapa de puntos “C14 “C14 " seleccione la opción Statistics” y luego “Spatial Correlation”. Seleccione la opción “Bidirectional “Statistics” “Spatial Correlation”. “Bidirectional”” y spacing” de 5000 metros. Especifique la dirección y tolerancia especifique un “Lag “Lag spacing” toleran cia width” desactivada. determinada en el ejercicio anterior, mantenga la opción “Band “Band width” Show” y Escriba “cebi5 “cebi5 " como nombre de la tabla de salida, active la opción ” Show” oprima el botón “OK “ OK”. ”. 24
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
El la tabla se debe visualizar el comportamiento de la semivarianza contra los clases de distancia para ambas direcciones: Graph”. En la caja En el menú “Options “Options”” de la tabla seleccione la opción “Show “Show Graph”. del diálogo seleccione la columna “Distance “Distance ” como eje X (izquierda) y la columna “SemiVar1” SemiVar1” como eje Y (derecha). Oprima el botón “OK “ OK”. ”. En la caja del diálogo Graph” seleccione “Point “Edit Graph” “Point”” como tipo de gráfico, cambie el tamaño del símbolo de “3 " a “ 5 ", ", seleccione el color rojo como “Color “Color”. ”. El rango del eje X no debe ser mayor a la mitad de la distancia di stancia máxima: especifique 45000 450 00 metros como máxima X y oprima el botón “OK “OK”. ”. Es posible desplegar los modelos en ventanas separadas pero por comodidad use la misma ventana: En el menú “Graph “Graph”” de la l a ventana del gráfico, seleccione la opción “Add “Add Graph from Columns”. Columns”. Seleccione las columnas “Distance “Distance ” contra “SemiVar2 “SemiVar2 ". ". Para el despliegue use símbolo de puntos, p untos, tamaño 5 y color c olor verde con c on la misma definició definiciónn del eje X. No cierre las ventanas de gráfico y tabla. 3.3 Modelar el semivariograma
El modelamiento del semivariograma es fundamental para el proceso del Kriging. Krigin g. Se basa en ajustar una función estadística estadísti ca al comportamiento calculado en la correlación correlación espacial. espacial. 3.3.1 Adicionar modelos Use el gráfico de la página 13 para definir cual o cuales modelos pueden ajustarse a los dos semivariogramas experimentales. En el menú “Graph “Graph”” de la ventana del gráfico, seleccione la opción “Add “Add Semivariogram Model”. Model”. Seleccione el semivariograma que más se ajuste y especifique las tres variables para definir el modelo. Oprima “OK “OK”” para adicionar el modelo a la ventana del gráfico. En caso que el modelo no se ajuste y se requiere edición de las l as variables, variables, entre al Management”. Seleccione el gráfico menú “Graph “Graph”” y seleccione seleccione la opción opci ón “Graph “Graph Management”. Graph”. Haga los ajustes necesarios y oprima el a editar y oprima el botón “Edit “Edit Graph”. Management”. botón “OK “OK”” para regresar a la ventana de “Graph “Graph Management”. Seleccione el siguiente gráfico a editar o oprima “OK “OK”” para cerrar la ventana. 25
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Anote los modelos para la semivarianza de la dirección perpendicular (SemiVar2, color verde): Modelo 1: ...................... Modelo2: ...................... Nugget: ...................... Nugget: ...................... Sill: ...................... Sill: ...................... Range: ...................... Range: ...................... Es más complejo definir el modelo para la dirección principal (SemiVar1 ( SemiVar1). ). Observe la figura abajo:
Figura 30: Modelamiento de semivariogramas
Como puede verse en la figura 30, la semivarianza correspondiente a las distancias 0, 5 y 15 km no se ajustan al modelo. 3.3.2 Anisotropía zonal Se puede notar que el comportamiento del modelo de semivariograma es distinto para la dirección principal y la dirección perpendicular. Esto significaría que se presenta anisotropía zonal en los datos y ILWIS ILW IS no soporta el procesamiento de anisotropía zonal. De todas formas, el cálculo cál culo no ha sido muy exacto, puesto que que para asegurar un resultado res ultado confidencial se requiere tener al menos 30 puntos por clases de distancia (por pixel) y solamente un rango pequeño cumple con esta condición (ver columnas “NrPairs1 “NrPairs1"" y “NrPairs2 "). "). Es posible que la diferencia en el número de puntos por pixel en el variograma representado como superficie cause la aparente anisotropía.
26
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
3.3.3 Correlación espacial espa cial omnidireccional Se recomienda en este caso hacer una correlación espacial omnidireccional y continuar con un proceso de Kriging descartando el “Anisotropic “Anisotropic Kriging”. En el menú “sensitivo “sensitivo al contexto ” del mapa de puntos “C14 “C14 " seleccione la opción Statistics” y luego “Spatial Correlation”. Seleccione la opción “Omnidirectional “Statistics” “Spatial Correlation”. “Omnidirectional”” spacing” de 5000 metros. Escriba el nombre “ceom5 y especifique un “Lag “Lag spacing” “ceom5 " como nombre de la tabla de salida, active la opción “Show “ Show”” y oprima “OK “ OK”. ”. Visualice el comportamiento de la semivarianza contra las clases de distancia: Graph”. En la caja En el menú “Options “Options”” de la tabla seleccione la opción “Show “Show Graph”. de diálogo seleccione la columna “Distance “Distance ” como eje X (izquierda) y la columna OK”. En la caja de diálogo “Edit “SemiVar ” como eje Y (derecha). Oprima el botón “OK”. “ Edit Graph” Graph” seleccione “Point “Point”” como tipo de gráfico, g ráfico, cambie cambie el tamaño del símbolo de “3 " a “5 " y seleccione el color rojo como “Color “ Color”. ”. Especifique 45000 metros como máximo X y oprima el botón “OK “OK”. ”. 3.3.4 Adicional modelos Se tiene que volver a modelar el semivariograma, puesto que su comportamiento es distinto al comportamiento calculado bi-direccionalmente: En el menú “Graph “Graph”” de la ventana del gráfico, seleccione la opción “Add “Add Semivariogram Model”. Model”. Seleccione el semivariograma que más se ajuste y especifique las tres variables para definir el modelo. Oprima “OK “OK”” para adicionar el modelo a la ventana del gráfico. Hay 4 modelos que se ajusten a la semivariograma experimental: Modelo 1: Nugget: Sill/Slope: Range/Power:
...................... ...................... ...................... ......................
Modelo 2: Nugget: Sill/Slope: Range/Power:
...................... ...................... ...................... ......................
Modelo 3: Nugget: Sill/Slope: Range/Power:
...................... ...................... ...................... ......................
Modelo 4: Nugget: Sill/Slope: Range/Power:
...................... ...................... ...................... ......................
27
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS 2
3.3.5 Verificación y selección de los modelos (Goodness of Fit R ) Es posible que existan varios modelos que se ajusten satisfactoriamente al semivariograma experimental. experimental. En este caso se puede hacer un cálculo cá lculo de “Goodness “Goodness of 2 Fit R ”. Los valores que resultan de esta operación permiten seleccionar el modelo más apropiada para el proceso. El “Goodness “Goodness of Fit R 2 ” se determina con la siguiente fórmula:
donde: ^ valor experimental de semivarianza semivarianza generado durante la operación de correlaci ón espacial o variograma de cruces. valor de semivarianza determinado por medio del modelamiento del semivariograma. N el número total de clases de distancias (Lags). Observe el resultado de analizar el “Goodness “Goodness of Fit R 2 " Se modeló 3 funciones estadísticas: estadísti cas: Power model (azul), Exponential model (magenta) y
28
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS 2
Spherical model (verde). Es necesario entender que no se debe calcular el R para toda la columna, puesto que esto únicamente tiene sentido para el intervalo de distancia que se utilizará como distancia limitante. En este caso se hizo 2 cálculos para dos intervalos de distancia: distanci a: La de abajo (de 2.5km hasta 37.5km) y la de arriba (de 2.5km 2 .5km hasta 12.5km y 17.5km hasta 37.5km). El último excluyó el valor correspondiente correspon diente a la distancia 15km y así así increm incrementa enta el valor de de 2 R (ejemplo: el Power model para el primer intervalo tuve un valor 0.63 y excluyendo la distancia 15km se incrementó a 0.92). Haciendo este cálculo cálcul o se puede determinar la distancia limitante necesaria en la mayoría mayoría de los procesos p rocesos de Kriging, Kri ging, y evaluar ev aluar cual modelo se ajuste mejor. mejor. Este Este proceso se hace en la tabla de la correlación espacial o la tabla de la variograma de cruces: Es necesario anotar para que rango de distancia los modelos se están ajustando: Modelo 1: ...................... Dist istanci anciaa mínima ima: ...................... Dist Distan anci ciaa máxim áxima: a: .... ...... .... .... .... .... .... .... .... .... ..
Modelo 2: ...................... Dist istancia ncia mí mínima: ima: ...................... Dist Distan anci ciaa máxim áxima: a: .... ...... .... .... ........ ...... .... .... ....
Modelo 3: ...................... Dist istanci anciaa mínima ima: ...................... Dist Distan anci ciaa máxim áxima: a: .... ...... .... .... .... .... .... .... .... .... ..
Modelo 4: ...................... Dist istancia ncia mí mínima: ima: ...................... Dist Distan anci ciaa máxim áxima: a: .... ...... .... .... ........ ...... .... .... ....
Cual modelo, según su criterio, se ajuste mejor a los valores experimentales? Cierre la ventana de los gráficos Entre al menú “Columns “Columns”” de la tabla y selecciona la opción “Semivariogram “Semivariogram”. ”. En la caja del diálogo seleccione el primer modelo (Power (Power model ), ), especifique 250 como “Nugget “Nugget”, ”, el valor 0.0001 como “Slope “Slope”” y el valor 1.5 como 1.5 como “Power “Power”. ”. Escriba “semipow como semipow como nombre de la columna colu mna y oprime “OK “OK”. ”. Acepta los valores que salen por defecto y oprime “OK “ OK”. ”. RESPUESTAS, Repite este proceso p roceso para los otros tres modelos modelos (ver capítulo 5. RESPUESTAS, punto 3.3.4 ), ), creando las columnas “semigau “semigau ” para el “Gaussian “Gaussian model ”, ”, “semirat ” para el “Rational “Rational Quadratic model ” y “semiwav “semiwav ” para el “Wave “Wave model ”. ”. No es útil determinar determinar el R 2 para toda la columna, como se indica en el manual de ILWIS, pues solamente estamos interesados en el comportamiento hasta la distancia limitante (entre 25 y 45km).
29
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Por tanto debemos crear una columna con dominio de tipo clase, c lase, e indicar indi car en ella las filas de la tabla a incluir en el cálculo: Column”. Especifique el Entre al menú “Columns “Columns”” y seleccione la opción “Add “Add Column”. nombre “usar “usar ” como nombre de la colu columna mna y oprima el botón pa para ra crear el dominio. domain” especifique el nombre En la caja caja de diálogo di álogo “Create domain” n ombre “usar “usar ” como nombre del dominio, asegurese que el tipo “Class “ Class”” está seleccionado y opri oprima ma “OK “OK”. ”. Luego en la ventana de dominio adicione el elemento “usar “usar ”: ”: entre al menú “Edit “ Edit”” y item”, luego especifique “usar seleccione la opción “Add “Add item”, “usar ” como nombre del elemento, no especifique ningún ningú n código y oprima “OK “OK”. ”. Cierre Cierre el dominio y oprima OK” en la caja de diálogo “Add Column”. “OK” “Add Column”. La columna se adiciona adicion a a la tabla y está llena lle na de interrogantes “?”. Es necesario necesario editar la columna: Haga clic en el primer campo de la columna “usar “ usar ”. ”. Haga clic sobre el botón que aparece en la parte pa rte derecha derecha del campo y seleccione selecci one el elemento “usar “usar ”. ”. Repita este proceso para la segunda seg unda y hasta la sexta fila. Confirme el último campo con la tecla y asegurese que los campos tienen un color blanco de fondo y no azul. Con este proceso se ha indicado que el intervalo de 0 a 25km de distancia será incluido en el cálculo. Más adelante, después de efectuar los cálculos para determinar el R 2 , editará esta columna para cambiar este intervalo y actualizará los R 2 para determinarlos para diferentes intervalos de distancia. Por medio de los siguientes cálculos obtendrá el R 2 . En la línea de comando escriba las siguientes form formulas ulas y acepte los valores de rango y precisión precisión por defecto d efecto (se estará estará model”, “gau ” para el utilizando “pow “pow ” para los cálculos relacionadas con el ”Power ”Power model”, model”, “rat ” para el “Rational model” y “wav ” para el “Wave “Gaussian model”, “Rational Quadratic model” “ Wave model”): model”): pow1=sq(semivar-semipow) gau1=sq(semivar-semigau) rat1=sq(semivar-semirat) wav1=sq(semivar-semiwav) Con estas formulas se generan 4 columnas con la diferencia entre la semivarianza experimental y la semivarianza de cada modelo en el cuadrado. Luego se calculará el promedio de la semivarianza experimental (columna “SemiVar “SemiVar ”) ”) usando la columna “usar “usar ”: ”: 30
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Entre al menú “Columns “Columns”” y seleccione la opción “Aggregation “Aggregation”. ”. Seleccione la columna “SemiVar “SemiVar ”, ”, asegurese de que la función “Average “Average”” está seleccionada, By” y asegurese de que la columna “usar active la opción “Group “Group By” “usar ” está seleccionada. seleccion ada. Escriba el nombre “AvgSV “AvgSV ” como nombre de salida y oprima el botón OK”. “OK”. El la nueva columna aparecen 2 valores: tanto para los elementos “usar “usar ” como para los interrogantes “? “? ” de la columna “usar “usar ” se determinó el promedio de la semivarianza experimental. Por medio del siguiente cálculo se determina la diferencia en el cuadrado entre la semivarianza experimental y el promedio de la semivarianza experimental: sqdif=sq(semivar-avgsv) Luego sumamos esta diferencia en el cuadrado usando nuevamente la columna “usar “usar ”: ”: Entre al menú “Columns “Columns”” y seleccione la opción “Aggregation “Aggregation”. ”. Seleccione la By” y columna “sqdif “sqdif ”, ”, cambie la función a “Sum “Sum”, ”, prenda la opción “Group “Group By” asegurese de que la columna “usar “usar ” está seleccionada. Escriba el nombre “SumSqDif ” como nombre de salida y oprima el botón “OK “ OK”. ”. Ahora se sumarán la primera columna que se generó en los cuatro modelos (“pow1 (“pow1", ", “gau1", gau1", “rat1 “rat1"" y “wav1 “wav1") ") usando en esta operación la columna “usar “usar ”: ”: Entre al menú “Columns “Columns”” y seleccione la opción “Aggregation “Aggregation”. ”. Seleccione la By” y columna “pow1 “pow1”, ”, cambie la función a “Sum “Sum”, ”, active la opción “Group “Group By” asegurese de que la columna “usar “usar ” está seleccionada. seleccionad a. Escriba el nombre “pow2 “pow2 ” como nombre de salida y oprima el botón “OK “OK”. ”. Repita este proceso para las columnas “gau1 “gau1", ", “rat1" rat1" y “wav1" wav1" generando las columnas “gau2 “gau2 ", ", “rat2 “rat2 " y “wav2 “wav2 ". ". Ya se puede determinar el R 2 para los cuatro modelos para el intervalo 0 a 25km: PowR2 = 1 - (pow2 / sumsqdif) GauR2 = 1 - (gau2 / sumsqdif) RatR2 = 1 - (rat2 / sumsqdif) WavR2 = 1 - (wav2 / sumsqdif) Se puede observar que los l os valores de los primeros primeros 6 registros fluctúan entre 0.5 y 0.6. Los valores mas cercanos a 1 representan un mejor modelo. También se puede notar que únicamente el “Wave “Wave model ” tiene valores positivo en los otros registros. 31
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Aunque el “Power model” se ajusta mejor en este intervalo, el “Wave model” representa el mejor modelo globalmente. Anote el R 2 para el intervalo de 0 a 25km de los 4 modelos en la siguiente tabla: Distancia
Modelo 1 (pow) Modelo 2 (gau)
Modelo 3 (rat)
Modelo 4 (wav)
0-25 km 0-35 km 0-40 km 0-45 km Promedio
Es posible editar la columna “usar” e incluir las distancia para obtener el segundo intervalo de 0-35km: Edite las filas 7 y 8 de la columna “usar “usar ” cambiando la interrogación “? “? ” por el . elemento “usar “usar ”. ”. Asegurese de confirmar el último cambio con la tecla . Para actualizar una columna se debe hacer doble clic sobre el nombre de la columna (empiece con la columna “PowR2 “PowR2 ") ") y en la caja del diálogo de las propiedades de la columna cambie el rango de valores de la columna a -1000 Date”. Los valores en la columna hasta 1 y oprima el botón “Make “Make Up to Date”. columna se actualizan. Repite este proceso para las columnas “GauR2 “GauR2 ", ", “RatR2 “RatR2 " y “WavR2 “WavR2 ". ". Es necesario hacer este proceso para los otros intervalos: 0-40km y 0-45km: Edite la fila 9 (intervalo 0-40km) de la columna “usar “ usar ” cambiando la interrogación “? ” por el elemento “usar “ usar ”. ”. Asegurese de confirmar el último cambio con la tecla . . Luego actualice las columnas y anote los R 2 en el intervalo de distancia correspondiente. Edite la fila 10 (intervalo (in tervalo 0-45km) de la columna “usar “usar ” cambiando la interrogación interrogación “? ” por el elemento “usar “usar ”. ”. Asegurese Asegurese de confirmar confirmar el último cambio cambio con la tecla 2 . . Luego actualice las columnas y anote los R en el intervalo de distancia correspondiente.
32
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Calcule también el promedio de los intervalos para cada modelo y anote los valores. Cual de los cuatro modelos se ajusta mejor al semivariograma experimental? Al incrementar el intervalo interv alo de distancia, dista ncia, el valor R 2 incrementa a 0.7667 para el intervalo 0-50km, luego se reduce. Así se puede concluir que la distancia limitante es de 50km. Cierre la tabla. 3.4 Interpolación por medio de Kriging
Ya que se ha establecido el modelo mas adecuado, se puede efectuar el proceso de Kriging. En este momento todavía se pueden usar tres tipos de Kriging: Simple Kriging Ordinary Kriging Universal Kriging Se puede descartar el “Simple “Simple Kriging ”, ”, puesto que el óptimo modelo debe usar la distancia limitante li mitante de 50 km y esto no es posible en “Simple “Simple Kriging Kriging ”. ”. El “Universal Kriging ” también también se puede puede excluir excluir puesto que qu e no se ha detectado ningún n ingún tipo de tendencia en los l os datos. De esta manera se puede concluir que se debe usar el “Ordinary “Ordinary Kriging ”: ”: En el menú “sensitivo “sensitivo al contexto ” del mapa de puntos “C14 “C14 " seleccione la opción Interpolation” y luego “Kriging “Interpolation” “Kriging”. ”. Asegurese de que el mapa “C14 “C14 " y la columna Kriging”. Especifique “C14 " están seleccionados y que el método es “Ordinary “Ordinary Kriging”. Especifique el el model” con el valor 250 como “Nugget “Wave model” “Nugget”, ”, el valor 650 como “Sill “Sill”” y el valor 50000, coloque el mínimo 7500 como “Range “Range”. ”. Especifique la l a distancia limitante limitante de 50000, 16. Escriba el nombre “OKwa50 número de puntos en 3 y deje el máximo en 16. “OKwa50 " como nombre, seleccione la georeferenciación “c14 “c14 " y active la opción “Show “Show ”. ”. Active Map” y escriba como nombre para este mapa “OKwa50e además la opción “Error “Error Map” “OKwa50e ”, ”, oprima el botón “OK “ OK”. ”. Este proceso es un poco demorado. Cuando aparezca la caja de diálogo de las opciones de despliegue haga clic en “OK “OK”” para desplegar el mapa. Despliegue también el mapa de errores “OKwa50e “OKwa50e ”. ”. Hagalo sobre el mapa de puntos “C14 “C14 " usando el símbolo “+”. Determine Determine el error mínimo y máximo del mapa: Error mínimo:
.................
Error máximo:
33
.................
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
3.5 Verificación del resultado (Confidence interval)
El mapa resultado del proceso proces o de Kriging contiene con tiene los estimat esti mativos ivos del del fenómeno fenómeno a estudiar y el mapa mapa de errores contiene una aproximación del error para el valor estimado. Es posible crear cre ar un mapa donde donde se indique con cierta certeza el rango en que se encuentran encu entran los valores del fenómeno estudiado. El rango se define de la siguiente manera: mínim ínimo: o: esti estim mativ ativoo Kri Krigi ging ng - (Er (Erro rorr * c ) máxim áximo: o: esti estim mativ ativoo Krig Krigin ingg + (Err (Error or * c ) “c ” representa el valor crítico y depende del intervalo de confiabilidad: intervalo de confiabilidad valor crítico “c ”
9 0%
9 5%
9 7 .5 %
98 %
99 %
9 9 .5 %
1.282
1 .6 4 5
1.960
2 .0 0 0
2 .3 2 6
2 .5 7 6
En la línea de comando del menú principal de ILWIS se pueden ejecutar las siguientes formulas para generar los mapas mapas de nivel nivel inferior infe rior y nivel superior s uperior para cada cad a pixel de C14. Es necesario corregir los valores en nuestro ejemplo para que queden entre 0 y 100%: Kinf:=iff(OKwa50-1.96*OKwa50e<0,0,OKwa50-1.96*OKwa50e) Ksup:=iff(OKwa50+1.96*OKwa50e>100,100,OKwa50+1.96*OKwa50e) En este ejemplo se está usando un valor crítico de 1.96 correspondiente al intervalo de confiabilidad de 97.5%. Cambie los valores a 0 hasta 100 con precisión 0.01 y oprima OK”. “OK”. Por medio de la siguiente formula se puede crear un mapa donde se combinen combinen los niveles nive les superiores e inferiores: conf975 := iff((%L mod 8) - (%C mod 8) <= 0,Ksup,Kinf) en donde: %L valo valorr de de núm númer eroo de de fila fila en el mapa apa de de sal salid idaa %C valor valor de númer númeroo de de colu colum mna en el mapa de salida salida mod mod valor residual residual de la división división enter enteraa (ej: (ej: 10 mod mod 3 = 1) El mapa que resulta tiene cuadros de 8x8 pixeles, pixeles, donde la parte inferior izquierda de cada cuadro representa el nivel inferior (mapa “Kinf “Kinf ”) ”) y la parte superior derecha representa el nivel superior (mapa “Ksup “Ksup ”). ”). (Ver figura 32)
34
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Figura 32: mapa de intervalo de confiabilidad
En el caso de que se cuente con un conjunto de datos relativamente grande se recomienda subdividirlo en dos partes con el fin de usar una parte para la verificación de los resultados del proceso de Kriging. Esto se puede hacer usando el operador operador de ILW IS “MAPVALUE”. Se debe abrir el mapa de puntos como tabla y escribir en la linea de comando de la tabla la siguiente formula: KRIG=MAPVALUE(OKWA50.MPR,COORD(X,Y)) para obtener el valor del resultado de Kriging en la tabla y: KERR=MAPVALUE(OKWA KERR=MAPVALUE(OKWA50E.MPR,COORD(X,Y)) 50E.MPR,COORD(X,Y)) para obtener el valor del mapa de errores en la tabla. Luego se puede adelantar el análisis del intervalo de confiabilidad en la tabla usando diferentes valores del valor crítico “c “ c ” y así evaluar eval uar si cada punto queda queda dentro del rango ran go de confiabilidad.
35
Ejercicio de Kriging
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
4. EXPOSICIÓN SOLAR En esta capítulo se impartirán los fundamentos que se requieren para adelantar un estudio de exposición exposición solar. Un proceso de modelamient modelamientoo de radiación solar incidente requiere datos auxiliares sobre la trayectoria solar, por lo que se presentarán las ecuaciones básicas que la definen. 4.1 La trayectoria solar
Los parámetros que definen la trayectoria del sol constituyen datos externos imprescindibles para llevar a cabo el modelamiento de la radiación solar. La localización aparente del sol s ol depende de d e una serie de parámetros de los cuales los más más importantes son: La latitud del punto a estudiar, que puede variar en un rango de -90° (polo Sur) hasta 90° (polo Norte). La declinación solar, solar, variable según la época del año en un rango de -23½° (invierno) y 23½° (verano) en el Hemisferio Norte. El ángulo horario, horario, dependiente de la hora hora del día y variable en e n un círculo de d e 360° centrado en el punto a analizar. La localización locali zación del sol se suele expresar en coordenadas coordenada s esféricas: esféricas: acimut y elevación angular sobre el horizonte , cuyas expresiones de cálculo son las l as siguientes: s e n α = ( s e n c o s ϕ =
D ⋅ s e n L)
( c o s L ⋅ s e n D)
+ (c o s
D ⋅ c o s L⋅ c o s H)
− ( c o s D ⋅ s e n L⋅ c o s
H)
c o s α
en donde: D es la Declinación solar L es la Latitud geográfica H es el ángulo Horario (Felicísimo, 1994 ) 4.1.1 Declinación solar La declinación decli nación solar sol ar para una determinada época del del año puede determinarse a partir de tablas o, aproximarse mediante expresiones empíricas. Una de las más simples es: D
= 2 3.5 ⋅ s e n [0 .9 8 6 ⋅ ( 2 8 4 + d ) ]
36
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
La siguiente figura muestra la relación entre el día del año y la declinación solar:
Figura 33: Declinación solar en función del día del año
En ILWIS para determinar D, se puede efectuar el siguiente cálculo en la línea de comando (se debe cambiar d por el día de interés): ? 23.5*SIN(DE 23.5 *SIN(DEGRAD(0.986*(284+ GRAD(0.986*(284+d ))) ))) Anexo se encuentra encuentra una tabla con la declinación solar para cada día del año. 4.1.2 Latitud geográfica La latitud influye en gran medida en la trayectoria del sol. Las influencias de las épocas del año no se notan cerca a la linea ecuatorial como es el caso de Colombia. Colombia está ubicada ubi cada entre en tre latitudes lati tudes ± 4°15' 4°15' Sur y ± 12°30' 12°30' Norte. Norte . 4.1.3 Ángulo horario En el hemisferio Norte (donde está ubicado la l a mayor parte de Colombia), el ángulo horario suele tomarse como cero en el mediodía, cuando cuand o el acimut solar es de 180° 1 80°y, por tanto, t anto, el sol está situado al Sur. Los ángulos son negativos hacia el Oeste y positivos hacia el Este, con intervalos de 15° por hora. Observe la siguiente tabla donde se relaciona el ángulo para cada hora entre las 6:00 am y las 6:00 pm.
37
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
Hora
Ángulo horario (H)
Hora
WINDOWS
Ángulo horario (H)
6:00 a.m.
90 °
1:00 p.m.
- 1 5°
7:00 a.m.
75 °
2:00 p.m.
- 3 0°
8:00 a.m.
60 °
3:00 p.m.
- 4 5°
9:00 a.m.
45 °
4:00 p.m.
- 6 0°
10:00 a.m.
30 °
5:00 p.m.
- 7 5°
11:00 a.m.
15 °
6:00 p.m.
- 9 0°
mediodía
0°
4.2 Ocultamiento topográfico
La existencia de zonas de sombra es una variable de gran interés en regiones montañosas, donde el relieve puede ser factor determinante del clima local. La sombra en un punto puede deberse a dos circunstancias: Autoocultamiento, Autoocultamiento, Se produce cuando el vector normal a la superficie forma un ángulo superior a los 90° con el vector solar, como sería el caso, por ejemplo, de una ladera orientada al Norte, con pendiente de 45°, cuando el sol ilumina desde el Sur, elevado solamente 30° sobre el horizonte. Ocultamiento por el relieve circundante, circundante, Se produce cuando la topografía interrumpe interrumpe la linea visual desde el sol hasta el punto analizado. (Felicísimo, 1994 ) 4.2.1 Definición de exposición solar Se define la exposición solar en un punto como el tiempo máximo que ese lugar puede estar sometido a la radiación solar directa en ausencia de nubosidad. (Felicísimo, (Felicísimo, 1994 ) En la mayoría de los estudios se asume este concepto como el tiempo teórico de incidencia solar directa sobre la superficie terrestre. Para crear un mapa de exposición solar, se deben usar diferentes intervalos de tiempo, empezando en el amanecer (6:00 am) y terminando en el ocaso (6:00 pm). En caso de usar intervalos de media hora, se puede generar un mapa de exposición solar donde la máxima exposición es de 12 horas con precisión de media hora.
38
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
4.2.2 Ángulo de incidencia
WINDOWS
Tanto para el cálculo de la exposición como para la radiación solar se requiere usar información acerca el terreno. El ángulo de incidencia es el ángulo entre el vector normal al terreno ( ) y el vector solar ( ). Ver formula abajo. cos( ϕ ) =
n ⋅s n
⋅s
Figura 34: ángulo de incidencia ( )
La determinación del vector normal al terreno ( ) se hace a partir del modelo digital de elevación. (Ver figura 35 y formula abajo):
∆ X 0 − ∆ Z x ∆ Y n = 0 ⋅ ∆ Y = − ∆ Z y ∆ X ∆ Z x ∆ Z y ∆ X ∆ Y Figura 35: vector normal al terreno ( )
La determinación del vector solar ( ) se hace tomando en cuenta el ángulo de la elevación del sol ( ) y el azimut del sol ( ). (Ver formula abajo y la figura 36):
s
c o s β s e n α = c o s β c o s β s e n β
Figura 36: vector solar ( )
4.2.3 Irradiación De esta manera se puede calcular la irradiación potencial (Meijerink (Meijerink et al., 1994): 1994): Ns a l id a
=
Ne n tr a d a 39
⋅ r⋅ c o s( s ( ϕ ) Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
En los anteriores gráficos y formulas se utilizaron las siguientes convenciones: N r
flujo de fotón reflectancia ángulo de incidencia vector normal al terreno vector solar ángulo de azimut del sol desde el Norte ángulo de elevación del sol tomado tomado desde la horizontal
Felicísimo (1994), utiliza la siguiente formula para determinar la irradiación (energía incidente) sobre una superficie perpendicular a la radiación solar I n : In
=
Io
⋅
pm
en donde: irradiación incidente I n constante solar (1367 (1367 W m-2) I o p la transmisividad media cenital m la masa óptica de aire En circunstancias reales, el valor teórico de I n se ve modificado por la pendiente y la orientación del terreno, lo cual obliga a considerar los ángulos de incidencia del vector solar sobre el mismo. La extensión del cálculo para periodos más largos obliga a considerar el cambio de valor de los parámetros que intervienen en la expresión anterior (ej.: el ángulo horario mínimo H1 y máximo H2). (Felicísimo, ( Felicísimo, 1994 ) H 2
In
=
Io
∫ p
m
⋅ c o s n∠ s ⋅ d H
H 1
4.3 Ejercicio práctico
Por la complejidad de los cálculos, cál culos, la exigencia exige ncia de variables varia bles involucradas involucradas para determinar la irradiación, restringiremos el ejercicio práctico a la generación de un mapa de exposición solar. 4.3.1 Operaciones de conectividad En cierta forma la exposición solar puede ser comparada con funciones funcion es de intervisibilidad intervisibi lidad,, tomando como punto de observación el sol. Las funciones funci ones de intervisibilidad intervisibilidad hacen parte de las operaciones de conectividad. Las operaciones de conectividad conecti vidad son funciones que acumulan valores sobre el área que se atraviesa. Se requiere que uno o más atributos sean evaluados y un recorrido acumulado esté incluido en el próximo paso. El recorrido acumulado puede ser cuantitativo, como la distancia acumulada viajando, o el tiempo acumulado viajando. 40
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
También puede ser cualitativo, como por ejemplo si un punto todavía es visible. (Bakker ( Bakker 2000 ) 4.3.2 Intervisibilidad Intervisibilidad Intervisibilidad ("viewshed modelling") es una operación acumulativa. Un "viewshed" es el área que se puede ver desde uno o más puntos de origen. Funciones de intervisibilidad intervi sibilidad son aplicadas para determinar la ubicación de torres de transmisión (aplicaciones en comunicaciones), determinar el impacto visual de la construcción cons trucción de torres de transmisión eléctrica (aplicaciones (apli caciones ambientales), determinar el área que el enemigo puede ver desde sus puntos de observación (aplicaciones (aplicaci ones militares). Las funciones de intervisibilidad intervisibilidad usan los modelos digitales digital es de elevación para pa ra definir la topografía de los alrededores. También se debería incluir elementos individuales (edificios, vegetación, etc). (Aronoff (Aronoff 1989 ) 4.3.3 Conectividad en ILWIS ILWIS soporta el procesamiento de pixeles en relación con su vecindad directa en un ambiente de 3x3 pixeles. Dado que las operaciones de conectividad exigen ambientes más grandes, ILWIS ofrece estos cálculos con iteraciones. Iteraciones son cálculos repetitivos donde el resultado de un paso es el mapa de entrada para el siguiente paso. El ambiente de 3x3 pixeles se maneja de la siguiente manera: 1
2
3
4
5
6
7
8
9
El pixel central (indicado con el número 5), obtiene el resultado de la operación de conectividad. En el cálculo de la exposición solar se debe tener en cuenta no solo el ángulo de elevación del sol, sino además su ángulo de azimut, puesto que éste ángulo define cuales pixeles pi xeles vecinos potencialmente potencia lmente pueden generar sombra sombra al pixel central (ver figura abajo):
Figura 37: Influencia de pixeles vecinos
41
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
La exposición solar del pixel central está directamente influenciada por los pixeles 4 y 7. En el ejemplo anterior (figura 37), el ángulo de azimut solar ( ) es de -60°. Este ángulo se determina a partir del Sur según se muestra en la siguiente figura figu ra (Meijerink (Meijerink et al.1994 ): ):
Figura 38: ángulos en función de la trayectoria solar
Como fue expresado anteriormente, no solamente se debe verificar los pixeles vecinos del pixel central, sino además los vecinos de sus vecinos y más alla. La siguiente figura muestra este efecto: El pixel p ixel indicado con h4 causa a los pixeles h5 y h6 sombra con el ángulo de elevación el evación . Puesto que esto no puede ser evaluado en un ambiente ambiente de 3x3 se asigna aall pixel h5 una nueva altura (h4 (h4--dz 45 ). En el siguiente paso esta nueva altura causa que el pixel h6 será identificado como sombra. Este proceso se debe repetir hasta que ninguna altura sea modificada. Luego se puede comparar comparar el nuevo Modelo Digital de Elevación (MDE) con el MDE modificado y donde se encuentre diferencia se asigna “ocultamiento ” y al resto “incidencia “ incidencia ”. ”. Figura 39: propagación de sombra
También También es posible posibl e asignar un valor val or expresado en horas según el intervalo de horario (ej: (ej: intervalo de ½hora se asigna el valor 0.5), y así sumar los mapas correspondientes a diferentes horas del día, obteniendo el mapa numérico de exposición solar. 42
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
La siguiente tabla relacione cuales pixeles tienen influencia sobre el pixel central dependiendo de diferentes ángulos de azimut solar y como calcular el valor p: p
m
n
0
< 45
| ta n ( ) |
9
8
45
< 90
| tan(90 - ) |
9
6
90
< 1 35
| tan(90 - ) |
3
6
135
< 1 80
| ta n ( ) |
3
2
180
< 2 25
| ta n ( ) |
1
2
225
< 2 70
| tan(90 - ) |
1
4
270
< 3 15
| tan(90 - ) |
7
4
315
< 3 60
| ta n ( ) |
7
8
Estos valores se usan para evaluar si el pixel central recibe luz solar directa o está ocultado por medio de la siguiente expresión: M D E5
en donde: MDE 5 p MDE m MDE n dx
< p ⋅ M D E m + (1 −
p) ⋅ M D E n
−
1+ p
2
⋅ d x ⋅ t a n α
La altura del pixel central La fracción según el ángulo de azimut solar La altura del pixel vecino m La altura del pixel vecino n El tamaño del pixel El ángulo de elevación solar
En caso que el pixel central obtenga luz directa, se mantiene el valor del pixel, en caso contrario se debe asignar el nuevo valor de pixel. Puesto que el proceso iterativo es un proceso demorado, se recomienda dejar el proceso por la noche, o usar una reducción de la precisión original. 4.3.4 Iniciar el proceso En este paso se reconocerá el modelo digital d e elevación con que se realizará el cálculo. cálcul o. Se trabajará con un MDE de gran parte de Cundinamarca (jurisdicción (jurisdicc ión de la l a CAR, Bogotá) Bogotá) remuestreado a un tamaño de pixel de 1000 metros (conjunto original tiene 30 metros/pixel). 43
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Cambie Cambie al sub-directorio “exposol ” del directorio “C:\ILWIS22\DATA\ “C:\ILWIS22\DATA\”. ”. Haga doble clic sobre el icono del mapa raster “mde “ mde ”. ”. En las opciones del despliegue, asegurase que la representación “fincur “fincur ” está seleccionada y haga clic en el botón OK”. “OK”. Determine los valores mínimos y máximos del MDE. Luego cierre la ventana del mapa. Como se puede observar el área es bastante montañosa montañosa y esto tiene como como efecto que habrá gran influencia de ocultamiento topográfico. Para este ejercicio se calculará la exposición solar durante el transcurso de la mañana (6:00 a.m. hasta 8:30 a.m.), cuando más se aprecian las diferencias en ocultamiento topográfico. Observe la siguiente tabla con los valores de los ángulos de elevación ( ) y azimut ( ) solar: Hora
elevación ( )
azimut ( )
p
m
n
6:00 a.m.
0 .4 °
9 0 .0 °
0 .0 0 0 0
3
6
6:30 a.m.
4 .1 °
8 3 .9 °
0 .1 0 6 9
9
6
7:00 a.m.
1 4 .3 °
74.2°
0 .2 8 3 0
9
6
7:30 a.m.
2 5 .0 °
64.5°
0 .4 7 7 0
9
6
8:00 a.m.
3 5 .8 °
54.3°
0 .7 1 8 6
9
6
8:30 a.m.
4 6 .5 °
42.2°
0 .9 0 6 7
9
8
La formula que se debe ingresar al proceso de iteraciones es la siguiente: max(mde,p*mde#[m]+(1-p)*mde#[n]-sqrt(1+p^2)*tpix*tan( )) donde se determina el máximo valor entre el MDE y los valores de los vecinos ajustados por el ángulo de elevación elevaci ón solar. A las 6:00 de la mañana (según la tabla en la página 43) p es igual a |tan(90-90.0)| 0.0000, m=3 y n=6. El tamaño del pixel (tpix) es de 1000 metros. Con estos valores se puede simplificar la formula de la siguiente manera: max(mde,mde#[6]-9.8732)
44
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Haga doble clic sobre el icono del programa “Iteration “ Iteration”” en el listado de Map” y escriba la formula operaciones. Seleccione el mapa “mde “mde ” como “Start “ Start Map” Criterium” está arriba en la caja de “Expression “Expression”. ”. Asegurese de que para el “Stop “Stop Criterium” changes” y desactive la opción “Propagation seleccionada la opción “Until “Until No changes” “Propagation”. ”. Escriba el nombre “mde60am “mde60am ” como mapa mapa de salida, seleccione se leccione el domino “value “value ”, ”, active la opción “Show “Show”” y oprima el botón “OK “ OK”. ”. El proceso es demorado y puede realizar al rededor de de 200 iteraciones, en caso de encontrar un elemento elemento ubicado en el oriente orie nte del área de estudio que cause un sombra (ocultamiento) a elementos en el extremo occidente. En este proceso iterativo se recalcula la formula infinitamente, por ello el usuario debe observar el número de cambios del paso anterior a nterior y detectar detecta r cuando este valor va lor sea constante cons tante (por unos 5 pasos). pasos). Cuando Cuando esto suceda se puede oprimir 1 vez el botón stop con el fin de que el sistema almacene el resultado (“Store Result”). Es importante notar que no se debe usar el dominio (“fincur (“fincur ”) ”) que corresponde al MDE. Aunque los valores que resultan estarán dentro el rango de este dominio, ILWIS internamente trabaja con valores de mayor precisión de la que este dominio puede manejar y usarlo daría como resultado un mapa con todos los valores indefinidos! Al terminar el proceso de iteración acepte las definiciones por defecto para desplegar el mapa. Despliegue el resultado al lado el mapa original “mde “mde ” con representación “pseudo “pseudo ”. ”. Explique las diferencias. Por medio de la siguiente formula se puede crear un mapa de exposición solar para las 6:00 de la mañana: sol60am:=iff(mde. En la caja del diálogo oprima . di álogo cree un u n nuevo dominio con nom nombre bre “sol “sol ”, ”, al cual NO se adiciona ningún elemento. Al oprimir el botón “OK “ OK”” en la caja de Definition” el sistema preguntará si desea adicionar los diálogo “Raster “Raster Map Definition” elementos “ocultamiento “ocultamiento ” y “incidencia “incidencia ” al dominio. En ambas casos elija “Si “ Si”. ”. Para crear un mapa que pueda ser usado de manera más eficiente en la generación del mapa de exposición solar para todo el día, se debe asignar el valor (en horas) correspondiente al intervalo que se utilizará: exp60am:=iff(mde
45
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Repita este proceso para los otros momentos de la mañana. Para el mapa “mde65am “mde65am ” (6:30 a.m.): max(mde,0.1069*mde#[9]+(0.8931)*mde#[6]-72.0893) Para el mapa “mde70am” “mde70am” (7:00 (7:00 a.m.): max(mde,0.2830*mde#[9]+(0.7170)*mde#[6]-264.9074) Para el mapa “mde75am “mde75am ” (7:30 a.m.): max(mde,0.4770*mde#[9]+(0.5230)*mde#[6]-516.6405) Para el mapa “mde80am” “mde80am” (8:00 (8:00 a.m.): max(mde,0.7186*mde#[9]+(0.2814)*mde#[6]-888.1254) Para el mapa “mde85am “mde85am ” (8:30 a.m.): max(mde,0.9067*mde#[9]+(0.0933)*mde#[8]-1422.4490) Cree también los mapas de exposición expresados en horas y nombrelos “exp65am “exp65am ”, ”, “exp70am ”, ”, “exp75am “exp75am ”, ”, “exp80am “exp80am ” y “exp85am “exp85am ”. ”. Despliegue los 6 mapas de exposición solar simultáneamente en la pantalla y compare los resultados. En que momento la diferencia entre dos cálculos es más grande? La suma de estos mapas generará un mapa de exposición solar para un intervalo de 3 horas. El máximo valor es 3 (correspondiente a 3 horas) y el mínimo es de 0 horas: exp60a85 := exp60am+exp65am+exp70am+exp75am+exp80am+exp85am Escribe esta formula formula en la linea de comando comando de la ventana principal de ILWIS y . Asegurese de que el rango de valores es de 0.0 a 3.0 con oprima . precisión de 0.1 y oprima el botón “OK “ OK”. ”. Despliegue el mapa usando la representación “Gray “Gray ” e interprete el resultado. Para su su mejor mejor interpretación se puede clasificar clasi ficar el resultado, asignar asig nar colores definido definid o por el usuario y desplegarlo en una vista tridimensional (ver figura 40):
46
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
La figura 40 muestra el resultado de esta operación, los colores más claros representan las zonas de mayor exposición solar. La posición relativa del sol es el la esquina superior derecha de la ventana. Figura 40: Exposición solar en 3D Los mapas “exp75am “exp75am ”, ”, “exp80am ” y “exp85am ” revelaron que casi ca si todo el terreno obtiene sol a partir de las la s 8:00 de la mañana. Al observar el terreno terren o (ver figura 41) se puede notar que el mayor ocultamiento topográfico se presenta por la mañana cuando el sol está ubicado al Oriente del área de estudio.
Figura 41: Corte transversal del terreno, trayectoria solar
De esta esta maner maneraa se puede asumir que todo el área áre a obtiene sol so l entre las 8:00 de la mañana y las 4:00 de la tarde. Se pueden hacer cálculos para exclusivamente este intervalo y sumarlos con el resultado “exp60a85 “exp60a85 " adicionando 7 horas para complementar el cálculo para todo el intervalo. 47
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
Hora
elevación ( )
azimut ( )
WINDOWS
p
m
n
4:00 p.m.
3 8 .6 °
3 0 5 .7 °
0.7186
7
4
4:30 p.m.
2 6 .2 °
2 9 5 .5 °
0.4770
7
4
5:00 p.m.
1 4 .7 °
2 8 5 .8 °
0.2830
7
4
5:30 p.m.
4 .3 °
276.1°
0 .1 0 6 9
7
4
6:00 p.m.
0 .5 °
270.0°
0 .0 0 0 0
7
4
Para el mapa “mde40pm “mde40pm ” (4:00 p.m.): max(mde,0.7186*mde#[7]+(0.2814)*mde#[4]-983.0267) Para el mapa “mde45pm “mde45pm ” (4:30 p.m.): max(mde,0.4770*mde#[7]+(0.5230)*mde#[4]-545.1736) Para el mapa “mde50pm “mde50pm ” (5:00 p.m.): max(mde,0.2830*mde#[7]+(0.7170)*mde#[4]-272.6482) Para el mapa “mde55pm “mde55pm ” (5:30 p.m.): max(mde,0.1069*mde#[7]+(0.8931)*mde#[4]-75.6188) Para el mapa “mde60pm “mde60pm ” (6:00 p.m.): max(mde,mde#[4]-8.7269) Cree los anteriores mapas. Después cree también los mapas de exposición solar para cada momento del día (“exp40pm (“exp40pm ”, ”, “exp45pm ”, ”, “exp50pm ”, ”, “exp55pm ”, ”, “exp60pm ”). ”). Sume los mapas para crear el mapa de exposición solar para todo el día: expdia := exp60a85+exp40pm+exp45pm+exp50pm+exp55pm+exp60pm+7 Puesto que estamos asignando ½hora a cada momento el intervalo real es de las 5:45 a.m. a.m. hasta las 6:15 p.m.. De esta manera el máximo valor de exposición expos ición solar sol ar será de 12½ horas. Observe las siguientes figuras (42 y 43) donde se muestra el resultado de esta operación.
48
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Los colores claros corresponden a las zonas con mayor incidencia solar, mientras mientras las más oscuras corresponden a las zonas de mayor ocultamiento ocultamiento topográfico.
Figura 42: Exposición solar en Cundinamarca
Las lineas en la parte cerca de los límites de Cundinamarca son debido a errores generados durante la interpolación de las curvas de nivel, por falta de información.
Figura 43: Vista tridimensional de la exposición solar en Cundinamarca
Se puede notar que los picos, crestas y planicies tienen ti enen mayor exposición solar, mientras las laderas, pozos y canales reciben menor cantidad de horas de sol. -- . -49
Exposición solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
5. RESPUESTAS En este capitulo se consignan las respuestas a los ejercicios expuestos 3.1.1 Despliegue y consulta de la información La distribución de los puntos en el mapa “C14" es: Aleatoria Medición de las distancias en los puntos: Distancia mínima: ± 150 metros Distancia máxima: ± 85000 metros Que tipo de dominio tiene el mapa (ver propiedades del mapa)? Dominio identificador Cuantos puntos tiene el mapa (ver propiedades del mapa)? 60 puntos Cuantos pares de puntos se pueden crear con estos puntos? n*(n-1)/2 = 60*59/2=1770 Cuales son los valores mínimos y máximos de la columna “C14 “ C14 "? "? Valor mínima: 0.35 Valor máxima: 81.00 La varianza de la columna “C14 “C14 " es: 581.901 3.1.2 Análisis de distribución espacial Con que distribución corresponde el comportamiento del gráfico? Aleatoria Este corresponde con lo que usted visualmente había concluido? Si 3.2.1 Variograma representado como superficie (Anisotropía) El ángulo principal de anisotropía anisotropía es: ± 20° 50
Respuestas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
El ángulo de tolerancia es: ± 15° 3.3.1 Adicionar modelos Los modelos que se ajustan al comportamiento perpendicular son: Modelo 1: Spherical model Modelo2: Nugget: 50 Nugget: Sill: 710 Sill: Range: 20000 Range:
del semivarianza de la dirección Circular model 50 690 19000
3.3.4 Adicionar modelos Los modelos que se ajustan al comportamiento del semivarianza omnidireccional son:
Modelo 1: Nugget: Slope: Power:
Power model 250 0.0001 1.5
Modelo 2: Nugget: Sill: Range:
Gaussian model 250 750 22000
Modelo 3: Nugget: Sill: Range:
Rational Quadratic 250 800 17500
Modelo 4: Nugget: Sill: Range:
Wave model 250 650 7500
51
Respuestas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS 2
3.3.5 Verificación y selección de los modelos (Goodness of Fit R ) Intervalo de distancia donde se ajusten los modelos: Modelo 1: Modelo 2: Power model Gaussian model Distancia mínima: 0 Distancia mínima: 0 Distancia máxima: 25km Distancia máxima: 35km Modelo 3: Rational Quadratic Modelo 4: Wave model Distancia mínima: 0 Distancia mínima: 0 Distancia máxima: 40km Distancia máxima: >45km Cual modelo, según su criterio, se ajusta mejor a los valores experimentales? Wave model (0-45km) Anote el R 2 para los 4 modelos para las 4 distancias (25, 35, 40 y 45 km): Distancia
Modelo 1 (pow) Modelo 2 (gau)
Modelo 3 (rat)
Modelo 4 (wav)
0-25 km
0.5637
0.5318
0.502
0 .5 5 4 9
0-35 km
0 .5 7 4 7
0 .7 2 8 2
0 .7 1 3 1
0.7436
0-40 km
0 .1 3 8 8
0 .7 4 6 5
0 .7 3 4 8
0.749
0-45 km
-0.7671
0 .7 4
0 .7 3 6 2
0.7647
0 .1 2 7 5
0 .6 8 6 6
0 .6 7 1 5
0.7031
Promedio
Cual de los cuatro modelos se ajusta mejor al semivariograma experimental? Wave model (tanto valor más alto y promedio más alto) 3.4 Interpolación por medio de Kriging
Los errores mínimos y máximos leídos en el mapa “OKwa50e “ OKwa50e ” son: Error mínimo: 16.704 Error máximo: 2 4 .1 7 6 Cerca de: pn t 3 8 Cerca de: pn t 6 0 a distancia: ½km a distancia: 9k m
52
Respuestas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Bibliografía Aronoff, S. 1989. Geographical Information Systems: A Management Perspective. WDL Publications Otawa. Canada. 294p. Bakker, X. and F. Pérez , 2000, Conceptos de análisis y modelamiento en el e l formato raster. Notas de clase. Especialización en Sistemas de Información Geográfica. CIAF Universidad de Tunja, Santa Fe de Bogotá, 57 p. Burrough, P.A. 1985. Principles Principle s of Geographical Information Systems for Land Resource Assessment. Clarendon Press, Oxford. 193p. Felicísimo, A.M. 1994. Modelos Digitales del Terreno. Introducción y aplicaciones en las ciencias ambientales. España. 118 p. ILWIS 1997, ILWIS 1997, ILWIS 2.1, User’s Guide. ITC, The Netherlands. ILWIS 1997, ILWIS 1997, ILWIS 2.1, Application Guide. ITC, The Netherlands. ILWIS 1998, ILWIS 1998, ILWIS 2.2, Installation & New Funcionality. ITC, The Netherlands. ILWIS 1999, ILWIS 1999, Ayuda en línea de ILWIS 2.23. ITC The Netherlands. Meijerink, A.M.J., H.A.M. de Brouwer, C.M. Mannaerts y C.R. Valenzuela , 1994. Introductionto the use of Geographic Information Information Systems for practical Hydrology. Hydrol ogy. ITC Publication number 23, Enschede, The Netherlands. 243p. Valenzuela, C, 1989, Análisis y modelamiento de datos, notas de clase.10p.
53
Respuestas
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
WINDOWS
Anexo: Declinación solar D F ec h a 1-Ene 2-Ene 3-Ene 4-Ene 5-Ene 6-Ene 7-Ene 8-Ene 9-Ene 10-Ene 11-Ene 12-Ene 13-Ene 14-Ene 15-Ene 16-Ene 17-Ene 18-Ene 19-Ene 20-Ene 21-Ene 22-Ene 23-Ene 24-Ene 25-Ene 26-Ene 27-Ene 28-Ene 29-Ene 30-Ene 31-Ene 1-Feb 2-Feb 3-Feb 4-Feb 5-Feb 6-Feb 7-Feb 8-Feb
d 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39
D -23.07 -22.99 -22.90 -22.81 -22.70 -22.60 -22.48 -22.36 -22.23 -22.10 -21.96 -21.81 -21.66 -21.50 -21.33 -21.16 -20.98 -20.79 -20.60 -20.40 -20.20 -19.99 -19.78 -19.55 -19.33 -19.09 -18.85 -18.61 -18.36 -18.11 -17.85 -17.58 -17.31 -17.03 -16.75 -16.47 -16.17 -15.88 -15.58
= 2 3.5 ⋅ s e n [0 .9 8 6 ⋅ ( 2 8 4 + d ) ] 9-Feb 10-Feb 11-Feb 12-Feb 13-Feb 14-Feb 15-Feb 16-Feb 17-Feb 18-Feb 19-Feb 20-Feb 21-Feb 22-Feb 23-Feb 24-Feb 25-Feb 26-Feb 27-Feb 28-Feb 1-Mar 2-Mar 3-Mar 4-Mar 5-Mar 6-Mar 7-Mar 8-Mar 9-Mar 10-Mar 11-Mar 12-Mar 13-Mar 14-Mar 15-Mar 16-Mar 17-Mar 18-Mar 19-Mar 20-Mar
40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 68 69 70 71 72 73 74 75 76 77 78 79 54
-15.27 -14.96 -14.65 -14.33 -14.01 -13.68 -13.35 -13.02 -12.68 -12.34 -11.99 -11.64 -11.29 -10.93 -10.57 -10.21 -9.84 -9.47 -9.10 -8.73 -8.35 -7.97 -7.59 -7.21 -6.82 -6.43 -6.04 -5.65 -5.26 -4.86 -4.47 -4.07 -3.67 -3.27 -2.87 -2.47 -2.06 -1.66 -1.26 -0.85
21-Mar 22-Mar 23-Mar 24-Mar 25-Mar 26-Mar 27-Mar 28-Mar 29-Mar 30-Mar 31-Mar 1-Abr 2-Abr 3-Abr 4-Abr 5-Abr 6-Abr 7-Abr 8-Abr 9-Abr 10-Abr 11-Abr 12-Abr 13-Abr 14-Abr 15-Abr 16-Abr 17-Abr 18-Abr 19-Abr 20-Abr 21-Abr 22-Abr 23-Abr 24-Abr 25-Abr 26-Abr 27-Abr 28-Abr 29-Abr
80 81 82 83 84 85 86 87 88 89 90 91 92 93 94 95 96 97 98 99 100 101 102 103 104 105 106 107 108 109 110 111 112 113 114 115 116 117 118 119
-0.45 -0.05 0.36 0.76 1.17 1.57 1.97 2.38 2.78 3.18 3.58 3.98 4.38 4.77 5.17 5.56 5.96 6.35 6.73 7.12 7.50 7.89 8.27 8.64 9.02 9.39 9.76 10.13 10.49 10.85 11.21 11.56 11.91 12.26 12.60 12.94 13.28 13.61 13.94 14.26
Anexo - Declinación solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
F ec h a 30-Abr 1-May 2-May 3-May 4-May 5-May 6-May 7-May 8-May 9-May 10-May 11-May 12-May 13-May 14-May 15-May 16-May 16-May 17-May 18-May 19-May 20-May 21-May 22-May 23-May 24-May 25-May 26-May 27-May 28-May 29-May 30-May 31-May 1-Jun 2-Jun 3-Jun 4-Jun 5-Jun 6-Jun 7-Jun 8-Jun 9-Jun 10-Jun
d 120 121 122 123 124 125 126 127 128 129 130 131 132 133 134 135 136 137 138 139 140 141 142 143 144 145 146 147 148 149 150 151 152 153 154 155 156 157 158 159 160 161
D 14.58 14.89 15.20 15.51 15.81 16.11 16.40 16.69 16.97 17.25 17.52 17.79 18.05 18.30 18.56 18.80 19.04 19.27 19.50 19.73 19.94 20.15 20.36 20.56 20.75 20.94 21.12 21.29 21.46 21.62 21.78 21.93 22.07 22.20 22.33 22.46 22.57 22.68 22.78 22.88 22.97 23.05
11-Jun 12-Jun 13-Jun 14-Jun 15-Jun 16-Jun 17-Jun 18-Jun 19-Jun 20-Jun 21-Jun 22-Jun 23-Jun 24-Jun 25-Jun 26-Jun 27-Jun 28-Jun 29-Jun 30-Jun 1-Jul 2-Jul 3-Jul 4-Jul 5-Jul 6-Jul 7-Jul 8-Jul 9-Jul 10-Jul 11-Jul 12-Jul 13-Jul 14-Jul 15-Jul 16-Jul 17-Jul 18-Jul 19-Jul 20-Jul 21-Jul 22-Jul
162 163 164 165 166 167 168 169 170 171 172 173 174 175 176 177 178 179 180 181 182 183 184 185 186 187 188 189 190 191 192 193 194 195 196 197 198 199 200 201 202 203 55
23.13 23.19 23.26 23.31 23.36 23.40 23.43 23.46 23.48 23.49 23.50 23.50 23.49 23.48 23.45 23.43 23.39 23.35 23.30 23.24 23.18 23.11 23.03 22.95 22.86 22.76 22.66 22.55 22.43 22.31 22.17 22.04 21.89 21.74 21.59 21.42 21.25 21.08 20.90 20.71 20.51 20.31
WINDOWS
23-Jul 24-Jul 25-Jul 26-Jul 27-Jul 28-Jul 29-Jul 30-Jul 31-Jul 1-Ago 2-Ago 3-Ago 4-Ago 5-Ago 6-Ago 7-Ago 8-Ag 8-Ago 9-Ago 10-Ago 11-Ago 12-Ago 13-Ago 14-Ago 15-Ago 16-Ago 17-Ago 18-Ago 19-Ago 20-Ago 21-Ago 22-Ago 23-Ago 24-Ago 25-Ago 26-Ago 27-Ago 28-Ago 29-Ago 30-Ago 31-Ago 1-Sep 2-Sep
204 205 206 207 208 209 210 211 212 213 214 215 216 217 218 219 220 221 222 223 224 225 226 227 228 229 230 231 232 233 234 235 236 237 238 239 240 241 242 243 244 245
20.11 19.90 19.68 19.45 19.22 18.99 18.75 18.50 18.25 17.99 17.73 17.46 17.19 16.91 16.63 16.34 16.04 15.75 15.44 15.14 14.82 14.51 14.19 13.86 13.54 13.20 12.87 12.53 12.18 11.83 11.48 11.13 10.77 10.41 10.05 9.68 9.31 8.94 8.56 8.18 7.80 7.42
Anexo - Declinación solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
F ec h a 3-Sep 4-Sep 5-Sep 6-Sep 7-Sep 8-Sep 9-Sep 10-Sep 11-Sep 12-Sep 13-Sep 14-Sep 15-Sep 16-Sep 17-Sep 18-Sep 19-Sep 20-Sep 21-Sep 22-Sep 23-Sep 24-Sep 25-Sep 26-Sep 27-Sep 28-Sep 29-Sep 30-Sep 1-Oct 2-Oct 3-Oct 4-Oct 5-Oct 6-Oct 7-Oct 8-Oct 9-Oct 10-Oct 11-Oct 12-Oct 13-Oct 14-Oct
d 246 247 248 249 250 251 252 253 254 255 256 257 258 259 260 261 262 263 264 265 266 267 268 269 270 271 272 273 274 275 276 277 278 279 280 281 282 283 284 285 286 287
D 7.04 6.65 6.26 5.87 5.48 5.08 4.69 4.29 3.89 3.49 3.09 2.69 2.29 1.89 1.48 1.08 0.67 0.27 -0.13 -0.54 -0.94 -1.35 -1.75 -2.15 -2.56 -2.96 -3.36 -3.76 -4.16 -4.55 -4.95 -5.34 -5.74 -6.13 -6.52 -6.91 -7.29 -7.67 -8.06 -8.43 -8.81 -9.18
15-Oct 16-Oct 17-Oct 18-Oct 19-Oct 20-Oct 21-Oct 22-Oct 23-Oct 24-Oct 25-Oct 26-Oct 27-Oct 28-Oct 29-Oct 30-Oct 31-Oct 1-Nov 2-Nov 3-Nov 4-Nov 5-Nov 6-Nov 7-Nov 8-Nov 9-Nov 10-Nov 11-Nov 12-Nov 13-Nov 14-Nov 15-Nov 16-Nov 17-Nov 18-Nov 19-Nov 20-Nov 21-Nov 22-Nov 23-Nov 24-Nov 25-Nov
288 289 290 291 292 293 294 295 296 297 298 299 300 301 302 303 304 305 306 307 308 309 310 311 312 313 314 315 316 317 318 319 320 321 322 323 324 325 326 327 328 329 56
-9.56 -9.92 -10.29 -10.65 -11.01 -11.37 -11.72 -12.07 -12.41 -12.75 -13.09 -13.42 -13.75 -14.08 -14.40 -14.72 -15.03 -15.34 -15.65 -15.94 -16.24 -16.53 -16.81 -17.09 -17.37 -17.64 -17.90 -18.16 -18.42 -18.67 -18.91 -19.15 -19.38 -19.60 -19.82 -20.04 -20.25 -20.45 -20.64 -20.83 -21.02 -21.20
WINDOWS
26-Nov 27-Nov 28-Nov 29-Nov 30-Nov 1-Dic 2-Dic 3-Dic 4-Dic 5-Dic 6-Dic 7-Dic 8-Dic 9-Dic 10-Dic 11-Dic 12-Dic 13-Dic 14-Dic 15-Dic 16-Dic 17-Dic 18-Dic 19-Dic 20-Dic 21-Dic 22-Dic 23-Dic 24-Dic 25-Dic 26-Dic 27-Dic 28-Dic 29-Dic 30-Dic 31-Dic
330 331 332 333 334 335 336 337 338 339 340 341 342 343 344 345 346 347 348 349 350 351 352 353 354 355 356 357 358 359 360 361 362 363 364 365
-21.37 -21.53 -21.69 -21.84 -21.99 -22.13 -22.26 -22.39 -22.51 -22.62 -22.73 -22.83 -22.92 -23.01 -23.08 -23.16 -23.22 -23.28 -23.33 -23.38 -23.41 -23.45 -23.47 -23.49 -23.50 -23.50 -23.50 -23.49 -23.47 -23.44 -23.41 -23.37 -23.33 -23.27 -23.22 -23.15
Anexo - Declinación solar
T HE I NT E GR A T E D L A N D A ND W A T ER I NF O RM A T I O N S YS T EM
for
INTERPRETACION DIGITAL Y SENSORES REMOTOS
57
WINDOWS
Anexo - Declinación solar