Airplane Noise Removal Demonstration: Difference between revisions

From Class Wiki
Jump to navigation Jump to search
(Created page with ' % Demonstration of LMS algorithm for noise cancellation. % Rich Kozick, Spring 1997 % Rob Frohne's modifications for Aircraft Noise Cancellation. % modifications for Macintos…')
 
No edit summary
Line 4: Line 4:
% modifications for Macintosh 2000. (commented out here)
% modifications for Macintosh 2000. (commented out here)
% Rob Frohne's modifications for Linux.
% Rob Frohne's modifications for Linux.

% Aircraft noise cancellation simulation.
% Aircraft noise cancellation simulation.

% Desired signal
% Desired signal

clear all
clear all
Totaltime=3;
Totaltime=3;
system("aoss espeak 'Hit a key and speak the signal.'");
system("aoss espeak 'Hit a key and speak the signal.'");

Fs=8000; %sample rate.
Fs=8000; %sample rate.
st=record(Totaltime,Fs);
st=record(Totaltime,Fs);
s = st';
s = st';
Ls = length(s);
Ls = length(s);

% Interference + random noise
% Interference + random noise
system("aoss espeak 'Hit a key and make the noise!'");
system("aoss espeak 'Hit a key and make the noise!'");
Line 35: Line 30:
nf = filter(bn,an,n);
nf = filter(bn,an,n);
%nf = [nf(Dn:(Ls)); zeros(Dn-1,1)] + Sign*randn(Ls,1);
%nf = [nf(Dn:(Ls)); zeros(Dn-1,1)] + Sign*randn(Ls,1);


y = s + nf;
y = s + nf;

%Sigx = 0.01;
%Sigx = 0.01;
%ax = [1]; % bx and ax are filtering on n to produce x
%ax = [1]; % bx and ax are filtering on n to produce x
Line 45: Line 37:
%x = filter(bx, ax, n);%x = [x(Dx:(Ls)); zeros(Dx-1,1)] + Sigx*randn(Ls,1);
%x = filter(bx, ax, n);%x = [x(Dx:(Ls)); zeros(Dx-1,1)] + Sigx*randn(Ls,1);
%x = x + Sigx*randn(1,Ls);
%x = x + Sigx*randn(1,Ls);

system("aoss espeak 'Here is the noisy signal.'");
system("aoss espeak 'Here is the noisy signal.'");
soundsc(y,Fs);
soundsc(y,Fs);

N = 20; % Length of adaptive filter
N = 20; % Length of adaptive filter

% LMS algorithm for adaptive noise cancellation
% LMS algorithm for adaptive noise cancellation

h = zeros(N,1);
h = zeros(N,1);
%mu = 1/(10*N*var(n));
%mu = 1/(10*N*var(n));
Line 62: Line 54:
h = h - mu*e(k)*xk'/(xk*xk'); % For the previous line.
h = h - mu*e(k)*xk'/(xk*xk'); % For the previous line.
end
end

% The signal estimate is in the vector e
% The signal estimate is in the vector e
%speak('Here is the cleaned signal.');
%speak('Here is the cleaned signal.');

Revision as of 10:39, 30 November 2010

% Demonstration of LMS algorithm for noise cancellation.
% Rich Kozick, Spring 1997
% Rob Frohne's modifications for Aircraft Noise Cancellation.
% modifications for Macintosh 2000. (commented out here)
% Rob Frohne's modifications for Linux. 
% Aircraft noise cancellation simulation.
% Desired signal
clear all
Totaltime=3;
system("aoss espeak 'Hit a key and speak the signal.'");
Fs=8000; %sample rate.
st=record(Totaltime,Fs);
s = st';
Ls = length(s);
% Interference + random noise
system("aoss espeak 'Hit a key and make the noise!'");
%pause();
nt = record(Totaltime, Fs);
n = nt';
%Sign = 0.01;
%Dn=20;  % Delay of the noise that appears in y.
%n = n(1:Ls) + Sign*randn(Ls,1);
%an = [0 .01 -.5 1 -.5 .1 .01 0];
%bn = 4*[0 0 0 0 0 0 0 0 0 0 0 0 0 0 .5 1 .5];
bn=rand(1,20);
an = [1];
%sys = tf(an,bn,1/Fs);
%bode(sys);
%figure(1);
nf = filter(bn,an,n);
%nf = [nf(Dn:(Ls)); zeros(Dn-1,1)] + Sign*randn(Ls,1);
y = s + nf;
%Sigx = 0.01;  
%ax = [1];       % bx and ax are filtering on n to produce x
%bx = [0 0 0 .1 1 .1];
%Dx = 1;        % Delay of n that appears in x
%x = filter(bx, ax, n);%x = [x(Dx:(Ls)); zeros(Dx-1,1)] + Sigx*randn(Ls,1);
%x = x + Sigx*randn(1,Ls);

system("aoss espeak 'Here is the noisy signal.'");
soundsc(y,Fs);

N = 20;         % Length of adaptive filter

% LMS algorithm for adaptive noise cancellation

h = zeros(N,1);
%mu = 1/(10*N*var(n));
mu = .05;
for k=N:Ls
   xk = n(k:-1:(k-N+1));
   nhat(k) = h'*xk';
   e(k) = - y(k) + nhat(k);
   h = h - mu*e(k)*xk'/(xk*xk'); % For the previous line.
end

% The signal estimate is in the vector e
%speak('Here is the cleaned signal.');
system("aoss espeak 'Here is the cleaned signal.'");
soundsc(e(1,6000:end),Fs);
%playaudio(e);  (For Linux without octave-forge.)