Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
Rectificacion Introduccion Las imágenes registradas por sensores montados en satélites son una representación de la superficie terrestre y presentan distorsiones. Entre las fuentes de distorsión geométricas pueden contarse las siguientes: • Rotación de la Tierra durante el tiempo de adquisición de los datos • Velocidad finita de barrido de algunos sensores • Ancho campo visual de algunos sensores • Curvatura de la Tierra • No linealidades de los sensores • Variaciones en altura, actitud y velocidad de la plataforma Siempre que se necesite procesar conjuntamente varias imágenes de la misma región geográfica, estas deben poder superponerse geométricamente previo al procesamiento. Por lo tanto en estos casos las correciones geométricas también pueden considerarse como una etapa de preprocesamiento. La rectificación es uno de los dos métodos que existen para corregir estos problemas y es un proceso mediante el cual la geometría de una imagen se hace planimétrica, es decir se la corrige con respecto a un mapa, lo que implica la selección de una dada proyección geográfica. En este proceso, la grilla de datos originales es proyectada a una nueva grilla y en ese caso es necesario hacer un remuestreo de los valores, lo que consiste en general en un proceso de interpolación para conocer los nuevos valores radiométricos de los pixeles a partir de los correspondientes a la grilla original (los valores radiométricos son aquellos que tiene que ver con el brillo de cada pixel y por lo tanto de toda la imagen). Þ Rectificacion es el proceso de transformar los datos desde un sistema de grilla a otro usando un orden polynomial n y un sistema de proyeccion de mapas dado. Esto ultimo, es cualquier sistema designado para representar la superficie de una esfera o esferoide -ej: la Tierra- en un plano. Debido a que los pixeles de la nueva grilla pueden no alinearse con los de la grilla original, los pixeles deben ser remuestreados. La rectificacion es necesaria en casos en donde la grilla de pixeles de una imagen debe ser cambiada para ajustarse a un sistema de proyeccion de mapas o a una imagen referencia. Usualmente, rectificacion es la conversion de las coordenadas del archivo de datos a alguna otra grilla y sistema de coordenadas, llamado “sistema de referencia”. Þ Resampling (remuestreo) es el proceso de extrapolar los valores de los datos desde la grilla original a la nueva grilla Þ Registracion es el proceso de hacer que una imagen conforme con otra imagen. En estos casos, un sistema de coordenada de mapas no esta necesariamente involucrado. Este proceso se aplica en casos en los cuales se dispone de imagenes de una misma area que fueron colectadas desde diferentes fuentes pero Página 1 de 10e
Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
que deben ser utilizadas juntas. Para poder comparar las diferentes imagenes pixel por pixel, la grilla de pixeles de cada imagen debe coincidir con las otras imagenes en la base de datos. Þ Georeferenciacion se refiere al proceso de asignar coordenadas de mapa a los datos de la imagen. La rectificacion, por definicion, involucra la georeferenciacion ya que todos los sistemas de proyeccion de mapas estan asociados con coordenadas de mapas. No obstante, la registracion imagen a imagen involucra georeferenciacion solamente si la imagen referencia esta ya georeferenciada. La georeferenciacion, por si misma, involucra cambiar solamente la informacion de coordenadas de mapa en el archivo imagen. La grilla de la imagen no cambia Þ Latitud/ Longitud es un sistema de coordenadas esferico que no esta asociado con una proyeccion de mapa. Lat./Long expresa locaciones en terminos de un esferoide, no de un plano. Aunque una imagen no es usualmente “rectificada” a Lat./Long, es posible convertir las imagenes a este sistema de coordenadas. Rectificar o Registrar una imagen incluye los siguientes pasos a)-Localizar puntos de control. Los GCP son pixeles que representan la misma posicion en la imagen a rectificar y en la imagen referencia (en la cual las coordenadas de esos puntos son conocidas). A partir de estos GCPs se extrapolan luego las coordenadas de rectificacion para todos los puntos restantes de la imagen. Por ello, cuanto mayor sea la precision en la seleccion de los puntos, mayor la cantidad que se tomen, y mayor la dispersion que presenten, mas fiable sera la rectificacion. El minimo numero de GCPs a tomar se obtiene con la siguiente ecuacion... (t+1)*(t+2)/2, donde t es el orden de transformacion a usar. b)-Computar y testear la matriz de transformacion. Esta matriz se crea a partir de los GCPs y consiste en coeficientes que son usados en ecuaciones polinomiales para convertir las coordenadas. El tamano de la matriz depende del orden de transformacion elegido. Lo ideal es utilizar la ecuacion polinomica para la cual exista la menor cantidad posible de error en el momento en que sea utilizada para transformar las coordenadas de los GCPs de la imagen referencia (denominada destination) en las coordenadas de la imagen a rectificar (denominada source). b1)-Evaluacion del error: Se analiza el error RMS (root mean square = raiz cuadrada media). El error para la totalidad de los GCPS existentes se obtiene automaticamente con cada computacion de la matriz de transformacion (ver debajo). El error RMS para cada punto indica la distancia entre la ubicacion del GCP en la imagen source y la ubicacion del mismo GCP en la imagen retransformada. En otras palabras, es la diferencia entre la coordenada deseada (o ideal) de salida para un GCP y la salida real obtenida cuando el punto es Página 2 de 10e
Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
retransformado con la matriz de transformacion. Luego de cada computaciond e la matriz de transformacion y calculo del error RMS hay 4 posibilidades si el error obtenido es mas alto que el esperado... * Descartar los GCPs con mas alto error RMS y volver a computar la matriz (siempre y cuando se cumpla con el numero minimo de puntos) * Tolerar un valor mas alto de error RMS * Incrementar el orden de transformacion (el problema es que esto crea alteraciones geometricas mas complejas en la imagen) * Seleccionar solo aquellos puntos para los cuales se tenga mayor confianza. b2)-Ordenes de Transformacion * Orden 1: Se aplica una ecuacion polinomica de grado 1, es decir que la transformacion es lineal. Se utiliza en rectificaciones de imagenes planas a proyecciones de mapas planas, proyecciones planas de maspas a otra proyeccion plana, rectificacion de imagenes que ocupan pequenas areas y datos que ya se encuentran proyectados dentro de un plano (ej: SPOT y Landsat1B estan ya transformadas a un plano pero pueden no estar rectificadas a una proyeccion de mapa deseada. * Orden 2: en adelante, son transformaciones no lineales, el proceso de correcion es conocido como rubber sheeting. En el caso del orden 2, puede usarse para convertir datos Lat/Long a una proyeccion plana, corregir datos que cubren un area grande, o datos distorsionados. * Orden 3: Puede usarse para corregir fotografias aereas distorsionadas e imagenes de radar El peligro de usar altos ordenes de transformacion es que cuanto mas complicada sea la ecuacion de transformacion menos regulares y predecibles seran los resultados. Por ello, es conveniente empezar con un orden 1 de transformacion a menos que se tenga certeza de que no funcionara. c)-Remuestreo de la imagen: Crear la imagen salida con la informacion de coordenadas nueva en el header. Los pixeles deben ser remuestreados para conformar la nueva grilla usando Resample con el metodo de remuestreo mas conveniente... * Nearest neighbor: Asigna a cada pixel de la imagen salida (rectificada) el valor del pixel mas cercano en la imagen source. Tiene como ventaja que transfiere los valores originales de los datos, sin promediarlos. Esto es importante cuando se discriminan tipos de vegetacion o se determinan distintos niveles de turbidez o temperaturas en el agua. Es bueno para utilizarlos antes de un proceso de clasificacion y util para imagenes tematicas (cualitativas). Las desventajas son que crea un efecto de “escalones” en los bordes cuando se remuestrea desde un tamano de grilla grande a uno chico. Ademas pueden perderse algunos valores (que no se asignan a ningun nuevo pixel) y otros duplicarse. Página 3 de 10e
Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
* Bilinear Interpolation: Considera el valor de los 4 pixeles mas cercanos en la imagen source para asignar el nuevo valor al pixel de la imagen salida. Las ventajas son que no existe el efecto de “escalones” en los bordes que puede aparece en nearest neighbor y ademas hay mejor exactitud espacial. Se usa a menudo cuando se cambia el tamano de las celdas en los datos. La desventaja es que como los pixeles son promediados tiene el efecto de una baja frecuencia de convolucion. Es decir, algunos extremos de los valores de los datos pueden perderse. * Cubic convolution: Es similar al anterior, solo que considera los 16 pixeles mas cercanos (una matriz de 4x4). La ventaja es que en la mayoria de los casos la media y varianza de los pixeles salida coincide con la media y varianza de los pixeles entrada mas que en cualquier otro metodo. Tambien puede tanto mejorar la forma de la imagen como remover ruido, aunque esto depende de los datos que se esten usando. Este metodo es recomendado cuando se esta cambiando drasticamente el tamano de las celdas en los datos (tal como una convergencia de LandsatTM y fotos aereas). La desventaja es que los datos son alterados a veces y que constituye el metodo mas lento de todos. Tamano de las celdas (pixeles) en la salida: Es el tamano que tendran los pixeles de la imagen salida, pero referido a las unidades de referencia que tenga la imagen destination, las cuales se indican en “map projection” de resample. Pueden plantearse 3 casos a groso modo... 1.Si la unidad de proyeccion de la imagen destination esta en metros, y el tamano de celda esta indicado en 5 por ejemplo, entonces el tamano del pixel en la salida sera de 5 metros por default, pero el usuario puede cambiarlo al valor -en este caso en metros- deseado (el tipo de transformacion se ubicara automaticamente en File to Map suponiendo que la imagen a remuestrear no tiene unidades de referencia previas). 2.Si no hay map projection en la imagen destination (o sea que tiene unidades de file, en pixeles), entonces el tipo de transformacion se ubicara en file to file y el tamano de las celdas en la salida no se podra modificar (aparece en 1); la imagen salida sera llevada automaticamente a la resolucion que tenga la imagen de referencia. 3.a) Si se quiere una resolucion distinta a la de la imagen referencia, hay que asignar un “map projection” trucho a la imagen referencia. Esto se hace desde el Image Info (ver Image Information), poniendo la proyeccion como “unknown” con todos los parametros por default, el tamano de los pixeles en la imagen referencia sera 1 y esto permitira modificar luego el tamano de las celdas en la imagen salida de acuerdo a lo que se necesite. Para que este cambio sea reconocido es necesario cerrar la imagen y volverla a desplegar. Ejemplo: Supongamos que se esta remuestreando una imagen Landsat TM con una resolucion de 30m tomando como referencia una imagen de Radar con resolucion de 12,5m. Página 4 de 10e
Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
Si la imagen Radar tiene un map projection “unknown”, el tamano X de las celdas en la salida de la TM remuestreada sera X*(unidad unknown de la imagen referencia)= X*1 para este ejemplo. Pero este 1 no es otra cosa que la resolucion original de la imagen referencia (12,5m para el Radar). Si quiero un tamano de celda de 30m en la salida, entonces mi X debera tener un valor tal que X*12,5 = 30, o sea que X=2,4. b) Que pasa si ya tome todos los puntos de control y despues recien me doy cuenta que deberia haber fijado el map projection de la imagen referencia en “unknown”? Si cambio el map projection en ese momento perdere todos los puntos de control (total que son 2 dias mas tomando puntitos!), entonces tengo que cambiar el map projection como se indico pero tengo que transformar mis GCPs al nuevo map projection trucho. 0
M pixeles
m-1
0 N lineas
1
M pixeles
M
N FILE (1) COORDINATES
n-1
N lineas
MAP "TRUCHO" (2) COORDINATES
1
Los diagramas indican como se ubican las coordenadas en un File y en un Map. En File las coordenadas empiezan en 0 y en Map en 1, ademas las filas se cuentan desde arriba hacia abajo para File y en Map es al reves. Entonces, los valores para cada X e Y en el pasaje de (1) a (2) se obtienen de acuerdo a las siguientes ecuaciones... Y map = N - Y file X map = 1 + X file Donde N = numero de filas (o height, altura) de la imagen destination. Estas ecuaciones son las que hay que aplicarle a todos los puntos Xdestination e Ydestination luego de haber fijado la proyeccion como “unknown”. Esto se hace seleccionando toda la columna de valores en el GCP editor y utilizando la opcion “formula”, a la cual se accede clickeando con el boton derecho del mouse mientras se esta sobre el encabezado de la columna. Una vez efectuado esto se puede calcular la matriz de transformacion y efectuar el remuestreo sin ningun problema. Rectificacion en Erdas 8.2 1)- Abrir en 2 viewers distintos la imagen a rectificar y la imagen referencia. 2)- Desde el viewer que contiene la imagen source (a rectificar), abrir GCP editor y, usando los menues del GCP editor, seleccionar bajo el modo “pairwise” la imagen referencia. “pairwise” permite registrar o rectificar una imagen en funcion de otra imagen, de una mesa digitalizadora o de coordenadas ingresadas desde el teclado... Página 5 de 10e
Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
Þ Raster Þ GCP Editor Þ Pairwise Þ Select Viewer (hacer click dentro de la imagen referencia) calcular matriz de transformacion calcular matriz automaticamente en cada paso anular select GCPs seleccionar GCPs bloquear/ desbloquear select GCPs
identificacion de cada GCP
buscar punto seleccionado en source o destination
coordenadas X-Y de cada GCP para source y destination
3)- Registrar los puntos de control (GCPs). Siempre desde los menues del GCP editor... Þ Source Þ Show Chip Extraction View Þ Destination Þ Show Chip Extraction View * Aparece para cada imagen un “mini” viewer que amplifica la parte de la imagen incluida dentro del cuadro de seleccion (este cuadro aparece sobre la imagen al seleccionar el chip extraction view wn cada imagen). Esto sirve de ayuda para ir seleccionando los puntos. Alinear ambas imagenes de manera que se visualice la misma area en ambas. Þ Hacer click sobre el icono para seleccionar un GCP en la imagen source. Luego volver a clickear sobre el icono para seleccionar el mismo GCP en la imagen destination. El primer par de GCPs estara seleccionado entonces. Þ Desplazar los cuadros de seleccion de los Chip Extraction View de ambas imagenes hasta donde se desee tomar un nuevo GCP. Repetir la operacion anterior para tomar un nuevo par de GCPs y asi sucesivamente hasta completar el total de puntos de control necesarios. 4)- Calculo de la Matriz de transformacion y error RMS Þ Hacer click sobre el icono que permite calcular la matriz, automaticamente apareceran 2 nuevas ventanas, una es el Transform Editor con el resultado de la matriz y otra con el valor de error RMS para el eje X, eje Y y total. Si el error obtenido es acorde con lo esperado pasar al paso 5). * Si el error RMS obtenido es demasiado alto proceder de acuerdo a lo indicado anteriormente (descartar puntos con alto error, aceptar un error total mayor o cambiar el orden de transformacion).
Página 6 de 10e
Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
recalcular matriz de transformacion resample
* Cambio del orden de transformacion: Desde el menu Edit Þ Transformation order, del Transform editor. Si el cambio de un orden 1 por 2 no reduce siginificativamente el error conviene remuestrear con orden 1 nomas. * Tipo de transformacion: En el transform editor Þ Edit Þ Transformation type. Puede haber 4 tipos: File to File, File to Map, Map to File y Map to Map. Erdas toma por defecto el tipo de acuerdo con las imagenes que intervienen en el proceso. Una imagen sera tomada como “map” siempre y cuando presente un sistema de coordenadas, aunque este sea “trucho”, es decir, colocado como “unknown” adrede por el usuario para algun fin determinado. Si la imagen no tiene sistema de coordenadas, sera considerada como “file”. El sistema de coordenadas puede verse y cambiarse desde el Image Info en el menu principal de Erdas o desde el menu Viewer Þ Image Info dentro del viewer. * Si se quiere anular ciertos GCPs para recalcular el error y ver que efcto tenian sobre el total de la transformacion, pero sin borrarlos todavia... usar el mouse + tecla “shift” para seleccionar los puntos que no quiero incluir en el calculo de la matriz, marcandolos desde el encabezado de fila del editor, luego siempre sobre el encabezado de las filas hacer click con el boton derecho y seleccionar “invert seleccion”. Luego ir al menu View Þ Use Current Seleccion dentro del GCP editor y recien entonces proceder a recalcular la matriz de transformacion. * Si se quiere borrar puntos definitivamente... seleccionar el punto desde el encabezado de fila y con el boton derecho del mouse desplegar el menu, seleccionando “delete”. El punto sera eliminado. *
Si se quiere cambiar el orden de transformacion, seleccionar desde los menues del Transform editor Þ Edit Þ Transformation order.
Página 7 de 10e
Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
* Si se desea que la matriz de transformacion y el error RMS se recalculen automaticamente con cada punto borrado o anulado... hacer click en el icono existente en el GCP editor. 5)- Remuestreo Desde el Transform Editor Þ Resample (icono) (Si tengo la matriz de transformacion guardada no hace falta desplegar las 2 imagenes y el GCP editor para remuestrear, se corre “resample” desde el menu Raster Þ Resample” del viewer que contiene la imagen a remuestrear. * Metodo de remuestreo: Ver arriba * Map Projection: Sera la que tenga la imagen que se tome por referencia (destination). * Output corners: Extremos de la imagen salida. Pueden modificarse si se quiere remuestrear un subset de la imagen original nada mas (ver abajo). indicar archivo de salida
metodo de remuestreo
fijar M ap Projection
tamano de las celdas en la salida
poner en "on" esta opcion
*
Tamano de las celdas (pixeles) en la salida: Es el tamano que tendran los pixeles de la imagen salida, pero referido a las unidades de referencia que tenga la imagen destination, las cuales se indican en “map projection” de resample (ver arriba)
Página 8 de 10e
Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
6)- Verificacion de la rectificacion Se visualizan al mismo tiempo la imagen rectificada y la imagen referencia (destination). Se puede hacer de 2 formas: (A)- Visualizando cada una de las imagenes en un viewer y utilizando la opcion link + el inquire cursor. Desde los menues de uno de los viewers... Þ Abrir ambas imagenes y utilizando Þ Session Þ Tile viewers, dentro del menu principal de Erdas, repartir la pantalla para ambas imagenes. Luego desde uno de los viewers... Þ View Þ Link/Unlink Þ hacer click dentro del otro viewer para vincular ambas imagenes. Luego, desde el mismo viewer... Þ View Þ Inquire cursor Þ moviendo el inquire cursor verificar si ambas imagenes concuerdan. (B)- Superponiendo ambas imagenes en el mismo viewer. Þ Abrir la imagen rectificada en un viewer y luego, desde el mismo viewer, clickear la opcion de “abrir archivo” nuevamente. Seleccionar la imagen destination y fijar la opcion “clear display” en off antes de dar el OK. Esto hace que la imagen destination quede superpuesta a la imagen rectificada, la que al quedar por debajo no se visualiza. Luego, desde los menues del viewer... Þ Utility Þ Swipe Aparece una cuadro de control, esta utilidad permite ir visualizando paso a paso la superposicion de la imagen que esta debajo sobre la primera. Guardado de los GCPs Si se quiere guardar los puntos de control, ya sea porque no se va a remuestrear en el mismo momento en que se tomaron, o simplemente como medida de seguridad, existen 2 maneras. 1- Desde el GCP editor Þ File Þ Save Source y luego repetir Þ Sfile Þ Save Destination. De esta manera, cuando se vuelvan a visualizar ambas imagenes al mismo tiempo y se seleccione la opcion “pairwise” desde el GCP editor, los GCPs apareceran en cada imagen. Se puede entonces continuar ingresando putos de control por ejemplo. 2- Con los puntos de control en el GCP editor Þ seleccionar desde el encabezado de columnas, con el mouse+ tecla “shift”, las 5 columnas correspondientes al # de GCP, valor de X, Y, X’ e Y’. Luego, siempre sobre el encabezado de columnas, seleccionar con el boton derecho del mouse la opcion “Export” y guardar la informacion en un archivo con extension .dat, que luego puede recuperarse con el procedimiento inverso usando la opcion “Import”. Guardado de la Matriz de Transformacion
Página 9 de 10e
Rectificación de imágenes satelitales, ejemplo usando ERDAS 8.3.1. y 8.2 José Bava (
[email protected])
Si tambien se desea se puede salvar la matriz de transformacion, esto se hace desde el Transform Editor Þ Save Þ se guarda el archivo con extension .cff.
Página 10 de 10e