Hi,
I understand that you are trying to perform noise cancellation using Wiener Filter. The provided code does not consider the correlation between the noise in the primary and the secondary sensers (v1 and v2). Consequently, the final signal does not match the desired signal.
Please refer to updated code below to perform noise cancellation using Wiener Filter.
clc;
close all;
%taking the recorded signal as d
[d, fs] = audioread('sound.mp3');
%taking noisy signal 2 (v2 in Monson hayes diagram) as g
n=length(d);
v2=0.5*randn(n,2);
p=20;
v2filt=filter(1, [1 -0.5] ,v2);
%%randn(1,n) is noisy signal 1 (v1) as per monson hayes diagram. we want to
% to approximate v1 using v2
%%%%%%%%%%%%%%
% v1=0.1*randn(n,2);
% Creating another random signal will make v1 and v2filt uncorrelated
% Generating v1 from the initial created random signal v2
v1=filter(1, [1 -0.1] ,v2);
%%%%%%%%%%%%%%
x=d+v1;
Rv1=covar(v2filt,p) ;
rxv=convm(x,p)'*convm(v2filt,p)/(n-1);
w=rxv(1,:)/Rv1;
v1hat=filter(w,1,v2filt);
appx=x-v1hat;
error=mse(d-appx);
%frequency domain mse
P = periodogram(d-appx,[],[],fs);
Hmss = dspdata.msspectrum(P,'Fs',fs,'spectrumtype','onesided');
%display
subplot(5,1,1),plot(d);
title('desired signal:');
subplot(5,1,2),plot(x);
title('desired + noise signal:');
subplot(5,1,3),plot(v1hat);
title('estimated signal:');
subplot(5,1,4),plot(appx);
title('Noise-removed signal');
subplot(5,1,5),plot(Hmss);
%title('frequency domain mean square error signal');
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%% Defining convm and covar functions
function y = convm(x, p)
% Generate a convolution matrix of signal x with given order p
N = length(x);
y = zeros(N-p, p+1);
for n = p+1:N
y(n-p, :) = x(n:-1:n-p);
end
end
function R = covar(x, p)
% Compute the covariance matrix of signal x using the given order p
N = length(x);
R = zeros(p+1, p+1);
for n = p+1:N
x_p = x(n-p:n);
R = R + x_p' * x_p;
end
R = R / (N - p);
end
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
The output of the code for a sample audio file has been shown below.
Hope this helps.
