In [119]:
#Le fichier doit être exécuté avec python 3.12 !
import importlib
import SamWavelet as sw 
importlib.reload(sw)
#le reste des libraries est importé dans SamWavelet.py
Out[119]:
<module 'SamWavelet' from 'c:\\Users\\darkz\\Desktop\\Cours\\TIPE\\1\\SamWavelet.py'>

Ceci est la première partie de mon TIPE. Elle constitue les premiers tests propres après plusieurs mois de programmation jusqu'à obtenir SamWavelet.py une bilbiothèque regroupant toutes les fonctions que j'ai écrites me permettant d'effectuer la transformation en ondelettes de signaux 1D (sons) et 2D (images). La librairie est disponible sur le dépôt github suivant : https://github.com/UseLess289/SamWavelet

Tehcniques de débruitage par lifting en ondelettes et seuillage dur avec ou sans invariance par translation¶

On utilise une transcription python du script de transformation en ondelettes fourni par Gabriel Peyre en 2009 pour MatLab

Le script est le suivant (en partie) :

In [ ]:
def perform_wavelet_transform(x, Jmin, direction, options=None):
    """
    Execute une transformation en ondelettes 1D ou 2D avec lissage symetrique
    et traitement des bords symetriques.
    
    Paramètres:
    -----------
    x : numpy.ndarray
        Signal ou image d'entrée
    Jmin : int
        Niveau de décomposition minimum
    direction : int
        1 pour transformée directe, -1 pour transformée inverse
    options : dict, optional
        liste des options:
        - filter: string ou array, coefficients du filtre d'ondelettes
          Peut être 'haar', 'linear'/'5-3', ou '9-7'
        - ti: bool, si vrai, calculer la transformée invariante de translation
        - separable: bool, si vrai, calculer la transformée 2D separable
        - ndims: int, force la dimension (1 ou 2) si elle est mal détectée
    
    Sortie:
    --------
    numpy.ndarray
        Coefficients d'ondelettes ou signal reconstruit
    
    """
        
    if options is None:
        options = {}
    
    # Options par défaut
    h = options.get('filter', '9-7')
    separable = options.get('separable', False)
    
   
    d = options.get('ndims', -1)  # Détection de la dimension
    if d < 0:
        if x.ndim == 1 or (x.ndim == 2 and (x.shape[0] == 1 or x.shape[1] == 1)):
            d = 1
        else:
            d = 2
    
    
    x = np.array(x, dtype=float)
    
    # Gérer l'orientation du vecteur 1D
    if d == 1 and x.ndim == 2 and x.shape[0] == 1:
        x = x.T
    
    # Convertir les noms de filtres en chaînes en coefficients
    if isinstance(h, str):
        h = h.lower()
        if h == 'haar':
            if d == 1 or not separable:
                return perform_haar_transform(x, Jmin, direction, options)
        elif h in ['linear', '5-3']:
            h = np.array([1/2, 1/4, np.sqrt(2)])
        elif h in ['9-7', '7-9']:
            h = np.array([1.586134342, -0.05298011854, -0.8829110762, 0.4435068522, 1.149604398])
        else:
            raise ValueError(f"Unknown filter: {h}")
    
    # Option translation invariant(clée pour les images)
    ti = options.get('ti', False)
    
    # Gérer la transformée 2D separable
    if d == 2 and separable:
        options_copy = options.copy()
        options_copy['ti'] = False
        if ti:
            print("Warning: Separable does not work for translation invariant transform")
        
        n = x.shape[0]
        if direction == 1:
            for i in range(n):
                x[:, i] = perform_wavelet_transform(x[:, i], Jmin, direction, options_copy)
            for i in range(n):
                x[i, :] = perform_wavelet_transform(x[i, :], Jmin, direction, options_copy)
        else:
            for i in range(n):
                x[i, :] = perform_wavelet_transform(x[i, :], Jmin, direction, options_copy)
            for i in range(n):
                x[:, i] = perform_wavelet_transform(x[:, i], Jmin, direction, options_copy)
        return x
    
    # Calculer les étapes de lifting
    n = x.shape[0]
    m = (len(h) - 1) // 2
    Jmax = int(math.log2(n)) - 1
    jlist = list(range(Jmax, Jmin-1, -1))
    if direction == -1:
        jlist = jlist[::-1]
    
    if not ti:
        # Transformée subsamplée
        for j in jlist:
            if d == 1:
                x[:2**(j+1)] = lifting_step(x[:2**(j+1)], h, direction)
            else:
                x[:2**(j+1), :2**(j+1)] = lifting_step(x[:2**(j+1), :2**(j+1)], h, direction)
                x[:2**(j+1), :2**(j+1)] = lifting_step(x[:2**(j+1), :2**(j+1)].T, h, direction).T
    else: #i.e if ti is True
        # Transformée invariante de translation
        nJ = Jmax - Jmin + 1
        if direction == 1:
            if d == 1: #pour un signal 1D,
                x = np.repeat(x[:, np.newaxis], nJ+1, axis=1) #on crée une matrice 3D avec nJ+1 colonnes
            elif d == 2: #pour une image 2D,
                x_3d = np.zeros((x.shape[0], x.shape[1], 3*nJ+1)) #on crée une matrice 3D avec 3*nJ+1 colonnes
                x_3d[:,:,0] = x #on copie la matrice d'entrée dans la première colonne de la matrice 3D
                x = x_3d #on affecte la matrice 3D à x
        
        for j in jlist:
            dist = 2**(Jmax-j) #distance entre les coefficients
            if d == 1: #pour un signal 1D,
                if direction == 1: #pour une transformée directe,
                    result = lifting_step_ti(x[:, 0], h, direction, dist) #on applique la transformée invariante de translation
                    x[:, 0] = result[:, 0] #on affecte la première colonne de la matrice 3D à x
                    x[:, j-Jmin+1] = result[:, 1] #on affecte la deuxième colonne de la matrice 3D à x
                else: #pour une transformée inverse,
                    x[:, 0] = lifting_step_ti(x[:, [0, j-Jmin+1]], h, direction, dist) #on applique la transformée inverse
            else:
                dj = 3*(j-Jmin)
                if direction == 1:
                    # Transformée horizontale
                    result = lifting_step_ti(x[:,:,0], h, direction, dist)
                    x[:,:,0] = result[:,:,0]
                    x[:,:,dj+1] = result[:,:,1]
                    
                    # Transformée verticale sur le signal original
                    result = lifting_step_ti(x[:,:,0].T, h, direction, dist)
                    x[:,:,0] = result[:,:,0].T
                    x[:,:,dj+2] = result[:,:,1].T
                    
                    # Transformée verticale sur les détails horizontaux
                    result = lifting_step_ti(x[:,:,dj+1].T, h, direction, dist)
                    x[:,:,dj+1] = result[:,:,0].T
                    x[:,:,dj+3] = result[:,:,1].T
                else:
                    # Inverse l'ordre des opérations
                    x[:,:,dj+1] = x[:,:,dj+1].T
                    x[:,:,dj+3] = x[:,:,dj+3].T
                    x[:,:,dj+1] = lifting_step_ti(x[:,:,[dj+1, dj+3]], h, direction, dist).T
                    
                    x[:,:,0] = x[:,:,0].T
                    x[:,:,dj+2] = x[:,:,dj+2].T
                    x[:,:,0] = lifting_step_ti(x[:,:,[0, dj+2]], h, direction, dist).T
                    
                    x[:,:,0] = lifting_step_ti(x[:,:,[0, dj+1]], h, direction, dist)
        
        if direction == -1:
            if d == 1:
                x = x[:, 0]
            else:
                x = x[:,:,0]
    
    return x

