bayuanshiliangzhen.m

矢量阵指向性程序

clear;
f0=2000;   %信号频率
fs=80*f0;   %采样频率 
Ns=2048;   %采样点数
c=1500;    %声速  
N1=8;      %阵元数
r1=1*c/f0;      %大子阵半径
r2=0.5*r1     %小子阵半径
fai=0;     %信号入射方向 
w=2*pi*f0;
ph1=-180;                   
ph2=180;                   
phn=1;                  
ax=(ph2-ph1)/phn+1;

%************产生入射信号****************%

for n=1:N1;
       a=2;
    if mod(n,a)~=0
       for k=1:Ns;
        p1(n,k)=exp(j*w*(k-1)/fs-j*w*r1*cos(fai*pi/180-(n-1)*2*pi/N1)/c);
        vx1(n,k)=p1(n,k)*cos(fai*pi/180);
        vy1(n,k)=p1(n,k)*sin(fai*pi/180);  
           x(k)=exp(j*w*(k-1)/fs);
       end
    else
       for k=1:Ns;
        p2(n,k)=exp(j*w*(k-1)/fs-j*w*r2*cos(fai*pi/180-(n-1)*2*pi/N1)/c);
        vx2(n,k)=p2(n,k)*cos(fai*pi/180);
        vy2(n,k)=p2(n,k)*sin(fai*pi/180);  
       end
    end
end


%**************对接收信号进行相移*******%
for i=1:ax;     
    num=ph1+(i-1)*phn;
    fain=num*pi/180;
     Xp1=0;
     Xp2=0;
    Xv1=0;
    Xv2=0;
      XP1=x(1);
    XV1=x(1)*(cos(fai*pi/180)+sin(fai*pi/180));
    for n=1:N1;
        a=2;
      if mod(n,a)~=0
        Xp01(n,1,i)=p1(n,1)*exp(j*w*r1*cos(fain-(n-1)*2*pi/N1)/c); 
        Xp1=Xp1+Xp01(n,1,i);
        Xv01(n,1,i)=Xp01(n,1,i)*(cos(fai*pi/180)*cos(fain)+sin(fai*pi/180)*sin(fain));
        Xv1=Xv1+Xv01(n,1,i);
       
        XP(n,1,i)=p1(n,1)*exp(j*w*r1*cos(fain-(n-1)*2*pi/N1)/c); 
        XP1=XP1+XP(n,1,i);
        XV(n,1,i)=XP(n,1,i)*(cos(fai*pi/180)*cos(fain)+sin(fai*pi/180)*sin(fain));
        XV1=XV1+XV(n,1,i);
      else
        Xp02(n,1,i)=p2(n,1)*exp(j*w*r2*cos(fain-(n-1)*2*pi/N1)/c); 
        Xp2=Xp2+Xp02(n,1,i);
        Xv02(n,1,i)=Xp02(n,1,i)*(cos(fai*pi/180)*cos(fain)+sin(fai*pi/180)*sin(fain));
        Xv2=Xv2+Xv02(n,1,i);
      end
    end
    Xp=Xp1+Xp2;
    Xv=Xv1+Xv2;
    Rxp(i)=abs(Xp);
    Rp(i)=abs(Xp1);
      RP(i)=abs(XP1*XV1)*XV1;
    Rp1(i)=abs(Xp*Xv);
    Rp2(i)=abs((Xp+Xv)*Xv);
    Rp3(i)=abs((Xp+Xv).^2);
%   Rp(i)=abs((Xp1+2*Xv1)*Xv1);
%   Rp(i)=abs((Xp1+2*Xv1).^2);
%   Rp(i)=abs((Xp1+2*Xv1)*(Xp1+Xv1));
end 
figure
subplot(321)
i=1:ax;
maxRxp=max(Rxp);
Rxpp=20*log10(Rxp(i)/maxRxp);
plot(ph1+(i-1)*phn,Rxpp);
grid on
hold on
axis([-180 180 -100 0])
title('八元声压阵')
xlabel('角度')
ylabel('dB')
subplot(322)
i=1:ax;
maxRp=max(Rp);
Rpp=20*log10(Rp(i)/maxRp);
plot(ph1+(i-1)*phn,Rpp);
grid on
hold on
axis([-180 180 -100 0])
title('四圆矢量阵(pv)')
xlabel('角度')
ylabel('dB')
subplot(323)
i=1:ax;
maxRp1=max(Rp1);
Rpp1=20*log10(Rp1(i)/maxRp1);
plot(ph1+(i-1)*phn,Rpp1);
grid on
hold on
axis([-180 180 -100 0])
title('八元矢量阵(pv)')
xlabel('角度')
ylabel('dB')
subplot(324)
i=1:ax;
maxRp2=max(Rp2);
Rpp2=20*log10(Rp2(i)/maxRp2);
plot(ph1+(i-1)*phn,Rpp2);
grid on
hold on
axis([-180 180 -100 0])
title('八元矢量阵((p+v)v)')
xlabel('角度')
ylabel('dB')
subplot(325)
i=1:ax;
maxRp3=max(Rp3);
Rpp3=20*log10(Rp3(i)/maxRp3);
plot(ph1+(i-1)*phn,Rpp3);
grid on
hold on
axis([-180 180 -100 0])
title('八元矢量阵((p+v)^2)')
xlabel('角度')
ylabel('dB')
subplot(326)
i=1:ax;
maxRP=max(RP);
RPP=20*log10(RP(i)/maxRP);
plot(ph1+(i-1)*phn,RPP);
grid on
hold on
axis([-180 180 -100 0])
title('五元矢量阵')
xlabel('角度')
ylabel('dB')