Applications of digital signal processing - Continuous and Discrete Time Signals and Systems - page 779

% 2D Bartlett window
>> x = bartlett(32);
>> for i = 1:32
zx(i,:) = x';
zy(:,i) = x;
>> end
>> bartlett2D = zx .* zy;
>>n=0;
% calculate power spectrum
>> P = zeros(256,256);
>> for (i = 1:16:320)
for (j = 1:16:288)
Isub = I(i:i+31,j:j+31).*bartlett2D;
P=P+fftshift(fft2(Isub,256,256));
n=n+1;
end
>> end
>> Pabs = (abs(P)/n).ˆ2;
>> mesh([-128:127]*2/256,[-128:127]*2/256,Pabs/
max(max(Pabs)));
The resulting power spectrum is shown in Fig. 17.20(b), where we see peaks at
frequencies [ Ω x , Ω y ] given by [0, 0], [0, 0.5 π ], and [ 0.5 π , 0]. The peak at
[0, 0] corresponds to the dc gain, whereas the remaining peaks are because of
the additive noise that has corrupted the image. We now attempt to eliminate
the noise with a lowpass filter.
Figure 17.21(a) shows the magnitude spectrum of the filter with the impulse
response specified in Eq. (17.27). We use the following M ATLAB code to plot
the magnitude spectrum:
>> h = 1/64*[12321;23432;34543;23432;1
2321];
>> H = fftshift(fft2(h,256,256));
% magnitude spectrum with 256-pt fft
% 2D mesh plot with frequency axis normalized to pi
>> mesh([-128: 127]*2*1/256, [-128:127]*2*1/256, abs(H));
Since the filter provides a higher gain at the lower frequencies and lower gain
at higher frequencies, it is clear that Fig. 17.21(a) corresponds to a lowpass
filter. Note that the gain at frequencies [0, 0.5 π ] and [ 0.5 π , 0] is zero,
therefore the lowpass filter would eliminate the additive noise. The filter2