Denoising via Filtering in MATLAB

In this post, i am going to show how to denoise a signal using a Gaussian filter in MATLAB

The MATLAB code below generates a noisy signal (i.e., un) from the function sech(x) (i.e., u) by adding a white noise via its frequency domain spectrum (ut is the fft of u and utn is the fft of un)

clear all; close all; clc;

T=30;
n=512;

t2=linspace(-T/2, T/2, n+1);

t=t2(1:n);

k=(2*pi/T)*[0:n/2-1 -n/2:-1];
ks=fftshift(k);

u=sech(t);

ut=fft(u);
noise=7;
utn=ut+noise*(randn(1, n)+i*randn(1, n));
un=ifft(utn);

subplot(2, 1, 1); plot(t, u, 'k', t, abs(un), 'm');
subplot(2, 1, 2); plot(ks, abs(fftshift(ut))/max(abs(fftshift(ut))), 'k', ks, abs(fftshift(utn))/max(abs(fftshift(utn))), 'm');

axis([-25 25 0 1]);

The figure below shows the resulting plots

The MATLAB code below shows the use of a gaussian filter which is multiplied with the fft of the noisy signal (i.e., utn) to obtained the fft of a filtered noisy signal (i.e., utnf) which is then converted back to the time domain filtered noisy signal (i.e., unf) that is denoised

clear all; close all; clc;

T=30;
n=512;

t2=linspace(-T/2, T/2, n+1);

t=t2(1:n);

k=(2*pi/T)*[0:n/2-1 -n/2:-1];
ks=fftshift(k);

u=sech(t);

ut=fft(u);
noise=7;
utn=ut+noise*(randn(1, n)+i*randn(1, n));
un=ifft(utn);

filter=exp(-k.^2); %guassian filter
utnf=filter.*utn;

unf=ifft(utnf);

subplot(2, 1, 1); plot(t, u, 'k', t, un, 'm', t, unf, 'g');
subplot(2, 1, 2); plot(ks, abs(fftshift(ut))/max(abs(fftshift(ut))), 'k', ks, abs(fftshift(utn))/max(abs(fftshift(utn))), 'm', ...
    ks, fftshift(filter), 'b', ...
    ks, abs(fftshift(utnf))/max(abs(fftshift(utnf))), 'g');

axis([-25 25 0 1]);

The figure below shows the resulting plots (where the green line shows the denoised signal in time and frequency domain

Note that in the above MATLAB code it assume that the user knows the frequency they want to look at (in this case it is centered at k=0), therefore the gaussian filter is center at that frequency range (i.e. that guassian is centered at k=0)

The following MATLAB code allows user to vary the width (i.e., sigma) and the location (i.e., the mu) of the Guassian signal

clear all; close all; clc;

T=30;
n=512;

t2=linspace(-T/2, T/2, n+1);

t=t2(1:n);

k=(2*pi/T)*[0:n/2-1 -n/2:-1];
ks=fftshift(k);

u=sech(t);

ut=fft(u);
noise=7;
utn=ut+noise*(randn(1, n)+i*randn(1, n));
un=ifft(utn);

sigma=1;
mu=0;

filter=1 / (sigma * (2*pi)^0.5) * exp(-(k-mu).^2/(2*sigma&2)); %guassian filter
utnf=filter.*utn;

unf=ifft(utnf);

subplot(2, 1, 1); plot(t, u, 'k', t, un, 'm', t, unf, 'g');
subplot(2, 1, 2); plot(ks, abs(fftshift(ut))/max(abs(fftshift(ut))), 'k', ks, abs(fftshift(utn))/max(abs(fftshift(utn))), 'm', ...
    ks, fftshift(filter), 'b', ...
    ks, abs(fftshift(utnf))/max(abs(fftshift(utnf))), 'g');

axis([-25 25 0 1]);