Dans le cadre du TIPE, on va préferer utiliser la transformée invariante de translation.

Principe :¶

Calculer et stocker la transformée en ondelette de notre signal, mais aussi de toutes ses translatées.

Avantages supposés :¶

  • réduit les artefacts proches des discontinuitées
  • préserve mieux les détails
  • plus efficace
  • complexité divisée par 2 par rapport à un FFT !

Désavantage :¶

  • elle est largement plus couteuse en mémoire

La méthode utilisée par Gabriel Peyre est particulièrement intéréssante de la part le moyennage des coefficients, la diminution des effets de bords et le traitement des bords eux mêmes. Toutes ces améliorations permettent d'obtenir un résultat vraiment propre et lisse. Mais au centre de cette fonction se trouve la fonction lifting :

In [121]:
def lifting_step_ti(x, h, direction, dist):
    """
    Effectue une étape de lifting pour la transformée en ondelettes invariante de translation.
    
    Paramètres:
    -----------
    x : numpy.ndarray
        Signal ou partie de la transformée
    h : numpy.ndarray
        Coefficients du filtre de lifting
    direction : int
        1 pour transformée directe, -1 pour transformée inverse
    dist : int
        Paramètre de distance pour la translation
    
    Returns:
    --------
    numpy.ndarray
        Résultat après l'étape de lifting
    """
    m = (len(h) - 1) // 2
    n = x.shape[0]
    
    # Calculer les indices décalés avec les conditions aux limites
    s1 = np.arange(n) + dist
    s2 = np.arange(n) - dist
    
    # Appliquer les conditions aux limites
    # Cette étape est cruciale car elle permet d'étendre le signal de manière symétrique aux bords
    # cela permet d'éviter les artefacts aux bords
    s1[s1 >= n] = 2*n - s1[s1 >= n] - 1
    s1[s1 < 0] = -s1[s1 < 0] - 1
    s2[s2 >= n] = 2*n - s2[s2 >= n] - 1
    s2[s2 < 0] = -s2[s2 < 0] - 1
    
    if direction == 1:
        # séparer les coefs
        d = x.copy()
        s = x.copy()
        
        # étapes de lifting
        for i in range(m):
            d = d - h[2*i] * (s[s1] + s[s2])
            s = s + h[2*i+1] * (d[s1] + d[s2])
        
        # scalling
        s = s * h[-1]
        d = d / h[-1]
        
        # Retourner les deux parties de fréquence basse et haute
        if x.ndim == 1:
            return np.column_stack([s, d])
        else:
            result = np.zeros((x.shape[0], x.shape[1], 2))
            result[:,:,0] = s
            result[:,:,1] = d
            return result
    else:
        # Pour la transformée inverse, x doit avoir deux composants
        if x.ndim == 2 and x.shape[1] == 2:
            s = x[:, 0].copy() / h[-1]
            d = x[:, 1].copy() * h[-1]
        else:
            s = x[:,:,0].copy() / h[-1]
            d = x[:,:,1].copy() * h[-1]
        
        # étapes de lifting inverses
        for i in range(m-1, -1, -1):
            s = s - h[2*i+1] * (d[s1] + d[s2])
            d = d + h[2*i] * (s[s1] + s[s2])
        
        # Fusionner
        return (s + d) / 2

