fangxiangjiaowucha.m

利用延时来估计单个声源位置

clear all;clc;
d=0.032;%麦克风间距；
c=340;%声速；
L=[1 2 3 4 5]'
n=5;
    R=0.2*L%声源与麦克风阵列的距离；
    x1=R/sqrt(4);y1=0.5*R*sqrt(3);%方向角为30°
    x2=0.12*L;y2=0.16*L;%方向角为36.78°
    x3=0.5*R*sqrt(2);y3=0.5*R*sqrt(2);%方向角为45°
    x4=0.16*L;y4=0.12*L;%方向角为53.13°
    x5=0.5*R*sqrt(3);y5=R/sqrt(4);%方向角为60°
    x6=R/4*(sqrt(6)+sqrt(2));y6=R/4*(sqrt(6)-sqrt(2));%方向角为75°
    
    m1=x1./sqrt(x1.^2+y1.^2);p1=asin(m1).*180/pi;
    m2=x2./sqrt(x2.^2+y2.^2);p2=asin(m2).*180/pi;
    m3=x3./sqrt(x3.^2+y3.^2);p3=asin(m3).*180/pi;
    m4=x4./sqrt(x4.^2+y4.^2);p4=asin(m4).*180/pi;
    m5=x5./sqrt(x5.^2+y5.^2);p5=asin(m5).*180/pi;
    m6=x6./sqrt(x6.^2+y6.^2);p6=asin(m6).*180/pi;
    p=[p1 p2 p3 p4 p5 p6];%理论上所有观测角度
    
for M=1:4  %4个麦克风
    t1(1:n,M)=((R-sqrt(R.^2-2.*(M-1).*d.*R.*m1+(M-1).^2.*(d.^2)))./c).*1000;%以第一个麦克风为麦克风阵列的中心，方向角30°，
    %各个距离下，各个麦克风相对于第一个麦克风的时间延迟
    %n1=(R.^2+(M-1).^2.*(d.^2)-(R-(c*t1)./1000).^2)./(2*(M-1)*d.*R)
    %P1(M)=asin(N1(M))*180/pi
    P1(1:n,M)=asin((c*t1(1:n,M))/((M-1)*d*1000))*180/pi;%以第一个麦克风为麦克风阵列的中心，方向角30°，各个距离下，各个麦克风相对于第一个麦克风估计出的方向角
    
    t2(1:n,M)=((R-sqrt(R.^2-2.*(M-1).*d.*R.*m2+(M-1).^2.*(d.^2)))./c).*1000;%以第一个麦克风为麦克风阵列的中心，方向角36.78°，各个距离下，各个麦克风相对于第一个麦克风的时间延迟
    %N2=(R.^2+(M-1).^2.*(d.^2)-(R-(c.*T2)./1000).^2)./(2.*(M-1).*d.*R)
    %P2(M)=asin(N2(M))*180/pi
    P2(1:n,M)=asin((c*t2(1:n,M))/((M-1)*d*1000))*180/pi;%以第一个麦克风为麦克风阵列的中心，方向角36.78°，各个距离下，各个麦克风相对于第一个麦克风估计出的方向角
    
    t3(1:n,M)=((R-sqrt(R.^2-2.*(M-1).*d.*R.*m3+(M-1).^2.*(d.^2)))./c).*1000;%以第一个麦克风为麦克风阵列的中心，方向角45°，各个距离下，各个麦克风相对于第一个麦克风的时间延迟
    %N3(M)=(R.^2+(M-1).^2.*(d.^2)-(R-(c.*T3(M))./1000).^2)./(2.*(M-1).*d.*R)
    %P3(M)=asin(N3(M))*180/pi
     P3(1:n,M)=asin((c*t3(1:n,M))/((M-1)*d*1000))*180/pi;%以第一个麦克风为麦克风阵列的中心，方向角45°，各个距离下，各个麦克风相对于第一个麦克风估计出的方向角
    
    t4(1:n,M)=((R-sqrt(R.^2-2.*(M-1).*d.*R.*m4+(M-1).^2.*(d.^2)))./c).*1000;%以第一个麦克风为麦克风阵列的中心，方向角53.13°，各个距离下，各个麦克风相对于第一个麦克风的时间延迟
    %N1(M)=(R.^2+(M-1).^2.*(d.^2)-(R-(c.*T4(M))./1000).^2)./(2.*(M-1).*d.*R)
    %P4(M)=asin(N4(M))*180/pi
     P4(1:n,M)=asin((c*t4(1:n,M))/((M-1)*d*1000))*180/pi;%以第一个麦克风为麦克风阵列的中心，方向角53.13°，各个距离下，各个麦克风相对于第一个麦克风估计出的方向角
    
    t5(1:n,M)=((R-sqrt(R.^2-2.*(M-1).*d.*R.*m5+(M-1).^2.*(d.^2)))./c).*1000;%以第一个麦克风为麦克风阵列的中心，方向角60°，各个距离下，各个麦克风相对于第一个麦克风的时间延迟
    %N5(M)=(R.^2+(M-1).^2.*(d.^2)-(R-(c.*T5(M))./1000).^2)./(2.*(M-1).*d.*R)
    %P2(M)=asin(N2(M))*180/pi
     P5(1:n,M)=asin((c*t5(1:n,M))/((M-1)*d*1000))*180/pi;%以第一个麦克风为麦克风阵列的中心，方向角60°，各个距离下，各个麦克风相对于第一个麦克风估计出的方向角
     
     t6(1:n,M)=((R-sqrt(R.^2-2.*(M-1).*d.*R.*m6+(M-1).^2.*(d.^2)))./c).*1000;%以第一个麦克风为麦克风阵列的中心，方向角75°，各个距离下，各个麦克风相对于第一个麦克风的时间延迟
    %N5(M)=(R.^2+(M-1).^2.*(d.^2)-(R-(c.*T5(M))./1000).^2)./(2.*(M-1).*d.*R)
    %P2(M)=asin(N2(M))*180/pi
     P6(1:n,M)=asin((c*t6(1:n,M))/((M-1)*d*1000))*180/pi;%以第一个麦克风为麦克风阵列的中心，方向角75°，各个距离下，各个麦克风相对于第一个麦克风估计出的方向角
