Travaux pratiques de traitement d’image

M.N

M 1

.

k=0

FreqIma[ m , n].e

2.j. pi .

n.i l.j

M.N

N 1

l=0

Ima est l'image et i,j sont les coordonn es du domaine spatial.

FreqIma est l'image et m,n sont les coordonn es du domaine fr quentiel.

M, N sont les dimensions de l'image.

3 . MANIPULATIONS :

3.1 . Transform e de Fourier 1D :

Dans cette premi re partie, on vous propose de revoir rapidement comment est construite la

Transform e de Fourier Discr te (FFT) sur un signal 2D :

Pour commencer, ouvrez l'image porte.bmp :

MaPorte=imread('porte.bmp');

Vous savez pr sent que limage est disponible sous la forme dune matrice de valeurs

cod es sur 8 bits. Pour certaines op rations, il est important de noter que le r sultat obtenu

est plus pr cis, voire tr s diff rents, si les calculs sont r alis s sur des valeurs r elles. Pour

convertir limage en flottant :

MaPorte = double(MaPorte) ;

Limage en niveau de gris est form e de Nlin lignes, Ncol colonnes et 1 seul plan. Ces

informations sont r cup rables facilement gr ce :

Pour afficher l'image :

= size(MaPorte) ;

figure(1); image(MaPorte);

matlab permet de calculer la FFT d'un signal 1D (un vecteur) gr ce la fonction :

Freq_MonSignal = fft(MaSignal);

Attention : votre image est un signal 2D r el et sa fft est complexe (la fonction abs() donne le module).

Essayez de regarder ce qu'il se passe si vous faites la fft sur la ligne Nlin/2 de votre image.

Tracez ce spectre.

Effectuez la m me op ration sur la premi re ligne de l'image. Conclusion ?

Sachant que le spectre d'un signal chantillonn est p riodique, effectuez les instructions

suivantes pour rendre votre espace des fr quences p riodique :

Freq_periodique=Freq_ligne;

for u=2:16

Freq_periodique= ;

end

Tracez les spectres Freq_ligne et Freq_periodique,

Calculez la fft inverse (fonction ifft(MonSpectre) )du nouveau spectre. Conclusion ?

Sur les spectres trac s, vous observez que l'axe des abscisses repr sente l'indice du point

2

Travaux pratiques de traitement dimage

P. Falgayrettes

calcul et non la fr quence. Vous pouvez aussi vous rendre compte que le trac repr sente

les fr quences entre 0 et Fech/2 sur les indices 0 Nlin/2 et les fr quences entre -Fech/2 et

0 pour les indices Nlin/2 Nlin.

Changez l' chelle de l'axe des abscisses.

plot (abs(Freq_periodique));

v=axis; % r cup re les valeurs des chelles sur x et y

axis( ); % modifie les valeurs sur

l'axe x uniquement !

Le spectre obtenue doit ressembler quelque chose que vous connaissez.

Sous matlab, il existe une fonction qui permet d'effectuer ce d calage d'axe. Essayez :

figure(3);

FFT_ligne=fft(Im1(Nlin/2,:));

FFT_ligne=fftshift(FFT_ligne); % recentrage de l'origine des

fr quences au centre de la fen tre

plot (abs(FFT_ligne)); % trace le module du spectre

3.2 . La Transform e de Fourier 2D :

Vous avez revue vos classiques sur la transform e 1D. Revenons aux images. Matlab

permet de calculer la fft d'une image (commande fft2, regardez l'aide)

Affichez la transform e de Fourier 2D de 'MaPorte.bmp'.

Sachant que les valeurs du module obtenues peuvent tre grandes, vous devrez ajuster les

amplitudes du spectre aux valeurs de couleurs disponibles.

Image_affiche=FFT_MaPorte.... % mise l' chelle des niveaux

de gris

image(Image_affiche);

colormap(pink(256)); % meilleure visibilit du spectre

axis('image');

Remarque : Vous pouvez utilisez la fonction 'imagesc()' de matlab la place de 'image()' qui permet de

recadrer automatiquement les valeurs dans l'intervalle des couleurs disponibles.

3.3 . Animez votre Transform e de Fourier 2D !

Vous ma trisez parfaitement les affichages dans le domaine fr quentiel.

On se propose de r aliser une petite application didactique :

a) Cr ez une image 512x512 contenant un cosinus aux fr quences fox=1/ Tox suivant x

et foy=1/Toy suivant y qui utilise tous les niveaux de gris disponibles. Affichez l'image

cr e.

b) Calculez la fft de cette image et affichez en le module

Faites une boucle qui reprend les tapes a) et b) en changeant la valeur de la fr quence

