Implementación y evaluación de técnicas de mejoras de la imagen

Equation Chapter 1 Section 1 Trabajo Fin de Grado Grado en Ingeniería de las Tecnologías de Telecomunicación Implementación y evaluación de técnicas de mejoras de la imagen Dpto. Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Autor: Juan Carlos Hidalgo Palomo Tutores: Carmen Serrano Gotarredona Begoña Acha Piñero Sevilla, 2023 iii Trabajo Fin de Grado Grado en Ingeniería de las Tecnologías de Telecomunicación Implementación y evaluación de técnicas de mejoras de la imagen Autor: Juan Carlos Hidalgo Palomo Tutores: Carmen Serrano Gotarredona Begoña Acha Piñero Dpto. de Teoría de la Señal y Comunicaciones Escuela Técnica Superior de Ingeniería Universidad de Sevilla Sevilla, 2023 v Trabajo Fin de Grado: Implementación y evaluación de técnicas de mejoras de la imagen Autor: Juan Carlos Hidalgo Palomo Tutores: Carmen Serrano Gotarredona Begoña Acha Piñero El tribunal nombrado para juzgar el Proyecto arriba indicado, compuesto por los siguientes miembros: Presidente: Vocales: Secretario: Acuerdan otorgarle la calificación de: Sevilla, 2023 El Secretario del Tribunal vii A mi familia A mis maestros ix Resumen En este proyecto se definen, en primer lugar, una serie de medidas estadísticas para la calidad de la imagen digital, denominadas referenciadas y no referenciadas, diferenciándose, fundamentalmente, en la comparación con la existencia de la imagen original o no. A través de estas medidas, se describen y aplican unos métodos, algunos sirven para eliminar imperfecciones de la imagen, o como se conoce en el ámbito del tratamiento digital de imágenes, el ruido, otros sirven para mejorar o visualizar mejor otras características de la imagen, como la croma, la intensidad, la nitidez de los bordes, etc. Finalmente, se implementaran y mostraran los resultados del tratamiento digital de las imágenes realizados en el software MATLAB. 1 INTRODUCCIÓN El objetivo técnico de este proyecto es evaluar los principales métodos o técnicas básicas de mejora y eliminación de ruido para una serie de imágenes médicas obtenidas de la base de datos pública International Skin Imaging Collaboration (ISIC). Para ello se emplearán diferentes medidas de calidad de la imagen, clasificadas en referenciadas y no referenciadas, distinguiendo un tipo del otro la existencia o no de la imagen original para su comparación. A las medidas referenciadas también se las conoce como medidas de fidelidad, y las mejoras o cambios que se le apliquen a la imagen están relacionados con la imagen original. Por otro lado, las medidas no referenciadas comparan propiedades y atributos de la propia imagen. Describiremos brevemente, qué es el ruido digital en la imagen, cómo afecta y qué principales técnicas usaremos para eliminarlo todo lo que sea posible. Finalmente, implementáramos los métodos de eliminación de ruido y de mejora de la imagen para ver los resultados de las medidas referenciadas y no referenciadas, a través de la aplicación MATLAB discutiendo los resultados. Implementación y evaluación de técnicas de mejoras de la imagen 12 12 2 MEDIDAS DE CALIDAD DE LA IMAGEN 2.1 Referenciadas En primer lugar, empezaremos definiendo qué son las medidas de calidad referenciadas, también conocidas como medidas de fidelidad. En este sentido, la fidelidad se conoce como una medida de referencia completa de la calidad de la imagen, ya que se requiere la referencia o la imagen original sin distorsiones para estimar dicha fidelidad. Por tanto, es solo un indicador de la calidad general de una imagen [1-2]. Las principales medidas de fidelidad que vamos a analizar en este proyecto son: 2.1.1 El error cuadrático medio, en inglés Mean Square Error (MSE) Esta medida es muy usada en la comparación de datos de forma general, mide el promedio de los errores al cuadrado, compara un valor predicho y un valor observado o conocido. En el caso de las imágenes vamos a comparar la imagen original, f, con la imagen a tratar g. En primer lugar, del espacio de color RGB, el cual está formado por 3 planos de cada color primario, R, G y B, vamos a pasar al espacio de color CIE L*a*b*, ya que éste permite separar la luminosidad o brillo (plano L) de las coordenadas cromáticas [a*= coordenadas rojo/verde (+a indica rojo, -a indica verde); b* = coordenadas amarillo/azul (+b indica amarillo, -b indica azul)]. Figura 1. Espacio de color CIE L*a*b*. El MSE para las imágenes se define como: 𝑀𝑆𝐸𝐶𝐼𝐸𝐿𝐴𝐵= 1 𝑁∑‖𝑓𝑛−𝑔𝑛‖2 1≤𝑛≤𝑁 (2–1) donde N es el número de píxeles de la imagen, || || representa la distancia Euclídea, 𝑓𝑛 y 𝑔𝑛 hace referencia a los píxeles en las imágenes en ese espacio de color acotado de 1 ≤ n ≤ N. 13 2.1.2 La proporción máxima o pico de la relación señal a ruido, del inglés Peak Signal-to- Noise Ratio, (PSNR) De forma general, este concepto mide la relación entre la máxima energía de una señal y el ruido que afecta a la misma. Debido a que muchas señales tienen un gran rango dinámico, el PSNR se expresa generalmente en escala logarítmica, utilizando como unidad el decibelio. Pero uno de sus usos principales es como medida de calidad y fidelidad en el ámbito del análisis, reconstrucción, compresión y estudio de las imágenes. Por tanto, para las imágenes se define como la relación entre los valores o píxeles máximos del rango y el error cuadrático medio (MSE). 𝑃𝑆𝑁𝑅𝐶𝐼𝐸𝐿𝐴𝐵= 10 log10 {2552 1 𝑁 ∑ ‖𝑓𝑛−𝑔𝑛‖2 1≤𝑛≤𝑁 } (2–2) 𝑃𝑆𝑁𝑅𝐶𝐼𝐸𝐿𝐴𝐵= 10 log10 { 2552 𝑀𝑆𝐸𝐶𝐼𝐸𝐿𝐴𝐵} (2–3) donde el rango dinámico de las componentes de color es de [0,255]. El error cuadrático medio guarda una relación inversamente proporcional al pico de la relación señal a ruido, es decir que, si para un mismo escenario, la primera aumenta, la segunda disminuye y viceversa. En lo que se refiere a la codificación para la compresión, el PSNR es mayor cuanto mejor es la codificación. Además el comité MPEG (Moving Picture Experts Group) definido como un grupo de trabajo de expertos que se formó por la ISO (Organización Internacional de Normalización) y la IEC (Comisión Electrotécnica Internacional) para establecer estándares para el audio y la transmisión de video, emplea un valor umbral informal de 0,5 dB en el incremento del PSNR para decidir si se incluye una determinada mejora en un algoritmo de codificación, ya que se considera que este aumento del PSNR es apreciable visualmente [3]. Estas medidas las hemos tomado en el plano L del espacio de color CIE Lab, ya que el rango está acotado a un solo plano, pero también las hemos obtenido del espacio de color RGB, con la única diferencia que el valor del MSE es la media aritmética de los tres planos de color R, G y B. Implementación y evaluación de técnicas de mejoras de la imagen 14 14 2.2 No referenciadas Se trata de medidas realizadas en una sola imagen, sin disponer de la imagen original con la que comparar, por tanto, todas las propiedades, atributos, transformaciones, etc., se realizan dentro de la misma imagen. En un estudio llevado a cabo por Sabakis AE, Stephen PE y Loui AC, llamado “Evaluation of image appeal in consumer photography”, se elaboró un ranking de los atributos de imagen más observados por los humanos y entre los más votados están tres medidas no referenciadas, como son: el contraste, el colorido (colorfulness) y la nitidez de los bordes de la imagen (sharpness). 2.2.1 El contraste Este término tiene varios significados, aunque de forma general la Real Academia Española (RAE) lo define como la oposición, contraposición o diferencia notable que existe entre personas o cosas. El contraste se usa en el diseño gráfico para diferenciar elementos a nivel de color, por forma o por disposición de los mismos. A nivel de color, podemos diferenciar siete tipos para contrastar: de color puro, de claro a oscuro, de temperatura, de colores complementarios, de saturación, de cantidad y simultáneo. Aunque en este proyecto nos vamos a centrar en la definición relacionada con el tratamiento de imágenes que, según la RAE, es la relación existente entre el brillo, luminancia o intensidad de cada uno de los puntos que forman la imagen. De forma general se representa como: 𝐷𝑖𝑓𝑒𝑟𝑒𝑛𝑐𝑖𝑎 𝑑𝑒 𝐿𝑢𝑚𝑖𝑛𝑎𝑛𝑐𝑖𝑎 𝐿𝑢𝑚𝑖𝑛𝑎𝑛𝑐𝑖𝑎 𝑃𝑟𝑜𝑚𝑒𝑑𝑖𝑜 (2–4) Y podemos medirlo gracias a la fórmula del contraste de Michelson [4]: 𝐼max − 𝐼𝑚𝑖𝑛 𝐼max + 𝐼𝑚𝑖𝑛 (2–5) donde Imax es la intensidad o luminancia máxima e Imin es la intensidad mínima. Y a la fórmula del contraste de Weber: 𝐼−𝐼𝑏 𝐼𝑏 (2–6) desarrollada por los científicos Ernst Heinrich Weber (1795–1878) y Gustav Theodor Fechner (1801–1887) en la ley de Weber-Fechner. Donde I es la luminancia de la características e Ib la luminancia del fondo. Usada en los casos donde existen una proporción pequeña de características dentro de un gran fondo uniforme, o lo que es lo mismo, la luminancia promedio es muy similar a la luminancia de fondo [5]. Esta última formula también aparece en muchos documentos de forma abreviada como: 𝛥𝐿 𝐿 (2–7) De esta forma, vamos a calcular el contraste localizado en un área limitada, realizando la diferencia entre su valor máximo y su valor mínimo. 15 2.2.2 El colorido (colorfulness) Se define como la presencia y la viveza del color en toda la imagen, da una sensación visual según la cual el color percibido de un área parece ser más o menos cromático. Según este concepto, cada color tiene dos atributos perceptivos [6]: 1. La saturación: es el colorido de un área o de toda la imagen con relación a su propio brillo. 2. La croma: es el colorido de un área relativo al brillo de otro color que aparece en blanco bajo condiciones similares de visualización. Podemos ver también algunos de estos atributos representados en el espacio de color CIE L*C*h, similar al espacio CIE L*a*b*, que se sirven del mismo diagrama, pero el primero usa coordenadas cilíndricas, mientras que el segundo coordenadas rectangulares. L (Lightness) es la luminosidad, C (Chroma) es la croma y h (Hue) es el matiz. Figura 2. Espacio de color CIE L*C*h. Para calcular el colorido de una imagen sumamos la media de todos los valores de la imagen más la variación estándar de la misma, la cual indica qué tan dispersos están los datos con respecto a la media. 2.2.3 La nitidez de la imagen (sharpness) De forma general, la nitidez es la claridad con la que se perciben los detalles. En las imágenes se suele medir el gradiente entre píxeles. El gradiente, en el caso de una imagen en dos dimensiones (caso discreto), se basa en las diferencias entre los niveles de grises de la imagen adyacentes en la misma fila y columna [7]: Implementación y evaluación de técnicas de mejoras de la imagen 16 16 ∇𝑓=[𝐺𝑥 𝐺𝑦] (2–8) Para su cálculo, realizaremos el promedio del módulo del gradiente. La nitidez especifica como los bordes de cada elemento existente en la imagen no se ven borrosos. 3 MÉTODOS DE ELIMINACIÓN DE RUIDO Y MEJORA EN IMÁGENES DIGITALES 3.1 El ruido en imágenes El ruido en las imágenes se define como una variación del color y la intensidad de manera aleatoria de algunos píxeles de la imagen, provoca una mala visualización de la calidad y resolución de las mismas, éste se debe a múltiples factores como: ruido de Poison debido a la naturaleza corpuscular de la luz, ruido térmico de los receptores (ejemplo: corriente de oscuridad en la CCD), ruido de cuantificación, ruido de compresión con pérdidas, etc. Surgen principalmente, en la fase de toma o adquisición y transmisión de la imagen [8]. Figura 3. Modelo del proceso de degradación/restauración de una imagen. Podemos modelarlos con dos tipos de ruidos principalmente, llamados: 3.1.1 Ruido Gaussiano Está asociado con la radiación electromagnética, muestra una densidad de probabilidad que responde a una distribución normal o distribución de Gauss, que siguen todos los píxeles de la imagen afectada. La función densidad de probabilidad (FDP), se define como la probabilidad relativa de que la variable aleatoria que tiene asociada tome un determinado valor [8], esta probabilidad viene dada por la formula: 𝑃[𝑎≤𝑋 ≤𝑏]= ∫𝑓𝑋(𝑥) 𝑑𝑥 𝑏 𝑎 (3 –1) La FDP de Gauss se representa con la expresión: 𝑝(𝑧)= 1 √2𝜋𝜎 𝑒−(𝑧− 𝑧) 2𝜎2 −∞<𝑧< ∞ (3 –2) Degradación + Filtro(s) de restauración f (x, y) g (x, y) 𝑓መ (x, y) Ruido Ƞ(x, y) 17 Figura 4. FDP Ruido gaussiano. Donde z representa la intensidad, 𝒛 es la media o valor esperado de la intensidad y 𝛔 es la desviación estándar. Estos dos últimos conceptos son muy importantes en el tratamiento y análisis del ruido, ya que aportan información de cómo están distribuidos los píxeles. También se usa el parámetro denominado varianza, que no es más que desviación estándar al cuadrado, 𝛔2 . Representación del ruido gaussiano en la Figura 5: Figura 5. Ruido Gaussiano. 3.1.2 Ruido Impulsivo (sal y pimienta) El valor del píxel no tiene relación con el valor original, sino que toma valores alrededor del máximo o alrededor del mínimo, en el caso de escala de grises, blanco (de ahí que se le llame sal) o negro (de ahí que se le llame pimienta), respectivamente. Se modela con una función de distribución no gaussiana o escalón. Vemos en la Figura 7 el efecto del ruido, donde algunos píxeles aleatoriamente toman unos valores y otros toman otros. Su función densidad de probabilidad (FDP) [8] se define como: 𝑝(𝑧)={𝑃𝑠 𝑧= 2𝑘−1 𝑃𝑝 𝑧=0 1−(𝑃𝑠+𝑃𝑝) 𝑧=𝑉 (3 –3) Implementación y evaluación de técnicas de mejoras de la imagen 18 18 donde el rango de intensidades, z, es [0, 2𝑘−1] y V es cualquier valor entero entre 0 < V < 2𝑘− 1 . Cuya representación gráfica se ve en la Figura 6: Figura 6. FDP Ruido sal y pimienta. A una imagen que ha sido degrada con ruido sal y pimienta, cuyos valores de intensidad sigue la PDF anterior, se le asigna un 0 (pimienta: negro) al ruido con una probabilidad 𝑃𝑝. De igual forma, el ruido toma un valor 2𝑘−1 (el valor máxino de intensidad, sal: blanco) con una probabilidad 𝑃𝑠. El resto valores permanecen igual con lo que resta de probabilidad, 1−(𝑃𝑠+𝑃𝑝) [8]. La probabilidad 𝑃, es la probabilidad con la que un píxel es corrompido por el ruido de sal y pimienta, 𝑃=𝑃𝑠+𝑃𝑝. A 𝑃 tambien se le conoce como densidad de ruido. Representación del ruido sal y pimienta en una imagen, Figura 7: Figura 7. Ruido Impulsivo. 19 3.2 Métodos de eliminación del ruido Existen métodos para solventar este problema y eliminarlo en gran parte, algunos funcionan mejor dependiendo del tipo y características del ruido. Nos vamos a centrar en dos de los principales. 3.2.1 Filtro de media: Para este método, consideremos un píxel rodeado por un grupo pequeño de píxeles dentro de la imagen, este grupo se denomina vecindario. Para poder aplicarlo es necesario los siguientes requisitos: - La imagen es relativamente constante en un vecindario, no hay apenas variación. - Las variaciones en los píxeles son debidas al ruido. Algunos de los vecindarios o agrupaciones de píxeles más utilizados son: un cuadrado 3x3 (donde el píxel está conectado con 8 vecinos), una ventana de 5 píxeles, donde el píxel central está conectado con 4 o un grupo de barras, donde el píxel solo está conectado arriba y abajo o a la derecha e izquierda [9]. Figura 8. Ejemplos de vecindarios. Para un vecindario 3x3, 8-conectado: Figura 9. Vecindario 3x3, 8 conectado de un píxel. 𝑔(𝑚,𝑛)=1 9 ∑ ∑𝑓(𝑚+𝑝,𝑛+𝑞) (3 –4) 1 𝑞=−1 1 𝑝=−1 donde f(m,n) es la imagen de entrada u original, p y q son las posiciones alrededor del píxel y g(m,n) es la imagen de salida tras el filtrado. Implementación y evaluación de técnicas de mejoras de la imagen 26 26 Donde 𝑘 especifica qué parte de la máscara es añadida, según el valor que tome: - Para k = 1 se conoce como enmascaramiento del desenfoque (USM). - Para k > 1 se trata del filtrado de refuerzo de frecuencias altas, debido a que se aumentan aún más las componentes de alta frecuencia, ya que se le está dando más peso a la imagen enmascarada. La máscara de enfoque no hace más que, en término de procesamiento de señales, amplificar las componentes de alta frecuencia. En ocasiones, aunque se mejora la nitidez de la imagen, puede suceder que se distorsionen elementos de la imagen. El enmascaramiento del desenfoque se puede controlar a través de tres parámetros: I. La cantidad, que controla la fuerza del efecto de la nitidez, se mide en porcentaje. A mayor porcentaje de cantidad, mejor será el contraste de los píxeles nítidos. II. El radio, el cual decide el área o la región alrededor de los bordes que se va a tratar. III. El umbral (threshold), que elige la región o regiones en la que se mejore la nitidez a partir de ese valor umbral. [14] Una popular máscara 3x3 usada en esta técnica es: [ −1 8−1 8−1 8 −1 8 2 −1 8 −1 8−1 8−1 8 ] (3−21) 3.3.3 Filtro de refuerzo de frecuencias altas Los filtros de altas frecuencia, debido a su naturaleza, son aquellos que refuerzan, tal y como su nombre dice, las altas frecuencias y eliminan o dan menos peso a las frecuencias bajas de la imagen. A nivel visual, resaltan los cambios bruscos de intensidad de la imagen digital, en cambio, reducen, minimizan o incluso eliminan en algunos casos los cambios suaves o degradados. Por tanto, son muy útiles para mejorar, resaltar y separar los bordes de las imágenes. Existen numerosos filtros de alta frecuencia, pero nos vamos a centrar en el filtro de refuerzo de frecuencias altas, más conocido popularmente por su traducción en inglés, High Boost Filter (HBF). Este filtro al ser de frecuencias altas enfatiza los valores de los píxeles de dichas frecuencias sin eliminar los componentes de bajas frecuencias, solo las suaviza. Por lo que son muy útiles en imágenes oscuras. 27 Agudiza los bordes de la imagen gracias a la amplificación de un factor que posee, obteniéndose una imagen más clara. Por lo que se logra una mejora de la nitidez en la imagen digital, sobre todo, en áreas de alto detalle, pero poco en áreas planas o suaves. Este filtrado se puede expresar con la siguiente formula [13]: 𝑓𝐻𝐵(𝑚,𝑛)= 𝐴∗𝑓(𝑚,𝑛)−𝑓𝐿𝑃𝐹(𝑚,𝑛) (3−22) donde 𝑓(𝑚,𝑛) es la imagen original, 𝐴 se denomina factor de amplificación que controla los pesos, 𝑓𝐿𝑃𝐹(𝑚,𝑛) es la imagen filtrada de paso bajo e 𝑓𝐻𝐵(𝑚,𝑛) es finalmente la imagen filtrada con HBF. Además, el valor del factor de amplificación tiene que ser: 𝐴 ≥1 Si ahora sumamos y restamos 𝑓(𝑚,𝑛), la formula queda: 𝑓𝐻𝐵(𝑚,𝑛)=(𝐴−1)∗𝑓(𝑚,𝑛)+𝑓(𝑚,𝑛)−𝑓𝐿𝑃𝐹(𝑚,𝑛) (3−23) Si a la imagen original le eliminamos las componentes de baja frecuencia, quedarían las de altas frecuencias: 𝑓(𝑚,𝑛)−𝑓𝐿𝑃𝐹(𝑚,𝑛)= 𝑓𝐻𝑃𝐹(𝑚,𝑛) (3−24) La fórmula queda entonces: 𝑓𝐻𝐵(𝑚,𝑛)=(𝐴−1)∗𝑓(𝑚,𝑛)+𝑓𝐻𝑃𝐹(𝑚,𝑛) (3−25) 𝑓𝐻𝐵(𝑚,𝑛)=(𝐴−1)∗𝑓(𝑚,𝑛)+𝑓(𝑚,𝑛)∗ℎ𝐻𝑃𝐹(𝑚,𝑛) (3−26) donde ℎ𝐻𝑃𝐹(𝑚,𝑛) es el filtro paso alto. La implementación del filtro de refuerzo de frecuencias altas en una imagen se logra aplicando una máscara adecuada y realizando la operación de convolución. Algunas de las máscaras típicas usada para esta técnica son [13]: [0 −1 0 −1 𝐴+4 −1 0 −1 0] [−1 −1 −1 −1 𝐴+8 −1 −1 −1 −1] (3−27) Estas dos máscaras las vamos a usar para ver cómo se comportan, tomamos el factor de amplificación como 𝐴=1, llamaremos a la primera mask1 y a la segunda mask2. Implementación y evaluación de técnicas de mejoras de la imagen 28 28 3.3.4 Comparativa de los métodos Tabla 3. Comparativa de técnicas de mejora de la imagen. Ventajas Desventajas Ecualización del histograma Consigue mejorar el contrate de la imagen con una distribución uniforme del histograma. Es un método sencillo y automático. Útil en imágenes médicas, como radiografías y TAC, ya que permite separar regiones óseas, musculares, tejidos y órganos. No funciona bien con imágenes con histogramas muy irregulares, ya que no se consigue hacerlo del todo uniforme. Puede producir efectos no deseados cuando la profundidad de color es pequeña. Unsharp masking El enmascaramiento del desenfoque es una técnica bastante útil para aumentar la nitidez de una imagen a través de la aplicación de una máscara. Posee parámetros para ajustar el método y hacerlo así más flexible. Puede crear efectos de borde llamativos no deseados. Puede aumentar el ruido de los píxeles de la imagen. Esta limitado por un factor, si se sobrepasa, se trataría entonces de un filtro High Boost. High Boost Filter Es una técnica para la mejora, sobre todo, de la nitidez de la imagen, al igual que el enmascaramiento del desenfoque. Produce un efecto más agresivo que el enmascaramiento del desenfoque debido a que el factor de amplificación es mayor ( ≥1 ). Funciona mejor en imágenes, regiones o áreas oscuras. El principal inconveniente es que la amplificación afecta tanto a los píxeles de información como los píxeles de ruido. Se aumenta el ruido de la imagen. 29 4 IMPLEMENTACIÓN EN MATLAB Se desarrolla ahora, en la aplicación MATLAB, los métodos usados para la eliminación del ruido y mejora de las imágenes. 4.1 Aplicación de los filtros de media y mediana en MATLAB En este apartado vamos a indicar el código, variables, métodos, etc, usados en MATLAB para eliminar el ruido gaussiano y de sal y pimienta, a través de los filtros de media y mediana. i1 = 'ISIC_0498990.jpg'; i2 = 'ISIC_0984817.jpg'; i4 = 'ISIC_4602463.jpg'; i5 = 'ISIC_6821316.jpg'; f=imread(i1); %Leemos la imagen que deseemos f_gray=rgb2gray(f); %Imagen original en escala de grises %Ruido Gaussiano con media 0 y varianza 0.01 fng_gray=imnoise(f_gray,"gaussian",0,0.01); fng=imnoise(f,"gaussian",0,0.01); %Añadimos ruido sal & pimienta, d=densidad (por defecto es 0.05) fni_gray=imnoise(f_gray,"salt & pepper"); fni=imnoise(f,"salt & pepper"); %Filtro de media sobre la imagen original con ruido gaussiano y ruido sal y pimienta. %Usamos ahora imfilter para aplicar el filtro h(creado por fspecial, el cual será de media) h=fspecial('average'); fmedia_gray=imfilter(fng_gray,h); fmedia=imfilter(fng,h); %Ahora con ruido sal y pimienta fmediaS_gray=imfilter(fni_gray,h); fmediaS=imfilter(fni,h); %% Aplicamos el filtro de mediana sobre la imagen con ruido gaussiano. fmedianaGauss_gray = medfilt2(fng_gray); %en escala de grises %En cada plano R, G y B y encadenamos al final fmedianaGauss_R = medfilt2(fng(:,:,1)); fmedianaGauss_G = medfilt2(fng(:,:,2)); fmedianaGauss_B = medfilt2(fng(:,:,3)); fmedianaGauss_color = cat(3,fmedianaGauss_R,fmedianaGauss_G,fmedianaGauss_B); %% Aplicamos el filtro de mediana sobre la imagen con ruido impulsivo. fmedianaSal_gray = medfilt2(fni_gray); %En cada plano R, G y B y encadenamos al final fmedianaSal_R = medfilt2(fni(:,:,1)); fmedianaSal_G = medfilt2(fni(:,:,2)); fmedianaSal_B = medfilt2(fni(:,:,3)); Implementación y evaluación de técnicas de mejoras de la imagen 30 30 fmedianaSal_color = cat(3,fmedianaSal_R,fmedianaSal_G,fmedianaSal_B); %% Comparativas en imágenes figure('Name','COMPARATIVA COLOR RUIDO GAUSSIANO') subplot(1,4,1),imshow(f); %title('Imagen original') subplot(1,4,2),imshow(fng); %title('Imagen con ruido gaussiano') subplot(1,4,3),imshow(fmediaG); %title('Filtro de media') subplot(1,4,4),imshow(fmedianaGauss_color); %title('Filtro de mediana') figure('Name','COMPARATIVA COLOR RUIDO SAL Y PIMIENTA') subplot(1,4,1),imshow(f); %title('Imagen original') subplot(1,4,2),imshow(fni); %title('Imagen con ruido sal y pimienta') subplot(1,4,3),imshow(fmediaS); %title('Filtro de media') subplot(1,4,4),imshow(fmedianaSal_color);%title('Filtro de mediana') figure('Name','COMPARATIVA RUIDO GAUSSIANO') subplot(1,4,1),imshow(f_gray); %title('Imagen original') subplot(1,4,2),imshow(fng_gray); %title('Imagen con ruido gaussiano') subplot(1,4,3),imshow(fmediaG_gray); %title('Filtro de media') subplot(1,4,4),imshow(fmedianaGauss_gray); %title('Filtro de mediana') figure('Name','COMPARATIVA RUIDO SAL&PIMIENTA') subplot(1,4,1),imshow(f_gray); %title('Imagen original') subplot(1,4,2),imshow(fni_gray); %title('Imagen con ruido sal y pimienta') subplot(1,4,3),imshow(fmediaS_gray); %title('Filtro de media') subplot(1,4,4),imshow(fmedianaSal_gray); %title('Filtro de mediana') • Funciones de MATLAB usadas: imread( filename ): Esta función lee la imagen del archivo especificado por el atributo de entrada filename, además tiene otros atributos de entrada que permiten especificar el formato del archivo (fmt), o incluso leer una imagen o imágenes de un grupo (idx). rgb2gray( RGB ): Este método convierte la imagen en color en el espacio RGB introducida por parámetro de entrada a una imagen en escala de grises. Por tanto, elimina la información de tono y saturación y conserva la de luminosidad. imnoise( I, ‘tipo de ruido’ ): Añade el tipo de ruido indicado como parámetro de entrada a la imagen I, normalizando los valores entre 0 y 1. Además según los tipos de ruido, se pueden especificar otros parámetros, señalamos algunos: - imnoise( I, ‘gausssian’, m, var_gauss): Introduce ruido gaussiano con media m y varianza var_gauss, por defecto la media vale 0 y la varianza 0.01. - imnoise( I, ‘salt & pepper’, d): Añade ruido de sal y pimienta, donde d equivale a la densidad de ruido, por defecto su valor es de 0.05, que quiere decir que afecta solo a 5% de los píxeles. - imnoise( I, ‘poisson’): Genera ruido de disparo a partir de los datos en lugar de añadir ruido artificial a los datos. fspecial ( type ): 31 Esta función crea un filtro en 2 dimensiones del tipo especificado como parámetro de entrada, el valor devuelto por fspecial actúa como núcleo de correlación, que es la forma adecuada para utilizar con otro método que describiremos a continuación, imfilter. Según el tipo de filtro, además admite unos atributos u otros: - fspecial (‘average’, hsize): Genera un filtro promediador de tamaño hsize, si no se especifica el tamaño, por defecto es 3 (3x3). - fspecial (‘disk’, radius): devuelve un filtro promediador circular situado dentro de la matriz de tamaño 2*radius+1. imfilter ( I, h ): Función que realiza el filtrado multidimensional de la imagen I a partir de un filtro indicado como parámetro de entrada, h. medfilt2 ( I ): Este método realiza el filtrado de mediana de la imagen I en 2 dimensiones. Por cada píxel de la imagen I, se obtiene un píxel resultante que contiene el valor de la mediana de un entorno 3x3 del píxel de entrada. cat ( dN, A1, A2, … , AN): Función que permite concatenar elementos A1, A2, … , AN en la dimensión dN. La concatenación se realiza desde el último elemento hacia el primero. Funciones para la representación gráfica de los resultados: figure ( Name, Value): Crea una ventana para representar una figura (puede ser una imagen, una gráfica, resultados, etc). Puede no tener parámetros de entrada, pero se puede modificar los valores de la figura, indicándolo en el parámetro Name y el valor de éste. Ejemplo: figure( ‘Color’, ‘red’). Estos parámetros establecerán el color del fondo como rojo. subplot ( n, m, p): Sirve para crear ejes en posiciones segmentadas, divide la figura en una cuadrícula de m por n y crea ejes en la posición especificada por p. imshow ( I ): Permite mostrar la imagen I. Utiliza el intervalo de visualización predeterminado para el tipo de datos de la imagen y optimiza la figura, los ejes y las propiedades del objeto de imagen para su visualización. Si se incluye [] como parámetro, representa la imagen entre los valores máximos y mínimos. 4.1.1 Comparativa en imágenes. Vamos a ver ahora, representado en la Tabla 4 y Tabla 5, el resultado del apartado anterior. Se observa cómo se comportan tanto el filtro de media cómo el de mediana en imágenes con ruidos gaussianos y ruidos impulsivos, tanto en escala de grises como en color. Usamos algunas de las imágenes del estudio, descargadas de la base de datos de imágenes médicas ISIC. Implementación y evaluación de técnicas de mejoras de la imagen 32 32 - Escala de grises (Plano intensidad): Tabla 4. Comparativa en imágenes, escala de grises. Imagen original Imagen con ruido Gaussiano/Sal y pimienta Filtro de media Filtro de mediana 33 A simple vista, observamos que para para la imagen distorsionada con ruido gaussiano funciona bastante mejor el filtro de media que el filtro de mediana. Y para la imagen distorsionada con ruido sal y pimienta, funciona mucho mejor el filtro de la mediana que el de media. En el apartado 4.6 de este capítulo, esta evidencia además se verá plasmada con las medidas de calidad de la imagen, MSE y PSNR. - En color: Tabla 5. Comparativa en imágenes a color. Imagen original Imagen con ruido Gaussiano/Sal y pimienta Filtro de media Filtro de mediana Implementación y evaluación de técnicas de mejoras de la imagen 34 34 En ambas tablas, y aunque viene indicado en el título de la columna, la primera fila de cada imagen corresponde al tratamiento del ruido gaussiano y la siguiente fila al tratamiento del ruido de sal y pimienta, así sucesivamente para cada imagen. De la misma forma, en color, vemos que el ruido sal y pimienta se elimina mejor con el filtro de mediana que con el filtro de media. Y para el ruido gaussiano, es el filtro de media el más efectivo. Analizaremos en el apartado 4.6 estas conclusiones con las medidas de fidelidad de la imagen, MSE y PSNR. 4.2 Binarización del histograma en MATLAB Vamos a mostrarlo en un ejemplo realizado en MATLAB para la imagen: ISIC_0498990.jpg, obtenida de la base de datos de ISIC. Figura 14. Binarización del histograma. Imágenes con sus correspondientes histogramas debajo. Donde se aprecia las regiones binarizadas a blanco y negro claramente en la imagen, y el número de valores de blanco y negro en el histograma. 35 Realizado con el siguiente código: i1 = 'ISIC_0498990.jpg'; f=imread(i1); %Leemos las imágenes fgray=rgb2gray(f); %Convertimos de RBG a escala de grises fB=imbinarize(fgray); %Binarizamos la imagen en escala de grises %Obtenemos los histogramas de la imagen h= imhist(f); hg= imhist(fgray); hb=imhist(fB); %Funciones que obtienen el valor del umbral, normalizado entre [0 1] To = otsuthresh(hg) Tg = graythresh(fgray) %Representación de los resultados figure('Name','Imagenes Binarizada') subplot(2,3,1),imshow(f);title('Imagen original') subplot(2,3,2),imshow(fgray);title('Imagen en grises') subplot(2,3,3),imshow(fB);title('Imagen binarizada') subplot(2,3,4),plot(h);title('Histrograma') subplot(2,3,5),plot(hg);title('Histrograma en grises') subplot(2,3,6),plot(hb);title('Histrograma binarizado') Descripciones de las funciones usadas aún no mencionadas: imhist ( I ): Calcula el histograma de la imagen I . Devuelve los recuentos del histograma. imbinarize ( I ): Crea una imagen binaria a partir de otra en 2 ó 3 dimensiones, seleccionando un valor umbral, de manera que todos los valores por debajo de ese umbral, le asignará un 0 (negro); y los que estén por encima de ese umbral le asigna un 1 (blanco). Por defecto utiliza el método de Otsu para obtener ese umbral, pero se puede indicar por parámetros un valor concreto. otsuthresh ( counts ): Calcula un umbral global a partir de los recuentos de un histograma, counts, usando el método de Otsu, el cual elige el que minimiza la varianza interclase de los píxeles blancos y negros pasados por el umbral. graythresh ( I ): Calcula un umbral global a partir de una imagen en escala de grises I empleando el método de Otsu. plot ( G ): Representa en una gráfica los nodos y los bordes. Implementación y evaluación de técnicas de mejoras de la imagen 42 42 4.6.1 Con filtro de media Mostraremos en este apartado dos tablas, Tabla 11 y Tabla 12, las cuales representan los resultados de las medidas referenciadas MSE y PSNR, la primera solo aplicando las transformaciones necesarias al plano L del espacio CIE, convirtiendo después al espacio RGB; la segunda aplicando la máscara a cada plano de color RGB. Vamos a comparar la imagen original con la misma imagen pero añadiéndole ruido gaussiano o de sal y pimienta. Y la imagen original con la imagen tras el filtrado de media, a cada caso. Los valores de la MSE están normalizado entre 0 y 1. Los valores de la PSNR están calculados según el rango dinámico de las imágenes, en este caso está entre 0 y 255. Tabla 11. MSE y PSNR. Filtro de media, plano L. Imagen con ruido gaussiano Imagen con filtro de media al ruido gaussiano Imagen con ruido de sal y pimienta Imagen con filtro de media al ruido sal y pimienta media = 0 ; varianza = 0,01 mse = 0,0285 psnr = 63,5842 mse = 0,0272 psnr = 63,7896 densidad = 0,05 mse = 0,0471 psnr = 61,3968 mse = 0,0337 psnr = 62,8591 media = 0,5 ; varianza = 0,1 mse = 0,3300 psnr = 52,9455 mse = 0,2829 psnr = 53,6139 densidad = 0,1 mse = 0,0929 psnr = 58,4497 mse = 0,0635 psnr = 60,1061 media = 0 ; varianza = 0,01 mse = 0,0290 psnr = 63,5117 mse = 0,0264 psnr = 63,9174 densidad = 0,05 mse = 0,0453 psnr = 61,5693 mse = 0,0323 psnr = 63,0403 43 media = 0,5 ; varianza = 0,1 mse = 0,2555 psnr = 54,0577 mse = 0,2094 psnr = 54,9207 densidad = 0,1 mse = 0,0913 psnr = 58,5260 mse = 0,0628 psnr = 60,1501 Tabla 12. MSE y PSNR. Filtro de media en RGB. Imagen con ruido gaussiano Imagen con filtro de media al ruido gaussiano Imagen con ruido de sal y pimienta Imagen con filtro de media al ruido sal y pimienta media = 0 ; varianza = 0,01 mse = 0,0285 psnr = 63,5842 mse = 0,0055 psnr = 70,7576 densidad = 0,05 mse = 0,0471 psnr = 61,3968 mse = 0,0080 psnr = 69,0740 media = 0,5 ; varianza = 0,1 mse = 0,3300 psnr = 52,9455 mse = 0,2829 psnr = 53,6139 densidad = 0,1 mse = 0,0929 psnr = 58,4497 mse = 0,0143 psnr = 66,5829 media = 0 ; varianza = 0,01 mse = 0,0290 psnr = 63,5117 mse = 0,0049 psnr = 71,2588 densidad = 0,05 mse = 0,0453 psnr = 61,5693 mse = 0,0071 psnr = 69,6386 Implementación y evaluación de técnicas de mejoras de la imagen 44 44 media = 0,5 ; varianza = 0,1 mse = 0,2555 psnr = 54,0577 mse = 0,2040 psnr = 55,0352 densidad = 0,1 mse = 0,0913 psnr = 58,5260 mse = 0,0133 psnr = 66,8998 Desarrollado con este código: i1 = 'ISIC_0498990.jpg'; i5 = 'ISIC_6821316.jpg'; f=imread(i1); f=im2double(f); %Normalizamos los valores entre 0 y 1 %Añadimos ruido blanco gaussiano de media y varianza indicada por parámetros fng=imnoise(f,"gaussian",0.5,0.1); %Añadimos ruido impulsivo o salt & pepper, d=densidad % fni=imnoise(f,"salt & pepper",0.1); %% División de la imagen con ruido gaussiano y ruido sal y pimienta en los canales RGB% fngR=fng(:,:,1); fngG=fng(:,:,2); fngB=fng(:,:,3); fniR=fni(:,:,1); fniG=fni(:,:,2); fniB=fni(:,:,3); %% Aplicamos el filtro (con la función fspecial) sobre la imagen original con ruido gaussiano y ruido sal y pimienta. % Usamos ahora imfilter para aplicar el filtro h(creado por fspecial, el cual será de media) h=fspecial('average'); fmediaG=imfilter(fng,h); fmediaS=imfilter(fni,h); %imfilter a cada plano: %Para el ruido gaussiano fmediaGaussR=imfilter(fng(:,:,1),h); fmediaGaussG=imfilter(fng(:,:,2),h); fmediaGaussB=imfilter(fng(:,:,3),h); fmediaGaussf=cat(3,fmediaGaussR,fmediaGaussG,fmediaGaussB); %encadenamos los planos %Para el ruido de sal y pimienta fmediaSalR=imfilter(fni(:,:,1),h); fmediaSalG=imfilter(fni(:,:,2),h); fmediaSalB=imfilter(fni(:,:,3),h); fmediaSalf=cat(3,fmediaSalR,fmediaSalG,fmediaSalB); %encadenamos los planos %Pintamos los resultados% figure('Name','Filtro de media con fspecial a img original') subplot(3,3,[1,4]),imshow(f); title('Imagen original') subplot(3,3,2),imshow(fng);title('Imagen con ruido gaussiano') subplot(3,3,3),imshow(fni);title('Imagen con ruido sal y pimienta') subplot(3,3,5),imshow(fmediaGaussf);title('Imagen filtrada por plano de ruido gaussiano') subplot(3,3,6),imshow(fmediaSalf);title('Imagen filtrada por plano de ruido sal y pimienta') 45 %% Aplicamos el filtro de media al plano L. % Convertimos la imagen de RGB a CIE. Imagen Original. fCIE=rgb2lab(f); %Al cambiar del espacio rgb a lab, cambia el rango dinámico de la imagen entre 0 y 100, por eso después lo representamos entre esos valores. fCIEL=fCIE(:,:,1); % Obtenemos el plano L (escala de grises). % Imagen con ruido GAUSSIANO fngCIE=rgb2lab(fng); fngCIEL=fngCIE(:,:,1); % Obtenemos el plano L (escala de grises). % Imagen con ruido SAL Y PIMIENTA fniCIE=rgb2lab(fni); fniCIEL=fniCIE(:,:,1); % Obtenemos el plano L (escala de grises). % Imagen filtrada de ruido gaussiano (filtro de media) al plano L fmediaCIEG=imfilter(fngCIEL,h); fmediaCIEGf=cat(3,fmediaCIEG,fngCIE(:,:,2),fngCIE(:,:,3)); %Encadenamos los planos. fmediaCIEG_RGB=lab2rgb(fmediaCIEGf); %Volvemos al espacio RGB fmediaCIEG_RGB(fmediaCIEG_RGB > 1.0 ) = 1;%Umbralizamos los valores max y min a 0 y 1 fmediaCIEG_RGB(fmediaCIEG_RGB < 0 ) = 0; % Imagen filtrada de ruido sal y pimienta (filtro de media) al plano L fmediaCIES=imfilter(fniCIEL,h); fmediaCIESf=cat(3,fmediaCIES,fniCIE(:,:,2),fniCIE(:,:,3)); %Encadenamos los planos fmediaCIES_RGB=lab2rgb(fmediaCIESf); %Al transformar otra vez a rgb se pierde el rango dinamico, ¿qué hacer? fmediaCIES_RGB(fmediaCIES_RGB > 1.0 ) = 1; fmediaCIES_RGB(fmediaCIES_RGB < 0 ) = 0; %% Representación en RGB tras media en el plano L de Lab. imshow(f,[]) cambia el rango dinámico entre los valores máximo y mínimo de la imagen. figure('Name','RGB To Lab 3 planos') subplot(3,3,[1,4]),imshow(f,[]); title('Imagen original') %Obtenemos el plano L subplot(3,3,2),imshow(fng,[]);title('Imagen con ruido gaussiano') subplot(3,3,3),imshow(fni,[]);title('Imagen con ruido sal y pimienta') subplot(3,3,5),imshow(fmediaCIEG_RGB,[]);title('Lab:Filtrado Gauss paso a RGB') subplot(3,3,6),imshow(fmediaCIES_RGB,[]);title('Lab:Filtrado S&P paso a RGB') %% ---------------------------- MEDIDADES DE FIDELIDAD ------------------------------ % Vamos calcular la MSE y la PSNR como medidas de fidelidad de las imágenes %Calculamos el error cuadrático en el espacio RGB tras solo aplicar la transformación al plano L squaredErrorImageG = double((f(:,:,1)-fng(:,:,1)) .^ 2+(f(:,:,2)-fng(:,:,2)) .^ 2+(f(:,:,3)-fng(:,:,3)) .^ 2); %Entre la imagen original y con ruido gaussiano. squaredErrorImageS = double((f(:,:,1)-fni(:,:,1)) .^ 2+(f(:,:,2)-fni(:,:,2)) .^ 2+(f(:,:,3)-fni(:,:,3)) .^ 2); %Entre la imagen original y con ruido sal y pimienta. squaredErrorImageGf = double((f(:,:,1)-fmediaCIEG_RGB(:,:,1)) .^ 2+(f(:,:,2)- fmediaCIEG_RGB(:,:,2)) .^ 2+(f(:,:,3)-fmediaCIEG_RGB(:,:,3)) .^ 2); %Entre la imagen original y filtrada del ruido gaussiano. squaredErrorImageSf = double((f(:,:,1)-fmediaCIES_RGB(:,:,1)) .^ 2+(f(:,:,2)- fmediaCIES_RGB(:,:,2)) .^ 2+(f(:,:,3)-fmediaCIES_RGB(:,:,3)) .^ 2); %Entre la imagen original y filtrada del ruido sal y pimienta. %Error cuadrático en el espacio RGB, con ruido gaussiano y de sal y pimienta squaredErrorImagecolorG = double((f(:,:,1)-fng(:,:,1)) .^ 2+(f(:,:,2)-fng(:,:,2)) .^ 2+(f(:,:,3)-fng(:,:,3)) .^ 2); %Entre la imagen original y con ruido gaussiano. squaredErrorImagecolorS = double((f(:,:,1)-fni(:,:,1)) .^ 2+(f(:,:,2)-fni(:,:,2)) .^ 2+(f(:,:,3)-fni(:,:,3)) .^ 2); %Entre la imagen original y con ruido sal y pimienta. Implementación y evaluación de técnicas de mejoras de la imagen 46 46 squaredErrorImagecolorGf = double((f(:,:,1)-fmediaGaussf(:,:,1)) .^ 2+(f(:,:,2)- fmediaGaussf(:,:,2)) .^ 2+(f(:,:,3)-fmediaGaussf(:,:,3)) .^ 2); %Entre la imagen original y filtrada del ruido gaussiano squaredErrorImagecolorSf = double((f(:,:,1)-fmediaSalf(:,:,1)) .^ 2+(f(:,:,2)- fmediaSalf(:,:,2)) .^ 2+(f(:,:,3)-fmediaSalf(:,:,3)) .^ 2); %Entre la imagen original y filtrada del ruido sal y pimienta % Obtenemos los píxeles de la imagen. [rows columns plane] = size(f); %Dividimos en filas, columnas y planos. %% Calculamos el MSE %Sumamos el error cuadrático de la imagen y dividimos por el número de elementos para obtener el error cuadrático medio. Será un escalar. %Tras aplicar al plano de intensidad mseG = sum(sum(squaredErrorImageG)) / (rows * columns); % MSE entre la imagen original y con ruido gaussiano mseS = sum(sum(squaredErrorImageS)) / (rows * columns); % MSE entre la imagen original y con ruido sal y pimienta mseGf = sum(sum(squaredErrorImageGf)) / (rows * columns); % MSE entre la imagen original y filtrada del ruido gaussiano mseSf = sum(sum(squaredErrorImageSf)) / (rows * columns); % MSE entre la imagen original y filtrada del ruido sal y pimienta %A color msecolorG = sum(sum(squaredErrorImagecolorG)) / (rows * columns); % MSE entre la imagen original en color y con ruido gaussiano msecolorS = sum(sum(squaredErrorImagecolorS)) / (rows * columns); % MSE entre la imagen original en color y con ruido sal y pimienta msecolorGf = sum(sum(squaredErrorImagecolorGf)) / (rows * columns); % MSE entre la imagen original en color y filtrada del ruido gaussiano msecolorSf = sum(sum(squaredErrorImagecolorSf)) / (rows * columns); % MSE entre la imagen original en color y filtrada del ruido sal y pimienta %% Calculamos la PSNR (Peak Signal to Noise Ratio) a partir del MSE de acuerdo con la fórmula. %Como estamos comparando la imagen original con la filtrada en el plano L, y convertimos en RGB, ponemos 255. PSNRg = 10 * log10( 255^2 / mseG); %De imagen original y ruido gaussiano PSNRs = 10 * log10( 255^2 / mseS); %De imagen original y ruido sal y pimienta PSNRgf = 10 * log10( 255^2 / mseGf);%De imagen original y filtrada de ruido gaussiano PSNRsf = 10 * log10( 255^2 / mseSf);%De imagen original y filtrada de ruido sal y pimienta %A color, en el espacio RGB si va de 0 a 256 PSNRcolorg = 10 * log10( 255^2 / msecolorG); % De imagen original en color y ruido gaussiano PSNRcolors = 10 * log10( 255^2 / msecolorS); % De imagen original en color y ruido sal y pimienta PSNRcolorgf = 10 * log10( 255^2 / msecolorGf); % De imagen original en color y filtrada de ruido gaussiano PSNRcolorsf = 10 * log10( 255^2 / msecolorSf); % De imagen original en color y filtrada de ruido sal y pimienta % Alertas con los resultados. %Tras aplicar al plano de intensidad 47 message1 = sprintf('\nMSE entre la imagen original y con ruido gaussiano = %.4f.\nThe PSNR = %.4f', mseG, PSNRg); message2 = sprintf('\nMSE entre la imagen original y con ruido sal y pimienta = %.4f.\nThe PSNR = %.4f', mseS, PSNRs); message3 = sprintf('\nMSE entre la imagen original y filtrada del ruido gaussiano = %.4f.\nThe PSNR = %.4f', mseGf, PSNRgf); message4 = sprintf('\nMSE entre la imagen original y filtrada del ruido sal y pimienta = %.4f.\nThe PSNR = %.4f', mseSf, PSNRsf); %A color message5 = sprintf('\nMSE entre la imagen original en color y con ruido gaussiano = %.4f.\nThe PSNR = %.4f', msecolorG, PSNRcolorg); message6 = sprintf('\nMSE entre la imagen original en color y con ruido sal y pimienta = %.4f.\nThe PSNR = %.4f', msecolorS, PSNRcolors); message7 = sprintf('\nMSE entre la imagen original en color y filtrada del ruido gaussiano = %.4f.\nThe PSNR = %.4f', msecolorGf, PSNRcolorgf); message8 = sprintf('\nMSE entre la imagen original en color y filtrada del ruido sal y pimienta = %.4f.\nThe PSNR = %.4f', msecolorSf, PSNRcolorsf); %Representación en cajas msgbox(message1,'Original/Ruido gaussiano'); msgbox(message2,'Original/Ruido sal y pimienta'); msgbox(message3,'Original/Fmedia gaussiano'); msgbox(message4,'Original/Fmedia sal y pimienta'); msgbox(message5,'Original/Ruido gaussiano espacio RGB'); msgbox(message6,'Original/Ruido sal y pimienta espacio RGB'); msgbox(message7,'Original/Fmedia gaussiano espacio RGB'); msgbox(message8,'Original/Fmedia sal y pimienta espacio RGB'); • Discusión de los resultados De forma general, observamos que tras convertir al espacio CIE Lab y aplicar sólo la máscara al plano L, y además dependiendo de la cantidad de píxeles con ruido que haya en la imagen y de sus características propias, el filtro de media elimina gran parte de ese ruido, mejor con ruido gaussiano. Por otro lado, al aplicar el filtro de media a cada plano R, G y B en lugar de uno solo, se reduce bastante más el error entre la imagen original y la imagen filtrada, por lo que la MSE es menor y la PSNR mayor. 4.6.2 Con filtro de mediana Al igual que en el apartado 4.6.1, vamos a representar ahora los resultados de las medidas de fidelidad MSE y PSNR, en la Tabla 13 y Tabla 14, entre la imagen original y la misma imagen pero con ruido gaussiano o de sal y pimienta, y entre la imagen original y la imagen filtrada con el filtro de mediana. Tanto para el plano de intensidad L como en el espacio RGB. Para la imagen en color aplicamos el filtro en cada plano, tal como anuncia el método de ordenamiento marginal. Implementación y evaluación de técnicas de mejoras de la imagen 48 48 Tabla 13. MSE y PSNR. Filtro de mediana, plano L. Imagen con ruido gaussiano Imagen con filtro de mediana al ruido gaussiano Imagen con ruido de sal y pimienta Imagen con filtro de mediana al ruido sal y pimienta media = 0 ; varianza = 0,01 mse = 0,0285 psnr = 63,5762 mse = 0,0260 psnr = 63,9787 densidad = 0,05 mse = 0,0475 psnr = 61,3683 mse = 0,0292 psnr = 63,4754 media = 0,5 ; varianza = 0,1 mse = 0,3312 psnr = 52,9295 mse = 0,3226 psnr = 53,0446 densidad = 0,1 mse = 0,0933 psnr = 58,4310 mse = 0,0564 psnr = 60,6195 media = 0 ; varianza = 0,01 mse = 0,0290 psnr = 63,5012 mse = 0,0256 psnr = 64,0451 densidad = 0,05 mse = 0,0458 psnr = 61,5216 mse = 0,0289 psnr = 63,5174 media = 0,5 ; varianza = 0,1 mse = 0,2555 psnr = 54,0569 mse = 0,2482 psnr = 54,1821 densidad = 0,1 mse = 0,0930 psnr = 58,4461 mse = 0,0580 psnr = 60,4930 49 Tabla 14. MSE y PSNR. Filtro de mediana, en RGB. Imagen con ruido gaussiano Imagen con filtro de media al ruido gaussiano Imagen con ruido de sal y pimienta Imagen con filtro de media al ruido sal y pimienta media = 0 ; varianza = 0,01 mse = 0,0285 psnr = 63,5762 mse = 0,0056 psnr = 70,6347 densidad = 0,05 mse = 0,0475 psnr = 61,3683 mse = 0,0005 psnr = 81,0822 media = 0,5 ; varianza = 0,1 mse = 0,3312 psnr = 52,9295 mse = 0,3579 psnr = 52,5929 densidad = 0,1 mse = 0,0933 psnr = 58,4310 mse = 0,0012 psnr = 77,5211 media = 0 ; varianza = 0,01 mse = 0,0290 psnr = 63,5012 mse = 0,0053 psnr = 70,8567 densidad = 0,05 mse = 0,0458 psnr = 61,5216 mse = 0,0002 psnr = 84,6893 media = 0,5 ; varianza = 0,1 mse = 0,2555 psnr = 54,0569 mse = 0,2855 psnr = 53,5742 densidad = 0,1 mse = 0,0930 psnr = 58,4461 mse = 0,0007 psnr = 79,9909 Implementación y evaluación de técnicas de mejoras de la imagen 50 50 Programado con: i1 = 'ISIC_0498990.jpg'; i5 = 'ISIC_6821316.jpg'; f=imread(i1); f=im2double(f); %Normalizamos los valores entre 0 y 1 %Añadimos ruido blanco gaussiano de media y varianza indicada por parámetros fng=imnoise(f,"gaussian",0.5,0.1); %Añadimos ruido impulsivo o salt & pepper, d=densidad % fni=imnoise(f,"salt & pepper",0.1); %Aplicamos el filtro de mediana a cada plano (se conoce como ordenamiento marginal) sobre la imagen con ruido gaussiano. fmedianaR_Gauss = medfilt2(fng(:,:,1)); fmedianaG_Gauss = medfilt2(fng(:,:,2)); fmedianaB_Gauss = medfilt2(fng(:,:,3)); fmedianaGauss = cat(3,fmedianaR_Gauss,fmedianaG_Gauss,fmedianaB_Gauss); %Aplicamos el filtro de mediana a cada plano sobre la imagen con ruido impulsivo. fmedianaR_Sal = medfilt2(fni(:,:,1)); fmedianaG_Sal = medfilt2(fni(:,:,2)); fmedianaB_Sal = medfilt2(fni(:,:,3)); fmedianaSal = cat(3,fmedianaR_Sal,fmedianaG_Sal,fmedianaB_Sal); %Pintamos los resultados del filtrado al ruido figure('Name','Filtro de mediana a imagen con ruido') subplot(3,3,[1,4]),imshow(f); title('Imagen original') subplot(3,3,2),imshow(fng);title('Imagen con ruido gaussiano') subplot(3,3,3),imshow(fni);title('Imagen con ruido sal y pimienta') subplot(3,3,5),imshow(fmedianaGauss);title('Imagen filtrada de ruido gaussiano') subplot(3,3,6),imshow(fmedianaSal);title('Imagen filtrada de ruido sal y pimienta') %% Convertimos la imagen de RGB a CIE. % Imagen Original fCIE=rgb2lab(f); % Al cambiar del espacio rgb a lab, cambia el rango dinámico de la imagen entre 0 y 100, por eso después lo representamos entre esos valores. fCIEL=fCIE(:,:,1); % Obtenemos el plano L (escala de grises) % Imagen con ruido GAUSSIANO fngCIE=rgb2lab(fng); fngCIEL=fngCIE(:,:,1); % Obtenemos el plano L (escala de grises) % Imagen con ruido SAL Y PIMIENTA fniCIE=rgb2lab(fni); fniCIEL=fniCIE(:,:,1); % Obtenemos el plano L (escala de grises) % Imagen filtrada de ruido gaussiano (filtro de mediana) al plano L fmedianaCIELG=medfilt2(fngCIEL); fmedianaCIEGf=cat(3,fmedianaCIELG,fngCIE(:,:,2),fngCIE(:,:,3));%Encadenamos los planos fmedianaCIEGf_RGB=lab2rgb(fmedianaCIEGf); %Pasamos a RGB fmedianaCIEGf_RGB(fmedianaCIEGf_RGB > 1.0 ) = 1; %Umbralizamos los valores max y min a 0 y 1 fmedianaCIEGf_RGB(fmedianaCIEGf_RGB < 0 ) = 0; % Imagen filtrada de ruido sal y pimienta (filtro de mediana) al plano L fmedianaCIELS=medfilt2(fniCIEL); fmedianaCIESf=cat(3,fmedianaCIELS,fniCIE(:,:,2),fniCIE(:,:,3));%Encadenamos los planos fmedianaCIESf_RGB=lab2rgb(fmedianaCIESf); %Pasamos a RGB fmedianaCIESf_RGB(fmedianaCIESf_RGB > 1.0 ) = 1; %Umbralizamos los valores max y min a 0 y 1 fmedianaCIESf_RGB(fmedianaCIESf_RGB < 0 ) = 0; 51 % Representación en RGB tras media en el plano L de Lab. figure('Name','RGB To Lab 3 planos') subplot(3,3,[1,4]),imshow(f,[]); title('Imagen original') %Obtenemos el plano L subplot(3,3,2),imshow(fng,[]);title('Imagen con ruido gaussiano') subplot(3,3,3),imshow(fni,[]);title('Imagen con ruido sal y pimienta') subplot(3,3,5),imshow(fmedianaCIEGf_RGB,[]);title('Lab:Filtrado Gauss paso a RGB') subplot(3,3,6),imshow(fmedianaCIESf_RGB,[]);title('Lab:Filtrado S&P paso a RGB') %----------------------- MEDIDADES DE FIDELIDAD ------------------------------------ % Vamos calcular la MSE y la PSNR como medidas de fidelidad de las imágenes %Calculamos el error cuadrático en el espacio RGB tras solo aplicar la transformacion al plano L squaredErrorImageG = double((f(:,:,1)-fng(:,:,1)) .^ 2+(f(:,:,2)-fng(:,:,2)) .^ 2+(f(:,:,3)-fng(:,:,3)) .^ 2); %Entre la imagen original y con ruido gaussiano squaredErrorImageS = double((f(:,:,1)-fni(:,:,1)) .^ 2+(f(:,:,2)-fni(:,:,2)) .^ 2+(f(:,:,3)-fni(:,:,3)) .^ 2); %Entre la imagen original y con ruido sal y pimienta squaredErrorImageGf = double((f(:,:,1)-fmedianaCIEGf_RGB(:,:,1)) .^ 2+(f(:,:,2)- fmedianaCIEGf_RGB(:,:,2)) .^ 2+(f(:,:,3)-fmedianaCIEGf_RGB(:,:,3)) .^ 2); %Entre la imagen original y filtrada del ruido gaussiano squaredErrorImageSf = double((f(:,:,1)-fmedianaCIESf_RGB(:,:,1)) .^ 2+(f(:,:,2)- fmedianaCIESf_RGB(:,:,2)) .^ 2+(f(:,:,3)-fmedianaCIESf_RGB(:,:,3)) .^ 2); %Entre la imagen original y filtrada del ruido sal y pimienta %Error cuadrático en el espacio RGB squaredErrorImagecolorG = double((f(:,:,1)-fng(:,:,1)) .^ 2+(f(:,:,2)-fng(:,:,2)) .^ 2+(f(:,:,3)-fng(:,:,3)) .^ 2); %Entre la imagen original y con ruido gaussiano. squaredErrorImagecolorS = double((f(:,:,1)-fni(:,:,1)) .^ 2+(f(:,:,2)-fni(:,:,2)) .^ 2+(f(:,:,3)-fni(:,:,3)) .^ 2); %Entre la imagen original y con ruido sal y pimienta. squaredErrorImagecolorGf = double((f(:,:,1)-fmedianaGauss(:,:,1)) .^ 2+(f(:,:,2)- fmedianaGauss(:,:,2)) .^ 2+(f(:,:,3)-fmedianaGauss(:,:,3)) .^ 2); %Entre la imagen original y filtrada del ruido gaussiano. squaredErrorImagecolorSf = double((f(:,:,1)-fmedianaSal(:,:,1)) .^ 2+(f(:,:,2)- fmedianaSal(:,:,2)) .^ 2+(f(:,:,3)-fmedianaSal(:,:,3)) .^ 2); %Entre la imagen original y filtrada del ruido sal y pimienta. % Obtenemos los píxeles de la imagen. [rows columns plane] = size(f); %Dividimos en filas, columnas y planos. %% Calculamos el MSE. Sumamos el error cuadrático de la imagen y dividimos por el número de elementos para obtener el error cuadrático medio. Será un escalar. %Tras aplicar al plano de intensidad mseG = sum(sum(squaredErrorImageG)) / (rows * columns); % MSE entre la imagen original y con ruido gaussiano. mseS = sum(sum(squaredErrorImageS)) / (rows * columns); % MSE entre la imagen original y con ruido sal y pimienta. mseGf = sum(sum(squaredErrorImageGf)) / (rows * columns); % MSE entre la imagen original y filtrada del ruido gaussiano mseSf = sum(sum(squaredErrorImageSf)) / (rows * columns); % MSE entre la imagen original y filtrada del ruido sal y pimienta %A color msecolorG = sum(sum(squaredErrorImagecolorG)) / (rows * columns); % MSE entre la imagen original en color y con ruido gaussiano msecolorS = sum(sum(squaredErrorImagecolorS)) / (rows * columns); % MSE entre la imagen original en color y con ruido sal y pimienta Implementación y evaluación de técnicas de mejoras de la imagen 58 58 • Discusión de los resultados Vemos que, en todos los casos, al aplicar la ecualización del histograma, el colorido empeora. Sin embargo, con los métodos USM y HBF hay una ligera mejora, dependiendo de las características de cada imagen. Tiene sentido ya que ninguno de estos métodos se centran en mejorar el colorido de una imagen. 4.7.3 Para la nitidez En el capítulo 2 describimos esta característica como el emborronado que poseen los bordes de los elementos de las imágenes digitales. Para ello usamos el gradiente de cada plano que la componen. Vamos a ver cómo se comporta cada método con esta medida. Desarrollo en MATLAB: i2 = 'ISIC_0984817.jpg'; i4 = 'ISIC_4602463.jpg'; i6 = 'ISIC_6948176.jpg'; f=imread(i6); fdouble= im2double(f); %Transformamos la imagen a doble y normalizada. %Separamos los planos RGB. r=fdouble(:,:,1); g=fdouble(:,:,2); b=fdouble(:,:,3); %% Obtenemos el valor de sharpness de la imagen inicial de cada plano para calcularle el gradiente. sharpness_originalR=estimate_sharpness(fdouble(:,:,1)); sharpness_originalG=estimate_sharpness(fdouble(:,:,2)); sharpness_originalB=estimate_sharpness(fdouble(:,:,3)); sharpness_original=sharpness_originalR + sharpness_originalG + sharpness_originalB %% Medidas de sharpness tras ecualización de histograma fHSV= rgb2hsv(fdouble); %Convertimos al espacio HSV feV=histeq(fHSV(:,:,3));%Calculamos la ecualización del histograma de la componente V fefinV=cat(3,fHSV(:,:,1),fHSV(:,:,2),feV); feHSV_RGB=hsv2rgb(fefinV); sharpness_ecualizaR=estimate_sharpness(feHSV_RGB(:,:,1)); sharpness_ecualizaG=estimate_sharpness(feHSV_RGB(:,:,2)); sharpness_ecualizaB=estimate_sharpness(feHSV_RGB(:,:,3)); sharpness_ecualiza=sharpness_ecualizaR + sharpness_ecualizaG + sharpness_ecualizaB %Máscara UNSHARP MASKING maskUSM=[-1/8 -1/8 -1/8;-1/8 2 -1/8;-1/8 -1/8 -1/8]; %Máscara HIGH BOOST maskHBF=[-1 -1 -1;-1 9 -1;-1 -1 -1]; %% Filtramos con el método USM a cada componente RGB. gUSM_R=filter2(maskUSM,r); gUSM_G=filter2(maskUSM,g); gUSM_B=filter2(maskUSM,b); gUSMfin=cat(3,gUSM_R,gUSM_G,gUSM_B); %Encadenamos las 3 componentes a la imagen. gUSMfin(gUSMfin > 1.0 ) = 1; gUSMfin(gUSMfin < 0 ) = 0; % Sharpness de la imagen tratada con USM. sharpness_usm3R=estimate_sharpness(gUSMfin(:,:,1)); sharpness_usm3G=estimate_sharpness(gUSMfin(:,:,2)); 59 sharpness_usm3B=estimate_sharpness(gUSMfin(:,:,3)); sharpness_usm3=sharpness_usm3R + sharpness_usm3G + sharpness_usm3B %% Filtramos con el método USM en el plano V gUSM_V=filter2(maskUSM,fHSV(:,:,3)); gUSM_Vfin=cat(3,fHSV(:,:,1),fHSV(:,:,2),gUSM_V); %Encadenamos las 3 componentes gUSM_HSV_RGB=hsv2rgb(gUSM_Vfin); gUSM_HSV_RGB(gUSM_HSV_RGB > 1.0 ) = 1; %Umbralizamos los valores max y min a 0 y 1 gUSM_HSV_RGB(gUSM_HSV_RGB < 0 ) = 0; % Sharpness de la imagen tratada con USM. sharpness_usmR=estimate_sharpness(gUSM_HSV_RGB(:,:,1)); sharpness_usmG=estimate_sharpness(gUSM_HSV_RGB(:,:,2)); sharpness_usmB=estimate_sharpness(gUSM_HSV_RGB(:,:,3)); sharpness_usm=sharpness_usmR + sharpness_usmG + sharpness_usmB %% Filtramos con el método HBF a cada componente RGB gHBF_R=filter2(maskHBF,r); gHBF_G=filter2(maskHBF,g); gHBF_B=filter2(maskHBF,b); gHBFfin=cat(3,gHBF_R,gHBF_G,gHBF_B); %Encadenamos las 3 componentes a la imagen gHBFfin(gHBFfin > 1.0 ) = 1; gHBFfin(gHBFfin < 0 ) = 0; % Sharpness de la imagen tratada con HBF sharpness_hbf3R=estimate_sharpness(gHBFfin(:,:,1)); sharpness_hbf3G=estimate_sharpness(gHBFfin(:,:,2)); sharpness_hbf3B=estimate_sharpness(gHBFfin(:,:,3)); sharpness_hbf3=sharpness_hbf3R + sharpness_hbf3G + sharpness_hbf3B %% Filtramos con el método HBF en el plano V gHBF_V=filter2(maskHBF,fHSV(:,:,3)); gHBF_Vfin=cat(3,fHSV(:,:,1),fHSV(:,:,2),gHBF_V);%Encadenamos las 3 componentes gHBF_HSV_RGB=hsv2rgb(gHBF_Vfin); gHBF_HSV_RGB(gHBF_HSV_RGB > 1.0 ) = 1; gHBF_HSV_RGB(gHBF_HSV_RGB < 0 ) = 0; % Sharpness de la imagen tratada con HBF sharpness_hbfR=estimate_sharpness(gHBF_HSV_RGB(:,:,1)); sharpness_hbfG=estimate_sharpness(gHBF_HSV_RGB(:,:,2)); sharpness_hbfB=estimate_sharpness(gHBF_HSV_RGB(:,:,3)); sharpness_hbf=sharpness_hbfR + sharpness_hbfG + sharpness_hbfB %% Representación gráfica figure('Name','Imágenes') subplot(3,2,1),imshow(f);title('Imagen Original') subplot(3,2,2),imshow(feHSV_RGB);title('Imagen ecualizada') subplot(3,2,3),imshow(gUSMfin,[]);title('Imagen filtrada USM en 3 planos') subplot(3,2,4),imshow(gHBFfin,[]);title('Imagen filtrada HBF en 3 planos') subplot(3,2,5),imshow(gUSM_HSV_RGB,[]);title('Imagen filtrada USM en V') subplot(3,2,6),imshow(gHBF_HSV_RGB,[]);title('Imagen filtrada HBF en V') %% Función que calcula la medida de sharpness a través del gradiente function [sharpness]=estimate_sharpness(G) [Gx, Gy]=gradient(G); S=sqrt(Gx.*Gx+Gy.*Gy); sharpness=sum(sum(S))./(numel(Gx)); %Promedio del módulo del gradiente. End Implementación y evaluación de técnicas de mejoras de la imagen 60 60 Mostramos los resultados en la Tabla 17. Tabla 17. Medidas de la nitidez. 61 • Discusión de los resultados Hemos aplicado las máscaras de USM y HBF en primer lugar a cada plano R, G y B, obteniendo las medidas sharpness_usm3 y sharpness_hbf; y en segundo lugar sólo al plano de intensidad V del espacio HSV. Tanto el método USM como HBF mejoran la nitidez de los bordes de la imagen, en menor medida el primero que el segundo, de forma general, al realzar las frecuencias altas. Por tanto, el método HBF es el que arroja mejores resultados en la mejora de la nitidez y por tanto el emborronado de los bordes. Implementación y evaluación de técnicas de mejoras de la imagen 62 62 5 CONCLUSIONES Para finalizar, podemos concluir que las medidas de calidad de la imagen, tanto referenciadas como no referenciadas, son de gran ayuda a la hora de comparar, analizar y discutir cuánto de “buena” puede ser una imagen, si existe mucho error entre los píxeles originales y los actuales, cómo es su contraste, el grado de intensidad que tiene en relación al blanco, su saturación, su nitidez en los bordes, etc. Vamos a mostrar en la Tabla 18 las conclusiones generales para cada medidad de calidad. Tabla 18. Conclusiones finales. MSE y PSNR El contraste La croma La nitidez Filtro de media Se mejora la imagen en cuanto al error, sobre todo si posee un ruido gaussiano. Filtro de mediana Funciona mucho mejor si la imagen tiene imperfecciones del tipo sal y pimienta. Ecualización del histograma Es sin duda el mejor método de mejora para esta característica, ya que medimos el parecido de los píxeles solo en el plano de intensidad. Se reduce al aplicar este método, ya que solo se actúa sobre el plano de intensidad para sumar la media y la desviación típica. Puede mejorar un poco, dependiendo de la imagen pero la nitidez de los bordes se resaltan mejor con los métodos USM y HBF. Unsharp masking Es el método que peor funciona para el contraste, no existe mucha mejora. No se ven apenas modificadas puesto que la mejora se aplica al plano intensidad. Este método junto con HBF son generalmente los mejores para esta característica. High Boost filtering Mejora bastante el contraste de la imagen si ésta posee frecuencias altas, ya que filtramos con una máscara que las enfatiza. Lo mismo sucede con este método, apenas hay variación. Es el mejor método para mejorar el emborronado de los bordes y por tanto la nitidez. 63 REFERENCIAS [1] Rangaraj M. Rangayyan, Begoña Acha, Carmen Serrano. Color image processing with biomedical applications, pages 117-118. [2] Winkler S. Visual fidelity and perceived quality: towards comprehensive metrics. In Rogowitz BE and Pappas TN, editors, Proceedings of SPIE, volume 4299, pages 114-125, 2001. doi:10.1117/12.429540. [3] A.N. Netravali and B.G. Haskell, Digital Pictures: Representation, Compression, and Standards (2nd Ed), Plenum Press, New York, NY (1995). [4] Michelson, A. (1927). Studies in Optics. U. of Chicago Press. [5] E. Peli (Oct 1990). "Contrast in Complex Images". Journal of the Optical Society of America A. 7 (10): 2032–2040. [6] Fedorovskaya EA, Ridder H, and Blommaert FJJ. Chroma variations and perceived quality of color images of natural scenes. Color Research and Application, 22(2): 96-110,1997. [7] Rafael C. Gonzalez, Richard E. Woods, Steven L. Eddins. Digital Image Processing Using MATLAB, pages 320-321. [8] Rafael C. Gonzalez, Richard E. Woods. Digital Image Processing, pages 320-323. [9] Rangaraj M. Rangayyan, Begoña Acha, Carmen Serrano. Color image processing with biomedical applications, pages 175-178. [10] Rangaraj M. Rangayyan, Begoña Acha, Carmen Serrano. Color image processing with biomedical applications, pages 215-216. [11] Rafael C. Gonzalez, Richard E. Woods. Digital Image Processing, pages 134-140. [12] Carl-Fredrik Westin Ron Kikinis, Hans Knutsson, in Handbook of Medical Imaging, 2000. Adaptive Image Filtering. [13] Kamal Sharma, Akhil Gupta, Bandana Sharma, Suman Lata Tripathi. Intelligent Communication and Automation Systems, pages 70-72. [14] Rangaraj M. Rangayyan, Begoña Acha, Carmen Serrano. Color image processing with biomedical applications, pages 229-233. Implementación y evaluación de técnicas de mejoras de la imagen 64 64