end

PH1=sum(P1(:,2:4)')'/3;%方向角为30°，各个距离下，每个麦克风估计出方向角的平均值
PH2=sum(P2(:,2:4)')'/3;%方向角为36.78°，各个距离下，每个麦克风估计出方向角的平均值
PH3=sum(P3(:,2:4)')'/3;%方向角为45°，各个距离下，每个麦克风估计出方向角的平均值
PH4=sum(P4(:,2:4)')'/3;%方向角为53.13°，各个距离下，每个麦克风估计出方向角的平均值
PH5=sum(P5(:,2:4)')'/3;%方向角为60°，各个距离下，每个麦克风估计出方向角的平均值
PH6=sum(P6(:,2:4)')'/3;%方向角为75°，各个距离下，每个麦克风估计出方向角的平均值
PH=[PH1 PH2 PH3 PH4 PH5 PH6];
err1=abs(p-PH);%实际测量方位角与理论方位角的误差(单位：°)'

%改进后
for i=1:6
    %p=[p1 p2 p3 p4 p5 p6]
    x=[x1 x2 x3 x4 x5 x6];
    y=[y1 y2 y3 y4 y5 y6];
    r21(1:n,i)=(sqrt(R./(sqrt((x(1:n,i)-d).^2+y(1:n,i).^2)))).^(1/8);
    r31(1:n,i)=(sqrt(R./(sqrt((x(1:n,i)-2*d).^2+y(1:n,i).^2)))).^(1/8);
    r41(1:n,i)=(sqrt(R./(sqrt((x(1:n,i)-3*d).^2+y(1:n,i).^2)))).^(1/8);
    
end
%改进后计算后得到的结果
P1_GA(1:n,2)=r21(1:n,1).*P1(1:n,2);
P2_GA(1:n,2)=r21(1:n,2).*P2(1:n,2);
P3_GA(1:n,2)=r21(1:n,3).*P3(1:n,2);
P4_GA(1:n,2)=r21(1:n,4).*P4(1:n,2);
P5_GA(1:n,2)=r21(1:n,5).*P5(1:n,2);
P6_GA(1:n,2)=r21(1:n,6).*P6(1:n,2);
PH1_GA=P1_GA(1:n,2);
PH2_GA=P2_GA(1:n,2);
PH3_GA=P3_GA(1:n,2);
PH4_GA=P4_GA(1:n,2);
PH5_GA=P5_GA(1:n,2);
PH6_GA=P6_GA(1:n,2);
PH_GA=[PH1_GA PH2_GA PH3_GA PH4_GA PH5_GA PH6_GA];%PH_GA:改进后计算后得到的结果
err_GA=abs(p-PH_GA);

%改进前计算后得到的结果
PH_BF=[P1(:,2) P2(:,2) P3(:,2) P4(:,2) P5(:,2) P6(:,2)];%PH_BF:改进前只用第一组角度进行计算后得到的结果
err4=abs(p-PH_BF);
figure(1);
plot(p(1,:),err4(1,:),'k*-',p(1,:),err4(2,:),'ks-', p(1,:),err4(3,:),'kv-',p(1,:),err4(4,:),'ko-',p(1,:),err4(5,:),'kd-');
xlabel('理论方位角/度');
ylabel('实际测量方位角与理论方位角的误差/度');
title('改进前近场球面前方向角误差曲线');
legend('R=0.2m','R=0.4m','R=0.6m','R=0.8m','R=1m');

figure(5);
plot(p(1,:),err_GA(1,:),'k*-',p(1,:),err_GA(2,:),'ks-', p(1,:),err_GA(3,:),'kv-',p(1,:),err_GA(4,:),'ko-',p(1,:),err_GA(5,:),'kd-');
xlabel('理论方位角(单位：°)');
ylabel('实际测量方位角与理论方位角的误差(单位：°)');
title('改进后近场球面前方向角误差曲线');
legend('R=0.2m','R=0.4m','R=0.6m','R=0.8m','R=1m');

%改进前后两幅曲线一起显示
figure(6);
plot(p(1,:),err4(1,:),'k*-',p(1,:),err4(2,:),'ks-', p(1,:),err4(3,:),'kv-',p(1,:),err4(4,:),'ko-',p(1,:),err4(5,:),'kd-',p(1,:),err_GA(1,:),'k*--',p(1,:),err_GA(2,:),'ks--', p(1,:),err_GA(3,:),'kv--',p(1,:),err_GA(4,:),'ko--',p(1,:),err_GA(5,:),'kd--');
xlabel('理论方位角/度');ylabel('实际测量方位角与理论方位角的误差/度');
title('近场球面前方向角误差曲线')
legend('R=0.2m','R=0.4m','R=0.6m','R=0.8m','R=1m');

%三个图放在一起
figure(7);
subplot(131);plot(p(1,:),err4(1,:),'r*-',p(1,:),err4(2,:),'bs-', p(1,:),err4(3,:),'gv-',p(1,:),err4(4,:),'mo-',p(1,:),err4(5,:),'kd-');
title('改进前(a)');xlabel('理论方位角(单位：°)');ylabel('实际测量方位角与理论方位角的误差(单位：°)');legend('R=0.2m','R=0.4m','R=0.6m','R=0.8m','R=1m');
subplot(132);plot(p(1,:),err_GA(1,:),'r*-',p(1,:),err_GA(2,:),'bs-', p(1,:),err_GA(3,:),'gv-',p(1,:),err_GA(4,:),'mo-',p(1,:),err_GA(5,:),'kd-');
title('改进后(b)');xlabel('理论方位角(单位：°)');ylabel('实际测量方位角与理论方位角的误差(单位：°)');legend('R=0.2m','R=0.4m','R=0.6m','R=0.8m','R=1m');
subplot(133);plot(p(1,:),err4(1,:),'r*-',p(1,:),err4(2,:),'bs-', p(1,:),err4(3,:),'gv-',p(1,:),err4(4,:),'mo-',p(1,:),err4(5,:),'kd-',p(1,:),err_GA(1,:),'r*-',p(1,:),err_GA(2,:),'bs-', p(1,:),err_GA(3,:),'gv-',p(1,:),err_GA(4,:),'mo-',p(1,:),err_GA(5,:),'kd-');
title('一起(c)');xlabel('理论方位角(单位：°)');ylabel('实际测量方位角与理论方位角的误差(单位：°)');legend('R=0.2m','R=0.4m','R=0.6m','R=0.8m','R=1m');