voir https://fr.wikipedia.org/wiki/Lifting_en_ondelettes#Exemple_de_la_transform%C3%A9e_lifting_en_ondelettes_de_Haar

Applications¶

Cas 1D¶

Signal sinusoïdal¶

In [122]:
#Signal pur
f0=2
nbpoints=2**12
t = np.arange(nbpoints)/ nbpoints
x = np.sin(2 * np.pi * f0 * t)
#x le signal sinusoidal pur
#bruit
sigma = 0.5
bruit = sigma * np.random.randn(len(x))
y=x+bruit 
#y le signal bruité

#affichage 
plt.figure(figsize=(12, 6))
plt.plot(y, 'r', alpha=0.7, label='Signal bruité')
plt.plot(x, 'b', linewidth=2, label='Signal original')
plt.legend()
plt.title('Signal original et signal bruité')
plt.grid(True)
plt.tight_layout()
#plt.savefig('signal_bruite.png', dpi=150)
plt.show()
No description has been provided for this image
In [123]:
#Mise en place des paramètres de débruitage
Jmin= 1
options = {'filter': '9-7'}

#on effectue la dwt des trois signaux
WTx = sw.perform_wavelet_transform(x, Jmin, 1, options)
WTbruit = sw.perform_wavelet_transform(bruit, Jmin, 1, options)
WTy = sw.perform_wavelet_transform(y, Jmin, 1, options)

