Tuesday, January 26, 2010

Coherent Noise Averaging




% James Eastham
% Member, IEEE
% Revision: R1
% Contact: james.eastham@ieee.org
%
% Description: This script shows the affect of coherent averaging to
% reduce random noise.

%% Affect of Additive Noise
close all; % close all plots
clear all; % clear all variables
fs = 500; % sample freq
Ts=1/fs; % sample period
n=0:1*fs; % time window (fs=samples per second), 1 second
tn=n*Ts; % sample index
x=cos(2*pi*50*tn); % 50Hz sinusoid
k=.5; %noise amplitude
n=randn(size(x)); % random noise
s=x+(k*n); %signal + random noise with amplitude k
% Time domain plot of x(t)=cos(2*pi*50*tn)
figure('Color',[1 1 1]);
plot(tn,s,'g','LineWidth',2);box off;hold on;
title('Sinusoid w/ Additive White Noise');
xlabel('Time(s)');
ylabel('Amplitude');
% Frequency Domain Plot
figure('Color',[1 1 1]);
[X,f]=ComputeSpectrum(s,fs,2^12);
orig_X=X;
h=plot(f,20*log(X));
box off;
grid on;
set(h,'Linewidth',1);
xlabel('Frequency(Hz)');
ylabel('Amplitude(dB)');
title(['Frequency Spectrum w/ Additive Noise k=.5']);
axis([0 250 -50 120]);

%% Example of Noise Averaging, 3 Coherent Averages
fs = 500; % sample freq
Ts=1/fs; % sample period
n=0:1*fs; % time window (fs=samples per second), 1 second
tn=n*Ts; % sample index
x=cos(2*pi*50*tn);
k=.2; %Noise amplitude
n1=randn(size(x));
n2=randn(size(x));
n3=randn(size(x));
s1=x+(k*n1);
s2=x+(k*n2);
s3=x+(k*n3);
average_s=((s1+s2+s3)/3);
% Time domain plot of x(t)=cos(2*pi*50*tn)
figure('Color',[1 1 1]);
h=plot(tn,average_s,'g','LineWidth',2);box off;hold on;
title('Sinusoid w/ Additive White Noise: 3 Averages');
xlabel('Time(s)');
ylabel('Amplitude');
% Frequency Domain Plot
figure('Color',[1 1 1]);
[X,f]=ComputeSpectrum(average_s,fs,2^12);
plot(f,20*log(X),f,20*log(orig_X),'--r');
box off;
grid on;
set(h,'Linewidth',1);
xlabel('Frequency(Hz)');
ylabel('Amplitude(dB)');
legend('3 Averages','Original');
title(['Frequency Spectrum w/ Noise Averaging: 3 Averages']);
axis([0 250 -50 120]);

%% Example of Noise Averaging, 10 Coherent Averages
fs = 500; % sample freq
Ts=1/fs; % sample period
n=0:1*fs; % time window (fs=samples per second), 1 second
tn=n*Ts; % sample index
x=cos(2*pi*50*tn);
k=.2; %Noise amplitude
n1=randn(size(x));
n2=randn(size(x));
n3=randn(size(x));
n4=randn(size(x));
n5=randn(size(x));
n6=randn(size(x));
n7=randn(size(x));
n8=randn(size(x));
n9=randn(size(x));
n10=randn(size(x));
s1=x+(k*n1);
s2=x+(k*n2);
s3=x+(k*n3);
s4=x+(k*n1);
s5=x+(k*n2);
s6=x+(k*n3);
s7=x+(k*n1);
s8=x+(k*n2);
s9=x+(k*n3);
s10=x+(k*n1);
average_s=((s1+s2+s3+s4+s5+s6+s7+s8+s9+s10)/10);
% Time domain plot of x(t)=cos(2*pi*50*tn)
figure('Color',[1 1 1]);
h=plot(tn,average_s,'g','LineWidth',2);box off;hold on;plot(tn,x,'--r');
title('Sinusoid w/ Additive White Noise: 10 Averages');
xlabel('Time(s)');
ylabel('Amplitude');
legend('composite signal','averaged');
% Frequency Domain Plot
figure('Color',[1 1 1]);
[X,f]=ComputeSpectrum(average_s,fs,2^12);
plot(f,20*log(X),f,20*log(orig_X),'--r');
box off;
grid on;
set(h,'Color',[0 0.8 1]);
set(h,'Linewidth',1);
xlabel('Frequency(Hz)');
ylabel('Amplitude(dB)');
title(['Frequency Spectrum w/ Noise Averaging: 10 Averages']);
legend('10 Averages','Original');
axis([0 250 -50 120]);

No comments: