lunes, 15 de abril de 2013

Convolución Bidimensional

En la entrada anterior me quedo pendiente profundizar algo en la base del código de la máscara de enfoque. La máscara de enfoque, y casi cualquier filtrado que se haga sobre una imagen, se implementa a partir de la convolución bidimensional discreta. En dicha convolución, el valor del píxel de salida se calcula mediante la suma ponderada de píxeles vecinos. En lo que se refiere al procesamiento de imágenes, la convolución se realiza entre la imagen y una matriz (los coeficientes del filtro) llamada máscara, o kernel. De manera matemática la convolución bidimensional se representa según:


Ya que lo más habitual es emplear una máscara de 3X3 elementos para realizar la convolución bidimensional, entonces la ecuación anterior se convierte en:


Como ejemplo para aclarar las áridas expresiones matématicas anteriores, tomamos una matriz de 4X4 que representa los valores de los píxeles de una imagen, y como máscara una matriz de 3X3

-Matriz imagen

 [15    20    101    100]
 [200  50      55        8]
 [10    11    230    202]
 [100  130  115    120]


-Máscara

 [7 8 3]
 [1 1 0]
 [1 2 1]


Los pasos para realizar la convolución consisten en:

1.- Rotar la máscara 180º a partir del elemento del centro.

 [1 2 1]
 [0 1 1]
 [3 8 7]


2.-Sobreponer el centro de la máscara sobre el elemento de interés.

 [15   1    20 2    101  1    100]
 [200 0    50 1      55  1        8]
 [10   3    11 8    230  7    202]
 [100       130     115        120]


3.-Multiplicar cada valor (peso) de la máscara rotada por el píxel de la matriz de imagen que se encuentra "bajo" la máscara.
4.-Sumar los productos individuales en el apartado anterior.

El resultado entonces es:

 [197   88]
 [181 190]


En el caso de trabajar en los extremos de la imagen lo que se hace es insertar ceros en los extremos, lo que se conoce como zero padding. El resultado del ejemplo es una matriz de 2X2 ya que no hemos trabajado en la zona de los extremos de la imagen.

El  código en Python, usando numpy y Scipy para realizar la convolución bidimensional es:

import numpy as np
import cv2

kernel=np.array([[1, 0, 0], [0, 1, 0], [0, 0, 1]], 'uint8')
signal=np.array([[1, 1, 1, 1], [1, 1, 1, 1], [1, 1, 1, 1], [1, 1, 1, 1]], 'uint8')

rows=signal.shape[0]-kernel.shape[0]+1
cols=signal.shape[1]-kernel.shape[1]+1
output=np.zeros((rows,cols), 'uint8')

kernel_reversed=np.rot90(np.rot90(kernel))

for i in range(0, output.shape[0]):
   for j in range(0, output.shape[1]):
       
#desplazando el kernel
        signal_patch=signal[i:i+len(kernel),j:j+len(kernel)]  
        output[i][j]=(kernel_reversed*signal_patch).sum()

print signal_patch
print output


Como referencias pueden valer los dos libros siguientes. El segundo es un clásico en el tratamiento digital de señales, viene a ser la continuación del, más clásico aún, Oppenheim.

-Gonzalez, R.C., and Woods, P.,
Digital Image Processing, Addison Wesley, 2002
-Proakis, John G., and Manolakis, Dimitris G.,
Digital Signal Processing: Principles, Algorithms and applications, Prentice-Hall Inc, 1996.

jueves, 11 de abril de 2013

Máscara de enfoque con Python


Código que implementa una máscara de enfoque para una imagen. Está escrito en Python con el apoyo de las librerías OpenCV, numpy y scipy. Para ejecutarlo es necesario emplear Python 2.7.4

En la próxima entrada comentaré el código en profundidad, explicando como se realiza la convolución bidimensional y el efecto de la máscara de enfoque sobre la imagen.


import numpy as np
import cv2
from scipy import ndimage, signal

            #filtro de enfoque
            # -*- coding: utf-8 -*-

