canshu.m

矢量水听器圆阵的常规波束形成仿真程序

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

%************产生入射信号****************%
for n=1:N1;
       a=2;
    if mod(n,a)~=0
       for k=1:Ns;
        p(n,k)=sin(w*(k-1)/fs-w*r1*cos(fai*pi/180-(n-1)*2*pi/N1)/c);
        
       end
       p1(n,:)=hilbert(p(n,:));
    else
       for k=1:Ns;
        p(n,k)=sin(w*(k-1)/fs-w*r2*cos(fai*pi/180-(n-1)*2*pi/N1)/c);
       end
        p1(n,:)=hilbert(p(n,:));
    end
end
% figure
% subplot(211)
% k=1:Ns;
% plot(k,p(1,:))
% subplot(212)
% p1=hilbert(p);
% plot(k,p1(2,:))
%**************对接收信号进行相移*******%
     Xp1=zeros(ax,Ns);
     Xvx1=zeros(ax,Ns);
     Xvy1=zeros(ax,Ns);
for i=1:ax;     
    num=ph1+(i-1)*phn;
    fain=num*pi/180;  
    for n=1:N1;
        a=2;
      if mod(n,a)~=0 
        Xp01(n,:,i)=p1(n,:)*exp(j*w*r1*cos(fain-(n-1)*2*pi/N1)/c); 
        Xp1(i,:)=Xp1(i,:)+Xp01(n,:,i);
        Xvx01(n,:,i)=Xp01(n,:,i)*cos(fai*pi/180);
        Xvy01(n,:,i)=Xp01(n,:,i)*sin(fai*pi/180);
         Xvx1(i,:)=Xvx1(i,:)+Xvx01(n,:,i);
          Xvy1(i,:)=Xvy1(i,:)+Xvy01(n,:,i);
          Xvc1(i,:)=Xvx1(i,:)*cos(fain)+Xvy1(i,:)*sin(fain);
      else
        Xp01(n,:,i)=p1(n,:)*exp(j*w*r2*cos(fain-(n-1)*2*pi/N1)/c); 
        Xp1(i,:)=Xp1(i,:)+Xp01(n,:,i);
       Xvx01(n,:,i)=Xp01(n,:,i)*cos(fai*pi/180-pi/4)*cos(pi/4)-Xp01(n,:,i)*sin(fai*pi/180-pi/4)*sin(pi/4);
        Xvy01(n,:,i)=Xp01(n,:,i)*cos(fai*pi/180-pi/4)*sin(pi/4)+Xp01(n,:,i)*sin(fai*pi/180-pi/4)*cos(pi/4);
        Xvx1(i,:)=Xvx1(i,:)+Xvx01(n,:,i);
        Xvy1(i,:)=Xvy1(i,:)+Xvy01(n,:,i);
         Xvc1(i,:)=Xvx1(i,:)*cos(fain)+Xvy1(i,:)*sin(fain);
      end
    end
 Xp=(sum(abs(Xp1').^2))/Ns;
Xvc=(sum(abs(Xvc1').^2))/Ns;
Rp(i)=(Xp(i)+Xvc(i))*Xvc(i);
end 
i=1:ax;
  num=ph1+(i-1)*phn;
  fain=num*pi/180;
maxRp=max(Rp);
Rpp=10*log10(Rp(i)./maxRp);
figure
plot(ph1+(i-1)*phn,Rpp);
legend('交叉矢量阵(p+vc)*vc',4);
axis([-180 180 -100 0])
grid on
hold on
xlabel('角度')
ylabel('dB')