Hello everyone, I'm trying to blur an image by removing the high frequencies from the DFT of that image. after reading the image i used the fft command, now i have the frequencies but can't figure out how to remove the high ones
any suggestion??
fft
Nx = 64; % Number of samples collected along first dimension
Ny = 32; % Number of samples collected along second dimension
dx = 1; % x increment (i.e., Spacing between each column)
dy = .1; % y increment (i.e., Spacing between each row)
x = 0 : dx : (Nx-1)*dx; % x-distance
y = 0 : dy : (Ny-1)*dy; % y-distance
data_spacedomain = randn(Ny,Nx); % random 2D matrix
Nyq_kx = 1/(2*dx); % Nyquist of data in first dimension
Nyq_ky = 1/(2*dy); % Nyquist of data in second dimension
dkx = 1/(Nx*dx); % x-Wavenumber increment
dky = 1/(Ny*dy); % y-Wavenumber increment
kx = -Nyq_kx : dkx : Nyq_kx-dkx; % x-wavenumber
ky = -Nyq_ky : dky : Nyq_ky-dky; % y-wavenumber
data_wavenumberdomain = zeros(size(data_spacedomain)); % initialize data
% Compute 2D Discrete Fourier Transform
for i1 = 1 : Nxfor j1 = 1 : Nyfor i2 = 1 : Nxfor j2 = 1 : Nydata_wavenumberdomain(j1,i1) = data_wavenumberdomain(j1,i1) + ...(2i*pi*ky(j1))*(2i*pi*kx(i1))* ...(data_spacedomain(j2,i2)*exp(-1i*(2*pi)*(kx(i1)*x(i2)+ky(j1)*y(j2))));endendendenddifferentiated_data = ifft2(ifftshift(data_wavenumberdomain));
Nx = 64; % Number of samples collected along first dimensionNy = 32; % Number of samples collected along second dimensiondx = 1; % x increment (i.e., Spacing between each column)dy = .1; % y increment (i.e., Spacing between each row)x = 0 : dx : (Nx-1)*dx; % x-distancey = 0 : dy : (Ny-1)*dy; % y-distancedata_spacedomain = randn(Ny,Nx); % random 2D matrixNyq_kx = 1/(2*dx); % Nyquist of data in first dimensionNyq_ky = 1/(2*dy); % Nyquist of data in second dimensiondkx = 1/(Nx*dx); % x-Wavenumber incrementdky = 1/(Ny*dy); % y-Wavenumber incrementkx = -Nyq_kx : dkx : Nyq_kx-dkx; % x-wavenumberky = -Nyq_ky : dky : Nyq_ky-dky; % y-wavenumber% Compute 2D FFT
data_wavenumberdomain = fft2(data_spacedomain);% Compute grid of wavenumbers
[KX, KY] = meshgrid(ifftshift(kx),ifftshift(ky));% Compute 2D derivative
data_wavenumberdomain_differentiated = (2i*pi)^2*KX.*KY.*data_wavenumberdomain; % Convert back to space domain
data_spacedomain_differentiated = ifft2(data_wavenumberdomain_differentiated );
d_MAT_x = 2i*pi*KX.*data_wavenumberdomain;d_MAT_y = 2i*pi*KY.*data_wavenumberdomain;
[KX,FQ]=meshgrid(ifftshift(kx),ifftshift(fq));
Nyq_fq = 1/(2*dt);Nyq_kx = 1/(2*dx);dfq = 1/(Nt*dt);dkx = 1/(Nx*dx);if mod(Nt,2) == 0 % Nt is even
fq = -Nyq_fq : dfq : Nyq_fq-dfq;else % Nt is odd
fq = [sort(-1*(dfq:dfq:Nyq_fq)) (0:dfq:Nyq_fq)];endif mod(Nx,2) == 0 % Nx is even
kx = -Nyq_kx : dkx : Nyq_kx-dkx;else % Nx is odd
kx = [sort(-1*(dkx:dkx:Nyq_kx)) (0:dkx:Nyq_kx)];end
M=8;N=16; N=8;M=16; dx=0.1;dy=0.2; f = randn(M,N); % Energy in time domain
energy_f = sum(sum(f.^2))*dx*dy dkx=1/(N*dx);dky=1/(M*dy); F = fftshift(fft2(f))*dx*dy; % Energy in wavenumber domain
energy_F = sum(sum(abs(F).^2))*dkx*dky % Energy using "half" the spectrum
energy_F2 = 2*sum(sum(abs(F(2:M/2, :).^2)))*dkx*dky + ... 2*sum(sum(abs(F(M/2+1,2:N/2).^2)))*dkx*dky + ... sum(sum(abs(F(M/2+1, 1).^2)))*dkx*dky + ... sum(sum(abs(F(M/2+1,N/2+1).^2)))*dkx*dky + ... 2*sum(sum(abs(F(1 ,2:N/2).^2)))*dkx*dky + ... sum(sum(abs(F(1 , 1).^2)))*dkx*dky + ... sum(sum(abs(F(1 ,N/2+1).^2)))*dkx*dky % Plot spectrum
Nyq_x = 1/2/dx; Nyq_y = 1/2/dy; kx = -Nyq_x : dkx : Nyq_x-dkx; ky = -Nyq_y : dky : Nyq_y-dky; imagesc(1:N,1:M,abs(F)) set(gca,'YTick',1:M,'YTickLabel',ky); set(gca,'XTick',1:N,'XTickLabel',kx);
Best Answer