def gauss2d(k, std):                                    #gaussiana bidimensional
        rows=2*k+1
        cols=2*k+1
        gaussian1d=signal.gaussian(cols,std)            #gaussiana unidimensional
        kernel_gauss=np.ndarray((rows,cols),"float")    #mascara/kernel
        for i in range(0,rows):
                for j in range(0,cols):
                        kernel_gauss[i,j]=gaussian1d[i]*gaussian1d[j]   
        kernel_gauss=kernel_gauss/kernel_gauss.sum()    #norm. de la gaussiana
        return kernel_gauss

def conv2d(img, krnl):                              #convolucion bidimensional
        rows=img.shape[0]-krnl.shape[0]+1
        cols=img.shape[1]-krnl.shape[1]+1
        output=np.ndarray((rows,cols), "float")
        kernel_reversed=np.rot90(np.rot90(kernel))
        for i in range(0, output.shape[0]):
                for j in range(0, output.shape[1]):
                        img_patch=img[i:i+len(krnl),j:j+len(krnl)]
                        y=max((kernel_reversed*img_patch).sum(),0)
                        z=min(y,255)
                        output[i,j]=z
        return output

def sharp(k,std):                  #creación de la mascara de enfoque           
        rows=2*k+1
        cols=2*k+1
        kernel_g=gauss2d(k,std)
        kernel=np.zeros((rows,cols),"float")
        kernel[k,k]=2
        kernel=kernel-kernel_g
        return kernel

def apply_filter(image, kernel):       #rutina para aplicarle el filtro a la imagen
        if len(img.shape) == 2:        #blanco y negro (1 canal y transparencia)
                img_filt = conv2d(img, kernel)
        else:                                 #color (3 canales y transparencia)
                img_filt = []
                for channel in range(img.shape[2]):
                    img_filt.append(conv2d(img[:,:,channel], kernel))
                img_filt = cv2.merge(img_filt)
        return img_filt

k=5                                  #radio de la gaussiana
std=3.0                           #desviacion estandar de la gaussiana

  #k y std son los parámetros que controlan la cantidad de enfoque que
  #se quieres aplicar a la imagen
img=cv2.imread('ruta de la imagen a leer')       #leer imagen
img=img.astype("float")
kernel=sharp(k,std)
img_filt=apply_filter(img, kernel)
cv2.imwrite('ruta de la imagen a escribir',img_filt)   #escribir imagen

lunes, 13 de agosto de 2012

Umbralización local o adaptativa.


En muchas aplicaciones, no se puede obtener un umbral global para un histograma, es decir, no se puede  obtener una buena segmentación con un único umbral  para toda la imagen. Esto ocurre cuando el fondo no es constante y el contraste de los objetos varía en la imagen, por lo que la umbralización daría buenos resultados en una parte de la imagen, pero para el resto de la imagen, la segmentación no sería la adecuada.

Si las variaciones del fondo de la imagen pueden ser descritas a través de funciones conocidas, dependientes de la posición en la imagen, se podría intentar corregir la segmentación utilizando la corrección de niveles de gris; para que, posteriormente, la aplicación de un único umbral a toda la imagen, diera buenos resultados en la segmentación.Otra solución sería el uso de la umbralización local o adaptativa.

Los umbrales locales pueden ser determinados dividiendo una imagen en subimágenes y calculando los umbrales de dichas subimágenes, o examinando las intensidades de la imagen en los alrededores de cada píxel. En el primer método, primero, se divide una imagen en subimágenes rectangulares solapadas, y se obtienen los histogramas de cada subimagen. Estas subimágenes deben ser los suficientemente grandes como para abarcar el fondo, y el objeto en cuestión. Si una subimagen tiene un histograma bimodal, el mínimo entre los dos picos determinará el umbral local. En caso que el histrograma sea unimodal, el umbral  se obtendrá por interpolación de los umbrales locales de las subimágenes adyacentes. Por último,  es necesario realizar una segunda interpolación para encontrar los umbrales adecuados a cada píxel.

En el último método, se puede seleccionar un umbral mediante el valor medio de la distribución local de intensidad. Para esto, se pueden usar otro tipo de medidas estadísticas, tales como la media de las desviaciones estándar, la media de los valores máximos y mínimos, o medidas basadas en los valores de los gradientes de intensidad. 

