用MATLAB實現頻域平滑濾波以及圖像去噪代碼

%%%%%%%%spatial frequency (SF) filtering by low pass filter%%%%%%%%

% the SF filter is unselective to orientation (doughnut-shaped in the SF

% domain).

[FileName,PathName,FilterIndex] = uigetfile ;

filename = fullfile(PathName, FileName) ;

[X map] = imread(filename, fmt); % read image

L = double(X); % transform to double

%%%%%%%%%%%%% need to add (-1)x+y to L

% calculate the number of points for FFT (power of 2)

fftsize = 2 .^ ceil(log2(size(L)));

% 2d fft

Y = fft2(X, fftsize(1), fftsize (2));

Y = fftshift(Y);

% obtain frequency (cycles/pixel)

f0 = floor([m n] / 2) + 1;

fy = ((m: -1: 1) - f0(1) + 1) / m;

fx = ((1: n) - f0(2)) / n;

[mfx mfy] = meshgrid(fx, fy);

% calculate radius

SF = sqrt(mfx .^ 2 + mfy .^ 2);

% SF-bandpass and orientation-unselective filter

filt = SF > k0;

A_filtered = filt .* A; % SF filtering

L_filtered = real(ifft2(ifftshift(A_filtered))); % IFFT

L_filtered = L_filtered(1: size(L, 1), 1: size(L, 2));

%%%%%%%%%%need to add (-1)x + y to L_filtered

% show

figure(1);

clf reset;

colormap gray;

% plot image

subplot(2, 2, 1);

imagesc(L);

colorbar;

axis square;

set(gca, 'TickDir', 'out');

title('original image');

xlabel('x');

ylabel('y');

imwrite(L, fullfile(FilePath, 'original image.bmp'), 'bmp') ;

% plot amplitude

A = abs(A);

A = log10(A);

% spectral amplitude

subplot(2, 2, 2);

imagesc(fx, fy, A);

axis xy;

axis square;

set(gca, 'TickDir', 'out');

title('amplitude spectrum');

xlabel('fx (cyc/pix)');

ylabel('fy (cyc/pix)');

imwrite(A, fullfile(FilePath, 'amplitude spectrum.bmp'), 'bmp') ;

% filter in the SF domain

subplot(2, 2, 3);

imagesc(fx, fy, filt);

axis xy;

axis square;

set(gca, 'TickDir', 'out');

title('filter in the SF domain');

xlabel('fx (cyc/pix)');

ylabel('fy (cyc/pix)');

imwrite(filt, fullfile(FilePath, 'filter in SF.bmp'), 'bmp') ;

% filtered image

subplot(2, 2, 4);

imagesc(L_filtered);

colorbar;

axis square;

set(gca, 'TickDir', 'out');

title('filtered image');

xlabel('x');

ylabel('y');

imwrite(filtered, fullfile(FilePath, 'filtered image.bmp'), 'bmp');

%%%%%%%%%%%%%%%%%median filter%%%%%%%%%%%%%%%%

[FileName,PathName,FilterIndex] = uigetfile ;

filename = fullfile(PathName, FileName) ;

[LNoise map] = imread(filename, fmt); % read image

L = medfilt2(LNoise, [3 3]); % remove the noise with 3*3 block

figure ;

imshow(LNoise) ;

title('image before fitlering') ;

figure

imshow(L)

title('filtered image') ;

imwrite(FilePath, 'filtered image.bmp', bmp)

上一篇:談SI住宅生態建築的模型建構？

下一篇:湖南軟件學院學費多少