|
查看: 1917|回复: 2
|
會在matlab里用cross correlation function的朋友請教下我
[复制链接]
|
|
|
本人在做著一份报告
关於noise reduction technique的
现在
我在matlab里輸入两种不一样程度的杂音
即左右耳各一个
我得把它們的frequency index分成二十一个subbands
然后左耳的第一个subband和右耳的第一个subband做cross colleration
左耳的第六个subband和右耳的第六个subband做cross colleration
以此类推
以下是我做的codes
但我每一次run都會有不同的結果
所以我应該做錯了
还有
如果在cross correlation我沒加上nlags
matlab會說"max lag must be integer"
什么意思?
要怎样做呢?
拜託了
很急
[ 本帖最后由 山羊座 于 8-11-2006 01:02 AM 编辑 ] |
|
|
|
|
|
|
|
|
|
|
|
楼主 |
发表于 7-11-2006 07:22 PM
|
显示全部楼层
n=512;
Fs=16000;
wavwrite(R_input,Fs,'c:\forgetnoiseR.wav');
wavwrite(L_input,Fs,'c:\forgetnoiseL.wav');
% Read the wave file in order to define the mixtured signal for both
% channels
x=wavread('C:\forgetnoiseR.wav'); % right channel
[x,fs,bits]=wavread('C:\forgetnoiseR.wav');
% wavplay (x,Fs)
y=wavread('C:\forgetnoiseL.wav'); % left channel
[y,fs,bits]=wavread('C:\forgetnoiseL.wav');
% wavplay (y,Fs)
% Buffering for 0.5 rate overlap-add, right channel
x_buffer=buffer(x,32,16,'nodelay'); %Window size for 16kHz is 32ms
length_a=length(x_buffer);
k=1;
for i=1:length_a,
for j=1:32,
row(k)=x_buffer(j,i);
k=k+1;
end;
end;
% Subsapmle
j=1;
for i=1:32:length(x( )
w(j)=x(i);
j=j+1;
end;
s=fft(w,n);
s=s*1000;
% Prepare the Hanning window
h = hanning(512, 'periodic');
% Power Density Spectrum
X=[];
for i=1:512;
X(i) = max(20 * log10(abs(1/n*s(i)* h)));
end;
% Calculation of the Power Sprectrum
X=abs(X).^2;
% Buffering for 0.5 rate overlap-add, left channel
y_buffer=buffer(y,32,16,'nodelay'); %Window size for 16kHz is 32ms
length_b=length(y_buffer);
k=1;
for i=1:length_b,
for j=1:32,
row(k)=y_buffer(j,i);
k=k+1;
end;
end;
% Subsapmle
j=1;
for i=1:32:length(y()
v(j)=y(i);
j=j+1;
end;
sy=fft(v,n);
sy=sy*1000;
% Prepare the Hanning window
h = hanning(512, 'periodic');
% Power Density Spectrum
Y=[];
for i=1:512;
Y(i) = max(10 * log10(abs(1/n*sy(i)* h))); % freq=1:n/2
end;
% Calculation of the Power Sprectrum
Y=abs(Y).^2;
% Table of the absolute threshold, available in Table D.1a
L = [
% Frequency, Hz(:,1) | Critical Band Rate, z(:,2) | Absolute Threshold, dB(:,3) |
31.25 0.309 58.23 ; 62.50 0.617 33.44 ;
93.75 0.925 24.17 ; 125.00 1.232 19.20 ;
156.25 1.538 16.05 ; 187.50 1.842 13.87 ;
218.75 2.145 12.26 ; 250.00 2.445 11.01 ;
281.25 2.742 10.01 ; 312.50 3.037 9.20 ;
343.75 3.329 8.52 ; 375.00 3.618 7.94 ;
406.25 3.903 7.44 ; 437.50 4.185 7.00 ;
468.75 4.463 6.64 ; 500.00 4.736 6.28 ;
531.25 5.006 5.97 ; 562.50 5.272 5.70 ;
593.75 5.533 5.44 ; 625.00 5.789 5.21 ;
656.25 6.041 5.00 ; 687.50 6.289 4.80 ;
718.75 6.532 4.62 ; 750.00 6.770 4.45 ;
812.50 7.004 4.29 ; 812.50 7.233 4.14 ;
843.75 7.457 4.00 ; 875.00 7.677 3.86 ;
906.25 7.892 3.73 ; 937.50 8.103 3.61 ;
968.75 8.309 3.49 ; 1000.00 8.511 3.37 ;
1031.25 8.708 3.26 ; 1062.50 8.901 3.15 ;
1093.75 9.090 3.04 ; 1125.00 9.275 2.93 ;
1156.25 9.456 2.83 ; 1187.50 9.632 2.73 ;
1218.75 9.805 2.63 ; 1250.00 9.974 2.53 ;
1281.25 10.139 2.42 ; 1312.50 10.301 2.32 ;
1343.75 10.459 2.22 ; 1375.00 10.614 2.12 ;
1406.25 10.765 2.02 ; 1437.50 10.913 1.92 ;
1468.75 11.058 1.81 ; 1500.00 11.199 1.71 ;
1562.50 11.474 1.49 ; 1625.00 11.736 1.27 ;
1687.50 11.988 1.04 ; 1750.00 12.230 0.80 ;
1812.50 12.461 0.55 ; 1875.00 12.684 0.29 ;
1937.50 12.898 0.02 ; 2000.00 13.104 -0.25 ;
2062.50 13.302 -0.54 ; 2125.00 13.493 -0.83 ;
2187.50 13.678 -1.12 ; 2250.00 13.855 -1.43 ;
2312.50 14.027 -1.73 ; 2375.00 14.193 -2.04 ;
2437.50 13.354 -2.34 ; 2500.00 14.509 -2.64 ;
2562.50 14.660 -2.93 ; 2625.00 14.807 -3.22 ;
2687.50 14.949 -2.49 ; 2750.00 15.087 -3.74 ;
2812.50 15.221 -3.98 ; 2875.00 15.351 -4.20 ;
2937.50 15.478 -4.40 ; 3000.00 15.602 -4.57 ;
3125.00 15.841 -4.82 ; 3250.00 16.069 -4.96 ;
3375.00 16.496 -4.98 ; 3500.00 16.496 -4.90 ;
3625.00 16.697 -4.70 ; 3750.00 16.891 -4.39 ;
3875.00 17.078 -3.99 ; 4000.00 17.259 -3.51 ;
4125.00 17.434 -2.99 ; 4250.00 17.605 -2.45 ;
4375.00 17.770 -1.90 ; 4500.00 17.932 -1.37 ;
4625.00 18.089 -0.86 ; 4750.00 18.242 -0.39 ;
4875.00 18.392 0.03 ; 5000.00 18.539 0.40 ;
5125.00 18.682 0.72 ; 5250.00 18.823 1.00 ;
5375.00 18.960 1.24 ; 5500.00 19.095 1.44 ;
5625.00 19.226 1.62 ; 5750.00 19.356 1.78 ;
5875.00 19.482 1.92 ; 6000.00 19.606 2.05 ;
6125.00 19.728 2.18 ; 6250.00 19.847 2.30 ;
6375.00 19.964 2.42 ; 6500.00 20.079 2.55 ;
6625.00 20.191 2.69 ; 6750.00 20.300 2.82 ;
6875.00 20.408 2.97 ; 7000.00 20.513 3.13 ;
7125.00 20.616 3.29 ; 7250.00 20.717 3.46 ;
7375.00 20.815 3.65 ; 7500.00 20.912 3.84 ;
]; |
|
|
|
|
|
|
|
|
|
|
|
楼主 |
发表于 7-11-2006 07:22 PM
|
显示全部楼层
LTm=[];
for n = 1 : 21 % 21 subbands
if (n == 1) % find minimum masking threshold for each subbands
for L = 1 : 3 % subband indexes
LTm(n)=(X(L));
end
elseif (n == 2)
for L = 4 : 7
LTm(n)=(X(L));
end
elseif (n == 3)
for L = 8 : 10
LTm(n)=(X(L));
end
elseif (n == 4)
for L = 11 : 13
LTm(n)=(X(L));
end
elseif (n == 5)
for L = 14 : 17
LTm(n)=(X(L));
end
elseif (n == 6)
for L = 18 : 21
LTm(n)=(X(L));
end
elseif (n == 7)
for L = 22 : 25
LTm(n)=(X(L));
end
elseif (n == 8)
for L = 26 : 30
LTm(n)=(X(L));
end
elseif (n == 9)
for L = 31 : 35
LTm(n)=(X(L));
end
elseif (n == 10)
for L = 36 : 40
LTm(n)=(X(L));
end
elseif (n == 11)
for L = 41 : 47
LTm(n)=(X(L));
end
elseif (n == 12)
for L = 48 : 51
LTm(n)=(X(L));
end
elseif (n == 13)
for L = 52 : 55
LTm(n)=(X(L));
end
elseif (n == 14)
for L = 56 : 61
LTm(n)=(X(L));
end
elseif (n == 15)
for L = 62 : 67
LTm(n)=(X(L));
end
elseif (n == 16)
for L = 68 : 74
LTm(n)=(X(L));
end
elseif (n == 17)
for L = 75 : 79
LTm(n)=(X(L));
end
elseif (n == 18)
for L = 80 : 84
LTm(n)=(X(L));
end
elseif (n == 19)
for L = 85 : 91
LTm(n)=(X(L));
end
elseif (n == 20)
for L = 92 : 99
LTm(n)=(X(L));
end
elseif (n == 21)
for L = 100 : 108
LTm(n)=(X(L));
end
end
end
LTn=[];
for n = 1 : 21 % 21 subbands
if (n == 1) % find minimum masking threshold for each subbands
for L = 1 : 3 % subband indexes
LTn(n)=(Y(L));
end
elseif (n == 2)
for L = 4 : 7
LTn(n)=(Y(L));
end
elseif (n == 3)
for L = 8 : 10
LTn(n)=(Y(L));
end
elseif (n == 4)
for L = 11 : 13
LTn(n)=(Y(L));
end
elseif (n == 5)
for L = 14 : 17
LTn(n)=(Y(L));
end
elseif (n == 6)
for L = 18 : 21
LTn(n)=(Y(L));
end
elseif (n == 7)
for L = 22 : 25
LTn(n)=(Y(L));
end
elseif (n == 8)
for L = 26 : 30
LTn(n)=(Y(L));
end
elseif (n == 9)
for L = 31 : 35
LTn(n)=(Y(L));
end
elseif (n == 10)
for L = 36 : 40
LTn(n)=(Y(L));
end
elseif (n == 11)
for L = 41 : 47
LTn(n)=(Y(L));
end
elseif (n == 12)
for L = 48 : 51
LTn(n)=(Y(L));
end
elseif (n == 13)
for L = 52 : 55
LTn(n)=(Y(L));
end
elseif (n == 14)
for L = 56 : 61
LTn(n)=(Y(L));
end
elseif (n == 15)
for L = 62 : 67
LTn(n)=(Y(L));
end
elseif (n == 16)
for L = 68 : 74
LTn(n)=(Y(L));
end
elseif (n == 17)
for L = 75 : 79
LTn(n)=(Y(L));
end
elseif (n == 18)
for L = 80 : 84
LTn(n)=(Y(L));
end
elseif (n == 19)
for L = 85 : 91
LTn(n)=(Y(L));
end
elseif (n == 20)
for L = 92 : 99
LTn(n)=(Y(L));
end
elseif (n == 21)
for L = 100 : 108
LTn(n)=(Y(L));
end
end
end
nlags=[];
for n = 1 : 21 % 21 subbands
if (n == 1) % find minimum masking threshold for each subbands
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 2)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 3)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 4)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 5)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 6)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 7)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 8)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 9)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 10)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 11)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 12)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 13)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 14)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 15)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 16)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 17)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 18)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 19)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 20)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
elseif (n == 21)
Sxy(n)=xcorr(LTm(n),LTn(n),nlags);
end
end
figure (1)
subplot(3,1,1);
plot(Z();
grid;
title('Target sound signal')
xlabel('nSamples')
ylabel('Magnitude')
subplot(3,1,2)
plot(R_input();
grid;
title('Mixtured sound signal of right channel')
xlabel('nSamples')
ylabel('Magnitude')
subplot(3,1,3)
plot(L_input();
grid;
title('Mixtured sound signal of left channel')
xlabel('nSamples')
ylabel('Magnitude')
figure (2)
subplot(2,1,1)
stairs(LTm()
grid;
title('Division of 21 subbands for right channel')
xlabel('Number of subbands')
ylabel('Magnitude (dB)')
subplot(2,1,2)
stairs(LTn()
grid;
title('Division of 21 subbands for left channel')
xlabel('Number of subbands')
ylabel('Magnitude (dB)')
figure (3)
stairs(Sxy()
figure (4)
subplot(2,1,1)
plot(X()
subplot(2,1,2)
plot(Y() |
|
|
|
|
|
|
|
|
|
| |
 本周最热论坛帖子
|
发表于 7-11-2006 07:22 PM
变色卡
千斤顶
1548 
237