La umbralización local es computacionalmente más costosa que la global. Es muy útil a la hora  de segmentar objetos en fondos no homogéneos, y para extraer regiones muy pequeñas y dispersas.

domingo, 12 de agosto de 2012

Umbralización global

Esta técnica de umbralización supone que la imagen posee un histograma bimodal, y por lo tanto, el objeto puede ser extraído del resto de la imagen mediante una simple operación que compare los valores de la imagen con el valor umbral T.

 Fig. 1: Histograma bimodal con umbral T.

Para una imagen f(x,y), con el histograma de la figura 1, los píxeles del objeto y del fondo tienen los niveles de gris agrupados en dos modos dominantes. La manera más obvia de extraer el objeto es seleccionar un umbral T que separe los dos modos.

La imagen g(x,y), resultante de aplicar la umbralización, viene definida por:

El resultado es una imagen binaria, donde los píxeles con valores de intensidad igual a 1 corresponden al objeto deseado; mientras que los píxeles con valor 0, corresponder al resto de la imagen.

La figura 2 muestra el resultado de la segmentación de una imagen mediante umbralización, donde la imagen original 2(a) contiene células blancas sobre un fondo negro, y las intensidades de los píxeles varían entre 0 y 255. El umbral T=127 se selecciona como el mínimo entre los dos modos del histograma 2(b), y el resultado de la segmentación se muestra en la figura 2(c), donde los píxeles con intensidades que superen el umbral (T=127) se representan en blanco. En última instancia, en la figura 2(d), los bordes de las células se obtienen aplicando el laplaciano de orden tres a la imagen umbralizada (figura 2c).

Fig. 2: (a) Imagen original, (b) histograma de la imagen, (c) resultado de la umbralización con T=127, y (d) bordes de las células.

Hay muchas maneras para seleccionar un umbral global, una de ellas se basa en un modelo de clasificación que minimice la probabilidad de error. Por ejemplo, para el caso de una imagen con un histograma bimodal, se puede calcular el error como el número total de píxeles del fondo que se clasifican como del objeto, y píxeles del objeto que se clasifican en el fondo. Una versión semiautomatizada de esta técnica, para medir volúmenes ventriculares de resonancia magnéticas en 3D; consiste en que un operador selecciona manualmente dos píxeles, uno del objeto y otro del fondo. Comparando la distribución de las intensidades de píxel en regiones circulares alrededor de los píxeles seleccionados, se calcula automáticamente el umbral y se corresponde con el mínimo número de píxeles mal clasificados entre las dos regiones. El resultado es de esta técnica se muestra como un mapa del contorno, superpuesto a la imagen original. Además, en caso que fuera necesario, el operador puede modificar manualmente cualquier parte del borde. Esta misma técnica también se aplicó para extraer linfomas de imágenes CT, pero se encontró muy dependiente de la selección de  los puntos interiores y exteriores por parte del usuario.

En muchas aplicaciones se obtiene una segmentación apropiada cuando el área, o perímetro, de los objetos, es sensible a pequeñas variaciones del nivel umbral seleccionado. En la figura 3(a) se muestra el perfil de intensidad de un objeto que es más brillante que el fondo, y tres umbrales (T1, T2 y T3) para llevar a cabo la segmentación. Una pequeña variación ΔT, en el nivel umbral más bajo, causaría cambios significativos en el área, o perímetro, del objeto segmentado; esto mismo se da, también, para el nivel umbral más alto. En cambio, una variación ΔT en el nivel umbral medio, tendría una mínima influencia sobre el área, o perímetro, del objeto.


Fig. 3: (a) Histograma con umbrales T1, T2 y T3, y (b) representación gráfica del área o el perímetro frente al umbral T.

El área A(T) y el perímetro P(T), son funciones que dependen del umbral T, y suelen presentar un aspecto similar al de la figura 3(b). Por lo tanto, el umbral que minimiza dA(T)/dT y dP(T)/dT, es una buena opción, especialmente ante la ausencia de un operador manual, y/o no se dispone de información previa del objeto.