#choix du seuil
seuil=np.max(np.abs(WTbruit))
print('le seuil est :', seuil)

#débruitage 
WTyHT=sw.ht1D(WTy,seuil)

#TO inverse
yHT=sw.perform_wavelet_transform(WTyHT, Jmin, -1, options)

#Visu
plt.figure(figsize=(12, 6))
plt.plot(y, 'r', alpha=0.7, label='Signal bruité')
plt.plot(x, 'b', linewidth=2, label='Signal original')
plt.plot(yHT, 'c', linewidth=2, label='Hard thresholding')
plt.legend()
plt.title('Débruitage par Hard Thresholding')
plt.grid(True)
plt.tight_layout()
#plt.savefig('debruitage_hard_thresholding.png', dpi=150)
plt.show()
le seuil est : 1.6885617494203897
No description has been provided for this image
In [124]:
#On va maintenant calculer l'erreur quadratique moyenne
mse_bruite = np.mean((y - x)**2)
mse_ht = np.mean((yHT - x)**2)

    
print(f"MSE du signal bruité: {mse_bruite:.6f}")
print(f"MSE après Hard Thresholding: {mse_ht:.6f}")
print(f"on a divisé la mse par un facteur de {mse_bruite/mse_ht:.2f}")
MSE du signal bruité: 0.248495
MSE après Hard Thresholding: 0.005681
on a divisé la mse par un facteur de 43.74

Le résultat est plutot acceptale malgré quelques irrégularitées ???

Signal Hölderien¶

In [125]:
#### on crée le signal
# Paramètres
alp = 0.5  # Régularité Höldérienne
n = 12     # Nombre de niveaux
Jmin = 1
    
# Créer la matrice des coefficients
matricepart = np.zeros((n, 2**(n-1)))
for j in range(n):
    matricepart[j, :2**j] = 2**(-j*(alp+1/2))
    
# Créer le vecteur des coefficients
WTpart = np.zeros(2**n)
for j in range(n):
    WTpart[2**j:2**(j+1)] = matricepart[j, :2**j]
    
# Reconstruction du signal
options = {'filter': '9-7'}
x = sw.perform_wavelet_transform(WTpart, Jmin, -1, options)
    
###on bruite le signal
y=x + 0.005 * np.random.randn(len(x))
bruit=y-x


#### Visualiser les signaux
plt.figure(figsize=(12, 6))
plt.plot(x, 'b', alpha =0.7, label='Signal original')
plt.plot(y, 'r', linestyle='--', linewidth=0.5, label='Signal bruité')
plt.title('Signal particulier (Höldérien) bruité')
plt.grid(True)
plt.legend()
plt.tight_layout()
#plt.savefig('signal_holderien.png', dpi=150)
plt.show()
No description has been provided for this image
In [126]:
#TO
WTy=sw.perform_wavelet_transform(y, Jmin, 1, options)
WTbruit=sw.perform_wavelet_transform(bruit, Jmin, 1, options)
#HT
seuil= np.max(np.abs(WTbruit))
WTyHT=sw.ht1D(WTy,seuil)
#TO inverse
yHT=sw.perform_wavelet_transform(WTyHT, Jmin, -1, options)

