实验报告三实验目的：实现常规波束形成及基于MUSIC 方法的方位估计。

实验内容：1）若干阵元的接收阵，信号频率为10KHz ，波束主轴12度，仿真给出常规波束形成的波束图。

2）16个阵元的均匀线列阵，信号频率为10KHz ，信号方位为12度，用MUSIC 方法完成目标定向，信噪比-5dB ，0dB ，5dB 。

i) 波束形成时的阵型设计为两种，一种是均匀线列阵，阵元16个；一种是均匀圆阵，阵元数为16个，比较这两种阵型的波束图。

ii ）比较不同信噪比下MUSIC 方法估计的性能（统计100次）。

实验原理：i)常规波束形成：如图所示，基阵的输出),(θt v 。

∑∑=*=*==Mm i i M m i iw t x t x w t v 11)()()()(),(θθθ采用向量符号则有，)()()()(),(H H θθθw x x w t t t v ==式中，x(t)和w(q )分别为观测数据向量和加权系数向量，),(θt v 图 1 波束形成器基本原理图T M 21])()()([)(t x t x t x t =x T M 21])()()([)(θθθθw w w =w基阵输出端的空间功率谱表示为：)()( )()]()([)( )]()()()([ )],(),([ ]),([)(H H H H H *2θθθθθθθθθθRw w w x x w w x x w =====t t E t t E t v t v E t v E P式中，R 为观测数据的协方差矩阵。

ii ）基于MUSIC 方法的方位估计：)()()()(1t n t s a t x i di +=∑=θT M 21])()()([)(t x t x t x t =x)()()()(t n t s A t x +=θ假设:(1 ) 信号源的数目d 是已知的, 且d < M ;(2 ) 各信号的方向矢量是相互独立的, 即)(θA 是一个列满秩矩阵;(3 ) 噪声)(t n 是空间平稳随机过程, 为具有各态历经性的均值为零、方差为σ2n 的高斯过程;(4 ) 噪声各取样间是统计独立的。

在上述假设条件下, 基阵输出的协方差矩阵可表示为：I A AR t x t x E R H s H 2])()([α+==其中, R s 为信号的协方差矩阵；I 为单位矩阵。

对R 进行特征分解, 并以特征值降值排列可得Hmm Md m m Hm m dm m ee e e R ∑∑+==+=11λλ信号子空间与噪声子空间正交。

若噪声子空间记为E N , 即∑+==Md m Hmm N eeE 1则有：0)(=i N a E θ因此，对矩阵R 进行特征分解后, 取噪声特征矢量, 就可得到信号源的到达角。

定义多重信号分类法的空间谱函数为：)()(1θθa E E a P H NN HMUSIC =实验过程及结论：1.常规波束形成，阵列为均匀线列阵：-100-80-60-40-20020*********-70-60-50-40-30-20-10均匀线列阵 12度主轴方位角/度归一化功率/d B2.常规波束形成，阵列为均匀圆阵：-200-150-100-50050100150200-80-70-60-50-40-30-20-100均匀圆阵 12度主轴方位角/度归一化功率/d B结论：由上面两个仿真结果可以看出，采用均匀线列阵和均匀圆阵都可以完成波束形成，都可以获得阵处理增益。

3.MUSIC算法方位估计：-100-80-60-40-20020406080100-60-50-40-30-20-10方位角 / 度归一化功率 / d B结论：从上图得出，信噪比越高，MUSIC方位估计越准确。

附录：MATLAB仿真程序代码1.常规波束形成%%clear all;close all;clc%%%%%%常规波束形成%信号参数f0=10e3;fs=10*f0;T=2e-3;t=1:1/fs:T-1/fs;%均匀线列阵阵列参数c=1500;lamda=c/f0;d=lamda/2;M=16;direct=12;%波束指向%阵列流行theta=-90:90;for i_theta=1:length(theta)fai=2*pi*d*sin(theta(i_theta)*pi/180)/lamda;for i_M=1:MA(i_M,i_theta)=exp(-sqrt(-1)*(i_M-1)*fai);endend%加权向量wfai_direct=2*pi*f0*d*sin(direct*pi/180)/c;for i_M=1:Mw(i_M)=exp(sqrt(-1)*(i_M-1)*fai_direct);endw=w';%无噪声输出P=w'*A;P=20*log10(abs(P)/max(abs(P)));figure(1)plot(theta,P);title('均匀线列阵 12度主轴');xlabel('方位角/度');ylabel('归一化功率/dB'); gridhold onclear all;clc%%%%%%常规波束形成%信号参数f0=10e3;fs=10*f0;T=2e-3;t=1:1/fs:T-1/fs;%元阵阵列参数c=1500;lamda=c/f0;d=lamda/2;M=16;i_M=0:M-1;fai=2*pi/M.*i_M;B=[d*cos(fai);d*sin(fai);d*zeros(1,M)];direct=12;%波束指向%阵列流行theta=-180:180;theta=theta*pi/180;for i_theta=1:length(theta)e=[cos(theta(i_theta));sin(theta(i_theta));0];tao=-e.'*B./c;for i_M=1:MA(i_M,i_theta)=exp(sqrt(-1)*2*pi*f0*tao(i_M));endend%加权向量wfai_direct=2*pi*f0*d*sin(direct*pi*180)/c;for i_M=1:Mw(i_M)=exp(sprt(-1)*2*pi*f0*d*cos(fai(i_M)-direct/180*pi)/c); endw=w';%无噪声输出P=w'*A;P=20*log10(abs(P)/max(abs(P)));figure(2)plot(theta./pi*180,P);title('均匀圆阵 12度主轴');xlabel('方位角/度');ylabel('归一化功率/dB'); grid2.MUSIC方位估计%%clc; clear; close all;M=16;c=1500;f0=10e3;lamda=c/f0;d=lamda/2;fs=10*f0;T=2;t=[0:1/fs:T];N=fs*T;signalnum=1;%%%信噪比SNR=-5;theta1=12;num=100;theta1=theta1*pi/180;tau0=d*sin(theta1)/c;for i_M=1:MS(i_M,:)=exp(sqrt(-1)*2*pi*f0*t)*exp(-sqrt(-1)*2*pi*f0*(i_M-1)*tau0); endP=zeros(1,181); %length(theta)for i=1:numX=awgn(S,SNR);R=X*X'/N; %协方差矩阵[V,I]=eig(R);[Order Index]=sort(diag(I),'descend'); %特征值降值排列VS=[]; %向量子空间P_music=[];for i_M=1:MVS=[VS,V(:,Index(i_M))]; %排序后的特征向量endE=VS(:,signalnum+1:M); %后M-signalnum向量构成向量子空间En=E*E'; %噪声子空间theta=-90:1:90;fai=2*pi*f0*d*sin(theta*pi/180)/c;for i_theta=1:length(theta)for i_M=1:Mw(i_theta,i_M)=exp(sqrt(-1)*(i_M-1)*fai(i_theta));%方向矢量endP_music(i,i_theta)=1./(w(i_theta,:)*En*En'*w(i_theta,:)');endendP=abs(P)/num;P=P+P_music(i,:);figure(3)plot(theta,10*log10(P/max(P)),'r');xlabel('方位角 / 度');ylabel('归一化功率 / dB');gridtitle('不同信噪比下MUSIC方法估计的性能（统计100次）')hold on%%%信噪比SNR=0;theta1=12;num=100;theta1=theta1*pi/180;tau0=d*sin(theta1)/c;for i_M=1:MS(i_M,:)=exp(sqrt(-1)*2*pi*f0*t)*exp(-sqrt(-1)*2*pi*f0*(i_M-1)*tau0); endP=zeros(1,181); %length(theta)for i=1:numX=awgn(S,SNR);R=X*X'/N; %协方差矩阵[V,I]=eig(R);[Order Index]=sort(diag(I),'descend'); %特征值降值排列VS=[]; %向量子空间P_music=[];for i_M=1:MVS=[VS,V(:,Index(i_M))]; %排序后的特征向量endE=VS(:,signalnum+1:M); %后M-signalnum向量构成向量子空间En=E*E'; %噪声子空间theta=-90:1:90;fai=2*pi*f0*d*sin(theta*pi/180)/c;for i_theta=1:length(theta)for i_M=1:Mw(i_theta,i_M)=exp(sqrt(-1)*(i_M-1)*fai(i_theta));%方向矢量endP_music(i,i_theta)=1./(w(i_theta,:)*En*En'*w(i_theta,:)');endendP=abs(P)/num;P=P+P_music(i,:);figure(3)plot(theta,10*log10(P/max(P)),'*');xlabel('方位角 / 度');ylabel('归一化功率 / dB');gridtitle('不同信噪比下MUSIC方法估计的性能（统计100次）')hold on%%%信噪比SNR=5;theta1=12;num=100;theta1=theta1*pi/180;tau0=d*sin(theta1)/c;for i_M=1:MS(i_M,:)=exp(sqrt(-1)*2*pi*f0*t)*exp(-sqrt(-1)*2*pi*f0*(i_M-1)*tau0); endP=zeros(1,181); %length(theta)for i=1:numX=awgn(S,SNR);R=X*X'/N; %协方差矩阵[V,I]=eig(R);[Order Index]=sort(diag(I),'descend'); %特征值降值排列VS=[]; %向量子空间P_music=[];for i_M=1:MVS=[VS,V(:,Index(i_M))]; %排序后的特征向量endE=VS(:,signalnum+1:M); %后M-signalnum向量构成向量子空间En=E*E'; %噪声子空间theta=-90:1:90;fai=2*pi*f0*d*sin(theta*pi/180)/c;for i_theta=1:length(theta)for i_M=1:Mw(i_theta,i_M)=exp(sqrt(-1)*(i_M-1)*fai(i_theta));%方向矢量endP_music(i,i_theta)=1./(w(i_theta,:)*En*En'*w(i_theta,:)');endendP=abs(P)/num;P=P+P_music(i,:);figure(3)plot(theta,10*log10(P/max(P)),'.');xlabel('方位角 / 度');ylabel('归一化功率 / dB');gridtitle('不同信噪比下MUSIC方法估计的性能（统计100次）') hold onlegend('-5dB','0dB','5dB');。