Una técnica relaciona que evalúa múltiples umbrales, está basada en la estimación del valor del gradiente en el contorno del objeto segmentado. El valor medio del gradiente medio viene dado por:

donde H(T) es la función de histograma. El umbral que maximiza el límite del gradiente medio es el que se selecciona.

       Si una imagen contiene más de dos tipos de regiones, es posible segmentarla usando varios umbrales simples, o usando la técnica de umbrales múltiples. Con el incremento del número de regiones, aumenta el número de modos a distinguir dentro del histograma; y por lo tanto, la selección del umbral se vuelve más difícil.

            El uso de umbrales globales es computacionalmente más simple y rápido. Funciona bien con imágenes que contienen objetos con una intensidad uniforme, y diferenciados del fondo. Sin embargo, falla si existe un bajo contraste entre el objeto y el fondo, si la imagen es ruidosa, o si la intensidad del fondo varía significativamente a lo largo de la imagen.

Dimensionalidad

La dimensionalidad hace referencia al dominio de la imagen, bidimensional o tridimensional, sobre el cual opera el método de segmentación. Los métodos que solamente se apoyan en las intensidades de la imagen son independientes del dominio de la imagen. De cualquier manera, ciertos métodos como los modelos deformables, los campos aleatorios de Markov (MRF – Markov random fields), y región creciente (region growing), incorporan información espacial y por lo tanto pueden operar de forma diferente dependiendo de la dimensionalidad de la imagen. Generalmente, los métodos 2D se aplican a imágenes 2D y los métodos 3D se aplican a imágenes 3D; pero en algunos casos, los métodos 2D se aplican de manera secuencial a los cortes de una imagen 3D. Esto se debe a razones prácticas como la facilidad de implementación, menor complejidad y coste computacional, y a que ciertas estructuras son más sencillas de definir mediante cortes bidimensionales.

miércoles, 8 de agosto de 2012

Pure Data-PD


Pure Data (PD) es un lenguaje de programación gráfico, en tiempo real, para audio, vídeo y  procesado de gráficos. Es la tercera rama más importante de la familia de lenguajes de  programación  patcher conocida como Max (Max / FTS, ISPW Max, Max / MSP,jMax, etc), y fue  originalmente desarrollado por Miller Puckette en el IRCAM. Aunque el núcleo de PD sigue siendo  desarrollado y mantenido por Miller Puckette, se incluye el trabajo de muchos desarrolladores, por  lo que todo el paquete supone un gran esfuerzo comunitario.

PD es un ejemplo de lenguaje de programación de "flujo de datos". Es decir, las funciones u  "objetos" son conectados, o "parcheados", unos con  otros en un ambiente gráfico que modela el  flujo de control y de los datos (audio, vídeo y/o gráficos). Por lo que PD es un lenguaje de programación especialmente orientado al procesado de señales

El entorno de programación posee una base modular  de código, en el que los objetos son utilizados como elementos de construcción en el desarrollo de programas en Pure Data. Además, el programa se hace arbitrariamente extensible a través de una API pública, por lo que alienta a otros desarrolladores a añadir sus propias rutinas, ya sea en el lenguaje de programación C o, con la ayuda de otros lenguajes externos. Estas unidades de código modulares y reusables, que resultan de la programación de algoritmos en PD, se suelen denominar "patches", y son usadas como programas independientes y compartidas libremente entre la comunidad de usuarios de PD. 

Para el desarrollo de los ejemplos gráficos que se presentarán en posteriores posts en este apartado se emplea la líbreria GEM (Graphics Environment for Multimedia),que está indicada para la generación de gráficos OpenGL en tiempo real. 

martes, 7 de agosto de 2012

Métodos de compresión de imágenes

Muchos formatos de almacenamiento de imágenes realizan una compresión de los datos para así reducir el espacio necesario para el almacenamiento de las imágenes, y aumentar la velocidad de transmisión de los archivos de imagen. La mayoría del software empleado en el procesado de imágenes permite al usuario especificar si desea compresión o no, y el método concreto para realizar la compresión. Para elegir el método de compresión, aparte de observar la compatibilidad del software a emplear, es necesario decantarse por un método con pérdida de información, o por uno que no represente pérdida de información. 