#Visu
plt.figure(figsize=(12, 6))
plt.plot(x, 'b', alpha =0.7, label='Signal original')
plt.plot(yHT, 'r', linestyle='--', linewidth=0.5, label='Signal débruité')
plt.legend()
plt.title('Débruitage du signal Höldérien')
plt.grid(True)
plt.tight_layout()
#plt.savefig('debruitage_holderien.png', dpi=150)
plt.show()
No description has been provided for this image
In [127]:
#On va maintenant calculer l'erreur quadratique moyenne
mse_bruite = np.mean((y - x)**2)
mse_ht = np.mean((yHT - x)**2)

print(f"MSE du signal bruité: {mse_bruite:.6f}")
print(f"MSE après Hard Thresholding: {mse_ht:.6f}")
print(f"on a divisé la mse par un facteur de {mse_bruite/mse_ht:.2f}")
MSE du signal bruité: 0.000025
MSE après Hard Thresholding: 0.000007
on a divisé la mse par un facteur de 3.61

Malgré quelques irrégularitées dans le résultat, il est assez satisfaisant.

Cas 2D, les images¶

Une image élémentaire¶

In [128]:
### Création d'une image élémentaire (matrice 2D)
image_size = 128  # Taille de l'image 128x128
x = np.zeros((image_size, image_size))  # Fond noir

# Dessiner un carré blanc au centre
square_size = 32  
start = (image_size - square_size) // 2
end = start + square_size
x[start:end, start:end] = 255  # Carré blanc

# Visualisation de l'image
plt.figure(figsize=(6, 6))
plt.imshow(x, cmap='gray', interpolation='nearest')
plt.title('Image élémentaire')
plt.axis('off')
plt.show()

### Maintenant on bruite l'image
var = 100
y = sw.bruit(x, var)
# visu
plt.figure(figsize=(6, 6))
plt.imshow(y, cmap='gray', interpolation='nearest')
plt.title('Image élémentaire bruitée')
plt.axis('off')
plt.show()

bruit= y-x
No description has been provided for this image
No description has been provided for this image
In [129]:
#TO
Jmin= 1
direction= 1
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 2}

WTy= sw.perform_wavelet_transform(y, Jmin, direction, options)
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 1}
WTbruit= sw.perform_wavelet_transform(bruit, Jmin, direction, options)


# Utilisation de la fonction
WTy = sw.hard_thresholding(WTy, seuil)

print(type(WTy))
#HT
seuil= 2*np.sqrt(var)
WTyHT= sw.hard_thresholding(WTy, seuil)
#TO inverse 
direction= -1
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 2}
yHT= sw.perform_wavelet_transform(WTyHT, Jmin, direction, options)

#visu
plt.figure(figsize=(6, 6))
plt.imshow(yHT, cmap='gray', interpolation='nearest')
plt.title('Image élémentaire débruitée')
plt.axis('off')
plt.show()
<class 'numpy.ndarray'>
No description has been provided for this image
In [130]:
# Calcul de l'erreur quadratique moyenne
mse_bruite = np.mean((y - x) ** 2)
mse_debruite = np.mean((yHT - x) ** 2)
print(f"MSE bruitée: {mse_bruite:.2f}")
print(f"MSE débruitée: {mse_debruite:.2f}")
MSE bruitée: 47.56
MSE débruitée: 15.25

Appréciations :¶

  • Le débruitage est clairement effectif
  • Néanmoins, le résultat présente des artefacts étranges qui pourraient être réduits par un lissage classique

Une image quelconque¶

In [131]:
x= cv.imread('juju.jpg', cv.IMREAD_GRAYSCALE)
x=sw.adjust_image_size(x)



#bruitage de l'image
var = 100
y = sw.bruit(x, var)
bruit = y-x

#affichage des images
plt.figure(figsize=(12, 12))
plt.imshow(x, cmap='gray')
plt.title("Image Originale")
plt.axis('off')
plt.show()
plt.close()

plt.figure(figsize=(12, 12))
plt.imshow(y, cmap='gray')
plt.title("Image bruitée")
plt.axis('off')
plt.show()
plt.close()