chaque passage dans la boucle (utilisez imshow(Image); drawnow; pour forcer

l'affichage).

Tox=0; % d finition de la p riode suivant x

Toy=0; % d finition de la p riode suivant y

while(Tox<=1024)

Tox=Tox+1;

3

Travaux pratiques de traitement dimage

P. Falgayrettes

Toy=Toy+1;

for i=1:128

for j=1:128

mon_cos(i,j)=double(128(1+cos(2pi*(i/Tox+j/Toy))));

end

end

figure (1);

imshow(((mon_cos-min(min(mon_cos)))*1/(max(max(mon_cos)))-

min(min(mon_cos))));

drawnow;

colormap(pink(256));

axis('image');

figure (2);

colormap(pink(256));

axis('image');

FFT_MonCos=fft2(mon_cos);

FFT_MonCos=fftshift(FFT_MonCos);

Image_affiche=abs(FFT_MonCos);

imshow(((Image_affiche-min(min(Image_affiche)))*256/

(max(max(Image_affiche)))-min(min(Image_affiche))));

drawnow;

end

4 . COMPRESSION :

Vous allez maintenant effectuer une compression JPEG. Cette compression utilise la

transform e en Cosinus Discret (DCT). Vous verrez que ce type de compression assume une

perte d'information lors de l'op ration.

4.1 . Principe :

Ce type de compression repose sur plusieurs tapes :

d couper l'image en blocs de taille identique (8x8 en pratique).

calculer la DCT de chaque bloc.

repr senter le r sultat de la DCT 2D en regroupant les pixels de tous les blocs par

ordre de coefficient.

supprimer les blocs de coefficients qui repr sentent les fr quences les plus lev es.

calculer les DCT inverse de chaque bloc.

regrouper les blocs dans une m me image.

4.2 . Manipulation :

Pour commencer, il est int ressant de regarder le r sultat obtenu en effectuant la

DCT de l'image compl te. Utilisez l'image 'Clown.bmp'.

MonClown=imread('clown.bmp');

figure (1);

imagesc(DCT);

colormap(gray(256));

DCT=dct2(MonClown);

4

Travaux pratiques de traitement dimage

P. Falgayrettes

Affichez le r sultat de cette DCT :

% affichage du r sultat

figure (2);

imagesc(DCT);

colormap(pink(256));

Remarquez que l'origine des fr quences se trouvent en haut gauche de l'image.

Ci-dessous, vous trouverez les codes qui vous permettent de r aliser ces diff rentes tapes.

Mettre z ro les pixels qui sont au-del de la bissectrice rouge figure ci-dessous.

DCT(0,0)

Pixels

conserver

Pixels mis 0

Calulez la DCT inverse et afficher l'image obtenue (image 'Ma_DCT_inv').

Estimez l'erreur induite par la mise z ro des pixels en faisant la diff rence entre

l'image originale et 'Ma_DCT_inv'.

D coupage des blocs & calcul DCT par bloc :

tbloc=8; % taille du bloc dans l'image (2,4,8,16,...)

%D coupage des blocs & calcul DCT par bloc

for i=1:round(Ncol/tbloc)

for j=1:round(Nlin/tbloc)

% D coupage des blocs

for k=1:tbloc

for l=1:tbloc

%i,j,k,l

bloc(k,l)=I((i-1)tbloc+k,(j-1)tbloc+l);

end

end

% calcul DCT d'un bloc

dct_image(i,j,:,:) = dct2(bloc);

end

end

Premier type d'affichage du r sultat :

Dans ce mode de repr sentation de la DCT2d, on conserve la localisation des blocs de

limage originale : observez lirr gularit spatiale de la r partition de linformation.

% reconstitution des blocs DCT en une image avec pr servation

de la position des blocs de limage originale

for i=1:round(Ncol/tbloc)

for j=1:round(Nlin/tbloc)

for k=1:tbloc

for l=1:tbloc

DCT1((i-1)tbloc+k,(j-1)tbloc+l) = dct_image(i,j,k,l);

end

end

end

end

% affichage du r sultat

figure (2);

5

Travaux pratiques de traitement dimage

P. Falgayrettes

imagesc(DCT1)

colormap(gray(256));

Deuxi me type d'affichage du r sultat :

Dans ce mode de repr sentation, on regroupe maintenant les pixels de tous les blocs par

coefficient (tous les premiers pixels de chaque bloc sont regroup s dans un bloc, tous les

second pixels sont regroup s dans un second bloc, etc.). Vous utiliserez ce mode d'affichage

pour valuer le niveau de d gradation acceptable de votre compression.

% mode de repr sentation de la DCT2d par regroupement des

pixels de tous les blocs par coefficient.

for k=1:tbloc

for l=1:tbloc

for i=1:Ncol/tbloc

for j=1:Nlin/tbloc

DCT2((k-1)Ncol/tbloc+i,(l-1)Nlin/tbloc+j) = dct_image(i,j,k,l);

end

end

end

end

% affichage du r sultat

figure (3);

imagesc(DCT2)

colormap(gray(256));

Reconstruction de l'image : d compression

%Reconstruction par DCT inverse de l'image

for i=1:round(Ncol/tbloc)

for j=1:round(Nlin/tbloc)

bloc=idct2(masque.*reshape(dct_image(i,j,:,:),tbloc,t

bloc));

for k=1:tbloc

for l=1:tbloc

ImReconstitue((i-1)tbloc+k,(j-1)tbloc+l) = bloc(k,l);

end

end

end

end

% affichage

figure (4);

imagesc(ImReconstitue)

colormap(gray(256));

Notez que dans cette derni re op ration appara t une variable nomm e

masque . Il s'agit d'une matrice de dimension (tbloc x tbloc) que vous

pouvez remplir de 0 et de 1. Un 0 signifie que vous perdrez cette partie de la

zone de fr quences, alors qu'un 1 signifie au contraire que vous souhaitez

conserver ces fr quences. Essayez diff rentes solutions.

Rappel : les basses fr quences sont en haut et gauche de votre image dans la repr sentation 2.

6