IHS
´ rdoba Universidad Nacional de Co ´ tica, Astronom´ıa y F´ısica Facultad de Matema
Filtros de ruido speckle y t´ ecnicas de clasificaci´ on para Im´ agenes SAR implementadas en Python
Autores: Mat´ıas Herranz Joaqu´ın Tita
Director: Dr. Oscar Bustos
“La omnipresencia de las c´ amaras insin´ ua de modo persuasivo que el tiempo consiste en acontecimientos interesantes, dignos de fotografiarse. Esto a su vez permite sentir f´ acilmente que a cualquier acontecimiento, una vez en marcha, y sea cual fuera su car´ acter moral, deber´ıa permit´ırsele concluir para que algo m´ as pueda a˜ nadirse al mundo, la fotograf´ıa. Una vez terminado el acontecimiento, la fotograf´ıa a´ un existir´ a, confiri´endole una especie de inmortalidad (e importancia) de la que jam´ as habr´ıa gozado de otra manera.”
Susan Sontag.
1
Agradecimientos Al Dr. Oscar Bustos, el profe, nuestro director. A todas las autoridades de FaMAF, sin las que nuestras aspiraciones acad´emicas y de investigaci´on no hubieran podido concretarse. A Pedro Mourelle, que colabor´o con buen humor y el dise˜ no del logo de PyRadar. A Meli Herranz, por sus cr´ıticas constructivas y revisiones generales. Mat´ıas Herranz y Joaqu´ın Tita. A mis padres, sin los cuales ni de casualidad hubiera llegado hasta ac´a. A mis hermanas, por entenderme siempre. A Paula, mi novia, por el incansable cari˜ no y apoyo, por estar siempre. A mis amigos y amigas, por los caminos transitados y todo lo aprendido juntos. Mat´ıas Herranz. A toda mi familia que siempre me apoy´o y me ense˜ n´o a nunca bajar los brazos. A mi novia Paula, a qui´en le debo mucho de todo lo que hoy soy. A cada uno de mis amigos, compa˜ neros y todo aquel que estir´o su mano para ayudarme cuando lo necesit´e. Joaqu´ın Tita.
2
´Indice 1. Definiciones y Contexto 1.1. ¿Por qu´e nos interes´ o el tema? . . . . . . . . . . 1.2. Repercusi´ on social del tema y de nuestro trabajo 1.3. Acerca de Python . . . . . . . . . . . . . . . . . . 1.3.1. ¿Por qu´e elegimos Python? . . . . . . . . 1.4. Sobre im´ agenes de Radar . . . . . . . . . . . . . 1.4.1. Prefacio hist´ orico . . . . . . . . . . . . . . 1.4.2. El sistema de im´agenes de Radar . . . . . 1.5. Algunas definiciones en torno a la teledetecci´on .
. . . . . . . .
. . . . . . . .
8 8 9 9 9 11 11 11 16
2. Introducci´ on e Historia 2.1. Historia y Enfoque . . . . . . . . . . . . . . . . . . . . . . . . . 2.1.1. Primeras fotograf´ıas e im´agenes a´ereas previas al avi´on . 2.1.2. Primeros usos del avi´on para la fotograf´ıa . . . . . . . . 2.1.3. Primera Guerra Mundial (1914–1918) . . . . . . . . . . 2.1.4. Per´ıodo de entre-guerras (1919-1939) . . . . . . . . . . . 2.1.5. Segunda Guerra Mundial (1939-1945) . . . . . . . . . . 2.1.6. Guerra Fr´ıa (1945–1991) . . . . . . . . . . . . . . . . . . 2.1.7. Investigaci´ on de Robert Colwell . . . . . . . . . . . . . . 2.1.8. Aplicaciones civiles de im´agenes a´ereas . . . . . . . . . . 2.1.9. Sat´elite de Teleobservaci´on . . . . . . . . . . . . . . . . 2.1.10. Teledetecci´ on Hiperespectral . . . . . . . . . . . . . . . 2.1.11. Informaci´ on Geoespacial . . . . . . . . . . . . . . . . . . 2.1.12. Teledetecci´ on P´ ublica . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . .
18 18 18 18 19 19 20 21 22 22 24 25 26 27
3. Definiciones algor´ıtmicas y referencias matem´ aticas 3.1. Algoritmos de Clustering . . . . . . . . . . . . . . . . 3.1.1. Introducci´ on a los Algoritmos de Clustering . . 3.1.2. Clasificaci´ on No Supervisada . . . . . . . . . . 3.1.3. Clasificaci´ on Supervisada . . . . . . . . . . . . 3.1.4. Definici´ on del Algoritmo Kmeans . . . . . . . . 3.1.5. Definici´ on del Algoritmo Isodata . . . . . . . . 3.2. Filtros de Ruido Speckle . . . . . . . . . . . . . . . . . 3.2.1. Introducci´ on a los Filtros de Ruido de Speckle 3.2.2. Modelo Matem´atico . . . . . . . . . . . . . . . 3.2.3. Filtro de Kuan . . . . . . . . . . . . . . . . . . 3.2.4. Filtro de Lee . . . . . . . . . . . . . . . . . . . 3.2.5. Filtro de Frost . . . . . . . . . . . . . . . . . . 3.2.6. Filtro de Lee Mejorado . . . . . . . . . . . . .
. . . . . . . . . . . . .
29 29 29 30 31 33 35 37 37 39 39 39 40 40
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . . . . . . .
. . . . . . . .
. . . . . . . . . . . . .
. . . . . . . .
. . . . . . . . . . . . .
. . . . . . . .
. . . . . . . . . . . . .
. . . . . . . . . . . . .
4. Implementaci´ on Computacional 41 4.1. Qu´e Hicimos y C´ omo lo Hicimos . . . . . . . . . . . . . . . . . . 41 4.2. Estructura de Archivos y Carpetas de PyRadar . . . . . . . . . . 42 4.3. Descripci´ on de los m´ odulos . . . . . . . . . . . . . . . . . . . . . 43 4.3.1. core, manejo de im´agenes SAR y algoritmos de ecualizaci´ on . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 43 4.3.2. classifiers, algoritmos de clasificaci´on supervisada y no supervisada . . . . . . . . . . . . . . . . . . . . . . . . . . . 47 3
4.3.3. comparator, es el m´odulo donde se encuentran rutinas para comparar im´agenes . . . . . . . . . . . . . . . . . . . . 4.3.4. filters, filtros de ruido . . . . . . . . . . . . . . . . . . . . 4.3.5. simulate, en este m´odulo se proveen algoritmos para simular im´ agenes por capas. . . . . . . . . . . . . . . . . . . 4.3.6. Im´ agenes de ejemplo generadas con el simulador de im´agenes 4.3.7. utils, es un m´ odulo en que se agrupan utilidades y herramientas generales . . . . . . . . . . . . . . . . . . . . . . . 4.3.8. examples, ejemplos de uso de PyRadar . . . . . . . . . . . 4.4. Apertura y disponibilidad del c´odigo. . . . . . . . . . . . . . . . . 4.4.1. Apertura del C´odigo . . . . . . . . . . . . . . . . . . . . . 4.4.2. Disponibilidad del C´odigo . . . . . . . . . . . . . . . . . . 5. Resultados 5.1. Filtros . . . . . . . . . . . . . . . . . . . . . . . . 5.1.1. Mean . . . . . . . . . . . . . . . . . . . . 5.1.2. Median . . . . . . . . . . . . . . . . . . . 5.1.3. Frost . . . . . . . . . . . . . . . . . . . . . 5.1.4. Lee . . . . . . . . . . . . . . . . . . . . . . 5.1.5. Lee Enhanced . . . . . . . . . . . . . . . . 5.1.6. Kuan . . . . . . . . . . . . . . . . . . . . 5.2. Algoritmos de Clasificaci´on . . . . . . . . . . . . 5.2.1. K-Means . . . . . . . . . . . . . . . . . . 5.2.2. Isodata . . . . . . . . . . . . . . . . . . . 5.3. Filtros y Algoritmos de Clasificaci´on Combinados 5.3.1. Lee Enhanced y K-Means . . . . . . . . . 5.3.2. Lee Enhanced e Isodata . . . . . . . . . . 5.3.3. Frost y K-Means . . . . . . . . . . . . . . 5.3.4. Frost e Isodata . . . . . . . . . . . . . . .
49 51 54 57 63 67 68 68 68
. . . . . . . . . . . . . . .
. . . . . . . . . . . . . . .
. . . . . . . . . . . . . . .
. . . . . . . . . . . . . . .
. . . . . . . . . . . . . . .
. . . . . . . . . . . . . . .
. . . . . . . . . . . . . . .
70 71 71 72 73 76 79 82 85 85 86 87 88 89 90 91
6. An´ alisis de los Resultados y Observaciones 6.1. Performance . . . . . . . . . . . . . . . . . . . . . . . . 6.2. Mejoras y Modificaciones de los Algoritmos Originales 6.3. Comparaciones con ENVI . . . . . . . . . . . . . . . . 6.3.1. PyRadar y ENVI: Frost Filter. . . . . . . . . . 6.3.2. PyRadar y ENVI: Lee Enhanced Filter. . . . . 6.3.3. PyRadar y ENVI: Lee. . . . . . . . . . . . . . . 6.3.4. PyRadar y ENVI: Kuan. . . . . . . . . . . . . 6.3.5. PyRadar y ENVI: Isodata. . . . . . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
. . . . . . . .
92 92 93 94 94 95 96 97 98
. . . . . . . . . . . . . . .
. . . . . . . . . . . . . . .
Referencias
99
A. Utilidades 101 A.1. timeit, mide el tiempo de ejecuci´on de una porci´on de c´odigo . . 101 A.2. Profiling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 103 B. Sistema Operativo, Hardware y Software B.1. Datos Espec´ıficos sobre Hardware Utilizado . . . . . . . . . . . . B.2. Datos Espec´ıficos sobre Software Utilizado . . . . . . . . . . . . . B.3. Herramientas de Software utilizadas . . . . . . . . . . . . . . . .
4
104 104 104 105
C. C´ omo ejecutar los ejemplos de c´ odigo 107 C.1. C´ odigo de los ejemplos . . . . . . . . . . . . . . . . . . . . . . . . 107 C.2. Ejecuci´ on de las scripts de ejemplo . . . . . . . . . . . . . . . . . 107 D. Instalaci´ on 108 D.1. Software B´ asico Requerido . . . . . . . . . . . . . . . . . . . . . . 108 D.2. Instrucciones de Instalaci´on . . . . . . . . . . . . . . . . . . . . . 108
´Indice de figuras 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.
El espectro electromagn´etico y las transmitancia de la atm´osfera en respecto de la distancia entre la tierra y el espacio exterior. . Diagrama de os elementos esenciales de un sistema de imagen por radar. [6] . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Un Lockheed U-2R/TR-1 en vuelo [15]. . . . . . . . . . . . . . . Sat´elite Meteorol´ ogico TIROS 1 [16]. . . . . . . . . . . . . . . . ERTS 1 (renombrado LANDSAT 1) luego de tests previos al lanzamiento en General Electric. [17] . . . . . . . . . . . . . . . . . . Imagen RapidEye [18] . . . . . . . . . . . . . . . . . . . . . . . . Imagen de Google Earth, a˜ no 2008. [20] . . . . . . . . . . . . . . Imagen cruda de entrada e imagen con la clasificaci´on resultante [10]. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Imagen SAR ERS con ruido speckle [19]. . . . . . . . . . . . . . . Simulaci´ on de imagen original con distribuci´on Gamma. . . . . . Simulaci´ on de imagen original con distribuci´on K. . . . . . . . . Imagen con distribuci´on Gamma con ruido de distribuci´on 2 . . . Imagen con distribuci´on K con ruido de distribuci´on 2 . . . . . . Plot de la imagen con distribuci´on . . . . . . . . . . . . . . . . . Plots de la imagen con distribuci´on y de la imagen con distribuci´ on K juntos. . . . . . . . . . . . . . . . . . . . . . . . . . . . Logotipo de PyRadar. . . . . . . . . . . . . . . . . . . . . . . . . Mean filter . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Median filter . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Frost: Im´ agenes A.1, A.2, A.3 y A.4. . . . . . . . . . . . . . . . . Frost: Im´ agenes B.1, B.2, B.3 y B.4. . . . . . . . . . . . . . . . . Frost: Im´ agenes C.1, C.2, C.3 y C.4. . . . . . . . . . . . . . . . . Lee: Im´ agenes A.1, A.2, A.3 y A.4. . . . . . . . . . . . . . . . . . Lee: Im´ agenes B.1, B.2, B.3 y B.4. . . . . . . . . . . . . . . . . . Lee: Im´ agenes C.1, C.2, C.3 y C.4. . . . . . . . . . . . . . . . . . Lee Enhanced: Im´ agenes A.1, A.2, A.3 y A.4. . . . . . . . . . . . Lee Enhanced: Im´ agenes B.1, B.2, B.3 y B.4. . . . . . . . . . . . Lee Enhanced: Im´ agenes C.1, C.2, C.3 y C.4. . . . . . . . . . . . Kuan: Im´ agenes A.1, A.2, A.3 y A.4. . . . . . . . . . . . . . . . . Kuan: Im´ agenes B.1, B.2, B.3 y B.4. . . . . . . . . . . . . . . . . Kuan: Im´ agenes C.1, C.2, C.3 y C.4. . . . . . . . . . . . . . . . . K-Means . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Isodata . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . Lee Enhanced y K-Means . . . . . . . . . . . . . . . . . . . . . . Lee Enhanced y K-Means . . . . . . . . . . . . . . . . . . . . . .
5
13 14 21 22 24 26 27 29 37 57 58 59 60 61 62 68 71 72 73 74 75 76 77 78 79 80 81 82 83 84 85 86 88 89
35. 36. 37. 38. 39. 40. 41. 42. 43. 44. 45.
Frost y K-Means . . . . . . . . . . . . . Frost e Isodata . . . . . . . . . . . . . . PyRadar y ENVI: Frost Filter. . . . . . PyRadar y ENVI: Lee Enhanced Filter. PyRadar y ENVI: Lee. . . . . . . . . . . PyRadar y ENVI: Kuan. . . . . . . . . . PyRadar y ENVI: Isodata. . . . . . . . . NINJA-IDE. . . . . . . . . . . . . . . . Texmaker. . . . . . . . . . . . . . . . . . GIT. . . . . . . . . . . . . . . . . . . . . Trac. . . . . . . . . . . . . . . . . . . . .
6
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
. . . . . . . . . . .
90 91 94 95 96 97 98 105 105 105 106
Resumen
En este trabajo se propone la construcci´on y documentaci´on de una librer´ıa en Python que incluya la implementaci´on de varios algoritmos para el filtrado del ruido speckle y dos de clasificaci´on de im´ agenes SAR (Synthetic Aperture Radar). Asimismo, se incluir´a en la librer´ıa una basta variedad de algoritmos, herramientas y utilidades para realizar diversas operaciones con im´agenes satelitales de Radar. Posteriormente se expondr´an resultados comparativos entre la librer´ıa desarrollada y los obtenidos con la aplicaci´on del software propietario ENVI. Finalmente, se analizar´a el rendimiento del producto obtenido al trabajar con im´agenes SAR reales obtenidas por gentileza de la CONAE (Comisi´on Nacional de Actividades Espaciales). Keywords: remote sensing, filters, speckle, kmeans, isodata, classification, supervised, unsupervised, Frost, Lee, Kuan, enhanced-Lee, mean, median.
7
1. 1.1.
Definiciones y Contexto ¿Por qu´ e nos interes´ o el tema?
Durante el curso de las materias optativas del Dr. Bustos sobre an´alisis estad´ıstico de im´ agenes satelitales tuvimos por primera vez contacto con el procesamiento de im´ agenes, los programas, algoritmos y m´etodos que se utilizan en esta amplia ´ area de investigaci´ on y aplicaci´on. Tambi´en fue durante el trancurso de estas optativas que nos topamos con el hecho no menor de que no exist´ıan alternativas libres a la hora de elegir un software con el cual trabajar. Y es aqu´ı donde empieza nuestra idea para esta tesis. Despu´es de los a˜ nos que la facultad lleva invertidos en nuestra educaci´on, una educaci´ on que nos abre y seguir´a abriendo puertas tanto acad´emicas como laborales, consideramos que a la hora de elegir un camino para recorrer nuestros u ´ltimos pasos en la carrera y realizar la tesis de grado, era casi nuestro deber retribuir a la sociedad, a la facultad o, al menos, a la comunidad acad´emica o del software de alg´ un modo. Pensando todo esto junto, fue que decidimos emprender como proyecto de tesis el desarrollo de un software que brindara una opci´on libre y gratuita para procesamiento de im´ agenes satelitales. Intentar cubrir todos los tipos de sensores y algoritmos hubiera sido una empresa poco realista y que rebasar´ıa el alcance de una tesis de grado. Por esto fue que decidimos enfocarnos en las im´agenes de radar y, a su vez, en dejar la infraestructura lista y armada para que otros programadores -como nosotros mismos, quiz´a- puedan unirse y continuar con el proyecto, haci´endolo crecer hacia cubrir otros dominios de procesamiento de im´ agenes satelitales.
8
1.2.
Repercusi´ on social del tema y de nuestro trabajo
A lo largo de los meses que trabajamos en nuestro proyecto de tesis, muchas veces nos topamos con un tipo de obst´aculo en particular que fue particularmente dif´ıcil de sortear: el de encontrar que nuestros algoritmos generaban resultados distintos que los softwares usados en el ´area (como ENVI, por ejemplo), y no poder comparar en c´ odigo o c´ omo se realizaba el procesamiento de la informaci´on en estos softwares, dado que no eran de c´odigo abierto. Y justamente este tipo de problemas, en alg´ un sentido, busc´abamos tener. Nuestro proyecto de tesis se trata, en buena medida, de un puntapi´e inicial, una base para que futuros investigadores, alumnos de carreras de grado y la comunidad en general pueda, partiendo desde hasta donde nosotros alcancemos a llegar, poder ir mucho m´ as lejos, avanzar m´as all´a. Pero poder hacer esto contando con un cuerpo, con un software amplio, extensible y de c´odigo abierto, que podr´ an modificar, extender, usar como base, e incluso mejorar. De este modo, intentamos devolver a la sociedad y a la comunidad un poco de todo el aprendizaje que nuestra facultad nos brind´o durante nuestra carrera de grado, brindando una opci´ on libre y extensible a la hora de elegir un software para procesamiento de im´ agenes satelitales.
1.3.
Acerca de Python
1.3.1.
¿Por qu´ e elegimos Python?
El lenguaje que elegimos para programar nuestra librer´ıa para procesamiento de im´ agenes de RADAR es Python. Las razones por las cuales elegimos este lenguaje y no otro son muchas, entre las cuales algunas de las m´as importantes son: Python es software libre: dado que nuestra tesis se enfocar´a en proveer una opci´ on libre a un abanico de opciones de softwares caros y privativos, elegir un lenguaje libre para implentarlo era algo indispensable. Python tiene una comunidad muy activa y cooperativa: La comunidad de programadores que programan con Python es muy grande (y tiende a seguir creciendo), y se caracteriza fuertemente por una muy buena aceptaci´ on de nuevos miembros y una gran predisposici´on a la ayuda mutua en medios como foros y listas de correos, como la lista de correo de PyAr (Python Argentina). Esto muestra que Python es un lenguaje vivo y en crecimiento. Python es un lenguaje con una sintaxis sencilla y es f´acil de entender y aprender: esto fue un factor muy importante, dado que el foco de nuestra tesis de grado es en buena medida colaborar haciendo una retribuci´on a la comunidad, el cual se pueda extender, modificar y mejorar por cualquier persona que lo necesite, tanto para usos acad´emicos, de aplicaci´on o los que fueran. Python es r´ apido: Python cuenta con una enorme cantidad de librer´ıas magn´ıficas implementadas en C (como NumPy o SciPy), C++ (como GDAL) que permiten realizar operaciones complejas en cuanto a procesamiento de manera sumamente r´apida y eficiente. Incluso si no existiera 9
una librer´ıa en particular en C para la funcionalidad que se desea, se puede implementar la porci´ on de c´odigo que realiza la funci´on cr´ıtica en C y utilizarla desde Python. Python tiene un gran soporte para computaci´on cient´ıfica: con librer´ıas s´ olidas, fuertemente mantenidas, bien documentadas y ´altamente eficientes, como SciPy, NumPy, matplotlib, VPython y Sage. Python es multiplataforma: el c´odigo Python de nuestro trabajo de tesis podr´ a ser utilizado pr´ acticamente en cualquier sistema operativo. Esto es una ventaja important´ısima y no limita al usuario a un sistema operativo o hardware en particular. Y esto nos parece un gran pro en pos de la libertad del usuario. [21]
10
1.4. 1.4.1.
Sobre im´ agenes de Radar Prefacio hist´ orico
La utilizaci´ on del radar para obtener im´agenes de la superficie terrestre y sus recursos no es reciente. Sistemas de microondas operaban en la d´ecada de 1960 en aeronaves, por delante de los sistemas ´opticos de imagen en las regiones visible e infrarroja del espectro. A la teledetecci´ on ´ optica se le dio un fuerte impulso con el lanzamiento de la primera de la serie de sat´elites Landsat a mediados de 1970. Aunque el sat´elite Seasat se puso en marcha en la misma ´epoca (1978) llevando un sensor de radar, funcion´ o s´ olo durante unos 12 meses y no hubo sistemas de microondas en la misma cantidad como s´ı en plataformas ´opticas en servicio durante la d´ecada de 1980. Como resultado, el control remoto de detecci´on de la comunidad a nivel mundial tiende a desarrollarse con fuerza en torno a im´agenes ´opticas hasta las misiones de transbordador Shuttle, a mediados de 1980 y sat´elites de imagen de vuelo libre a mediados de 1990, junto con varias sofisticadas plataformas en aviones. Desde entonces, y en particular con las capacidades u ´nicas y la flexibilidad de im´ agenes de radar, se ha producido un enorme aumento de inter´es en la tecnolog´ıa de im´ agenes de microondas. A diferencia de imagen ´optica, la comprensi´ on de los fundamentos te´oricos de im´agenes de radar puede ser un reto, sobre todo cuando se es nuevo en el campo. 1.4.2.
El sistema de im´ agenes de Radar
Por qu´ e utilizar microondas Lograr un entendimiento de la obtenci´on remota de im´agenes con radar puede ser m´ as dificultosa que en el contexto de un sistema ´optimo dado que la tecnolog´ıa misma es m´ as complicada y los datos de recolectados en las im´agenes captadas son m´ as variados. En primer lugar, intentaremos establecer las razones del inter´es en los sistemas de imagen de radar como modalidad de muestreo remoto. Una respuesta sencilla puede encontrarse al examinar la longitud de onda de la radiaci´on captada, comparada con la de las radiaciones visibles e infrarrojas empleadas en sistemas ´ opticos. Los sistemas de imagen ´optica operan en longitudes de onda del orden de alrededor 1µ m- es decir, una millon´esima parte de un metro. Las im´ agenes de radar, por otro lado, se basan en microondas cuya longitud de onda ronda los 10cm- aproximadamente 100 mil veces mayores que las ´opticas. Dada esta gran disparidad en las longitudes de onda, es de esperarse que las caracter´ısticas de la superficie terrestre se muestren de distinto modo en im´agenes provenientes de muestreo con microondas y medios ´opticos. Este es, en efecto, el caso. Esto no hace un tipo de im´agenes mejores que otro: los hace distintos, y dado lo antes expuesto, es que en m´ ultiples aplicaciones se acude a una combinaci´ on de ambos m´etodos, lo cual deriva en un modelo de espectro m´as amplio y sistemas geogr´ aficos con informaci´on m´as rica y valiosa de las zonas muestreadas. Otra diferencia importante surge del nivel de penetraci´on a trav´es de medios tales como agua y follaje. Las longitudes de onda mayores de las ondas de los sistemas de radar permiten cierto grado de penetraci´on, en tanto que esto es
11
mucho m´ as limitado con medios ´opticos. Es por esto que mientras las im´agenes opticas generalmente representan los elementos de la superficie del paisaje, las ´ im´ agenes de radar cuentan con informaci´on de mayor nivel de complejidad, incluyendo datos volum´etricos y de sub-superficie. Por u ´ltimo, con el uso de radar, es posible tener el control sobre las propiedades de la energ´ıa incidente. Esto permite recolectar una vasta variedad de tipos de datos y da lugar a aplicaciones innovadoras como mapeos topogr´aficos, detecci´ on de cambios en el paisaje y, para un entorno acotado de situaciones, incluso modelado tridimensional de los vol´ umenes de objetos superficiales. Obtenci´ on de im´ agenes con microondas En la b´ usqueda de formar una imagen, cualquier sea la tecnolog´ıa que se emplee, la primer consideraci´ on a tener en cuenta es de d´onde provendr´a la energ´ıa que permitir´ a observar el terreno. En el caso de im´agenes ´opticas, se trata de de luz solar visible e infrarroja, o energ´ıa termal de la Tierra misma. A pesar de que existe una peque˜ na cantidad de energ´ıa de microondas provenientes de la Tierra y el sol, esta es tan ´ınfima que generalmente es necesario proveer una fuente de radiaci´ on a nuestros prop´ositos. Un transmisor de microondas se lleva en la plataforma de sensores remotos y sirve para iluminar la tierra. La energ´ıa dispersada de vuelta a la plataforma se recibe y se usa para crear una imagen del paisaje. Podr´ıa haber dos plataformas, uno que lleva la fuente de energ´ıa y el otro (o incluso varios) que reciben la energ´ıa dispersada. Hasta la fecha la mayor´ıa de los sistemas de teledetecci´on radar han utilizado la misma plataforma y se llaman monoest´atico. Cuando se utilizan dos plataformas el sistema de radar se denomina bi-est´atico, que ahora est´a surgiendo como una importante modalidad de percepci´on remota. La energ´ıa de microondas es s´olo una forma de radiaci´on electromagn´etica, es parte de un espectro continuo como se muestra en la figura 1. El espectro tambi´en incluye la energ´ıa visible e infrarroja que es la base de teledetecci´on optica. La diferencia m´ ´ as significativa en las propiedades es la longitud de onda. En principio, podr´ıa preverse el uso de cualquier longitud de onda para obtener im´ agenes de la superficie de la Tierra, el u ´nico l´ımite real es el nivel de energ´ıa disponible en la superficie. Pero, se describ´ıa reci´en, es necesario considerar tambi´en una plataforma que lleva su propia fuente de energ´ıa. Entonces necesitamos preguntarnos si existe alguna limitaci´ on fundamental para el uso de cualquier rango de longitud de onda particular para fines de detecci´on a distancia. Existe una en particular: la atm´ osfera de la Tierra no es transparente en todas las longitudes de onda. Eso es muy afortunado ya que la absorci´on de una proporci´on significativa de la radiaci´ on ultravioleta del sol es importante para nuestro bienestar. Existe tambi´en la absorci´ on sustancial atmosf´erica en el infrarrojo lejano. Las caracter´ısticas de absorci´ on de la atm´ osfera son bastante complejas debido a su composici´on molecular. Figura 1 muestra la transmitancia atmosf´erica como una funci´on de longitud de onda, que cubre el intervalo de longitud de onda hasta ultravioleta para el espectro de ondas de radio. Varios aspectos son dignos de menci´on. La mayor parte del espectro de los componentes de la atm´ osfera de vapor de agua, ox´ıgeno y di´oxido de carbono bloquean selectivamente la transmisi´on de energ´ıa electromagn´etica a trav´es de la atm´ osfera. Regiones en las que hay poca absorci´on se refieren a menudo como 12
Figura 1: El espectro electromagn´etico y las transmitancia de la atm´osfera en respecto de la distancia entre la tierra y el espacio exterior. ventanas atmosf´ericas, el ser m´as importante en la regi´on visible e infrarrojo cercano (⇠ 0,3 1,3µm), el infrarrojo medio (⇠ 1,5 1,8µm, ⇠ 2,0a2,6µm, ⇠ 3,0 3,6µm, ⇠ 4,2 5µm) y el infrarrojo t´ermico (⇠ 7, 0 15µm). Por debajo de la gama visible el contenido de ozono de la atm´osfera superior bloquea la radiaci´ on solar. La atm´ osfera tambi´en est´a pr´acticamente cerrada a la radiaci´on de longitudes de onda m´ as all´a del infrarrojo t´ermico hasta que nos encontramos con la parte de las ondas de radio del espectro. Para longitudes de onda m´ as all´ a de unos 3 cm la atm´ osfera se considera transparente. Para aplicaciones terrestres esto se aplica de forma indefinida. Sin embargo, en los caminos desde el espacio a la superficie de la tierra, o viceversa, la regi´on de la atm´osfera denominada ion´ osfera, que consiste en un d´ebil pero significativamente ionizado conjunto de capas, refleja la energ´ıa electromagn´etica con longitudes de onda m´ as largas que alrededor de unos 10 metros, tanto radiada hacia arriba de la tierra, o hacia abajo desde un veh´ıculo espacial. Para fines de detecci´on a distancia por lo tanto, consideran que la atm´osfera potencialmente es un problema en esas longitudes de onda m´ as largas. Componentes de un sistema de imagen por radar La figura 2 resume la tecnolog´ıa de im´agenes de radar; ilustra los componentes esenciales del sistema que necesitan ser entendidos en el desarrollo de una valoraci´ on global de la esfera. La primera consideraci´on es ser capaz de resolver el campo de inter´es en c´elulas resoluci´on, o pixels. Principios diferentes se utilizan para crear la resoluci´on en la direcci´on paralela al movimiento de la plataforma (a lo largo de la pista o azimut), y ortogonal a la misma (a trav´es de la pista o rango). La irradiaci´on del paisaje utiliza pulsos de energ´ıa, el tiempo que toman de la transmisi´ on del paisaje y de vuelta al radar determina la distancia a la parte del paisaje. Innovadoras t´ecnicas de procesado de se˜ nal se emplean para dar la mayor resoluci´on espacial posible en esta dimensi´on. 13
En la direcci´ on de la pista a lo largo del movimiento de la plataforma en relaci´ on con el paisaje se ver´ a para dar un cambio Doppler en la frecuencia de la radiaci´ on que se utiliza para iluminar el paisaje. Al mantener un registro del efecto Doppler en la plataforma pasa a las regiones de inter´es que veremos, una vez m´ as, los m´etodos de procesamiento de se˜ nales se pueden utilizar para alcanzar resoluciones espaciales muy altas en azimut incluso desde un punto de vista del espacio transmitidas - aqu´ı es donde el concepto de “de apertura sint´etica” entra en juego. El siguiente aspecto importante para comprender es c´omo una franja de datos de imagen se puede establecer, a partir de las im´agenes seleccionadas. A pesar de que es, en principio, s´olo una propiedad de la antena llevada en la plataforma que establece el ancho de la franja, nos encontraremos con que la franja se ve limitada por la aparici´on de se˜ nales ambiguas si es demasiado ancha. Como consecuencia vamos a considerar la aplicaci´on de un principio que ahora se conoce como ScanSAR para lograr una cobertura muy amplia en el lado de la plataforma. Para utilizar los datos de manera significativa es necesario entender las distorsiones que pueden haber sido introducidos en la imagen grabada. Ellas pueden ser muy graves y se requiere orientaci´on y gu´ıa experimentada sobre c´omo las im´ agenes de radar debe ser configuradas para reducir al m´ınimo los tipos de distorsiones que afectan a aplicaciones particulares.
Figura 2: Diagrama de os elementos esenciales de un sistema de imagen por radar. [6] Una vez que llegamos a esta etapa vamos a tener una comprensi´on de c´omo funciona el radar de apertura sint´etica. A continuaci´on necesitamos entender c´ omo se dispersa la radiaci´ on incidente desde el paisaje ya que la energ´ıa devuelta contiene informaci´ on sobre las propiedades de la parte de la superficie de la Tierra que se est´ a estudiando. Esa es una consideraci´on importante, ya que constituye la base de la percepci´on remota con microondas. Es un estudio significativo en derecho propio y su tratamiento es de gran inter´es de estudio. No
14
s´ olo es importante para comprender las propiedades de dispersi´on de materiales de la superficie de la tierra, sino que es deseable ser capaz de modelar, ya que puede ser un paso importante en la interpretaci´on de im´agenes por radar. Es u ´til introducir terminolog´ıa en este punto. En la teledetecci´on, por lo general, se resuelve la escena de inter´es en p´ıxeles. Lo mismo es cierto en el radar, sin embargo a menudo se utiliza “c´elula de resoluci´on” como sin´onimo de p´ıxel. Adem´ as, a menudo se llama al pixel “objetivo”. Esto se produce debido a la herencia de la tecnolog´ıa de radar como un medio para detectar objetos discretos. A veces, nuestros p´ıxeles se ver´an como objetivos concretos, como cuando est´ an dominados por un objeto de dispersi´on simple como un gran ´arbol o un edificio, o, posiblemente, un barco en la superficie del oc´eano. Se agrega a´ un m´ as complejidad por el hecho de que la tierra responde de manera diferente para diferentes polarizaciones y longitudes de onda de la energ´ıa incidente, del mismo modo que para diferentes ´angulos con los que el paisaje se ve. Una generalizaci´ on importante subyacente emerge - las propiedades de radar de una regi´ on sobre la superficie terrestre son dependientes de su geometr´ıa y su contenido de humedad. Una caracter´ıstica especial, y particularmente indeseable, de la teledetecci´on utilizando t´ecnicas consistentes en el uso de radiaciones electromagn´eticas es la aparici´ on de ruido de speckle. Veremos que ese es el resultado de la interferencia de la energ´ıa reflejada por los m´ ultiples y variados dispersores elementales que se producen dentro de una celda de resoluci´on (p´ıxeles). El ruido de Speckle necesita entenderse correctamente, como as´ı tambi´en los medios para reducir su impacto, que son tema constitutivo de este trabajo. Adem´ as, debido a la la naturaleza b´asicamente pura (o al menos coherente) de la energ´ıa utilizada en la obtenci´on de im´agenes de microondas, es posible desarrollar t´ecnicas intereferom´etricas partir de las cuales las caracter´ısticas topogr´ aficas del paisaje se pueden obtener con una resoluci´on espacial muy alta y con la que los cambios espaciales, ya sea a corto o de largo plazo, pueden ser detectados y mapeados con una precisi´on muy alta. Una t´ecnica de imagen asociado utiliza una muy peque˜ na cantidad de energ´ıa de microondas emitida naturalmente por la tierra como un medio para la formaci´ on de im´ agenes. Aunque no es una t´ecnica de radar, tiene especial relevancia para las aplicaciones oceanogr´ aficas y la humedad del suelo. Para concluir, una comprensi´on profunda de los sistemas de imagen por radar involucra: saber de d´ onde proviene la energ´ıa, y sus propiedades; saber c´ omo resolver la escena en p´ıxeles; apreciaci´ on de las propiedades de dispersi´on de caracter´ısticas de la superficie de la tierra; la comprensi´ on de la dependencia de la dispersi´on sobre los par´ametros del sistema tales como la longitud de onda de la radiaci´on, su polarizaci´on y el ´ angulo con el que la radiaci´on intersea la superficie de la tierra; saber c´ omo se forma una imagen; saber interpretar las im´ agenes grabadas;
15
comprender c´ omo la cartograf´ıa tem´atica puede llevarse a cabo a partir de los datos de radar, y apreciaci´ on de las aplicaciones especiales de im´agenes de microondas como interferometr´ıa.
1.5.
Algunas definiciones en torno a la teledetecci´ on
El concepto de teledetecci´ on ha sido redefinido muchas veces a lo largo del tiempo: Seg´ un Fischer et al. (1986 p. 34), la teledetecci´on es el arte o ciencia de decir algo sobre un objecto pero sin tocarlo. Seg´ un Lintz y Simonett (1976, p. 1), la teledetecci´on es la adquisici´on f´ısica de datos sobre un objeto sin tocarlo ni tener contacto con ´este. Seg´ un la Sociedad Americana de Fotograf´ıa, las im´agenes son adquiridas utilizando un sensor, a diferencia de otros m´etodos como el uso de las c´ amaras convencionales. El mismo, graba la escena utilizando escaneo electr´ onico, pudiendo tambi´en utilizar radiaci´on fuera del rango visual normal de las c´ amaras o las filmaciones. Ejemplos de ´esto, son las microondas, radares, im´ agenes t´ermicas, infrarrojas, ultravioletas, as´ı como tambi´en, multiespectrales. T´ecnicas especiales son aplicadas para procesar e interpretar las im´ agenes de teledetecci´on, con el fin de aportar datos a distintos sectores como agricultura, arqueolog´ıa, informaci´on forestal, geograf´ıa, geolog´ıa y otras ´areas. La teledetecci´ on es la observaci´on de un objetivo utilizando alg´ un dispositivo alejado por alguna distancia (Barrett y Curtis, 1976 p. 3) El t´ermino “teledetecci´ on” en un sentido m´as amplio, significa “reconocimiento a distancia” ( Colwell, 1966, p. 71) La teledetecci´ on, a pesar de no estar precisamente definida, incluye todos los m´etodos para obtener im´agenes u otras formas de grabaciones electromagn´eticas de la superficie terrestre, desde la distancia, junto con el tratamiento y procesamiento de estos datos. Luego, la teledetecci´on en el sentido amplio, se refiere a la detecci´on y grabaci´on de radiaci´on electromagn´etica desde el punto de vista del instrumento sensor. Esta radiaci´on pudo haber sido originada directamente desde el sensor o desde una fuente ajena, como la energ´ıa solar reflejada. (White, 1977, p. 1–2). Teledetecci´ on, es un t´ermino usado por muchos cient´ıficos para denominar el estudio de objetos remotos desde grandes distancias (superficies como la de la Tierra, la Luna, la Atm´ostfera, incluso fen´omenos estelares y gal´ acticos). En t´erminos generales, la teledetecci´on denota la utilizaci´on de modernos sensores, equipamento de procesamiento de datos, teor´ıa de la informaci´ on y metodolog´ıa de procesamiento, teor´ıa de la comunicaci´on y dispositivos, veh´ıculos espaciales y terrestres, sumados a sistemas te´oricos y pr´ acticos, con el prop´osito de llevar a cabo reconocimientos a´ereos o en el espacio de la superficie de la tierra (National Academy of Sciences, 1970, p. 1). 16
La teledetecci´ on es la ciencia de derivar informaci´on de un objeto utilizando mediciones hechas a distancia y sin tomar contacto con ´el. La medida comunmente usada para sistemas de teledetecci´on, es la energ´ıa electromagn´etica que emana el objeto de inter´es, aunque existen otras posibilidades como ondas s´ısmicas, ondas sonoras y tambi´en la fuerza gravitacional (D. A. Landgrebe, quoted in Swain and Davis, 1978, p. 1). La teledetecci´ on es la pr´ actica de obtener informaci´on de la superficie de la Tierra usando im´ agenes adquiridas desde una perspectiva a´erea utilizando para ello radiaci´ on electromagn´etica reflejada o emitida desde la superficie misma en una o m´ as regiones del espectro electromagn´etico (Introduction to Remote Rensing, Fifth Edition 2011, Campbell, Wynne).
17
2. 2.1.
Introducci´ on e Historia Historia y Enfoque
Las im´ agenes en general envuelven mucha informaci´on sobre posiciones, tama˜ nos, relaciones entre elementos e incluso hasta ideas. Por su naturaleza, pueden retratar informaci´on compleja sobre elementos que reconocemos como objetos y estos objetos poseer un profundo nivel de significado. Por otro lado, las personas poseen una alta capacidad de deducci´on de informaci´ on a partir de im´ agenes pero, a su vez, encuentran dificultad para interpretar aquellas que son visualmente complejas. El ser humano es tan h´ abil en dicha tarea, que cuando se quiere replicar esta capacidad en programas de computadoras, se toma conciencia de cuan poderosa es esta habilidad para derivar este tipo de informaci´on intrincada. Las im´ agenes pueden ser utilizadas en beneficio del ser humano, es por ello que se las materializa en muchas formas distintas: cartograf´ıa, fotograf´ıas a´eras, im´ agenes de radar, im´ agenes m´edicas, videos, esquemas, incluso hasta conceptos e ideas. El potencial conceptual que poseen las im´agenes es infinito. 2.1.1.
Primeras fotograf´ıas e im´ agenes a´ ereas previas al avi´ on
Dado que la teledetecci´ on, en su parte pr´actica, se centra en el an´alisis de im´ agenes de la superficie de la Tierra, su origen naturalmente comienza en la fotograf´ıa. Los primeros intentos para crear fotograf´ıas datan en 1800, cuando un grupo de cient´ıficos condujeron experimentos con qu´ımicos fotosensitivos, lamentablemente, olvidados por la historia. En 1839, Louis Daguerre, public´o los resultados de sus experimentos con qu´ımicos fotosensitivos. Es por esto que se toma a ´esta como la primera fecha arbitraria para definir el nacimiento de la fotograf´ıa. La primera captura de una fotograf´ıa a´erea es atribu´ıda a Gaspard F´elix Tournachon, tambi´en conocido por su seud´onimo “Nadar”. En 1858, Nadar, tom´ o la primera fotograf´ıa a´erea desde un globo aerost´atico de aire caliente sobre la ciudad de Par´ıs en Francia. En los a˜ nos siguientes, se realizaron numerosas mejoras en la tecnolog´ıa fotogr´afica y en los m´etodos de adquisici´on de las fotograf´ıas que eran tomadas desde globos y las cometas. Estas im´agenes a´ereas de la Tierra est´ an entre las primeras que se ajustan a la definici´on de teledetecci´ on. 2.1.2.
Primeros usos del avi´ on para la fotograf´ıa
En 1909, Wilbur Wright pilote´o el avi´on que captur´o durante su vuelo las primeras im´ agenes del paisaje italiano cerca de la zona de Centocelli; esas se dicen que fueron las primeras fotograf´ıas a´ereas tomadas desde un avi´on. La maniobrabilidad del avi´ on facilit´o el control de la velocidad, altitud y direcci´on necesaria para el manejo de la c´amara. Si bien hubo muchos intentos de combinar la c´ amara con el avi´ on, los instrumentos de la ´epoca no estaban adaptados para su uso con otros.
18
2.1.3.
Primera Guerra Mundial (1914–1918)
La primera Guerra Mundial marc´o el comienzo de la adquisici´on de fotograf´ıas a´ereas como rutina b´ asica en operaciones militares. Si bien las c´amaras utilizadas para tomar las fotograf´ıas durante el conflicto fueron dise˜ nadas espec´ıficamente para ser usadas desde los aviones, la adaptaci´on entre ambos segu´ıa siendo rudimentaria para los est´andares de la ´epoca. A medida que la guerra continuaba su curso, la fotograf´ıa a´erea para el reconocimiento y vigilancia militar se fue volviendo cada vez m´as valiosa, y as´ı las aplicaciones se convirtieron cada vez m´as sofisticadas. Sobre el final de la guerra, el rol de la fotograf´ıa a´erea en operaciones militares fue reconocido como un gran ´exito, a pesar de no tener programas especializados, estructuras organizativas ni doctrinas operacionales. 2.1.4.
Per´ıodo de entre-guerras (1919-1939)
Durante estos a˜ nos, se dieron numerosas mejoras. Los dise˜ nos de las c´amaras fueron mejorados y adaptados espec´ıficamente para ser utilizadas en aviones militares. La ciencia de la fotogrametr´ıa —procedimiento por el cual se realizan mediciones precisas por medio de fotograf´ıas— se aplic´o a la fotograf´ıa a´erea por primera vez, junto con el desarrollo de instrumentos espec´ıficos dise˜ nados para el an´ alisis de fotograf´ıas a´ereas. A pesar de que los fundamentos de la fotogrametr´ıa fueron definidos muchos a˜ nos antes, el campo se termin´o de desarrollar a fines de los a˜ nos veinte con la aplicaci´ on especializada de instrumentos fotogram´etricos. Desde estos or´ıgenes, un nuevo hito fue establecido en el campo: la aplicaci´on de la fotograf´ıa a´erea en los programas del gobierno. En un principio en el ´area de la cartogr´ afica topogr´ afica, y luego para relevamiento de suelo, cartograf´ıa geol´ ogica, relevamiento forestal e incluso tambi´en agr´ıcola. Durante este per´ıodo, muchas de las innovaciones eran conducidas por pioneros visionarios quienes se establec´ıan en la industria privada para desarrollar aplicaciones civiles de cartograf´ıa a´erea. Sherman Fairchild fund´o numerosas compa˜ n´ıas, incluyendo Fairchild Surveys y Fairchild Camera and Instrument, que se convirtieron en l´ıderes de aviaci´on y dise˜ no de c´amaras a´ereas. Talbert Abrams condujo muchas innovaciones en reconocimiento a´ereo, redise˜ no de c´ amaras y operaciones comerciales a nivel global. Mientras este per´ıodo se desarrollaba, el libro de Lee, “The Face of the Earth as Seen from the Air”, expuso aplicaciones para fotograf´ıas a´ereas en muchas disciplinas desde la perspectiva de aquella ´epoca. A pesar de que Lee expuso muchas aplicaciones, la mayor´ıa eran alcanzadas a largo plazo y a peque˜ nos pasos, aunque fue suficiente para capturar el inter´es del gobierno para asegurar el desarrollo del ´area. Sin embargo, el uso de la fotograf´ıa a´erea, tuvo un poco de resistencia entre los tradicionalistas, por la imperfecci´on de las t´ecnicas y los equipos; adem´as de incertidumbre con respecto al rol de la fotograf´ıa a´erea en la investigaci´on cient´ıfica y las aplicaciones pr´ acticas de ella. La crisis econ´ omica mundial (1929-1939), no solo impact´o a nivel econ´omico y financiero, sino tambi´en en el medio ambiente. Las preocupaciones acerca de los impactos sociales y econ´ omicos del desarrollo econ´omico rural, sumado a la erosi´ on generalizada del suelo y la escasez de suministros de agua, condujeron a las primeras aplicaciones del gobierno en el reconocimiento a´ereo con el objetivo
19
de registrar y monitorear el desarrollo econ´omico rural. En Estados Unidos, el departamento de Agricultura y el Tennessee Valley Authority — departamento federal creado en 1933 para proveer navegaci´on, electricidad, control de inundaciones y desarrollo econ´ omico en el ´area de Tennessee durante la crisis— llev´o a cabo planes para aplicar la fotograf´ıa a´erea para guiar la planificaci´on ambiental y desarrollo econ´ omico. Dichos esfuerzos, fueron de gran importancia a la institucionalizaci´ on del uso de la fotograf´ıa a´erea por parte del gobierno y la creaci´ on de una base experimental en la aplicaci´on de la fotograf´ıa a´erea. 2.1.5.
Segunda Guerra Mundial (1939-1945)
Durante los a˜ nos de guerra, el uso del aspecto electromagn´etico fue extendido desde el uso pr´ acticamente exclusivo del espectro visible electromagn´etico hacia otras regiones, en particular a la regi´on del infrarrojo y microondas (lejos del rango de la visi´ on humana). El conocimiento de estas regiones del espectro ha sido desarrollado tanto en las ciencias b´asicas como en las ciencias aplicadas durante los anteriores 150 a˜ nos. Sin embargo, durante los a˜ nos de guerra, la aplicaci´ on y el desarrollo de este conocimiento se aceler´o, como as´ı tambi´en los medios para aplicarla. Aunque los investigadores comprend´ıan el potencial del espectro no visible, el equipamiento, los materiales y la experiencia necesarios para aplicarla a problemas pr´acticos en la ´epoca no estaban disponibles. La investigaci´ on durante la guerra y la experiencia operacional obtenida, proporcionaron tanto el conocimiento te´orico como la pr´actica necesaria para el uso diario del espectro no visible en la teledetecci´on. Adem´ as, la experiencia y el entrenamiento de un gran n´ umero de pilotos, operadores de c´ amaras e int´erpretes fotogr´aficos, cre´o un conjunto de personal experimentado que una vez terminada la guerra ser´ıan capaces de transmitir sus habilidades y experiencias a aplicaciones civiles. Muchos de estos cient´ıficos tomaron posiciones importantes en la industria cient´ıfica y programas gubernamentales para aplicar fotograf´ıas a´ereas y teledetecci´on a un gran rango de problemas. Durante la Segunda Guerra Mundial, la expansi´on del reconocimiento a´ereo desde la t´ actica militar, con nacimiento en la Primera Guerra Mundial en torno a la inteligencia dentro del territorio enemigo para comprender la infraestructura industrial y del transporte, adem´as de se˜ nalar planes y objetivos a largo plazo, aument´ o considerablemente la importancia de la inteligencia a´erea. Tal vez, la mejor publicidad de estas capacidades fueron vistas con el papel de los int´erpretes fotogr´ aficos Brit´ anicos, quienes aportaron componentes claves de inteligencia militar para detectar los misiles alemanes V-1 y V-2, mucho antes de su lanzamiento. De esta manera, eliminando cualquier sorpresa e ideando un plan r´ apido para contra arrestarlos (Babington-Smith, 1957). La Segunda Guerra Mundial tambi´en gener´o una redefinici´on en cuanto al enfoque tem´ atico del reconocimiento a´ereo. Considerando que los int´erpretes fotogr´ aficos de la Primera Guerra Mundial, estaban enfocados en la identificaci´ on y examinaci´ on de equipamiento y fortificaciones militares, sus hom´ologos en la Segunda Guerra Mundial, se encargaron adem´as de examinar la topograf´ıa, vegetaci´ on, transitabilidad y otras cualidades del terreno. De esta forma, expandiendo el enfoque inicial y conocimiento b´asico y pr´actico de la foto interpretaci´ on.
20
2.1.6.
Guerra Fr´ıa (1945–1991)
El ´exito estrat´egico que tuvo la interpretaci´on fotogr´afica durante la Segunda Guerra Mundial, propuls´ o el inter´es en la vigilancia a´erea durante la Guerra Fr´ıa. Las tendencias tecnol´ ogicas establecidas durante la Segunda Guerra Mundial fueron continuadas y mejoradas. Sin embargo, a medida que el conflicto avanzaba, la interpretaci´ on fotogr´afica como medio estrat´egico, fue uno de los pocos medios confiables que hab´ıa para obtener informaci´on. Dentro de este contexto, se desarroll´ o el avi´ on U-2 (Figura 3) para la vigilancia a una gran altitud y adem´ as, nuevos sistemas de c´amaras para extender la aviaci´on y los sistemas ´ opticos mucho m´ as lejos de los l´ımites esperados.
Figura 3: Un Lockheed U-2R/TR-1 en vuelo [15]. En 1960, se construy´ o el sat´elite de reconocimiento estrat´egico CORONA, que proporcion´ o la recolecci´ on de im´agenes desde el espacio. La mayor contribuci´ on de parte del reconocimiento fotogr´afico durante la Guerra Fr´ıa, se present´ o en “la Crisis de los misiles en Cuba”. En 1962, int´erpretes fotogr´aficos de Estados Unidos pudieron detectar con certeza las etapas iniciales de la introducci´ on sovi´etica de misiles en Cuba, mucho antes de lo que los estrategas sovi´eticos lo hubieran anticipado, logrando as´ı establecer el escenario para mitigar uno de los m´ as serios incidentes de la Guerra Fr´ıa (Brugioni, 1991). La era de la Guerra Fr´ıa dio lugar a una bifurcaci´on de las capacidades, con una amplia separaci´ on entre las aplicaciones dentro de la parte civil y de la defensa y los establecimientos de seguridad. Los inicios de la era de la Guerra Fr´ıa entre los bloques occidental-capitalista y oriental-comunista crearon el ambiente para el desarrollo de t´ecnicas avanzadas de reconocimiento, que generalmente eran guardadas como secretos de defensa y no siempre estaban disponibles para aplicaciones civiles. A medida que nuevos y m´as sofisticados instrumentos eran desarrollados, las tecnolog´ıas eran reemplazadas y las sustituidas se liberaban para aplicaciones civiles.
21
2.1.7.
Investigaci´ on de Robert Colwell
Entre los desarrollos m´ as importantes para el ´area civil de la teledetecci´on, est´ a el trabajo de Robert Colwell (1956), que aplic´o pel´ıcula infrarroja en color (mas conocida como “pel´ıcula para detecci´on de camuflaje” desarrollada para la Segunda Guerra Mundial), a problemas en los cuales se quer´ıa identificar peque˜ nos cultivos de cereales, las enfermedades pose´ıan y otros problemas en la investigaci´ on de plantas. Aunque la mayor´ıa de los principios b´asicos de sus desarrollos fueron establecidos mucho antes, su investigaci´on sistem´atica de dimensiones pr´ acticas constituye un hito claro en el desarrollo del campo de la teledetecci´ on. Incluso en esta ´epoca, Colwell deline´o las bases para sensores remotos modernos y anticip´ o muchas dificultades y ventajas en este campo de investigaci´ on. 2.1.8.
Aplicaciones civiles de im´ agenes a´ ereas
Para fines de los a˜ nos 50, las im´agenes a´ereas se institucionalizaron tanto en aplicaciones gubernamentales como en el ´area civil como fuente de informaci´on cartogr´ afica. En 1960, ocurrieron una serie de desarrollos importantes de forma muy veloz. El primer sat´elite meteorol´ogico, TIROS-1, fue lanzado por la NASA el 1o de abril de 1960.
Figura 4: Sat´elite Meteorol´ogico TIROS 1 [16]. Este sat´elite fue dise˜ nado para las observaciones climatol´ogicas y meteorol´ ogicas, pero tambi´en sirvi´ o de base para el posterior desarrollo de sat´elites
22
para la observaci´ on de la superficie terrestre. Durante este per´ıodo, algunos de los instrumentos de teledetecci´on desarrollados para el reconocimiento militar, y clasificados como “secretos de defensa” fueron liberados para uso civil a medida que los dise˜ nos m´ as avanzados estuvieron disponibles para uso militar. Estos instrumentos extendieron el alcance de la observaci´on fuera del espectro electromagn´etico visible, hacia las zonas del infrarrojo y microondas, que antes no eran alcanzadas.
23
2.1.9.
Sat´ elite de Teleobservaci´ on
El 23 de Julio de 1972, se lanz´o el Landsat 1, el primero de muchos sat´elites con ´ orbita terrestre dise˜ nado para la observaci´on de la superficie terrestre, marcando as´ı otro hecho significativo en la historia aeroespacial. El sat´elite Landsat 1, prove´ıa por primera vez, observaci´on repetitiva sistem´atica de la superficie de la Tierra. Cada imagen Landsat 1 representaba grandes ´areas con un nivel bajo de detalle en distintas regiones del espectro electromagn´etico pero a la vez siendo suficiente para las mayor´ıa de las aplicaciones pr´acticas en distintos campos. La importancia de este sat´elite no era totalmente apreciada pero fue posible identificar tres contribuciones fundamentales a la teledetecci´on. Primero, la posibilidad de obtener informaci´ on multiespectral de amplias ´areas de la superficie terrestre, expandi´ o enormemente la experiencia de las personas y el inter´es en el an´ alisis de informaci´ on multiespectral. La informaci´on multiespectral ya estaba disponible a˜ nos antes pero solo para laboratorios de investigaci´on. Esta informaci´ on que repentinamente empezaron a proveer los sat´elites Landsat, expandi´o la cantidad de cient´ıficos interesados en an´alisis de informaci´on multiespectral. La segunda contribuci´ on fue crear un incentivo para la r´apida expansi´on y ampliar los usos del an´ alisis digital aplicado a teledetecci´on. Antes de poseer informaci´on Landsat, el an´ alisis de las im´ agenes era completado examinando impresiones y transparencias de im´ agenes a´ereas.
Figura 5: ERTS 1 (renombrado LANDSAT 1) luego de tests previos al lanzamiento en General Electric. [17] Los an´ alisis de im´ agenes digitales por computadora fueron posibles principalmente en las instituciones de investigaci´on; las computadoras personales, y la variedad de programas para an´alisis de im´agenes que ahora son comunes, no exist´ıan. Aunque los datos provenientes del Landsat se utilizaron inicialmente
24
como impresiones o transparencias, posteriormente se entregaron tambi´en en formato digital. La disponibilidad de los datos digitales en un formato est´andar, cre´ o el contexto que permiti´ o popularizar el an´alisis digital y sentar las bases para el desarrollo de software para an´alisis de im´agenes, que hoy en d´ıa es com´ un. Durante esta etapa, el proceso fotogram´etrico que originalmente era llevado a cabo utilizando instrumentos mec´anicos fue redefinido como an´alisis digital, llevando a mejoras en la precisi´on y en la serializaci´on de adquisici´on, procesamiento, producci´ on y distribuci´on de los datos obtenidos por teledetecci´ on. Una tercera contribuci´ on del programa Landsat, fue su papel como modelo para el desarrollo de otros sat´elites de observaci´on, dise˜ nados y operados por diversas organizaciones en todo el mundo. A principios de los a˜ nos 80, una segunda generaci´on de instrumentos para recolecci´ on de im´ agenes satelitales aportaron una resoluci´on espacial m´as precisa con niveles de detalle de 30m, 20m y 10m. Ya para la d´ecada de los 90s las im´ agenes alcanzaron el metro y tama˜ nos menores de resoluci´on. A finales de los a˜ nos 90s, surgi´ o la posibilidad de explotar el desarrollo comercial de la adquisici´ on de im´ agenes satelitales, con compa˜ n´ıas como Geoeye (sat´elite IKONOS), RapidEye (sat´elite RapidEye AG). Este surgimiento, abri´o las puertas para el desarrollo de nuevas aplicaciones civiles que anteriormente estaba disponible s´ olo a trav´es del uso de la fotograf´ıa a´erea. Es importante se˜ nalar que estos progresos en el campo de la teledetecci´on avanzada, sumado a su vez al avance de los sistemas de informaci´on geogr´aficos (GIS), proveyeron la posibilidad de llevar los datos de teledetecci´on y otros datos geoespaciales en un marco anal´ıtico com´ un, de esta forma mejorando los productos y la apertura de nuevos mercados, como por ejemplo, infraestructura urbana, apoyo a la agricultura de alta precisi´on, y el apoyo de la cartograf´ıa de inundaci´ on. 2.1.10.
Teledetecci´ on Hiperespectral
Durante los a˜ nos 80, cient´ıficos en el “Jet Propulsion Laboratories” en California, Estados Unidos, comenzaron a desarrollar con acci´on conjunta de la NASA, instrumentos que pudieran crear im´agenes de la Tierra a nivel de detalle espectral sin precedentes. Los sensores multiespectrales hasta ahora solo recolectaban datos en unas pocas bandas espectrales distintas, con estos nuevos instrumentos se pod´ıa recolectar 200 o incluso m´as bandas espectrales de forma precisa y con gran nivel de detalle. Estos nuevos instrumentos crearon el campo de la teledetecci´ on hiperespectral, que incluso todav´ıa hoy est´a en v´ıas de desarrollo como campo de investigaci´on. La teledetecci´on hiperespectral aport´o un avance significativo al an´ alisis de la teledetecci´on, llev´andola a nuevos niveles y formando una base para una comprensi´on m´as profunda de la mejor forma de desarrollar las futuras capacidades de los sensores remotos.
25
Figura 6: Imagen RapidEye [18] 2.1.11.
Informaci´ on Geoespacial
A comienzos de la d´ecada de 1980 pero sin madurar completamente sino hasta mediados del a˜ no 2000, varias tecnolog´ıas comenzaron a converger, creando sistemas con una sinergia u ´nica, mejorando y reforzando el valor de los otros, para crear informaci´ on geoespacial, un t´ermino aplicado colectivamente a varias tecnolog´ıas, principalmente relevamiento remoto, GIS, y GPS. Muchos de estos sistemas son implementados en el plano digital, reemplazando as´ı aplicaciones manuales y mec´ anicas desarrolladas muchas d´ecadas antes. A pesar de que muchas de estas tecnolog´ıas fueron desarrolladas como “tecnolog´ıas interrelacionadas”, para principios de la d´ecada del 2000 llegaron a formar sistemas complejos que pod´ıan adquirir im´agenes de alta precisi´on de posicionamiento y que permit´ıan la integraci´on de otros tipos de informaci´on e im´agenes. As´ı, recolectar y analizar informaci´on geoespacial se transform´o de depender de varias tecnolog´ıas diferentes, vagamente conectadas, a la formaci´on de instrumentos completamente integrados y sin´ergicos que reforzaban el valor de los dem´ as. Estos incrementos en calidad y flexibilidad de la informaci´on geoespacial, vinculados tambi´en a la descentralizaci´on de los servicios de tecnolog´ıa de la informaci´ on y la adquisici´ on de im´agenes a´ereas, incrementaron significativamente la cantidad de emprendedores a desarrollar productos innovadores y adaptados a mercados espec´ıficos. Dentro de las preocupaciones principales se inclu´ıan las necesidades para la planificaci´on y la construcci´on, las necesidades civiles en la agricultura y la silvicultura, la hidrolog´ıa y los servicios de consumo y productos (aunque, debe tenerse en cuenta, a menudo en nichos muy espec´ıficos dentro de cada una de estas ´areas). La amplia disponibilidad de datos de gran resoluci´ on provistos por sat´elites hab´ıa abierto nuevos mercados en los me-
26
dios de comunicaci´ on, organizaciones sin fines de lucro y no gubernamentales, adem´ as de continuar los mercados ya existentes en las comunidades militares y de seguridad nacional. Despu´es de los eventos asociados al ataque del 11 de Septiembre y el Hurac´ an Katrina, el valor de las im´agenes de sat´elite para la seguridad nacional fue reconocido m´as formalmente e integrado en la respuesta de emergencias y planificaci´ on de la seguridad nacional en los planos nacional, regional y local. 2.1.12.
Teledetecci´ on P´ ublica
Durante la primera d´ecada del siglo veintiuno, el constante crecimiento de Internet comenz´ o a potenciar el acceso p´ ublico de im´agenes de teledetecci´on, en parte fue a trav´es del dise˜ no de productos para im´agenes con una veloz distribuci´ on v´ıa Internet, y por otra parte, a trav´es del dise˜ no de productos de consumo y servicios, con base en im´agenes de teledetecci´on, presentadas en formato de mapa. Mientras la historia previa de la teledetecci´on puede ser vista como el trabajo de especialistas para generar productos especializados para el uso y consumo de otros especialistas, estos desarrollos fomentaron el dise˜ no de productos para el uso del p´ ublico en general. Google Earth, lanzado al p´ ublico en el a˜ no 2005, forma una representaci´on virtual de la superficie de la Tierra utilizando para ello muchas im´agenes digitales compuestas, aplicando conceptos b´asicos desarrollados por Keyhole, Inc., quien fue posteriormente adquirida por Google Inc. en el 2004. Google Earth fue dise˜ nado para comunicarse con una audiencia espec´ıfica, el p´ ublico que no posee conocimientos especializados en im´agenes de teledetecci´on, en contraste a a˜ nos anteriores donde ´este se requer´ıa como base inapelable. Este nuevo producto solamente necesitaba como premisa que el usuario est´e familiarizado con Internet para poder manejarlo.
Figura 7: Imagen de Google Earth, a˜ no 2008. [20]
27
Una innovaci´ on b´ asica de Google Earth es el reconocimiento del valor de un amplio rango de im´ agenes producidas por georeferencias precisas de im´agenes adquiridas en varias fechas, escalas, y resoluciones. Dichas composiciones pueden ser vistas usando una intuitiva interface para la navegaci´on, cambio de escala, orientaci´ on, y controlar el nivel de detalle. Google Earth ha a˜ nadido herramientas especializadas para adaptar el software a las necesidades de usuarios espec´ıficos, mejorando las opciones de visualizaci´on, y la integraci´on de otros datos dentro del producto Google Earth. El producto se basa en la idea de que Google Earth podr´ıa atraer a una poblaci´on de personas no especializadas, mientras que al mismo tiempo otorgue capacidades especializadas a los usuarios en la defensa, seguridad, respuesta a emergencias, y el comercio que requieran acceso inmediato y simult´ aneo a un cuerpo com´ un de datos geoespaciales. Google Earth y muchos servicios similares online, representan una nueva clase de aplicaciones de im´ agenes de teledetecci´on que contrastan fuertemente con las aplicaciones de las ´epocas anteriores. En 2010, Campbell, Hardy, y Bernard, delinearon el contexto para el desarrollo de muchas de estas nuevas aplicaciones cartogr´ aficas de datos de teledetecci´on. Entre ellas, se encuentran: una pol´ıtica p´ ublica que mantiene restricciones laxas sobre la adquisici´on de datos satelitales de alta resoluci´on. pol´ıticas de privacidad personal que favorecen la recolecci´on de im´agenes. disponibilidad y cantidad de dispositivos confiables personales de navegaci´ on. incremento de aplicaciones portadas a dispositivos m´oviles. Dichos desarrollos originaron una clase de productos cartogr´aficos derivados de datos obtenidos por teledetecci´on, incluyendo por ejemplo, aquellos creados por MapQuest y softwares relacionados con navegaci´on, cuya base es la red de carreteras que es continuamente actualizada con el an´alisis de im´agenes de teledetecci´ on. [10]
28
3. 3.1. 3.1.1.
Definiciones algor´ıtmicas y referencias matem´ aticas Algoritmos de Clustering Introducci´ on a los Algoritmos de Clustering
El an´ alisis de agrupamiento o clustering es la tarea de asignar un conjunto de objetos en grupos (llamados cl´ uster) de modo que los objetos en el mismo grupo son m´ as similares (en alg´ un sentido) entre s´ı que a los de otros grupos. La denominaci´ on “an´ alisis de cl´ uster” o “clustering” no implica el uso de un determinado m´etodo estad´ıstico o modelo, como s´ı lo har´ıa un an´alisis discriminante, an´ alisis factorial y la regresi´on. A menudo no se tiene que hacer ninguna suposici´ on sobre la distribuci´on subyacente de los datos. Utilizando el an´ alisis de cl´ uster, tambi´en se pueden formar grupos de variables relacionadas, de forma similar a lo que se har´ıa en el an´alisis factorial. Hay muchas maneras de ordenar los casos en cl´ uster o grupos. La elecci´on de un m´etodo depende, entre otras cosas, del tama˜ no de los datos de entrada. Los m´etodos utilizados para conjuntos peque˜ nos de datos no son pr´acticos para los archivos de datos con miles de casos. El clustering es una tarea principal de la miner´ıa de datos de exploraci´on, y una t´ecnica com´ un para el an´alisis de datos estad´ısticos utilizados en muchos campos, incluyendo machine learning, reconocimiento de patrones, an´alisis de im´ agenes, recuperaci´ on de informaci´on, y la bioinform´atica. Para realizar un agrupamiento por clases o clustering, se comienza con una serie de casos que luego se desea subdividir en grupos homog´eneos. En primer lugar, se deben elegir las variables en las que se desea que los grupos sean similares. A continuaci´ on, se debe decidir si estandarizar las variables de alguna manera para que los datos contribuyan por igual a la distancia o similitud entre los casos. Por u ´ltimo, se debe decidir el procedimiento a utilizar, bas´andose en el n´ umero de casos y tipos de variables que se desee utilizar para la formaci´on de los grupos.
Figura 8: Imagen cruda de entrada e imagen con la clasificaci´on resultante [10]. El an´ alisis de conglomerados en s´ı no es un algoritmo espec´ıfico, sino la tarea general que hay que resolver. Esto se puede lograr por medio de varios algoritmos que difieren significativamente en su noci´on de lo que constituye un cluster 29
y en c´ omo encontrarlos eficientemente. Las nociones populares de grupos incluyen los grupos con bajas distancias entre los miembros del cl´ uster, las zonas densas del espacio de datos, intervalos o particulares distribuciones estad´ısticas. La agrupaci´ on por lo tanto se puede formular como un problema de optimizaci´ on multi-objetivo. Los correspondientes valores de agrupaci´on de algoritmos y par´ ametros (incluidos los valores tales como la funci´on de la distancia a usar, un umbral de densidad o el n´ umero de grupos previstos) depender´a de los datos individuales establecidos y el uso previsto de los resultados. El an´alisis de cl´ uster, como tal, no es una tarea autom´atica, sino un proceso iterativo de descubrimiento de conocimiento o interactivo de optimizaci´on multi-objetivo que implica ensayo y error. A menudo ser´a necesario modificar los par´ametros de procesamiento hasta que el resultado alcance las propiedades deseadas. Adem´ as del t´ermino “clustering”, hay una serie de t´erminos con significados similares, incluida la clasificaci´on autom´atica, la taxonom´ıa num´erica y el an´ alisis tipol´ ogico. Las sutiles diferencias son a menudo en el uso de los resultados: mientras que en la miner´ıa de datos, los grupos resultantes son la materia de inter´es, en la clasificaci´ on autom´atica su poder discriminativo es de inter´es. Esto conduce a menudo a malentendidos entre los investigadores procedentes de los campos de miner´ıa de datos y de machine learning, ya que utilizan los mismos t´erminos y con frecuencia los mismos algoritmos, pero tienen diferentes objetivos. [13, 14] 3.1.2.
Clasificaci´ on No Supervisada
La clasificaci´ on no supervisada puede ser definida como la identificaci´on de grupos naturales o estructuras dentro de informaci´on multiespectral. La noci´on de la existencia de agrupaciones naturales, inherentes a los valores espectrales dentro de una escena, puede no ser intuitivamente obvia, pero se puede demostrar que las im´ agenes de sensores remotos se componen generalmente de las clases espectrales que son razonablemente uniforme internamente con respecto a la intensidad en varios canales espectrales. La clasificaci´on no supervisada es la definici´ on, identificaci´ on, etiquetado, y la asignaci´on de estas clases naturales. Ventajas No se requiere conocimiento previo de la regi´on de inter´es. M´as precisamente, la naturaleza de los conocimientos necesarios para la clasificaci´on no supervisada difiere de la exigida para la clasificaci´on supervisada. Para llevar a cabo una clasificaci´on supervisada, es un requisito tener conocimiento detallado de la zona a examinar para poder seleccionar muestras representativas de cada clase a asignarse. Para llevar a cabo una clasificaci´ on no supervisada, no se necesita ning´ un conocimiento previo, pero el conocimiento de la regi´ on es necesario para posteriormente poder interpretar el significado de los resultados producidos por el proceso de clasificaci´on. El error humano por la interacci´on humana se minimiza. Para realizar una clasificaci´ on no supervisada, el usuario tal vez s´olo especifique el n´ umero de categor´ıas que desee (o, posiblemente, el m´ınimo y m´aximo en el n´ umero de categor´ıas) y, a veces las limitaciones que rige la distinci´on y la homogeneidad de los grupos. Muchas de las decisiones requeridas para la
30
clasificaci´ on supervisada no son necesarias para la clasificaci´on no supervisada, por lo que al analista se presenta con menos probabilidad de error. Si el analista tiene ideas imprecisas con respecto a la regi´on, esas ideas tendr´ an pocas posibilidades de incidir negativamente en el resultado de la clasificaci´ on. Las clases definidas por la clasificaci´on no supervisada son a menudo mucho m´ as uniforme con respecto a la composici´on espectral que las generadas por la clasificaci´on supervisada. Las clases son reconocidas como unidades distintas. Estas clases, quiz´as algunas ´ areas muy peque˜ nas, puedan permanecer no reconocidas en el proceso de clasificaci´ on supervisada e inadvertidamente se las podr´ıa incorporar en otras clases, generando un error de clasificaci´on y provocando confusi´ on al momento de interpretaci´on de la clase completa. Desventajas y Limitaciones Las desventajas y limitaciones de la clasificaci´on no supervisada se deben principalmente a la dependencia de los grupos “naturales” y las dificultades para hacer coincidir estos grupos con las categor´ıas que se consideran de inter´es para el analista. La clasificaci´ on no supervisada identifica clases espectralmente homog´eneas dentro de los datos pero que no necesariamente corresponden a las categor´ıas que son de inter´es para el analista. Como consecuencia, el analista se enfrenta con el problema de hacer coincidir las clases espectrales generadas por la clasificaci´on con las categor´ıas que son requeridas por el usuario final de la informaci´on. Rara vez existe una relaci´on simple entre los dos conjuntos. El analista tiene un control muy limitado sobre las clases y sus identidades espec´ıficas. Si es necesario generar un grupo espec´ıfico de clases de informaci´ on (por ejemplo, para que coincida con otras clasificaciones para otras fechas o regiones adyacentes), el uso de la clasificaci´on no supervisada puede resultar inadecuado. Las propiedades espectrales de las clases seguramente cambien con el paso del tiempo (en un base estacional, as´ı como durante los a˜ nos). Como resultado, las relaciones entre las clases naturales y las clases espectrales no sean constantes, y las relaciones definidas para una imagen puedan no extenderse a las dem´ as. 3.1.3.
Clasificaci´ on Supervisada
Se puede definir de manera informal como el proceso de utilizaci´on de muestras de identidad conocida (es decir, los pixeles ya tienen una asignaci´on a algunas clases) para clasificar los pixeles de identidad desconocida (es decir, para asignar pixeles sin clasificar a una de muchas clases disponibles). Las muestras de identidad conocida son los pixeles situados dentro de las ´areas de formaci´on, o de los campos de entrenamiento. El analista define las ´areas de formaci´on mediante la identificaci´ on de regiones en la imagen que pueden coincidir claramente con las ´ areas de identidad conocida en la imagen. Tales ´areas deben representar propiedades espectrales de las categor´ıas y, por supuesto, deben ser homog´eneas 31
(o en la mayor proporci´ on posible) en lo que respecta a la categor´ıa de informaci´ on que debe clasificarse. Es decir, las ´areas de formaci´on no deber´ıan incluir las regiones inusuales, ni sobrepasar fronteras entre las categor´ıas. El tama˜ no, forma posici´ on deber´ an favorecer la identificaci´on tanto en la imagen como en el suelo. Los p´ıxeles situados dentro de estas ´areas forman las muestras de entrenamiento que se utilizan para guiar el algoritmo de clasificaci´on donde se asignan valores espec´ıficos del espectro de las clases. Es evidente que la selecci´on de estos datos de entrenamiento es un paso clave en la clasificaci´on supervisada. Ventajas Las ventajas de la clasificaci´on supervisada, en relaci´on a la clasificaci´on no supervisada, se pueden enumerar de la siguiente manera: El analista tiene el control total de la selecci´on de las clases con un prop´osito espec´ıficas y regi´ on geogr´afica. Esta propiedad puede ser muy importante si es necesario generar una clasificaci´on para un prop´osito espec´ıfico de comparaci´ on con otras clasificaciones de la misma zona en fechas diferentes, o si la clasificaci´on debe ser compatible con las de las regiones vecinas. Bajo tales circunstancias, la asignaci´on imprevisible de clases (es decir, con respecto al n´ umero, la identidad, el tama˜ no, y patr´on) puede afectar la calidad de las categor´ıas generadas por clasificaci´on no supervisada resultando inconveniente o inadecuado. La clasificaci´ on supervisada est´a vinculada a ´areas espec´ıficas con identidad conocida, determinadas a trav´es del proceso de selecci´on de ´areas de entrenamiento. El analista, con la clasificaci´on supervisada, no se enfrenta con el problema de hacer coincidir las categor´ıas espectrales en el mapa final con las categor´ıas de inter´es (esta tarea se realiza durante el proceso de selecci´on de datos de entrenamiento). El operador puede ser capaz de detectar errores en la clasificaci´on mediante la examinaci´ on de los datos de entrenamiento. Datos de entrenamiento inexactos tienden a generar serios problemas en la clasificaci´on. Pero a pesar de la correcta clasificaci´on de los datos de entrenamiento no siempre concluye en una clasificaci´on correcta de los datos. Desventajas y Limitaciones Las desventajas de la clasificaci´on supervisada son numerosas. El analista, en efecto, impone una estructura de clasificaci´on de los datos. Estas clases definidas por el operador pueden no coincidir con las clases naturales que existen dentro de los datos y por lo tanto pueden no ser distinguibles o bien definidas en el espacio de datos multidimensional. Los datos de entrenamiento se definen principalmente con referencia a las categor´ıas deseadas y s´ olo en segundo lugar con referencia a las propiedades espectrales. Un ´ area que es “100 % bosque” puede ser precisa con
32
respecto a “forestaci´ on” , pero todav´ıa puede ser confusa en cuanto a densidad, la edad, de la sombra, y similares caracter´ısticas, y por lo tanto formar un ´ area de escasa informaci´on. Los datos de entrenamiento seleccionados por el analista pueden no ser representativos con los encontrados dentro de la imagen. Esto puede ocurrir a pesar de los esfuerzos del analista, especialmente si el ´area que clasifica es grande, compleja o de dif´ıcil acceso. La selecci´ on consciente de los datos de entrenamiento puede ser una tarea que demande mucho tiempo, sea costosa y tediosa. El analista puede experimentar problemas para adaptar las ´areas posibles de entrenamiento tal como se encuentran en los mapas a´ereos y fotograf´ıas para ser clasificados. La clasificaci´ on supervisada puede no ser capaz de reconocer y representar a las categor´ıas especiales o u ´nicas que no est´an representadas en los datos de entrenamiento, posiblemente porque no son de conocimiento para el analista o porque ocupan ´areas muy peque˜ nas en la imagen.
A continuaci´ on se presentan dos algoritmos de clustering de amplia aplicaci´ on y uso en im´ agenes SAR: los algoritmos K-Means e Isodata. 3.1.4.
Definici´ on del Algoritmo Kmeans
Dado un conjunto de datos {xi , i = 1, 2, . . . , N } que se desea clasificar y sean K e IT ER M AX par´ ametros prefijados. El algoritmo se describe a continuaci´on: (0)
(0)
(0)
(0)
Paso 1: Elegir arbitrariamente K centros de clase iniciales m1 , m2 , m3 , . . . , mK , tomados del conjunto de datos (el super´ındice indica el n´ umero de iteraci´ on). Paso 2: Asignar cada uno de los elementos del conjunto de datos {xi , i = 1, 2, . . . , N } a un cluster de acuerdo a la distancia que se encuentre de cada centro de clase: (l)
(l)
x 2 !j if DL (x, mj ) = min{DL (x, mi ), i = 1, 2, 3, . . . , K}
(1)
(l)
donde !i denota el i-´esimo cluster cuyo centro de clase es mi de la l-´esima iteraci´ on y DL es la funci´on que calcula distancia eucl´ıdea. (l)
(l)
(l)
Paso 3: Actualizar los centros de clase m1 , m2 , . . . , mK : (l+1)
mj
=
1 X x, (j = 1, 2, . . . , K) Kj x2!
(2)
j
(l)
donde Kj ´esima.
(l)
es la cantidad de elementos dentro de !j
33
en la iteraci´on l-
Paso 4: Terminar el algoritmo si los centros de clases no son alterados (o se alteran en proporci´ on) o si se excede el n´ umero de iteraciones m´aximas permitidas IT ER M AX. Sino, incrementar el n´ umero de iteraci´on l e ir al Paso 2[8].
34
3.1.5.
Definici´ on del Algoritmo Isodata
Dado un conjunto de datos {xi , i = 1, 2, . . . , N } que se desea clasificar y sean los siguiente par´ ametros prefijados: K: n´ umero de cl´ uster esperado. I: m´ aximo n´ umero de iteraciones. ⇥M : umbral para el m´ınimo n´ umero de elementos que puede contener un cluster. ⇥S : umbral para la desviaci´on est´andard (para hacer “splits” de cl´ uster) . ⇥C : umbral para la distancia entre cl´ uster (para hacer “merge” de cl´ uster). P : cantidad m´ axima de pares para hacer “merge” de cl´ uster. El algoritmo se describe a continuaci´on: Paso 1: Elegir arbitrariamente k (no necesariamente igual a K) centros de cl´ uster iniciales: m1 , m2 , m3 , . . . , mk , tomados de los datos {xi , i = 1, 2, . . . , N }. Paso 2: Asignar cada uno de los elementos al centro de cluster m´as cercano: x 2 !j if DL (x, mj ) = min{DL (x, mi ), i = 1, 2, 3, . . . , k}
(3)
donde !i denota el i-´esimo cluster cuyo centro de clase es mi y DL es la funci´ on que calcula distancia eucl´ıdea. Paso 3: Eliminar los cl´ usters que contengan pocos elementos, es decir, los cuales cumplan que Kj < ⇥M , donde Kj es la cantidad de elementos que pertenecen a !j . Luego descartar !j y decrementar k. Paso 4: Actualizar cada centro de clase: 1 X mj = x, (j = 1, 2, . . . , k) Kj x2!
(4)
j
Paso 5: Computar la distancia promedio Dj de los elementos de cada cluster !j hacia sus respectivos centros clase mj : 1 X Dj = DL (x, mj ), (j = 1, 2, . . . , k) (5) Kj x2! j
Paso 6: Computar la distancia total de los elementos hacia los centros de clase: D = Pk
1
j=1
o equivalentemente, ´ D=
k X j=1
Kj
k X
Kj D j ,
(6)
j=1
Kj Pj Dj , donde Pj = Pk j=1
35
Kj
(7)
K Paso 7: Si k , es decir, hay pocos cl´ uster, ir al Paso 8. 2 Si k > 2K, es decir, tenemos demasiados cl´ uster, ir al paso Paso 11. Sino ir al Paso 13. Paso 8 (split): Calcular la desviaci´on est´andard de cada cluster: s 1 X = (x mj )2 , (j = 1, 2, . . . , k) j Kj x2!
(8)
j
donde mj y j , son el centro del cluster y la desviaci´on est´andard del cluster !j y Kj la cantidad de elementos en !j . Paso 9: Buscar la m´ axima desviaci´on est´andard max
Paso 10: Si para max
max
= max{
1,
max :
2, . . . ,
k}
(9)
se cumplen los siguientes:
> ⇥S ,
Dj > D, Kj > 2⇥M entonces dividir el cluster !j (cuyo centro es mj ) en dos nuevos centros de clase, mj y m+ j . Para ello, se adiciona ± a mj . Luego eliminar mj del conjunto de centros de cl´ uster y k = k + 1. Sino, ir al Paso 13. Paso 11 (merge): Computar la distancia entre cl´ uster de a pares Dij : Dij = DL (mi , mj ), (1 i, j k, i > j)
(10)
y ordernarlos de forma ascendente, Di1 j1 Di2 j2 . . . Dik jk Paso 12: Encontrar las P distancias que cumplan Dij < ⇥C y unir de a pares los centros de clase mi con mj como sigue debajo. Para l = 1, 2, . . . , P : m⇤ =
1 (Kil mil + Kjl mjl ) K il + K j l
(11)
Luego, eliminar mil y mjl . Adem´as agregar m⇤ al conjunto de centros de clase y por u ´ltimo k = k 1. Paso 13: Terminar el algoritmo si el n´ umero de iteraciones I es alcanzado. Sino, ir al Paso 2[8].
36
3.2. 3.2.1.
Filtros de Ruido Speckle Introducci´ on a los Filtros de Ruido de Speckle
Las im´ agenes de Radar de apertura sint´etica (SAR), se han vuelto una importante fuente de informaci´ on sobre la superficie de la Tierra y han sido ampliamente usadas en muchos campos, como la ecolog´ıa, hidrolog´ıa, geolog´ıa y oceanograf´ıa. Sin embargo, las im´agenes SAR, por su naturaleza, siempre exhiben una alteraci´ on en su calidad fotogr´afica. Esta modificaci´on, se le denomina “ruido speckle”. El ruido speckle se manifiesta en im´agenes con iluminaci´on coherente, como las de scanner ultras´onico, sonar y Radar de apertura sint´etica (SAR). Este fen´ omeno se desv´ıa del modelo cl´asico, en el que se supone que el ruido es Gaussiano e independiente de la se˜ nal y adicionado al verdadero valor (modelo aditivo de ruido gaussiano). El speckle, en cambio, es multiplicativo y no Gaussiano (en los formatos intensidad y amplitud), y dificulta la interpretaci´ on de las im´ agenes porque el “efecto de sal y pimienta” degenera la informaci´ on [11].
Figura 9: Imagen SAR ERS con ruido speckle [19]. El speckle es creado por las dispersiones de la iluminaci´on del Radar por ser demasiado peque˜ nas con respecto a longitud de onda como para ser capturadas individualmente. Dado que la se˜ nal de radar se dice que es coherente (es transmitida en un rango estrecho de onda), la energ´ıa dispersada tiende a ser reforzada (interferencia constructiva) o a ser suprimida (interferencia destructiva). Cuando picos altos de amplitud de onda coinciden, se crea un retorno con brillo; contrariamente, cuando picos de amplitud alto coinciden con los valores de picos de baja amplitud tienden a cancelarse entre s´ı, creando un retorno oscuro. Debido a que el speckle constituye una forma de ruido (porque no transmite informaci´ on u ´til), generalmente se procesan bien mediante procesos MultiLook
37
que iluminan la escena con frecuencias ligeramente diferentes que producen retornos independientes que luego pueden ser promediadas para reducir el efecto de la speckle. Estas dos estrategias de filtrado se caracterizan ya sea como no adaptativas o adaptativas. Los filtros no adaptativos requieren menos c´omputo ya que se aplica un solo filtro a toda la imagen de manera uniforme, mientras que los filtros adaptativos complejos se ajustan a las propiedades locales del entorno, con enfoque en la preservaci´on de los bordes y los l´ımites. Todas las estrategias, tienen la intenci´ on de extraer estimaciones precisas de la dispersi´on, pero siempre se corre el riesgo de eliminar la informaci´on original de la escena, as´ı que siempre hay que equilibrar los beneficios y las p´erdidas dentro de cada imagen en particular [10]. Por este motivo, para tener un mayor control del filtrado, se aborda utilizar filtros adaptativos. Siempre es deseable reducir o filtrar el ruido previamente a una interpretaci´ on para luego realizar el an´alisis de la misma. La mayor´ıa de los filtros para reducci´ on de speckle com´ unmente utilizados, tienen buenos resultados para suavizar el ruido. Aunque, las im´agenes resultantes siempre est´an sujetas a una previsible degradaci´ on de la resoluci´on tanto espacial como radiom´etrica. Este deterioro concluye en una p´erdida de informaci´on de la imagen en mayor o menor medida, dependiendo del comportamiento del filtro utilizado. La cantidad deseada de reducci´ on de speckle debe estar balanceada con la cantidad de detalle requerido para la escala espacial y para la naturaleza de la aplicaci´ on en cuesti´ on. Para interpretaciones a gran escala, los detalles minuciosos puede ser en muchos casos ignorados. As´ı, un filtrado amplio del ruido speckle y, consecuentemente, con p´erdida de detalles de la imagen puede ser aceptable o apropiado. Mientras que para aplicaciones donde los detalles finos y de alta resoluci´ on son necesarios, el filtro de speckle debe enfocarse en la preservaci´on de detalles. Este aspecto puede ser tan importante como a su vez la reducci´on del ruido speckle [9]. A continuaci´on se presentan una serie de filtros conocidos para la reducci´ on de ruido speckle.
38
3.2.2.
Modelo Matem´ atico
Se supone el modelo multiplicativo, I(t) = R(t) · u(t)
(12)
donde: t = (x, y), representan las coordenadas espaciales de la imagen, I(t), es el valor de intensidad observado en la imagen, R(t), es la reflectividad correspondiente al terreno, u(t), es ruido speckle multiplicativo estad´ısticamente independiente de R(t), con media u ¯ y varianza u2 . 3.2.3.
Filtro de Kuan
En el filtro desarrollado por Kuan, el modelo multiplicativo es transformado a un modelo aditivo dependiente de la se˜ nal. Luego se aplica una minimizaci´on al criterio de error cuadr´ atico medio (MMSE), obteniendo como resultado la misma ecuaci´ on que Lee, pero diferente funci´on de ponderaci´on. ˆ = I(t) · W (t) + I(t) ¯ · (1 R(t)
W (t))
(13)
La ecuaci´ on refleja una suma ponderada de los valores observados I(t) y la ¯ media local I(t). La funci´ on de ponderaci´on utilizada es: W (t) =
1
Cu2 /Cl2 (t) 1 + Cu2
(14)
donde Cu = u /¯ u es el coeficiente de variaci´on del ruido, est´ andard y la media del ruido respectivamente [4, 5]. 3.2.4.
u
yu ¯ son la desviaci´on
Filtro de Lee
En el filtro formulado por Lee, se aproxima el modelo multiplicativo con un modelo lineal. Luego tambi´en se minimiza el criterio de error cuadr´atico medio (MMSE) aplicado al modelo lineal dando como resultado: ˆ = I(t) · W (t) + I(t) ¯ · (1 R(t)
W (t))
(15)
donde la funci´ on de ponderaci´on es: Cu2 /Cl2 (t)
W (t) = 1
(16)
y Cu =
¯ I(t)
u
, Cl = ¯ u ¯ I(t)
(17)
son los coeficientes de variaci´ on del ruido speckle u(t) y de la imagen I(t) respectivamente [1, 2, 5].
39
3.2.5.
Filtro de Frost
El filtro de Frost difiere de la definici´on de Lee y de Kuan con respecto a que la reflectividad de la escena es estimada realizando una convoluci´on entre la imagen observada y la respuesta a un impulso de un sistema SAR. La respuesta a un impulso del sistema SAR, se obtiene realizando una minimizaci´on del error cuadr´ atico medio (MMSE) entre la imagen observada y el modelo de reflectividad de la escena que se asume como un proceso autoregresivo. El filtro resultante luego de algunas simplificaciones, queda como: KCl (t0 )|t|
m(t) = e
(18)
donde: K, es una constante que controla el factor de amortiguamiento o damping factor de la funci´ on de respuesta a un impulso, t0 , denota el pixel a ser filtrado, Cl (t0 ), denota el coeficiente de variaci´on dentro de la ventana con centro t0 . Se puede ver que si Cl (t0 ) tiene un valor peque˜ no, el filtro se comporta como un filtro de paso bajo, filtrando el ruido speckle, pero cuando Cl (t0 ) tiene un valor alto, el filtro tiene tendencia a preservar el valor original observado de la imagen [1, 3]. 3.2.6.
Filtro de Lee Mejorado
A la definici´ on inicial del filtro de Lee, solo se le agregan condiciones extras sobre los coeficientes de variaci´on [5]. En general, la imagen se la puede dividir en ´ areas de tres clases distintas: la primera clase, corresponde a las ´areas homog´eneas de la imagen donde s´ olo se debe aplicar un filtro de paso bajo para eliminar el speckle existente, la segunda clase, corresponde a las ´areas heterog´eneas de la imagen donde se debe reducir el speckle pero a su vez conservando la textura de la imagen, la u ´ltima clase, es donde se encuentran objetivos aislados donde se debe preservar el valor observado. 8 ¯ > : I(t0 ),
W (t0 )),
si CI (t0 ) Cu si Cu CI (t0 ) < Cmax si CI (t0 ) Ccmax
(19)
donde:
W (t0 ) = exp[ K(CI (t0 )
Cu )/(Cmax
CI (t0 ))]
K, es el factor de amortiguamiento, Cmax , es el m´ aximo coeficiente de variaci´on del ruido [1, 2, 5].
40
(20)
4. 4.1.
Implementaci´ on Computacional Qu´ e Hicimos y C´ omo lo Hicimos
La tarea de implementaci´ on que desarrollamos consisti´o en la construcci´on de una librer´ıa, escrita en Python, para el procesamiento de im´agenes SAR. Bautizamos a nuestra librer´ıa PyRadar. En la siguiente p´ agina se incluye un diagrama de ´arbol que muestra la estructura de archivos, carpetas y m´odulos de la librer´ıa PyRadar. Seguido del diagrama de ´ arbol, iremos describiendo, comentado y dando ejemplos de cada uno de los m´odulos.
41
4.2.
Estructura de Archivos y Carpetas de PyRadar
pyradar/ init .py classifiers init .py isodata.py kmeans.py comparator init .py comparator utils.py image comparator.py core init .py equalizers.py sar.py examples filters init .py frost.py kuan.py lee.py lee enhanced.py mean.py median.py utils.py simulate ExampleImages DemoSet gamma img layer.png gamma noisy img.png k img layer.png k noisy img.png demo K.png demo gamma.png init .py image simulator.py utils init .py sar debugger.py statutils.py system info.py timeutils.py
42
4.3.
Descripci´ on de los m´ odulos
4.3.1.
core, manejo de im´ agenes SAR y algoritmos de ecualizaci´ on
Este m´ odulo permite el manejo de im´agenes SAR y adem´as posee algoritmos para su ecualizaci´ on.
Las funciones que contiene este m´odulo son: create_dataset_from_path(image_path) Desde el path de una imagen SAR, crea una estructura de datos, un dataset, para manejar su informaci´on. Esta estructura de datos, proviene de una librer´ıa externa llamada Gdal 1 . get_band_from_dataset(dataset) Obtiene la u ´nica banda de utilidad a nuestros fines del “dataset” pasado como argumento. Es importante notar que las im´agenes que utilizamos poseen s´ olo una banda, dado que son im´agenes de radar, en blanco y negro. get_band_min_max(band) Retorna una tupla de Python con el m´aximo y m´ınimo de la banda “band ”. read_image_from_band(band, xoff, yoff, win_xsize, win_ysize) Lee la banda convirti´endola en un arreglo bidimensional(o matriz ) de numpy, facilitando as´ı el manejo de la informaci´on en un formato simple y m´ as natural de manejar. Los par´ ametros significan lo siguiente: • “xo↵ ”, “yo↵ ”, indican con qu´e o↵set sobre el eje x y el eje y se deben leer los datos desde la imagen, • “win xsize”, “win ysize”, indican el tama˜ no de la ventana a leer de la imagen en largo y ancho. get_geoinfo(dataset, cast_to_int) Esta funci´ on extrae la informaci´on de georreferenciaci´on de un “dataset”, con “cast to int” indicando si los datos de georreferenciaci´on se encuentran como strings o como n´ umeros crudos. De momento, no es utilizada por ninguno de los algoritmos, pero si en un futuro se quisiese extender PyRadar con funcionalidades que requieran el uso de georreferenciaci´ on, la base est´a preparada. save_image(img_dest_dir, filename, img) 1 Gdal
es una librer´ıa para leer y escribir datos provenientes de rasters geoespaciales, y est´ a liberada bajo licencia MIT por la Open Source Geospatial Foundation.
43
Esta funci´ on se encarga de guardar una imagen “img” representada por un arreglo de numpy en le directorio “img dest dir’ ’ con nombre de archivo “filename”. Es importante destacar que “img” debe poseer valores en el rango [0 : 255] para que el procedimiento resulte exitoso. Si la imagen llegase a poseer valores fuera de este rango, se deber´an normalizar los valores previamente utilizando los algoritmos de equalizaci´on provistos por este m´odulo. equalization_using_histogram(img) Esta funci´ on normaliza los valores de “img” al rango [0 : 255], utilizando como procedimiento intermedio “equalize histogram”. Dicho algoritmo est´ a basado en [7]. equalize_histogram(img, histogram, cfs) Dado el histograma y la funci´on de distribuci´on acumulada de “img”, este procedimiento normaliza los valores al rango [0 : 255]. Cabe notar que esta funci´ on es utilizada internamente por “equalization using histogram”. naive_equalize_image(img, input_range, output_range) Esta funci´ on es una implementaci´on sencilla y sin optimizaciones que sirve para normalizar im´ agenes desde el rango “input range” al rango “output range”.
44
Ejemplo de uso: 1 2 3 4 5 6
# -*- coding: utf-8 *-* from pyradar.core.sar import from pyradar.core.sar import from pyradar.core.sar import from pyradar.core.sar import from pyradar.core.sar import
create_dataset_from_path get_band_from_dataset get_geoinfo read_image_from_band save_image
7 8
from pyradar.core.equalizers import equalization_using_histogram
9 10 11 12
IMAGE_PATH = "./img_sar/DAT_01.001" IMG_DEST_DIR = "."
13 14 15 16 17 18 19
# create dataset dataset = create_dataset_from_path(IMAGE_PATH) # get band from dataset band = get_band_from_dataset(dataset) # get geo info from dataset geoinfo = get_geoinfo(dataset, cast_to_int=True)
20 21 22 23
#usually both values are zero xoff = geoinfo[’xoff’] yoff = geoinfo[’yoff’]
24 25 26 27 28
# window size in coord x win_xsize = 512 # window size in coord y win_ysize = 512
29 30
image = read_image_from_band(band, xoff, yoff, win_xsize, win_ysize)
31 32 33 34 35
#equalize img to 0:255 image_eq = equalization_using_histogram(image) # save img in current directory save_image(IMG_DEST_DIR, "image_sar", image_eq) Este es el procedimiento est´ andar para abrir, leer y guardar una imagen SAR utilizando PyRadar. A continuaci´on, en los ejemplos de los dem´as m´odulos, omitiremos estos pasos y haremos referencia a los mismos como “pasos b´ asicos de lectura” desde los imports hasta (inclusive) la llamada del procedimiento “read image from band ”.
45
El siguiente ejemplo ilustra como utilizar “naive equalize image”. Para ello se deben seguir los “pasos b´ asicos de lectura”, y luego agregar el siguiente c´odigo:
32 33 34 35 36 37 38 39
# get actual range input_range = image.min(), image.max() # set new range output_range = 0, 255 # equalize image image_eq = naive_equalize_image(image, input_range, output_range) # save image in current directory save_image(IMG_DEST_DIR, "image_sar", image_eq)
46
4.3.2.
classifiers, algoritmos de clasificaci´ on supervisada y no supervisada
El prop´ osito de este m´ odulo es ofrecer algoritmos de clasificaci´on. Si bien los algoritmos desarrollados en este trabajo son Isodata y Kmeans, el m´odulo puede ser extendido con algoritmos nuevos o incluso versiones modificadas de los existentes. Las definiciones matem´ aticas de Isodata y Kmeans se encuentran en el cap´ıtulo 3.1 junto con sus ventajas y desventajas.
El m´ odulo cuenta con las siguientes funciones: kmeans_classification(img, k=K_CLASSES, iter_max=ITER_MAX) Computa el algoritmo Kmeans con “k ” clases fijas y con n´ umero m´aximo de iteraciones “iter max ”. De no ser especificados, se tomar´an los siguientes valores por defecto k = 5 y iter max = 100. Esta versi´ on de Kmeans utiliza la distancia eucl´ıdea definida como sigue: p DL ((x1 , x2 , · · · , xn ), (y1 , y2 , · · · , yn )) = (x1 y1 )2 + (x2 y2 )2 + · · · + (xn isodata_classification(img, parameters=None)
Calcula el algoritmo Isodata con los siguientes par´ametros por defecto si los mismos no fuesen especificados “params” en un diccionario Python: • K = 5, la cantidad deseada de cl´ uster,
• I = 100, n´ umero m´ aximo de iteraciones que realiza el algoritmo, • P = 4, m´ axima cantidad de pares de cl´ uster a “unir ”,
• T HET A S = 1, valor de la desviaci´on est´andard entre cl´ usters para decidir si se debe “dividir ” un cl´ uster, en dos nuevos cl´ uster, • T HET A M = 10, es el n´ umero m´ınimo de elementos por cluster; de tener el cluster menos elementos que T HET A M = 10 ser´a eliminado por el algoritmo, • T HET A C = 20, distancia entre dos cl´ uster para decidir si “unirlos”, es decir, el algoritmo eliminar´a los dos cl´ uster y formar´a un cluster nuevo, • T HET A O = 0,05, porcentaje de cambio de los centros de los cl´ uster entre cada iteraci´ on. Si los centros de los cl´ uster presentan un cambio menor a T HET A O entonces el algoritmo se detiene. • k = 5, cantidad de cl´ uster iniciales con que arranca el algoritmo, tomando como centros “k ” valores aleatorios dentro de la imagen.
47
yn ) 2
Ejemplo de uso de Kmeans: Para correr el algoritmo Kmeans, primero se deben seguir los “pasos b´ asicos de lectura” en la secci´ on 4.3.1 de la p´agina 45 para poder as´ı obtener la imagen a clasificar en la variable “image”.
32 33
# this should be placed at the top with all the imports from pyradar.classifiers.kmeans import kmeans_classification
34 35 36 37 38 39 40
# number of clusters k= 4 # max number of iterations iter_max = 1000 # run K-Means class_image = kmeans_classification(image, k, iter_max)
41 42 43 44 45 46 47 48 49
# equalize class image to 0:255 class_image_eq = equalization_using_histogram(class_image) # save it save_image(IMG_DEST_DIR, "class_image_eq", class_image_eq) # also save original image image_eq = equalization_using_histogram(image) # save it save_image(IMG_DEST_DIR, "image_eq", image_eq)
Ejemplo de uso de Isodata: Para correr el algoritmo Isodata, al igual que Kmeans se necesitan ejecutar los “pasos b´ asicos de lectura” previamente.
32 33
# this should be placed at the top with all the imports from pyradar.classifiers.isodata import isodata_classification
34 35 36
params = {"K": 15, "I" : 100, "P" : 2, "THETA_M" : 10, "THETA_S" : 0.1, "THETA_C" : 2, "THETA_O" : 0.01}
37 38 39
# run Isodata class_image = isodata_classification(img, parameters=params)
40 41 42 43 44 45 46 47 48
# equalize class image to 0:255 class_image_eq = equalization_using_histogram(class_image) # save it save_image(IMG_DEST_DIR, "class_image_eq", class_image_eq) # also save original image image_eq = equalization_using_histogram(image) # save it save_image(IMG_DEST_DIR, "image_eq", image_eq)
48
4.3.3.
comparator, es el m´ odulo donde se encuentran rutinas para comparar im´ agenes
Este m´ odulo permite realizar comparaciones entre im´agenes y la obtenci´on de m´etricas para medir su simitud.
Funciones y Clases del M´ odulo class ComparatorException(Exception) Custom excepction class definida para manejar las situaciones de excepci´ on del m´ odulo de comparaci´on de im´agenes. class BaseImageComparator(object) Clase base del comparador de im´agenes que agrupa los m´etodos y atributos del comparador de im´ agenes, el m´etodo de inicializaci´on y el de validaci´on de argumentos de inicializaci´on. Posee los siguientes m´etodos: • def __init__(self, image_1, image_2)
Funci´ on de inicializaci´on. Internamente, antes de inicializar los atributos, llama a la funci´on de validaci´on de argumentos.
• def validate_images_are_comparable(self, image1, image2)
Funci´ on que realiza m´ ultiples validaciones de las dos im´agenes pasadas como par´ ametros. Retorna True si las dos im´agenes son comparables; eleva una excepci´on si alguna de las im´agenes no est´a definida o si por alguna raz´on (ejemplo: no tienen las mismas dimensiones) no son comparables.
class ImageComparator(BaseImageComparator) Clase que hereda de BaseImageComparator y que la extiende con la implementaci´ on de las funciones de comparaci´on de im´agenes. Posee los siguientes m´etodos: • def compare_by(self, strategy, params)
Entry point del las funciones de comparaci´on. Toma como argumento: strategy : Estrategia o algoritmo con el cual realizar la comparaci´ on. Las opciones son: ⇧ rmse1 : Una de las dos implemantaciones que realizamos del c´ alculo de Root Mean Square Error(RMSE). ⇧ rmse2 : Segunda de las dos implementaciones que realizamos del c´ alculo de Root Mean Square Error(RMSE). ⇧ mae: Computa el Mean Square Error(MAE). ⇧ pearson: Calcula la correlaci´on de Pearson. params: Diccionario de par´ametros necesarios para el c´alculo de cada algoritmo de comparaci´on.
49
Ejemplo de uso: 1 2 3
from pyradar.comparator.image_comparator import ImageComparator from pyradar.examples.misc_resources.sample_matrixes import (numpy_image, numpy_image1)
4 5
im = ImageComparator(numpy_image, numpy_image1)
6 7 8 9 10
print print print print
’rmse1: ’, im.compare_by(’rmse1’, None) ’rmse2: ’, im.compare_by(’rmse2’, None) ’mae: ’, im.compare_by(’mae’, None) ’pearson: ’, im.compare_by(’pearson’, None)
50
4.3.4.
filters, filtros de ruido
En este m´ odulo se encuentran los siguientes filtros de ruido speckle: Frost Kuan Lee Lee Mejorado Los mismos siguen las definiciones matem´aticas de la secci´on 3.2. Adem´as, se encuentran tambi´en las implementaciones de los filtros cl´asicos de media y de mediana. Como PyRadar tiene entre sus objetivos seguir expandi´endose y creciendo, este m´ odulo puede ser extendido con nuevos filtros. Como complemento a estos algoritmos, se desarrollar´o una serie funciones que verifican la consistencia de los algoritmos en tiempo de ejecuci´on. La finalidad de esto es que, al extender el m´odulo con nuevos filtros, el desarrollador no necesite escribir nuevo c´ odigo para verificar estas condiciones en sus algoritmos para verificar consistencia. A continuaci´ on se detallan las funciones del m´odulo:
frost_filter(img, damping_factor, win_size) Esta funci´ on implementa el filtro de Frost sobre una imagen img, tomando como argumentos el damping factor y el tama˜ no de ventana win size. Si estos argumentos no fueran especificados, se toma por defecto: damping factor=2.0 y win size=3. kuan_filter(img, win_size, cu) Esta funci´ on aplica el filtro de Kuan sobre una imagen img tomando como argumentos el tama˜ no de ventana win size y coeficiente de variaci´on del ruido cu. De no ser especificados los siguientes valores se toman por defecto: win size=3 y cu=0.25 lee_filter(img, win_size, cu) Esta funci´ on aplica el filtro de Lee sobre una imagen img, tomando como argumentos el tama˜ no de ventana win size y coeficiente de variaci´on del ruido cu. De no ser especificados, se toman los mismos valores por defecto que en el filtro de Kuan.
51
lee_enhanced_filter(img, win_size, k, cu, cmax) Esta funci´ on aplica el filtro de Lee Mejorado con los siguientes argumentos: • la imagen img,
• el tama˜ no win size (por defecto 3),
• k el factor de amortiguamiento (por defecto 1.0),
• coeficiente de variaci´on del ruido cu (por defecto 0.523),
• coeficiente de variaci´on m´aximo en la imagen cmax (por defecto 1.73). mean_filter(img, win_size) Esta funci´ on ejecuta un filtro de paso bajo cl´ asico, como lo es el filtro de media. Los argumentos que toma son la imagen img y el tama˜ no de la ventana win size. Por defecto win size toma el valor de 3. median_filter(img, win_size) Esta funci´ on ejecuta el segundo filtro de paso bajo cl´ asico que contiene el m´ odulo: el filtro de mediana. Los argumentos de este filtro son la imagen img y el tama˜ no de ventana win size, tomando por defecto 3 para este u ´ltimo. Funciones desarrolladas para verificar consistencia de los filtros en tiempo de ejecuci´ on assert_window_size(win_size) Verifica que el tama˜ no de ventana sea m´ ultiplo de 3 y positivo. De no cumplirse la condici´ on, levanta una excepci´on de Python. assert_indices_in_range(width, height, xleft, xright, yup, ydown) Verifica que los ´ındices de la ventana deslizante se encuentren dentro de los valores normales. Es decir, que siempre se cumpla lo siguiente: (0 xlef t ^ xright width ^ 0 yup ^ ydown height) De no ser cierta la expresi´ on booleana anterior, se levanta una excepci´on de Python.
52
Ejemplo de uso de los filtros: Para correr los algoritmos de los filtros antes mencionados se necesitan ejecutar los “pasos b´ asicos de lectura”, para tener as´ı la imagen img a usar en la variable “image”.
32 33 34 35 36 37
from from from from from from
pyradar.filters.frost import frost_filter pyradar.filters.kuan import kuan_filter pyradar.filters.lee import lee_filter pyradar.filters.lee_enhanced_filter import lee_enhanced_filter pyradar.filters.median import median_filter pyradar.filters.mean import mean_filter
38 39 40 41 42 43 44 45 46 47 48 49 50 51
# filters parameters # window size winsize = 9 # damping factor for frost k_value1 = 2.0 # damping factor for lee enhanced k_value2 = 1.0 # coefficient of variation of noise cu_value = 0.25 # coefficient of variation for lee enhanced of noise cu_lee_enchanced = 0.523 # max coefficient of variation for lee enhanced cmax_value = 1.73
52 53 54 55 56 57 58 59 60 61 62 63 64 65
# frost filter image_frost = frost_filter(image, damping_factor=k_value1, win_size=winsize) # kuan filter image_kuan = kuan_filter(image, win_size=winsize, cu=cu_value) # lee filter image_lee = lee_filter(image, win_size=winsize, cu=cu_value) # lee enhanced filter image_lee_enhanced = lee_enhanced_filter(image, win_size=winsize, k=k_value2, cu=cu_lee_enchanced, cmax=cmax_value) # mean filter image_mean = mean_filter(image, win_size=winsize) # median filter image_median = median_filter(image, win_size=winsize)
53
4.3.5.
simulate, en este m´ odulo se proveen algoritmos para simular im´ agenes por capas.
En este m´ odulo se proveen algoritmos para simular im´agenes por capas: Capa de imagen original: • Imagen con distribuci´on (Gamma): Recomendable para simular zonas m´ as bien homog´eneas, como pasturas. • Imagen con distribuci´on K: Recomendable para simular zonas m´as bien inhomog´eneas, como ciudades. Capa de ruido, con distribuci´on
2
.
Capa de imagen con ruido: Combinaci´on multiplicativa de las primeras dos capas. Funciones y Clases del M´ odulo class SimulatorException(Exception) Custom excepction class definida para manejar las situaciones de excepci´ on del m´ odulo de simulaci´on de im´agenes. • def __init__(self, width, height)
Inicializa las matrices numpy para las capas de ruido, imagen origina e imagen con ruido, con las dimensiones indicadas. ‘width’ refiere a las columnas, ‘height’ a las filas.
• def generate_image_layer(self, distribution, params)
Genera la capa de imagen original(simulada), sin ruido. Almacena el resultado en self.image layer. Si se utiliza una distribuci´on ‘gamma’, params debe ser un diccionario con las claves ‘shape’ and ‘scale’, con valores ambos mayores que cero. En el caso de la distribuci´on k, se debe proveer en el diccionario de par´ ametros la clave ‘df’(por degrees of freedom, grados de libertad), con valor tambi´en mayor que cero. La funci´ on
• generate image layer() utiliza internamente las siguientes dos funciones: def _generate_img_with_gama_dist(self, params) def _generate_img_with_k_dist(self, params) • def generate_noise_layer(self)
Genera la capa de ruido(simulada). Almacena el resultado en self.noisy image. Para esto, combina la capa de imagen sin ruido y la capa de ruido, previamente generadas.
54
• def generate_noisy_layer(self)
Genera la capa de imagen con ruido(simulada). Almacena el resultado en self.noise layer. La funci´ on generate noisy layer() utiliza internamente la funci´on generate image with chisquare dist().
• def export_image_layer(self, layer_name, filename, path_to) Exporta un archivo de imagen para la capa ‘layer name’ al path indicado por ‘path to’ con el nombre ‘filename’. Internamente, se utiliza la funci´ on get layer securely() para obtener la capa solicitada con seguridad, arrojando(haciendo raise) una excepci´on que informe detalladamente al usuario si el nombre de capa es incorrecto.
• def plot_layer_histogram(self, layer_name, filename)
Genera un plot de la capa ‘layer name’ y lo guarda en un archivo llamado ‘filename’ en el directorio actual. Al igual que en la funci´on export image layer(), internamente, se utiliza la funci´on get layer securely() para obtener la capa solicitada con seguridad, arrojando(haciendo raise) una excepci´ on que informe detalladamente al usuario si el nombre de capa es incorrecto.
55
Ejemplo de uso: 1 2 3
from pyradar.simulate.image_simulator import ImageSimulator from pyradar.utils.timeutils import Timer pylab.close()
4 5 6
timer = Timer() width, height = 2000, 2000
7 8 9 10
gamma_ims = ImageSimulator(width, height) k_ims = ImageSimulator(width, height) noise_layer_ims = ImageSimulator(width, height)
11 12 13 14
gamma_params = {’scale’: 2.0, ’shape’: 3.0} k_params = {’mean’: 2.0, ’shape’: 2.0} noise_layer_params = {’df’: 3}
15 16 17 18
gamma_ims.generate_image_layer(distribution=’gamma’, params=gamma_params) k_ims.generate_image_layer(distribution=’k’, params=k_params) noise_layer_ims.generate_noise_layer(distribution=’chisquare’, params=noise_layer_params)
19 20 21 22 23 24
# Make some noise! gamma_ims.noise_layer = noise_layer_ims.noise_layer k_ims.noise_layer = noise_layer_ims.noise_layer gamma_ims.generate_noisy_layer() k_ims.generate_noisy_layer()
25 26 27 28 29 30 31 32 33 34 35 36
timer.calculate_time_elapsed(print_value=True) # Export the files: gamma_ims.export_image_layer(layer_name=’image_layer’, filename=’gamma_img_layer’, path_to=’.’) k_ims.export_image_layer(layer_name=’image_layer’, filename=’k_img_layer’, path_to=’.’) gamma_ims.export_image_layer(layer_name=’noisy_image’, filename=’gamma_noisy_img’, path_to=’.’) k_ims.export_image_layer(layer_name=’noisy_image’, filename=’k_noisy_img’, path_to=’.’) timer.calculate_time_elapsed(print_value=True)
37 38 39 40 41 42
# Make a plot: print ’Making a plot to "plot_img.png":’ pylab.close() gamma_ims.plot_layer_histogram(layer_name=’image_layer’, filename=’plot_gamma_img’) k_ims.plot_layer_histogram(layer_name=’image_layer’, filename=’plot_k_img’)
43 44 45
timer.stop_timer() timer.calculate_time_elapsed(print_value=True)
56
4.3.6.
Im´ agenes de ejemplo generadas con el simulador de im´ agenes
Capa de imagen original, distribuci´ on Par´ ametros: scale: 2.0 shape: 3.0
Figura 10: Simulaci´ on de imagen original con distribuci´on Gamma.
57
Capa de imagen original, distribuci´ on K Par´ ametros: mean: 2.0 shape: 2.0
Figura 11: Simulaci´ on de imagen original con distribuci´on K.
58
Imagen
, con ruido de distribuci´ on
2
Par´ ametros: df: 3
Figura 12: Imagen con distribuci´on Gamma con ruido de distribuci´on
59
2
.
Imagen K, con ruido de distribuci´ on
2
Par´ ametros: df: 3
Figura 13: Imagen con distribuci´on K con ruido de distribuci´on
60
2
.
Plot de la imagen con distribuci´ on Par´ ametros: scale: 2.0 shape: 3.0
Figura 14: Plot de la imagen con distribuci´on .
61
Plots de la imagen con distribuci´ on K juntos
y de la imagen con distribuci´ on
Par´ ametros distribuci´ on : scale: 2.0 shape: 3.0 Par´ ametros distribuci´ on K: df: 3
Figura 15: Plots de la imagen con distribuci´on K juntos.
62
y de la imagen con distribuci´on
4.3.7.
utils, es un m´ odulo en que se agrupan utilidades y herramientas generales
timer, un peque˜ no timer para cronometrar el tiempo de ejecuci´ on de porciones de c´ odigo Python Se desarroll´ o este peque˜ no m´ odulo con el fin de cronometrar el “tiempo de pared” de ejecuci´ on de algunas funciones. El tiempo de pared de ejecuci´on es el tiempo total que un c´ alculo permanece en ejecuci´on. Se le llama de “pared” porque dentro del sistema operativo la ejecuci´on de un proceso tambi´en acarrea otra operaciones b´asicas adem´as de la algoritmia programada. Operaciones como cambios de contexto del sistema operativo, carga y descarga de librer´ıas, volcados de datos a disco y otras operaciones agregan tiempo extra a la medici´ on. Si bien la medici´on no es de alta precisi´on, su valor es servir como medici´ on de referencia de tiempo de ejecuci´on para realizar optimizaciones. Ejemplos de uso: A diferencia de las dem´ as utilidades antes mencionadas, esta utilidad la utilizamos dentro del c´ odigo mismo.
1
from pyradar.utils.timeutils import Timer
2 3 4 5 6 7 8 9 10
# crea y arranca el timer simple_timer = Timer() # procedimiento que queremos medir result = function(arg1, arg2) # paramos el timer simple_timer.stop_timer() #imprimimos los resultados y los guardamos en diff diff = simple_timer.calculate_time_elapsed(print_value=True)
Sar Debugger El prop´ osito de este m´ odulo de PyRadar es agrupar funcionalidades y herramientas para realizar tareas de debugging sobre algoritmos que manipulen im´agenes SAR. De esta forma, el m´ odulo ir´a satisfaciendo las necesidades de la comunidad con nuevos features. De momento s´olo posee una funci´on, take snapshot(). def take_snapshot(img) es una funci´ on que toma una fotograf´ıa instant´ anea de la imagen img y la guarda en el disco. El prop´osito de esta funci´on, es poder exportar capturas de los algoritmos de clasificaci´on a medida que van evolucionando en el c´ omputo.
63
Ejemplo de uso: 1
from pyradar.utils.sar_debugger import take_snapshot
2 3
MAX_ITER = 1000
4 5 6 7
for iter in xrange(0, MAX_ITER): image = some_algorithm(image) take_snapshot(image, iteration_step=iter)
system info, obtiene informaci´ on del Sistema Operativo, Hardware y Software Se desarroll´ o un m´ odulo que permite obtener informaci´on del Sistema Operativo de manera simple y detallada. El objetivo de este m´odulo es que cuando alg´ un usuario o desarrollador tenga alg´ un problema con la librer´ıa pyradar, ´este pueda comparar la informaci´ on espec´ıfica sobre su Sistema Operativo, Hardware y Software, con el fin de descartar(o confirmar) la posibilidad de que su problema provenga de estas fuentes. Ejemplo de uso: 1 2
from pyradar.utils.system_info import get_system_info from pyradar.utils.system_info import print_info
3 4 5
info = get_system_info() print_info(info)
64
statutils, utilidades estad´ısticas Este m´ odulo provee algunas funciones estad´ısticas que pueden llegar a ser de utilidad para la comunidad de procesamiento de im´agenes, incluso m´as all´a del contexto espec´ıfico de la librer´ıa. Al igual que el m´odulo Sar Debugger, su objetivo tambi´en es ser extendido por la comunidad a medida que la necesidad lo demande. def compute_cfs() Recibiendo como argumento un histograma generado con la librer´ıa numpy, esta funci´ on genera una tabla con todas las frecuencias acumuladas. def calculate_pdf_for_pixel() Calcula la probabilidad de que un valor en particular aparezca en la imagen, donde la probabilidad est´a dada como: la cantidad de ocurrencias efectivas del valorcantidad total de elementos. def calculate_cdf_for_pixel() Calcula el valor de un pixel en la funci´on distribuci´on acumulada. def compute_cdfs() Esta funci´ on computa todas la probabilidades de la distribuci´on de frecuencia acumulada de cada valor en la imagen.
65
Ejemplos de uso: 1 2 3 4 5 6 7 8 9 10 11
import numpy as np from pyradar.utils.statutils import from pyradar.utils.statutils import from pyradar.utils.statutils import from pyradar.utils.statutils import arr = np.array([31, 49, 19, 62, 24, 54, 26, 57, 37, 43, 39, 52, 35, 51, 63, max_value = arr.max() min_value = arr.min() start, stop, step = int(min_value),
compute_cfs calculate_pdf_for_pixel calculate_cdf_for_pixel compute_cdfs 45, 23, 51, 55, 60, 40, 35, 65, 18, 41, 50, 56, 4, 54, 42])
int(max_value + 2), 1
12 13 14 15 16 17 18 19
histogram, bin_edge = np.histogram(arr, xrange(start, stop, step) compute_cfs(histogram) >>> array([ 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 2, 3, 3, 3, 3, 4, 5, 5, 6, 6, 6, 6, 6, 7, 7, 7, 7, 9, 9, 10, 10, 11, 12, 13, 14, 15, 15, 16, 16, 16, 16, 17, 18, 20, 21, 21, 23, 24, 25, 26, 26, 26, 27, 27, 28, 29, 29, 30])
20 21 22 23 24 25 26
calculate_pdf_for_pixel(arr, histogram, bin_edge, 54) >>> 0.066666666666666666 calculate_pdf_for_pixel(arr, histogram, bin_edge, 20) >>> 0.0 calculate_pdf_for_pixel(arr, histogram, bin_edge, 18) >>> 0.033333333333333333
27 28 29 30 31
calculate_cdf_for_pixel(arr, histogram, bin_edge, 4) >>> 0.033333333333333333 calculate_cdf_for_pixel(arr, histogram, bin_edge, 50) >>> 0.59999999999999998
32 33 34 35 36 37 38 39 40 41 42 43 44 45 46
compute_cdfs(arr, histogram, bin_edge) >>> array([ 0.03333333, 0.03333333, 0.03333333, 0.03333333, 0.03333333, 0.03333333, 0.03333333, 0.03333333, 0.03333333, 0.1 , 0.1 , 0.1 , 0.16666667, 0.16666667, 0.2 , 0.2 , 0.2 , 0.23333333, 0.23333333, 0.3 , 0.3 , 0.36666667, 0.4 , 0.43333333, 0.5 , 0.53333333, 0.53333333, 0.56666667, 0.6 , 0.66666667, 0.76666667, 0.8 , 0.83333333, 0.86666667, 0.9 , 0.9 , 0.96666667, 1. ])
66
0.03333333, 0.03333333, 0.03333333, 0.1 , 0.2 , 0.23333333, 0.33333333, 0.46666667, 0.53333333, 0.7 , 0.86666667, 0.93333333,
0.03333333, 0.03333333, 0.06666667, 0.13333333, 0.2 , 0.23333333, 0.33333333, 0.5 , 0.53333333, 0.7 , 0.86666667, 0.96666667,
4.3.8.
examples, ejemplos de uso de PyRadar
Este m´ odulo aloja todos los ejemplos mencionados en este trabajo de Tesis. pyradar/ init .py classifiers comparator core examples init .py example comparator.py example core.py example core2.py example filtros.py example isodata.py example kmeans.py example sar debugger.py example simulate.py example statutils.py example system info.py example timer.py misc resources compare images.py examples.py histogram.py kmean example.py sample matrixes.py filters simulate utils
67
4.4. 4.4.1.
Apertura y disponibilidad del c´ odigo. Apertura del C´ odigo
Como mencionamos en el cap´ıtulo 1.2, buena parte del prop´osito de nuestro proyecto de tesis es el de devolver algo a la sociedad en base a la educaci´on que hemos recibido en nuestros a˜ nos de carrera universitaria. Es por esto que la librer´ıa para procesamiento de im´agenes SAR y todo el c´odigo que escribimos ser´ a liberado con licencia LGPL 2 . 4.4.2.
Disponibilidad del C´ odigo
Repositorio Oficial
Figura 16: Logotipo de PyRadar. El c´ odigo estar´ a p´ ublicamente disponible en el siguiente repositorio GIT: https://github.com/PyRadar/pyradar All´ı mismo se podr´ an reportar bugs, comentar sobre nuevas features que fueran de inter´es de la comunidad para implementar y tambi´en otros desarrolladores podr´ an colaborar aportando c´ odigo o uni´endose al proyecto como desarrolladores.
2 Se adjuntan los textos completos de las licencias GPLv3 y LGPL, tanto en formato LATEX como compilados en formato PDF.
68
Python Package Index(PyPi) Para facilitar la instalaci´ o de PyRadar a los usuarios y para proceder de la manera que se acostumbra para con paquetes Python, cuando liberemos PyRadar tambi´en lo empaquetaremos y pondremos a disposici´on p´ ublica en el Python Package Index(PyPi)3 [22]. De este modo, se podr´a instalar PyRadar f´acilmente con alguno de los siguientes comandos: easy_install-2.7 pyradar o bien: pip install pyradar
3 PyPI - the Python Package Index: El Python Package Index es un repositorio de software Python. Actualmente cuenta con 20770 paquetes p´ ublicamente disponibles[22].
69
5.
Resultados
En todos los casos a continuaci´on expuestos se utiliz´o como base un raster ERS SAR, gentileza de CONAE. La imagen original tiene 8000 por 8200 p´ıxeles, por lo cual utilizamos una ventana y no la imagen para los siguientes resultados. La ventana utilizada fue la siguiente: Tama˜ no: 512x512 p´ıxeles. xo↵set: 4000. yo↵set: 4000. Las tiras de im´ agenes de resultados que se muestran, por razones evidentes de tama˜ no, se muestran reducidas. Asimismo, se adjuntan y ponen a disposici´on todas las im´ agenes y tiras de im´agenes en tama˜ no real que obtuvimos.
70
5.1.
Filtros
5.1.1.
Mean
Imagen 1 • window size: 3
• Run took: 0 mins, 3 secs, 359130 microsecs. Imagen 2 • window size 5
• Run took: 0 mins, 3 secs, 515199 microsecs. Imagen 3 • window size 7
• Run took: 0 mins, 3 secs, 610979 microsecs. Imagen 4 • window size 11
• Run took: 0 mins, 3 secs, 547104 microsecs.
Figura 17: Mean filter
71
5.1.2.
Median
Imagen 1 • window size: 3
• Run took: 0 mins, 6 secs, 163640 microsecs. Imagen 2 • window size 5
• Run took: 0 mins, 6 secs, 541002 microsecs. Imagen 3 • window size 7
• Run took: 0 mins, 6 secs, 431973 microsecs. Imagen 4 • window size 11
• Run took: 0 mins, 6 secs, 847116 microsecs.
Figura 18: Median filter
72
5.1.3.
Frost
damping factor: 0.1 • Imagen A.1
window size: 3 Run took: 0 mins, 22 secs, 162560 microsecs.
• Imagen A.2
window size 5 Run took: 0 mins, 22 secs, 545563 microsecs.
• Imagen A.3
window size 7 Run took: 0 mins, 22 secs, 874566 microsecs.
• Imagen A.4
window size 11 Run took: 0 mins, 23 secs, 284587 microsecs.
Figura 19: Frost: Im´agenes A.1, A.2, A.3 y A.4.
73
damping factor: 1 • Imagen B.1
window size: 3 Run took: 0 mins, 22 secs, 281812 microsecs.
• Imagen B.2
window size 5 Run took: 0 mins, 22 secs, 421805 microsecs.
• Imagen B.3
window size 7 Run took: 0 mins, 23 secs, 681723 microsecs.
• Imagen B.4
window size 11 Run took: 0 mins, 27 secs, 670154 microsecs.
Figura 20: Frost: Im´agenes B.1, B.2, B.3 y B.4.
74
damping factor: 10 • Imagen C.1
window size: 3 Run took: 0 mins, 22 secs, 791972 microsecs.
• Imagen C.2
window size 5 Run took: 0 mins, 22 secs, 571320 microsecs.
• Imagen C.3
window size 7 Run took: 0 mins, 23 secs, 106602 microsecs.
• Imagen C.4
window size 11 Run took: 0 mins, 24 secs, 74359 microsecs.
Figura 21: Frost: Im´agenes C.1, C.2, C.3 y C.4.
75
5.1.4.
Lee
cu = 0.10 • Imagen A.1
window size: 3 Run took: 0 mins, 11 secs, 764738 microsecs.
• Imagen A.2
window size 5 Run took: 0 mins, 12 secs, 954831 microsecs.
• Imagen A.3
window size 7 Run took: 0 mins, 12 secs, 882499 microsecs.
• Imagen A.4
window size 11 Run took: 0 mins, 12 secs, 731251 microsecs.
Figura 22: Lee: Im´agenes A.1, A.2, A.3 y A.4.
76
cu=0.20 • Imagen B.1
window size: 3 Run took: 0 mins, 11 secs, 986142 microsecs.
• Imagen B.2
window size 5 Run took: 0 mins, 11 secs, 923538 microsecs.
• Imagen B.3
window size 7 Run took: 0 mins, 12 secs, 78370 microsecs.
• Imagen B.4
window size 11 Run took: 0 mins, 12 secs, 869182 microsecs.
Figura 23: Lee: Im´agenes B.1, B.2, B.3 y B.4.
77
cu=0.30 • Imagen C.1
window size: 3 Run took: 0 mins, 11 secs, 986142 microsecs.
• Imagen C.2
window size 5 Run took: 0 mins, 11 secs, 923538 microsecs.
• Imagen C.3
window size 7 Run took: 0 mins, 12 secs, 78370 microsecs.
• Imagen C.4
window size 11 Run took: 0 mins, 12 secs, 869182 microsecs.
Figura 24: Lee: Im´agenes C.1, C.2, C.3 y C.4.
78
5.1.5.
Lee Enhanced
Para todos los casos: K DEF AU LT = 1,0 CM AX DEF AU LT = 1,73 cu = 0.10 • Imagen A.1
window size: 3 Run took: 0 mins, 13 secs, 938038 microsecs.
• Imagen A.2
window size 5 Run took: 0 mins, 15 secs, 219922 microsecs.
• Imagen A.3
window size 7 Run took: 0 mins, 14 secs, 884585 microsecs.
• Imagen A.4
window size 11 Run took: 0 mins, 14 secs, 943042 microsecs.
Figura 25: Lee Enhanced: Im´agenes A.1, A.2, A.3 y A.4.
79
cu=0.20 • Imagen B.1
window size: 3 Run took: 0 mins, 13 secs, 415071 microsecs.
• Imagen B.2
window size 5 Run took: 0 mins, 13 secs, 913231 microsecs.
• Imagen B.3
window size 7 Run took: 0 mins, 14 secs, 309572 microsecs.
• Imagen B.4
window size 11 Run took: 0 mins, 14 secs, 691692 microsecs.
Figura 26: Lee Enhanced: Im´agenes B.1, B.2, B.3 y B.4.
80
cu=0.30 • Imagen C.1
window size: 3 Run took: 0 mins, 13 secs, 312316 microsecs.
• Imagen C.2
window size 5 Run took: 0 mins, 13 secs, 466702 microsecs.
• Imagen C.3
window size 7 Run took: 0 mins, 14 secs, 20081 microsecs.
• Imagen C.4
window size 11 Run took: 0 mins, 13 secs, 547806 microsecs.
Figura 27: Lee Enhanced: Im´agenes C.1, C.2, C.3 y C.4.
81
5.1.6.
Kuan
cu = 0.10 • Imagen A.1
window size: 3 Run took: 0 mins, 13 secs, 821943 microsecs.
• Imagen A.2
window size 5 Run took: 0 mins, 14 secs, 317794 microsecs.
• Imagen A.3
window size 7 Run took: 0 mins, 14 secs, 646845 microsecs.
• Imagen A.4
window size 11 Run took: 0 mins, 14 secs, 722029 microsecs.
Figura 28: Kuan: Im´agenes A.1, A.2, A.3 y A.4.
82
cu=0.20 • Imagen B.1
window size: 3 Run took: 0 mins, 13 secs, 559700 microsecs.
• Imagen B.2
window size 5 Run took: 0 mins, 13 secs, 827549 microsecs.
• Imagen B.3
window size 7 Run took: 0 mins, 14 secs, 116431 microsecs.
• Imagen B.4
window size 11 Run took: 0 mins, 14 secs, 869885 microsecs.
Figura 29: Kuan: Im´agenes B.1, B.2, B.3 y B.4.
83
cu=0.30 • Imagen C.1
window size: 3 Run took: 0 mins, 13 secs, 484662 microsecs.
• Imagen C.2
window size 5 Run took: 0 mins, 13 secs, 561343 microsecs.
• Imagen C.3
window size 7 Run took: 0 mins, 14 secs, 128048 microsecs.
• Imagen C.4
window size 11 Run took: 0 mins, 13 secs, 936172 microsecs.
Figura 30: Kuan: Im´agenes C.1, C.2, C.3 y C.4.
84
5.2.
Algoritmos de Clasificaci´ on
5.2.1.
K-Means
Imagen 1 • N´ umero de clases: 3
• Iteraciones realizadas / Iteraciones m´aximas: 20/1000 • Run took: 0 mins, 0 secs, 287946 microsecs. Imagen 2 • N´ umero de clases: 5
• Iteraciones realizadas / Iteraciones m´aximas: 34/1000 • Run took: 0 mins, 0 secs, 801951 microsecs. Imagen 3 • N´ umero de clases: 8
• Iteraciones realizadas / Iteraciones m´aximas: 82/1000 • Run took: 0 mins, 2 secs, 748033 microsecs.
Figura 31: K-Means
85
5.2.2.
Isodata
Para todos los casos: P=4 THETA M = 10 THETA S = 1 THETA C = 20 THETA O = 0.05 Resultados: Imagen 1 • Clases obtenidas / deseadas: 3 / 3
• Iteraciones realizadas / Iteraciones m´aximas: 13/1000 • Run took: 0 mins, 0 secs, 311133 microsecs. Imagen 2 • Clases obtenidas / deseadas: 4 / 5
• Iteraciones realizadas / Iteraciones m´aximas: 9 / 1000 • Run took: 0 mins, 0 secs, 302457 microsecs. Imagen 3 • Clases obtenidas / deseadas: 5 / 8
• Iteraciones realizadas / Iteraciones m´aximas: 10/1000 • Run took: 0 mins, 0 secs, 447103 microsecs.
Figura 32: Isodata
86
5.3.
Filtros y Algoritmos de Clasificaci´ on Combinados
Para todos los casos: iteraciones 100 win size 7 cmax 1.73
87
5.3.1.
Lee Enhanced y K-Means
Imagen 1 • cu: 0.20 • clases 3
• iteraciones 19/100
• Run took: 0 mins, 14 secs, 892865 microsecs. Imagen 2 • cu: 0.20 • clases 5
• iteraciones 40/100
• Run took: 0 mins, 15 secs, 235338 microsecs. Imagen 3 • cu: 0.30 • clases 3
• iteraciones 18/100
• Run took: 0 mins, 13 secs, 459135 microsecs. Imagen 4 • cu: 0.30 • clases 5
• iteraciones 69/100
• Run took: 0 mins, 15 secs, 475884 microsecs.
Figura 33: Lee Enhanced y K-Means
88
5.3.2.
Lee Enhanced e Isodata
Imagen 1 • cu: 0.20
• clases usadas/deseadas 3/3 • iteraciones 11/100
• Run took: 0 mins, 14 secs, 749540 microsecs. Imagen 2 • cu: 0.20
• clases usadas/deseadas 4/5 • iteraciones 13/100
• Run took: 0 mins, 14 secs, 804483 microsecs. Imagen 3 • cu: 0.30
• clases usadas/deseadas 3/3 • iteraciones 11/100
• Run took: 0 mins, 13 secs, 860447 microsecs. Imagen 4 • cu: 0.30
• clases usadas/deseadas 4/5 • iteraciones 17/100
• Run took: 0 mins, 14 secs, 314145 microsecs.
Figura 34: Lee Enhanced y K-Means
89
5.3.3.
Frost y K-Means
Imagen 1 • damp=0.10 • clases 3
• iteraciones 20/100
• Run took: 0 mins, 25 secs, 985338 microsecs. Imagen 2 • damp=0.10 • clases 5
• iteraciones 34/100
• Run took: 0 mins, 26 secs, 761625 microsecs. Imagen 3 • damp=1.0 • clases 3
• iteraciones 20/100
• Run took: 0 mins, 26 secs, 985227 microsecs. Imagen 4 • damp=1.0 • clases 5
• iteraciones 34/100
• Run took: 0 mins, 26 secs, 987382 microsecs.
Figura 35: Frost y K-Means
90
5.3.4.
Frost e Isodata
Imagen 1 • damp=0.10 • clases 3/3
• iteraciones 11/100
• Run took: 0 mins, 26 secs, 834296 microsecs. Imagen 2 • damp=0.10 • clases 5/5
• iteraciones 11/100
• Run took: 0 mins, 26 secs, 184869 microsecs. Imagen 3 • damp=1.0 • clases 3/3
• iteraciones 10/ 100
• Run took: 0 mins, 26 secs, 523221 microsecs. Imagen 4 • damp=1.0 • clases 4/5
• iteraciones 15/1000
• Run took: 0 mins, 26 secs, 395668 microsecs.
Figura 36: Frost e Isodata
91
6. 6.1.
An´ alisis de los Resultados y Observaciones Performance
Las implementaciones de los diversos algoritmos de este trabajo de Tesis, las realizamos de manera progresiva: comenzamos con implementaci´on simple, verificamos la correcci´ on del funcionamiento del algoritmo. Luego identificamos las porciones de c´ odigo que m´as tiempo de procesamiento insumen, hacemos profiling, las reimplementamos y as´ı sucesivamente. Si bien nuestra experiencia en el campo de algoritmos de alta performance en Python era inicial, r´ apidamente aprendimos nuevos m´etodos com´ unmente utilizados y recomendados por la comunidad para aplicarlos a nuestro software. La librer´ıa que proveemos puede mejorarse incluso a´ un mas en t´erminos de performace. Una alternativa es aplicar el uso de threads o hilos que provee Python. Los algoritmos de clasificaci´ on pueden paralelizarse aplicando t´ecnicas cl´asicas de paralelizaci´ on. Los algoritmos de filtrado tambi´en pueden paralelizarse pero de forma no trivial, ya que el filtrado de un pixel depende de sus vecinos, lo cual genera una fuerte dependencia a la hora de paralelizar. Otra alternativa ampliamente utilizada en la actualidad es migrar hacia el paradigma de la computaci´ on cient´ıfica de alta performance aunque aunque alej´andonos del lenguaje Python. En cuyo caso perder´ıamos la simplicidad y expresividad del lenguaje que elegimos ya que los lenguaje de programaci´on que se suelen utilizar para programar GPUs son lenguajes de bajo nivel a contraposici´on con Python4 . El uso de unidades de procesamiento gr´afico (GPUs o placas aceleradoras gr´ aficas) para aplicar estos algoritmos utilizando paralelizaci´on incrementar´ıa la performance muchos ordenes de magnitud en tiempo de procesamiento. Aunque como desventaja posee que las GPUs suelen tener un elevado costo de adquisici´ on.
4 Dependiendo del fabricante pueden existir wrappers para utilizar otros lenguajes de programaci´ on para desarrollar algoritmos en GPUs.
92
6.2.
Mejoras y Modificaciones de los Algoritmos Originales
En el algoritmo de Lee, la funci´on de peso es la siguiente: w(t) = 1
Cu2 /Cl (t)2
Donde: Cl (t) es el coeficiente de variaci´on de la ventana centrada en el p´ıxel t. Cu es el coeficiente de variaci´on del ruido en la imagen(en la imagen completa, no en la ventana). Al implementar el algoritmo y comenzar realizar pruebas, observamos que el algoritmo, tal cual sedefine en [1], al filtrar la imagen, asignaba valores negativos a algunos p´ıxeles. Investigando este mal comportamiento llegamos a la conclusi´on de que el problema se suscitaba cuando Cl (t) era mucho menor que Cu , caso en el cual el cociente Cu2 /Cl (t)2 se hace mayor que 1 y la funci´on de peso w(t) retorna vaˆ lores negativos. Esto, a su vez, generaba que la funci´on R(t) asignara valores negativos a los p´ıxeles filtrados. Si bien este problema no es inherente a la implementaci´on sino que al algoritmo mismo, nos pareci´ o que el mejor camino a seguir era buscar una soluci´on que estabilice el comportamiento del algoritmo en estos casos. La soluci´on que implementamos fue la de, en los casos en que el algoritmo genere un valor negativo para un p´ıxel, ignorar este valor y asignarle en vez la media de la ventana. La mejor´ıa que esta soluci´ on report´o no s´olo evit´o por completo la aparici´on de valores negativos en la imagen filtrada sino que tambi´en mejor´o notablemente las im´ agenes resultantes en cuanto a su apariencia visual. En el algoritmo de Kuan se observ´o exactamente el mismo problema y se solucion´ o del mismo modo y con similares resultados.
93
6.3.
Comparaciones con ENVI
Comparaci´ on de los resultados obtenidos en PyRadar y en ENVI, con el mismo input y con los mismos par´ ametros para los algoritmos. Cabe mencionar que la diferencia de tono que se percibe entre los colores de l´as im´ agenes obtenidas con PyRadar y las generadas con ENVI se deben enteramente a diferencias entre los respectivos algoritmos de equalizaci´on empleados. Los resultados de los algoritmos, tanto de filtrado como de clasificaci´on, son pr´ acticamente id´enticos. Dado que el software ENVI es propietario y de c´odigo cerrado, no podemos saber de qu´e modo ni con qu´e procedimientos equaliza las im´agenes. 6.3.1.
PyRadar y ENVI: Frost Filter.
Par´ ametros: window size: 7 damping factor: 1.0
Figura 37: PyRadar y ENVI: Frost Filter.
94
6.3.2.
PyRadar y ENVI: Lee Enhanced Filter.
Par´ ametros: Cµ : 7 damping factor: 1.0 Cmax : 1.73
Figura 38: PyRadar y ENVI: Lee Enhanced Filter.
95
6.3.3.
PyRadar y ENVI: Lee.
Par´ ametros: Cµ : 7 damping factor: 1.0 Cmax : 1.73
Figura 39: PyRadar y ENVI: Lee.
96
6.3.4.
PyRadar y ENVI: Kuan.
Par´ ametros: window size: 7 Cµ : 0.3
Figura 40: PyRadar y ENVI: Kuan.
97
6.3.5.
PyRadar y ENVI: Isodata.
Par´ ametros: clases: 8 iteraciones m´ aximas: 1000 P=4 THETA M = 10 THETA S = 1 THETA C = 20 THETA O = 0.05
Figura 41: PyRadar y ENVI: Isodata.
98
Referencias [1] Armand Lopes, Ridha Touzi, Edmond Nezry, “Adaptive Speckle Filters and Scene Heterogeneity”, IEEE Transactions on Geoscience and Remote Sensing, Vol 28, 6 November 1990. [2] Jon-Sen Lee, “Digital Image Enhancement And Noise Filtering by Use of Local Statistics”, IEEE Transactions on Pattern Analysis and Machine Intelligence, Vol Pami-2, 2 March 1980. [3] Victor S. Frost, Josephine Abott Stiles, K. S. Shanmugan, Julian C. Holtzman, “A Model for Radar Images and Its Application to Adaptive Digital Filtering of Multiplicative Noise”, IEEE Transactions on Pattern Analysis and Machine Intelligence, Vol Pami-4, 2 March 1982. [4] Darwin Kuan, Alexander A. Sawchuk, Timothy C. Strand, Pierre Chavel, “Adaptive Restoration of Images with Speckle”, IEEE Transactions on Acoustics, Speech and Signal Processing, Vol. ASSP-35, No. 3, March 1987. [5] Zhenghao Shi, Ko B. Fung, “A Comparison of Digital Speckle Filters”, Canada Center for Remote Sensing, 588 Booth Street, Ottawa, Ontario, Canada K1A0Y7. [6] John A. Richards, “Remote Sensing with Imaging Radar”, Springer, ANU College of Engineering and Computer Science The Australian National University, Canberra, ACT, 0200, Australia. [7] Ji-Hee Han, Sejung Yang, Byung-Uk Lee, “A Novel 3-D Color Histogram Equalization Method with Uniform 1-D Gray Scale Histogram”, IEEE Trans. on Image Processing, Vol. 20, No. 2, pp. 506-512, Feb. 2011. [8] Rafael C. Gonzalez, Richard E. Woods, “Digital Image Processing”, Third Edition, Prentice Hall, 2008. [9] Jingfeng Xiao, Jing Li, A. Moody, “A detail-preserving and flexible adaptive filter for speckle suppression in SAR imagery”, Department of Geography, University of North Carolina at Chapel Hill, Chapel Hill, NC 27599-3220, USA, Institute of Remote Sensing & GIS, Peking University, Beijing, 100871, P. R. China. [10] James B. Campbell, Randolph H. Wynne, “Introduction to Remote Sensing”, Fifth Edition, The Guilford Press 2011, 72 Spring Street, New York, NY 10012. [11] O. H. Bustos, M. G. Palacio, A. C. Frery, “Filtros interactivos reductores de ruido speckle en im´ agenes”, Revista Espa˜ nola de Teledetecci´on. 2002, 17: 61-70, Asociaci´ on Espa˜ nola de Teledetecci´on [12] Marija Norusis, “SPSS 17.0 Statistical Procedures Companion”, ISBN10: 0321621417, ISBN-13: 9780321621412, Publisher: Pearson, Copyright: 2010, Published: 12/28/2008.
99
[13] R. Sibson, “SLINK: an optimally efficient algorithm for the single-link cluster method”, The Computer Journal (British Computer Society), 16 (1): 30–34, 1973. [14] Martin Ester, Hans-Peter Kriegel, Jorg Sander, Xiaowei Xu, “A densitybased algorithm for discovering clusters in large spatial databases with noise”, In Evangelos Simoudis, Jiawei Han, Usama M. Fayyad, Proceedings of the Second International Conference on Knowledge Discovery and Data Mining, (KDD-96). AAAI Press. pp. 226–231., ISBN 1-57735-004-9, 1996. [15] US Air Force, http://www.af.mil/, “The Official Web site of the United States Air Force”, Last checked: 30/04/2012. [16] NOAA Satellite Information System (NOAASIS), http://noaasis. noaa.gov/, “NOAA Satellite Information System Home Page - Office of Satellite and Product, Operations.”, Last checked: 30/04/2012. [17] NASA History Program Office, http://history.nasa.gov/SP-4402/ ch2.htm, “SP-4402 Origins of NASA Names”, Last checked: 30/04/2012. [18] RapidEye AG, http://www.rapideye.net/gallery/index.htm, “RapidEye Image Gallery”, Last checked: 30/04/2012. [19] Centre for Remote Imaging, Sensing and Processing, CRISP, http:// www.crisp.nus.edu.sg/~research/tutorial/sar_int.htm, “Interpreting SAR Images”, Last checked: 30/04/2012. [20] Google Inc., http://earth.google.com/, “Google Earth”, Last checked: 30/04/2012. [21] Python Software Foundation, http://www.python.org/, “Python Programming Language – Official Website”, Last checked: 30/04/2012. [22] Python Software Foundation, Python Package Index, http://pypi. python.org/pypi, “PyPI - the Python Package Index”, Last checked: 05/05/2012.
100
A. A.1.
Utilidades timeit, mide el tiempo de ejecuci´ on de una porci´ on de c´ odigo
Este m´ odulo provee una manera simple de cronometrar el tiempo de ejecuci´on de peque˜ nas porciones de c´ odigo Python. Se puede utilizar desde la l´ınea de comandos como as´ı tambi´en desde otras interfaces. Timeit utiliza una funci´ on de tiempo espec´ıfica por plataforma para proporcionar el c´ alculo del tiempo de la manera m´as exacta posible. Es uno de los m´as utilizados ya que evita un n´ umero importante de trampas comunes al momento de medir tiempos de ejecuci´ on. Utilizaci´ on desde la l´ınea de comandos: python -m timeit [-n N] [-r N] [-s S] [-t] [-c] [-h] [statement ...] Donde los par´ ametros representan lo siguiente: n N, --number=N cuantas veces ejecutar “statement” -r N, --repeat=N cuantas veces repetir el timer (por defecto 3) -s S, --setup=S la porci´ on de c´ odigo que debe ser ejecutada solo inicialmente (por defecto: pass) -t, --time utiliza time.time() (por defecto en la mayor´ıa de las plataformas excepto Windows) -c, --clock utiliza time.clock() (por defecto en Windows) -v, --verbose muestra los resultados de la cronometraci´on en crudo -h, --help muestra una peque˜ na ayuda
101
Ejemplos de uso 5 : Desde l´ınea de comandos: python -m timeit -s ’’d={}’’ ’’for i in range(1000):’’ ’’
d[str(i)] = i’’
1000 loops, best of 3: 281 usec per loop
python -m timeit -s ’’import numpy’’ ’’rndNums = numpy.random.randint(0,100,1600000)’’ 10 loops, best of 3: 27.2 msec per loop
5 Algunos ejemplos fueron inspirados de http://docs.python.org/library/timeit.html y el libro Dive Into Python de Mark Pilgrim.
102
A.2.
Profiling
¿Qu´ e es un Profiler? Un profiler es un programa que describe la performance o desempe˜ no en tiempo de ejecuci´ on de un programa, proporcionando un gran n´ umero de estad´ısticos u ´tiles. La librer´ıa est´ andard de Python provee tres profilers distintos: cProfile, es una extensi´ on de Python hecha en C con una sobrecarga razonable para realizar profiling de programas de larga duraci´on. Est´a basado en “lsprof”, un profiler desarrollado por Brett Rosen y Ted Czotter. profile, es un m´ odulo Python cuya interface imita cProfile. Este profiler agrega una sobrecarga importante para realizar profiling. hotspot, era un m´ odulo experimental escrito en C, cuya meta era reducir la sobrecarga del profiling a expensas de los largos tiempos de postprocesamiento. Este profiler ya no es mantenido y no se incluye en las u ´ltimas versiones de Python. De los anteriores, nosotros utilizamos cProfile por ser el m´as recomendado por la comunidad. Ejemplo de uso: 1
import cProfile
2 3 4 5 6 7 8
def add(x, y): while x != 0: x = x - 1 y = y + 1 return y
9 10 11 12
def f(): return add(10 ** 7, 10 ** 7)
13 14
cProfile.run(’f()’) 4 function calls in 1.468 CPU seconds Ordered by: standard name ncalls 1 1 1 1
tottime 0.000 0.000 1.468 0.000
percall 0.000 0.000 1.468 0.000
cumtime 1.468 1.468 1.468 0.000
103
percall filename:lineno(function) 1.468 :1() 1.468 example.py:11(f) 1.468 example.py:4(add) 0.000 {method ’disable’ of ’_lsprof.Profiler’ objects}
B.
Sistema Operativo, Hardware y Software
La mayor parte de los algoritmos fueron desarrollados en el Sistema Operativo Mac OS X Lion de Apple. Al haber elegido Python como lenguaje de programaci´ on permite utilizarlos en distintas plataformas, como Linux y Windows, sin problemas.
B.1.
Datos Espec´ıficos sobre Hardware Utilizado
Computadora 1: Procesador: 2.7 GHz Intel Core i7 Memoria Ram: 4 GB 1333 Mhz DDR3 Tipo de Arquitectura: Darwin-11.3.0-x86 64-i386-64bit Computadora 2: Procesador: 2,66 GHz Intel Core 2 Duo Memoria Ram: 4 GB 1067 MHz DDR3 Tipo de Arquitectura: Darwin-11.3.0-x86 64-i386-64bit
B.2.
Datos Espec´ıficos sobre Software Utilizado
(Iguales en ambas computadoras) Python 2.6.7 (r267:88850, Jul 31 2011, 19:30:54) PIL Version 1.1.7 GCC 4.2.1 (Based on Apple Inc. build 5658) (LLVM build 2335.15.00) numpy Version 1.5.1 matplotlib Version 1.0.1 ipython Version 0.11 ipdb Version 0.6 readline Version 6.2.1 git version 1.7.4.4
104
B.3.
Herramientas de Software utilizadas
Durante la realizaci´ on de nuestro trabajo de tesis, usamos varias herramientas de software libre que nos fueron indispensables y de una utilidad crucial. Sirva esta menci´ on tambi´en como se˜ na de nuestro sincero agradecimiento a la comunidad por compartir su esfuerzo. Algunas de estas herramientas fueron: NINJA-IDE: IDE que utilizamos para escribir todo el c´odigo Python de PyRadar6 .
Figura 42: NINJA-IDE. TexMaker: IDE que utilizamos para escribir todo el c´odigo LATEX de la Tesis7 .
Figura 43: Texmaker. GIT: Sistema de control de versiones que usamos para todo el c´odigo absolutamente todo el contenido de nuestra tesi como tambi´en de PyRadar y la documentaci´ on. 8 .
Figura 44: GIT. 6 http://ninja-ide.org/ 7 http://www.xm1math.net/texmaker/ 8 http://git-scm.com/
105
TRAC: Este Issue tracker su sistema de Wiki, integrado con nuestro repositorio GIT fue la herramienta de base en torno a la cual todo nuestro trabajo se organiz´ o. 9 .
Figura 45: Trac.
9 http://trac.edgewall.org/
106
C. C.1.
C´ omo ejecutar los ejemplos de c´ odigo C´ odigo de los ejemplos
El c´ odigo de los ejemplos citados se encuentra dentro de la carpeta ‘pyradar/examples’. Cada archivo contiene un ejemplo autocontenido de uso de la librer´ıa PyRadar.
C.2.
Ejecuci´ on de las scripts de ejemplo
Las scripts de ejemplo son autocontenidas y, para ejecutar cualquiera de ellas, basta con los siguientes comandos:
1 2
cd /pyradar/examples python Un ejemplo particular podr´ıa ser el siguiente:
1 2
cd /pyradar/examples python example_kmeans.py
107
D. D.1.
Instalaci´ on Software B´ asico Requerido
Para poder instalar PyRadar asumimos instalados, configurados y funcionando lo siguiente: Python 2.7. pip Para instalar pip se puede utilizar el siguiente comando: easy_install-2.7 pip GDAL(Geospatial Data Abstraction Library): La instalaci´on de este paquete variar´ a entre distintos los sistemas operativos. Para instrucciones espec´ıficas, referirse al sitio oficial: http://www.gdal.org/
D.2.
Instrucciones de Instalaci´ on
Luego, el primer paso ser´ a obtener una copia del c´odigo. El repositorio GIT, que es p´ ublico acceso, est´a en: https://github.com/pyradar/pyradar/ La documentaci´ on, tickets, recursos adicionales, referencias y dem´as se encuentran en: https://github.com/pyradar/pyradar/wiki https://github.com/pyradar/pyradar/issues Para hacer checkout del repositorio GIT, ejecutar el siguiente comando:
1
git clone
[email protected]:PyRadar/pyradar.git/ Procedimiento recomendado de instalaci´on(en un sistema Linux, Mac o similar). Desde una consola, ejecutar los siguientes comandos:
1 2 3 4 5 6
cd ~ mkdir Dev cd Dev mkdir Code cd Code git clone
[email protected]:PyRadar/pyradar.git Una vez clonado el repositorio, el resto de los requerimientos se instalar´an f´acilmente utilizando el comando “pip”. 108
1 2
cd ~/Dev/Code/ pip install -r proj/setup/pip-requirements.txt Tambi´en creamos un pip-bundle, en caso de que por alg´ un motivo no se pudiera realizar la instalaci´ on directamente desde PyPi con el comando anterior. Para instalar las dependencias desde el pip-bundle, utilizar el siguiente procedimiento: Una vez clonado el repositorio, el resto de los requerimientos se instalar´an f´acilmente utilizando el comando “pip”.
1 2
cd ~/Dev/Code/ pip install proj/setup/PyRadar_Bundle.pybundle En el archivo “proj/setup/pip-requirements.txt” est´an los pip-requirements y las instrucciones de instalaci´ on. Tambi´en hay un setup recomendado para Bash.
Comentarios La librer´ıa fue desarrollada integramente en Mac OS X. Sin embargo, deber´ıa funcionar sin modificaciones en Linux, en tanto cuenten con los requerimientos mencionados al principio de esta secci´on.
109
C´ ordoba, Argentina, 2012
Mat´ıas Herranz Joaqu´ın Tita
110