#mesure de l'erreur quadratique de l'image bruitée :
mse_bruite = np.mean((y - x) ** 2)
print(f"MSE bruitée: {mse_bruite:.2f}")
print(f"SNR de l'image bruitée : {sw.SNR(x, y):.2f} dB")
No description has been provided for this image
No description has been provided for this image
MSE bruitée: 33.14
SNR de l'image bruitée : 30.85 dB
In [132]:
####sans invariance par translation
#TO 
Jmin= 1
direction= 1
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 2}

WTy= sw.perform_wavelet_transform(y, Jmin, direction, options)
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 1}
WTbruit= sw.perform_wavelet_transform(bruit, Jmin, direction, options)

#HT
seuil= 2*np.sqrt(var)
WTyHT= sw.hard_thresholding(WTy, seuil)
#TO inverse
direction= -1
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 2}
yHT= sw.perform_wavelet_transform(WTyHT, Jmin, direction, options)

#visu
plt.figure(figsize=(12, 12))
plt.imshow(yHT, cmap='gray', interpolation='nearest')
plt.title('Image débruitée sans transformée invariante par translation')
plt.axis('off')
plt.show()


###avec invariance par translation
#TO
Jmin= 1
direction= 1
options= {'filter': 'haar', 'ti': True, 'separable': False, 'ndims': 2}

WTy= sw.perform_wavelet_transform(y, Jmin, direction, options)
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 1}
WTbruit= sw.perform_wavelet_transform(bruit, Jmin, direction, options)

#HT
seuil= 2*np.sqrt(var)
WTyHT= sw.hard_thresholding(WTy, seuil)
#TO inverse
direction= -1
options= {'filter': 'haar', 'ti': True, 'separable': False, 'ndims': 2}
yHT= sw.perform_wavelet_transform(WTyHT, Jmin, direction, options)

#visu
plt.figure(figsize=(12, 12))
plt.imshow(yHT, cmap='gray', interpolation='nearest')
plt.title('Image débruitée avec transformée invariante par translation')
plt.axis('off')
plt.show()
No description has been provided for this image
No description has been provided for this image
In [133]:
#Calcul de l'erreur quadratique moyenne et des rapports signaux/bruits
mse_bruite = np.mean((y - x) ** 2)
mse_debruite = np.mean((yHT - x) ** 2)
snr= sw.SNR(x, yHT)
print(f"MSE bruitée: {mse_bruite:.2f}")
print(f"MSE débruitée: {mse_debruite:.2f}")
print(f"SNR de l'image débruitée : {snr:.2f} dB")
print(f"Gain en dB par débruitage : {snr - sw.SNR(x, y):.2f} dB")
MSE bruitée: 33.14
MSE débruitée: 25.54
SNR de l'image débruitée : 34.06 dB
Gain en dB par débruitage : 3.21 dB

Appréciations :¶

  • Le débruitage sans invariance par translation donne un résultat étrange (voir carré uniforme)
  • Le débruitage avec invariance par translation donne un résultat uniforme
  • Les deux débruitages présentent des points gris aléatoirement distribués dans le fond noir et les bords de Jupiter sont nettement dégradés
  • On peut supposer que le fort contraste entre objet et fond en est à l'origine
  • Les images sont moins bonnes visuellement après débruitage malgré un gain de SNR la réduction de l'erreur quadratique moyenne

Origines des problèmes :¶

Les problèmes rencontrés avec l'image de Jupiter sont liés à la nature même de l'image. L'image de Jupiter a un contraste extrême, des bords nets et des détails particulièrements fins ainsi qu'un fond absolument noir ce qui a pour conséquence de révéler la moinde imperfection.

Un traitement standard comme utilisé n'est donc pas adapté à ce type d'image

Un seconde exemple plus adapté¶

In [134]:
x= cv.imread('dm.jpg', cv.IMREAD_GRAYSCALE)
x=sw.adjust_image_size(x)



#bruitage de l'image
var = 100
y = sw.bruit(x, var)
bruit = y-x

