2006-2007
Digitale Image Processing Labo 1 : Verwijderen van periodische ruis uit beelden
Bart Vanrumste
Alexander Alderweireldt
1 Maak gebruik van “periodic.m” en maak een som van 5 verschillende sinussen. Kies de frequenties (u,v) arbitrair. Bereken de FFT van dit beeld. Source m-file : % opgave 1 : % ========== M=256 N=256 % row vectors c(olum) en r(ow) r = 0:M-1; c = 0:N-1; % c en r -> array C en R [C,R] = meshgrid(c,r) % 5 frequenties (u,v) u0=2*pi/256*10 v0=2*pi/256*4 u1=2*pi/256*6 v1=2*pi/256*7 u2=2*pi/256*9 v2=2*pi/256*4 u3=2*pi/256*15 v3=2*pi/256*10 u4=2*pi/256*18 v4=2*pi/256*22 % 5 verschillende sinussen s0 = sin(u0*R + v0*C); s1 = sin(u1*R + v1*C); s2 = sin(u2*R + v2*C); s3 = sin(u3*R + v3*C); s4 = sin(u4*R + v4*C); % som van de 5 sinussen s = s0+s1+s2+s3+s4; % resultaat weergeven figure, imshow(mat2gray(s),[]) % FFT berekenen van het resultaat en weergeven [a,b]=size(s); S = fft2(s,2*a,2*b); Sc = 1+abs(fftshift(S)); figure, imshow(Sc,[]); % inverse van de vorige afbeelding % is veel duidelijker (en inkt-besparend!!) S_invers = imcomplement(1+abs(S)); Sc_invers = imcomplement(Sc); figure, imshow(S_invers,[]); figure, imshow(Sc_invers,[]);
Figuren :
Afbeelding 1: Som van 5 verschillende sinussen
Afbeelding 2: De FFT van de som van de 5 sinussen (het negatief)
Afbeelding 3: De gecentreerde FFT van de som van de 5 sinussen (het negatief)
Commentaar : Door middel van de FFT in een afbeelding weer te geven, kunnen we het spectrum van de afbeelding makkelijk visueel analyseren. Omwille van het periodisch karakter worden de resultaten in de vier hoeken van de afbeelding weergegeven (zie afbeelding 2). Daarom herleggen we de oorsprong van de transformatie naar het middelpunt van de afbeelding met behulp van de matlab functie 'fftshift' (zie afbeelding 3). De FFT wordt berekend voor twee maal de grootte van het originele beeld om zo ook de cyclische convolutie in rekening te brengen. De afbeelding illustreerd ook duidelijk het symmetrisch karakter van de FFT. Elke sinus wordt weergeven in de FFT door middel van twee punten. (5 x 2 punten) die symmetrisch ten opzichte van de oorsprong van elkaar liggen.
2 Bereken de FFT van het beeld (Fig5.05(a).jpg). Pas een logaritmische schaling toe. Wat merk je op. Source m-file : % opgave 2 : % ========== % de afbeelding inladen en weergeven f=double(imread('cb_periodic.tif')); figure, imshow(mat2gray(f),[]) % de afbeelding grootte ophalen [a,b]=size(f); % FFT en shift F=fft2(f,2*a,2*b); F=fftshift(F); S = log(1+abs(F)); figure, imshow(S,[]); figure, imshow(S,[0.55*max(S(:)) max(S(:))]); % inverse van de vorige afbeelding % is veel duidelijker (en inkt-besparend!!) S_invers = imcomplement(S); figure, imshow(S_invers,[]); figure, imshow(S_invers,[min(S_invers(:)) 0.55*min(S_invers(:))]);
Figuren :
Afbeelding 4: Afbeelding van een printplaat die ruis bevat
Afbeelding 5: De FFT van de afbeelding (het negatief)
Afbeelding 6: De FFT van de afbeelding met zijn waarden begrensd tussen 55% van de maximale waarde en de maximale waarde (het negatief)
Commentaar : Als we afbeelding 6 bestuderen merken we allereerst 2 punten op, die weer symmetrisch tegenover elkaar liggen ten opzichte van de oorsprong. De verbindingslijn tussen deze twee punten geeft de richting van de periodieke storing weer die we in afbeelding 4 konden waarnemen. Uit de eerste opgave hebben we geleerd dat deze 2 punten overeenkomen met een bepaalde frequentie. De punten rond de oorsprong representeren de gemiddelde grijswaarden van de afbeelding.
3 Maak gebruik van een notch filter om de periodisch ruis te verwijderen. Met n =2 en uo en vo de te verwijderen frequenties, te bepalen uit de FFT van opgave 2. Maak voor alle periodische signalen samen één notch filter. Gebruik ook “example4_8.m”
1
H (u , v) =
n
D02 1+ D1 (u , v) D2 (u , v) D1 (u , v) = (u − u 0 ) 2 + (v − v0 ) 2 1/ 2
[ D (u , v) = [(u + u ) D = [ (u ) + (v ) ]
2
] + (v + v ) ]
2
0
2
0
0
2 1/ 2
0
2 1/ 2
0
Formule 1 : De Notch Filter Source m-file : % opgave 3 : % ========== % de afbeelding inladen en weergeven f=double(imread('cb_periodic.tif')); figure, imshow(mat2gray(f),[]) FFT_opg2 = fft2(f,2*a,2*b); u = -a+1:a; v = -b+1:b; [V,U]=meshgrid(v,u); % coördinaten (x=v, y=u) % middelpt : x = 465 y = 449 % pt linksboven : x = 385 y = 249 % U0, V0 geven de afstand van 1 punt tot de oorsprong weer U0 = 449-249; V0 = 465-385; % Notch filter (eps = 2^(-52)) % Afhankelijk van n, krijgen we een groter of kleiner filtergebied n = 2; D0 =(U0.^2+V0.^2).^0.5; D1 =((U-U0).^2+(V-V0).^2).^0.5; D2 =((U+U0).^2+(V+V0).^2).^0.5; H = 1./(1+(D0^2./((D1.*D2)+eps)).^n); % geef de Notch filter weer figure, imshow(log(1+H),[]); % resulting spectrum G : H = fftshift(H); G = H.*FFT_opg2; % -> inverse fft -> beeld G = real(ifft2(G)); figure, imshow(G(1:a,1:b),[]); % de FFT na het filteren G = H.*Ftt_opg2; figure, imshow(log(1+abs(fftshift(G))),[]);
Figuren :
Afbeelding 7: Afbeelding van een printplaat die ruis bevat
Afbeelding 8: weergave van de coördinaten van een ruispunt bekomen door de FFT op de originele afbeelding
Afbeelding 9: De Notch Filter
Afbeelding 10: De afbeelding na een Notch filtering
Afbeelding 11: De FFT na de Notch filtering
Commentaar : De bedoeling van deze opgave is de periodieke ruis die aanwezig is bij de oorspronkelijke afbeelding (zie afbeelding 7) te verwijderen door middel van een Notch filter (zie afbeelding 9). In opgave 1 hebben we duidelijk gezien dat periodieke ruis overeenkomt met twee symmetrische punten. We trachten nu met de Notch filter deze twee ruispunten weg te filteren. De filter wordt gerealiseerd aan de hand van de volgende formule :
1
H (u , v) =
n
D02 1+ D1 (u , v) D2 (u , v) D1 (u , v) = (u − u0 ) 2 + (v − v0 ) 2 1/ 2
[ D (u, v) = [(u + u ) D = [(u ) + (v ) ]
2
2
0
2
0
0
] + (v + v ) ]
2 1/ 2
0
2 1/ 2
0
u = x-richting v = y-richting H = de Notch filter. Deze wordt berekend aan de hand van de afstand tussen het ruispunt en de oorsprong van de FFT. N = afhankelijk van n krijgen we een groter of een kleiner frequentiegebied.
Aan afbeelding 10 kunnen we zien dat de Notch filter de periodieke ruis verwijdert heeft. Als we afbeelding 11 verder analyseren merken we dat de ruispunten die verantwoordelijk waren voor de periodieke ruis uit de FFT zijn gefilterd.