Matlab Low Pass Filter Using fft

Question

I implemented a simple low-pass filter in Matlab using forward and reverse fft. It works in principle, but the minimum and maximum values ​​are different from the original.

signal = data; %% fourier spectrum % number of elements in fft NFFT = 1024; % fft of data Y = fft(signal,NFFT)/L; % plot(freq_spectrum) %% apply filter fullw = zeros(1, numel(Y)); fullw( 1 : 20 ) = 1; filteredData = Y.*fullw; %% invers fft iY = ifft(filteredData,NFFT); % amplitude is in abs part fY = abs(iY); % use only the length of the original data fY = fY(1:numel(signal)); filteredSignal = fY * NFFT; % correct maximum clf; hold on; plot(signal, 'g-') plot(filteredSignal ,'b-') hold off; 

the resulting image is as follows:

What am I doing wrong? If I normalize both data, the filtered signal looks correct.

Answer 1

Just to remind yourself how MATLAB stores frequency content for Y = fft(y,N) :

• Y(1) - constant offset
• Y(2:N/2 + 1) - many positive frequencies
• Y(N/2 + 2:end) is a set of negative frequencies ... (usually we will build it to the left of the vertical axis)

To create a true low pass filter, we must keep both low positive frequencies and low negative frequencies.

Here is an example of this with a multiplier rectangle filter in the frequency domain, as you did:

 % make our noisy function t = linspace(1,10,1024); x = -(t-5).^2 + 2; y = awgn(x,0.5); Y = fft(y,1024); r = 20; % range of frequencies we want to preserve rectangle = zeros(size(Y)); rectangle(1:r+1) = 1; % preserve low +ve frequencies y_half = ifft(Y.*rectangle,1024); % +ve low-pass filtered signal rectangle(end-r+1:end) = 1; % preserve low -ve frequencies y_rect = ifft(Y.*rectangle,1024); % full low-pass filtered signal hold on; plot(t,y,'g--'); plot(t,x,'k','LineWidth',2); plot(t,y_half,'b','LineWidth',2); plot(t,y_rect,'r','LineWidth',2); legend('noisy signal','true signal','+ve low-pass','full low-pass','Location','southwest') 

A full low-frequency fitler does a better job, but you will notice that the reconstruction is a bit “wavy”. This is because the multiplication with the rectangle function in the frequency domain coincides with the convolution with the sinc function in the temporary domain . The convolution with the sinc function replaces each point with a very uneven average value of its neighbors, therefore, the "wave" effect.

The Gaussian filter has more pleasant properties of the low-pass filter, because the Fourier transform of the Gaussian is Gaussian . Gaussian decays to zero, so distant neighbors are not included in the weighted average during convolution. Here is an example with a Gaussian filter preserving positive and negative frequencies:

 gauss = zeros(size(Y)); sigma = 8; % just a guess for a range of ~20 gauss(1:r+1) = exp(-(1:r+1).^ 2 / (2 * sigma ^ 2)); % +ve frequencies gauss(end-r+1:end) = fliplr(gauss(2:r+1)); % -ve frequencies y_gauss = ifft(Y.*gauss,1024); hold on; plot(t,x,'k','LineWidth',2); plot(t,y_rect,'r','LineWidth',2); plot(t,y_gauss,'c','LineWidth',2); legend('true signal','full low-pass','gaussian','Location','southwest') 

As you can see, the reconstruction is much better.