#affichage des images
plt.figure(figsize=(12, 12))
plt.imshow(x, cmap='gray')
plt.title("Image Originale")
plt.axis('off')
plt.show()
plt.close()

plt.figure(figsize=(12, 12))
plt.imshow(y, cmap='gray')
plt.title("Image bruitée")
plt.axis('off')
plt.show()
plt.close()

#mesure de l'erreur quadratique de l'image bruitée :
mse_bruite = np.mean((y - x) ** 2)
print(f"MSE bruitée: {mse_bruite:.2f}")
print(f"SNR de l'image bruitée : {sw.SNR(x, y):.2f} dB")
No description has been provided for this image
No description has been provided for this image
MSE bruitée: 35.21
SNR de l'image bruitée : 30.47 dB
In [135]:
####sans invariance par translation
#TO 
Jmin= 1
direction= 1
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 2}

WTy= sw.perform_wavelet_transform(y, Jmin, direction, options)
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 1}
WTbruit= sw.perform_wavelet_transform(bruit, Jmin, direction, options)

#HT
seuil= 2*np.sqrt(var)
WTyHT= sw.hard_thresholding(WTy, seuil)
#TO inverse
direction= -1
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 2}
yHT= sw.perform_wavelet_transform(WTyHT, Jmin, direction, options)

#visu
plt.figure(figsize=(12, 12))
plt.imshow(yHT, cmap='gray', interpolation='nearest')
plt.title('Image débruitée sans transformée invariante par translation')
plt.axis('off')
plt.show()

mse_debruite_nti = np.mean((yHT - x) ** 2)


###avec invariance par translation
#TO
Jmin= 1
direction= 1
options= {'filter': 'haar', 'ti': True, 'separable': False, 'ndims': 2}

WTy= sw.perform_wavelet_transform(y, Jmin, direction, options)
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 1}
WTbruit= sw.perform_wavelet_transform(bruit, Jmin, direction, options)

#HT
seuil= 2*np.sqrt(var)
WTyHT= sw.hard_thresholding(WTy, seuil)
#TO inverse
direction= -1
options= {'filter': 'haar', 'ti': True, 'separable': False, 'ndims': 2}
yHT= sw.perform_wavelet_transform(WTyHT, Jmin, direction, options)

#visu
plt.figure(figsize=(12, 12))
plt.imshow(yHT, cmap='gray', interpolation='nearest')
plt.title('Image débruitée avec transformée invariante par translation')
plt.axis('off')
plt.show()
No description has been provided for this image
No description has been provided for this image
In [137]:
#Calcul de l'erreur quadratique moyenne et des rapports signaux/bruits
mse_bruite = np.mean((y - x) ** 2)
mse_debruite = np.mean((yHT - x) ** 2)
snr= sw.SNR(x, yHT)
print(f"MSE bruitée: {mse_bruite:.2f}")
print(f"MSE débruitée avec transfo invariante: {mse_debruite:.2f}")
print(f"MSE débruitée sans transfo invariante: {mse_debruite_nti:.2f}")
print(f"SNR de l'image débruitée : {snr:.2f} dB")

print(f"Gain en dB par débruitage : {snr - sw.SNR(x, y):.2f} dB")
MSE bruitée: 35.21
MSE débruitée avec transfo invariante: 59.70
MSE débruitée sans transfo invariante: 44.39
SNR de l'image débruitée : 30.37 dB
Gain en dB par débruitage : -0.10 dB

Appréciations :¶

  • Le débruitage par transformation invariante est un peu grisé mais satisfaisant
  • Celui avec est moins bien, apparition d'un grain

Origine des problèmes :¶

  • Pour ce qui est de la transformation invariante par translation, cela s'explique par les motifs répétitifs de la page ainsi que par les forts contrastes présents entre page vide et lignes ou écritures

Un paysage¶

In [142]:
x= cv.imread('img4.jpg', cv.IMREAD_GRAYSCALE)
x=sw.adjust_image_size(x)



