Deblurring-numpy

La presente repository ha puramente scopi educativi. Consiste in una semplice implementazione in Python3 di un algoritmo di deconvoluzione per recuperare un'immagine corrotta da una sfocatura gaussiana ("deblurring"). Siete liberi di utilizzare il codice come meglio credete, citate l'autore originale (io) in caso di pubblicazione e nel caso vi fosse utile fatemelo sapere.

Prerequisiti

• Python 3.11.1
• Numpy 1.23.0
• Matplotlib 3.8.4

Ho cercato di limitare al minimo le librerie esterne, solo Numpy è necessario (semplicità di gestione di array 2D e presenza di DFT), ma non escludo un giorno di creare uno script totalmente indipendente. Numpy può essere installato con pip3 install numpy. Se volete salvare dei grafici anche Matplotlib deve essere installato con pip3 install matplotlib.

Cenni teorici

La trasformata di Fourier discreta (DFT) di un segnale 2D $f(x,y)$ (un'immagine, per esempio) con dimensioni $N \times M$ è:

$$ F(u,v) = \mathcal{F} [ f(x,y) ] = \sum_{x=0}^{N-1} \sum_{y=0}^{M-1} f(x,y) exp[-j 2 \pi (\frac{ux}{M} + \frac{vy}{N}) ] $$

La convoluzione discreta 2D è:

$$ g(x,y)= \omega(x,y) * f(x,y) = \sum_{dx=-\infty}^\infty \sum_{dy=-\infty}^\infty \omega (dx,dy)f(x-dx,y-dy) $$

Un filtro gaussiano consiste in una convoluzione 2D tra un'immagine $f(x,y)$ e una matrice $\omega(x,y)$ che rappresenta una curva normale discreta (il "kernel"). Il kernel 3x3, per esempio, è definito come

$$ \frac{1}{16} \begin{bmatrix} \ \ 1 &\ \ 2 &\ \ 1 \\ \ \ 2 &\ \ 4 &\ \ 2 \\ \ \ 1 &\ \ 2 &\ \ 1 \end{bmatrix} $$

Nel dominio delle frequenze:

$$ G(u,v) = \mathcal{F} [ g(x,y) ] = \mathcal{F} [ f(x,y) * \omega(x,y) ] = F(u,v) \Omega(u,v) $$

dove $\Omega(u,v) = \mathcal{F} [ \omega(x,y) ]$ è la DFT del kernel. Di conseguenza è possibile recuperare il segnale originale calcolando

$$ f(x, y) = \mathcal{F}^{-1} [ F(u,v) ] = \mathcal{F}^{-1} \left[ \frac{G(u,v)}{\Omega(u,v)} \right] $$

dove $\mathcal{F}^{-1}$ è la trasformata di Fourier discreta inversa (IDFT).

Note

I precedenti calcoli non tengono in considerazione il rumore introdotto ad ogni passaggio. In quel caso bisognerebbe stimare $\hat{f}(x,y)$ che minimizza l'errore quadratico medio $\mathbb{E} \left| f(x,y) - \hat{f}(x,y) \right|^2$ (v. deconvoluzione di Wiener).

Esperimento 1

Descrizione

Come input è stata usata un'immagine RGB con estensione .ppm. Il file viene convertito in scala di grigi e salvato come .pgm. Nella repository sono già presenti varie immagini di esempio $512 \times 512$, tra cui la classica Lenna¹ e la meno classica Sabrina Salerno². Altre immagini .ppm possono essere utilizzate, basta aggiungerle alla medesima cartella.

L'immagine in bianco e nero viene sfocata tramite convoluzione. L'energia presente nell'immagine sfocata è confrontata con quella dell'immagine originale e viene calcolato l'SNR. Viene calcolato l'inverso della trasformata di Fourier del kernel scelto (a sua volta un'altra curva gaussiana) che viene poi moltiplicato per lo spettro dell'immagine sfocata. Del risultato viene conservata solo la parte reale, dato che la parte immaginaria conterrà solo rumore. I valori superiori a 255 e quelli inferiori a 0 vengono tagliati e il tutto viene convertito in byte. Viene calcolato anche l'SNR dell'immagine risultante e confrontato con quello dell'immagine sfocata.

Problemi

Problema 1

I bordi...

Problema 2

La Symmetrization dell'input è fondamentale per evitare di introdurre alte frequenze e invalidare ogni misura.

Problema 3

Le alte frequenze del kernel sono nulle o prossime allo zero e invertirle dà origine a valori molto elevati. Moltiplicare questi valori per le alte frequenze dell'immagine sfocata porterà all'introduzione di rumore. Di conseguenza è stata definita una soglia sul valore assoluto oltre la quale le frequenze non sono state invertite, ma portate a 1, al fine di non modificare le alte frequenze dell'immagine, comunque già di per sè poco rilevanti in un'immagine naturale. Un processo iterativo a portato a definire 0,07 come valore di soglia medio grazie al quale è stato possibile recuperare circa 2 dB di SNR dalle immagini usate nel test. Ogni immagine avrà, naturalmente, un preciso valore in corrispondenza del quale si potrà recuperare il massimo di dB di SNR.

Risultati

Di seguito si riportano i risultati ottenuti con Lena (qui le immagini sono convertite in .png).

L'immagine è stata sfocata con una convoluzione (3x3, 5x5, 7x7):

La de-convoluzione (soglia: 0,07) ha restituito le seguenti immagini:

File	SNRs [dB]
	3x3			5x5			7x7
	Blurred	De-blurred	Difference	Blurred	De-blurred	Difference	Blurred	De-blurred	Difference
sabrina.ppm	9.42	12.18	2.77	6.99	10.26	3.27	6.92	10.03	3.11
lena.ppm	6.41	9.70	3.29	4.25	7.76	3.50	4.31	7.91	3.60
sara.ppm	5.17	10.55	5.38	2.20	7.07	4.88	2.02	6.52	4.50
laura.ppm	5.89	9.07	3.18	3.84	6.88	3.04	3.84	6.95	3.11

Esperimento 2

Che succede se si prova a de-sfocare un'immagine con l'inverso di un kernel diverso da quello usato per sfocarla?

Risulta evidente che la conoscenza del filtro di sfocatura è fondamentale per recuperare al meglio l'immagine originale. Tuttavia il de-blurring con dimensioni inferiori a quelle della sfocatura ha comunque permesso di recuperare qualche dB. Inoltre i filtri 5x5 e 7x7 paiono molto simili.

Footnotes

1. Playboy, vol. 19, n. 11, novembre 1972, Playboy Enterprises. ↩

2. Playmen, anno XXII, n. 9, settembre 1988, Tattilo Editrice. ↩