Los métodos de compresión con pérdidas descartan información que se considera irrelevante a la percepción humana; es decir, en el caso de imágenes, se explotan las limitaciones del ojo humano y, entonces, es posible llegar a lograr  una reducción del 80-90% en el tamaño del  archivo. 

En cambio, los métodos de compresión sin pérdidas, reducen el tamaño de los datos mediante  procesos totalmente reversibles, que sólo eliminarán los datos redundantes. Por lo tanto, el  contenido de las imágenes comprimidas mediante métodos sin pérdidas serán idénticas, en el contenido de información, a sus homólogas sin comprimir. 

Los tipos de redundancia que se pueden dar en un archivo de imagen son: 

1. Redundancia de código: 
Este tipo de redundancia existe cuando el método de codificación tiene más precisión que la  que realmente es necesaria para una imagen en particular. Se puede minimizar usando una  profundidad de bit menor y mediante el empleo de tablas de búsqueda. 

2. Redundancia espacial: 
Se da cuando hay grandes regiones con valores de píxeles idénticos, o por lo menos muy  similares. Por ejemplo, el fondo negro de una imagen de rayos-X. Esta redundancia se puede  reducir mediante un método que codifique una descripción de esas regiones homogéneas. 

3. Redundancia psico-visual: 
Esta es la información que no puede ser percibida; por ejemplo, los seres humanos no responden bien a los pequeños cambios en la intensidad y el color. Esta redundancia puede ser eliminada homogeneizando las pequeñas variaciones de intensidad y color. 

Un formato de archivo de mapa de bits almacena una  imagen  m x n como una matriz  rectangular de  m x n intensidades de píxeles. Este formato es, en general, una manera muy costosa de almacenar una imagen. En primer lugar, la mayoría de las imágenes tienen extensiones considerables, con valores de intensidades idénticas o casi idénticas (Redundancia espacial). En una  imagen médica, por ejemplo, la mayoría de los fondos suelen ser negros, conteniendo únicamente ruido.  La segunda limitación radica en el hecho que a menudo hay muchos menos valores de intensidad de píxeles presentes en la imagen que los que realmente se pueden codificar con la profundidad de bits nominal (Redundancia de código). 

Los métodos de compresión de imagen se aprovechan de las redundancias espacial, de código y psicovisual para llegar a reducir drásticamente la cantidad de espacio necesaria para el almacenamiento de imágenes y; por lo tanto, reducir el tiempo necesario para la transmisión de la imagen. Como ejemplo, si las primeras 100 filas de  una matriz de imagen de  m x n representan un fondo negro, entonces en vez de usar 100 x n x 8 bits, todos a cero; para almacenar esta información podemos simplemente utilizar un código que diga que los píxeles n de 1 a 100 tienen un valor cero. Esta codificación no sólo ahorra una gran cantidad  de espacio, sino que no conlleva pérdida de información en la imagen. 

Por otra parte, si se asume la pérdida de algo de  información que se considera irrelevante, una alternativa es ajustar valores de píxeles muy similares a un valor idéntico y, por lo tanto, reducir el número total de diferentes intensidades que habría que codificar. Cuando la imagen contiene menos valores discretos de intensidad que la profundidad de bits nominal puede codificar, entonces se puede ahorrar espacio codificando las intensidades en una tabla de búsqueda. 

El análisis estadístico de los datos de imagen puede dar lugar a nuevas mejoras en la compresión. Si las intensidades de los píxeles no representan información significativa de la imagen, entonces se pueden omitir de la tabla de búsqueda cambiándolos por el valor cercano más comun. Del mismo modo, se puede decidir que los píxeles individuales o pequeños grupos de píxeles, que no se ajustan a algún patrón observado en su vecindad y, que no aportan información; se pueden sustituir por sus valores originales en la codificación almacenada. Todo esto hace que se ahorre más espacio de almacenamiento, pero se pierde más información de la imagen original.