#bruitage de l'image
var = 100
y = sw.bruit(x, var)
bruit = y-x

#affichage des images
plt.figure(figsize=(12, 12))
plt.imshow(x, cmap='gray')
plt.title("Image Originale")
plt.axis('off')
plt.show()
plt.close()

plt.figure(figsize=(12, 12))
plt.imshow(y, cmap='gray')
plt.title("Image bruitée")
plt.axis('off')
plt.show()
plt.close()

#mesure de l'erreur quadratique de l'image bruitée :
mse_bruite = np.mean((y - x) ** 2)
print(f"MSE bruitée: {mse_bruite:.2f}")
print(f"SNR de l'image bruitée : {sw.SNR(x, y):.2f} dB")
No description has been provided for this image
No description has been provided for this image
MSE bruitée: 60.80
SNR de l'image bruitée : 28.12 dB
In [143]:
####sans invariance par translation
#TO 
Jmin= 1
direction= 1
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 2}

WTy= sw.perform_wavelet_transform(y, Jmin, direction, options)
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 1}
WTbruit= sw.perform_wavelet_transform(bruit, Jmin, direction, options)

#HT
seuil= 2*np.sqrt(var)
WTyHT= sw.hard_thresholding(WTy, seuil)
#TO inverse
direction= -1
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 2}
yHT= sw.perform_wavelet_transform(WTyHT, Jmin, direction, options)

#visu
plt.figure(figsize=(12, 12))
plt.imshow(yHT, cmap='gray', interpolation='nearest')
plt.title('Image débruitée sans transformée invariante par translation')
plt.axis('off')
plt.show()

mse_debruite_nti = np.mean((yHT - x) ** 2)


###avec invariance par translation
#TO
Jmin= 1
direction= 1
options= {'filter': 'haar', 'ti': True, 'separable': False, 'ndims': 2}

WTy= sw.perform_wavelet_transform(y, Jmin, direction, options)
options= {'filter': 'haar', 'ti': False, 'separable': False, 'ndims': 1}
WTbruit= sw.perform_wavelet_transform(bruit, Jmin, direction, options)

#HT
seuil= 2*np.sqrt(var)
WTyHT= sw.hard_thresholding(WTy, seuil)
#TO inverse
direction= -1
options= {'filter': 'haar', 'ti': True, 'separable': False, 'ndims': 2}
yHT= sw.perform_wavelet_transform(WTyHT, Jmin, direction, options)

#visu
plt.figure(figsize=(12, 12))
plt.imshow(yHT, cmap='gray', interpolation='nearest')
plt.title('Image débruitée avec transformée invariante par translation')
plt.axis('off')
plt.show()
No description has been provided for this image
No description has been provided for this image

Conclusion générale¶

On remarque que bien que les résultats puissent être satisfaisants, il est nécéssaire d'adapter la fonction à la nature de l'image.

En effet, les logiciels de traitement d'image conventionels utilisent des algorithmes adaptatifs pour traiter le bruit. L'outil Adobe Sensei analyse le contenu des images pour identifier les zones de détail. Topaz DeNoise AI utilise des réseaux de neurones entrainés sur des paires d'images bruitées/propres pour aider à séparer le bruit du détail.

Cependant, nous avons développé une méthode basique qui se révèle fonctionnelle dans certains cas adaptés. Seulement il faut remarquer qu'ici, le débruitage passe par un seuillage dur des coefficients. Pour cela voir, voir mon MCOT. Donc la valeur du seuil est cruciale. Sauf que dans tous les exemples qui précèdent, on crée un signal x, un bruit et on définit y=x+bruit. On choisit ensuite le seuil à partir de la variable bruit et non y. Alors que dans un cas d'usage du programme, on ne dispose que de l'image y (bruitée donc) et pas des caractéristiques du bruit. La prochaine étape consiste de donc à visiter plusieurs méthode pour estimer la valeur de seuil la plus performante.