ESPRIT是借助旋转不变技术估计信号参数(estimating signal parameter via rotational invariance techniques)方法的英文缩写,最早由Roy等人[233]于1986年提出。ESPRIT方法现已成为现代信号处理中一种代表性方法,并得到了广泛的应用。
和MUSIC方法存在阵元空间和波束空间两类算法一样,ESPRIT方法也有阵元空间E-SPRIT算法和波束空间ESPRIT算法两种类型,并且还有专门针对复观测数据的酉ESPRIT方法。

算法一、基本ESPRIT算法步骤

步骤1、利用已知的观测数据x(1),…..,x(N)估计自协方差函数.

步骤2、由估计的自协方差函数构造MxM自协方差矩阵Rxx和MxM互协方差Rxy。

步骤3、求Rxx的特征值分解。对于M>p,最小特征值的平均值作为噪声方差的估计。

步骤4、利用计算Cxx=Rxx-I和Cxy=Rxy-Z。

步骤5、求矩阵束{Cxx,Cxy}的广义特征值分解,得到位于单位圆上的p个广义特征值,它们直接给出谐波频率wi的估计。然后,利用估计波
达方向θ1,…,θp。

Roy和Kailath指出[232],最小二乘算子会导致在求解广义特征值问题的某些潜在的数值困难。现在已广泛认识到,奇异值分解(SVD)和总体最小二乘(TLS)的应用可以将一个较大维数(MxM)病态广义特征问题转化为一个较小维数(p×p)的无病态广义特征问题。

算法二、TLS-ESPRIT算法步骤

TLS-ESPRIT算法matlab仿真实验

TLS-ESPRIT函数程序如下:

function estimate = tls_esprit(dd,cr,Le)twpi =2.0*pi;derad = pi / 180.0;radeg = 180.0 / pi;%对接收信号协方差矩阵进行特征值分解[K,KK] = size(cr);[V,D]=eig(cr);EVA = real(diag(D)');[EVA,I] = sort(EVA);EVA=fliplr(EVA);EV=fliplr(V(:,I));%构造E_{xy}和E_xys=E_{xy}HE_{xy}Exy =[EV(1:K-1,1:Le) EV(2:K,1:Le)];E_xys = Exy'*Exy;%对E_xys进行特征值分解[V,D]=eig(E_xys);EVA_xys=real(diag(D)');[EVA_xys,I] =sort(EVA_xys);EVA_xys=fliplr(EVA_xys);EV_xys=fliplr(V(:,I));%将EV_xys分解Gx = EV_xys(1:Le,Le+1:Le*2);Gy=EV_xys(Le+1:Le*2,Le+1:Le*2);%计算Psi=-Gx[Gy]{-1}Psi = -Gx/Gy;%对Psi进行特征值分解[V,D]=eig(Psi);EGS = diag(D).';[EGS,I] = sort(EGS);EGS=fliplr(EGS);EVS=fliplr(V(:,I));%估计DOAephi = atan2(imag(EGS), real(EGS));ange = -asin( ephi / twpi / dd ) * radeg;estimate(1,:)=ange;%功率估计T=inv(EVS);powe = T*diag(EVA(1:Le)-EVA(K))*T';powe = abs(diag(powe).')/K;estimate(2,:)=powe;

主函数程序如下:

clear allclose allderad = pi/180;%角度->弧度radeg = 180/pi;%弧度->角度twpi=2*pi;kelm = 8;%阵元数dd=0.5;%阵元间距d=0:dd:(kelm-1)*dd;iwave = 3; %信源数theta =[10 20 30]; %波达方向snr = 10;%信噪比n=500; %采样数(快拍)A=exp(-1i*twpi*d.'*sin(theta*derad));%方向矢量S=randn(iwave,n ); %信源信号snr0=0:3:100 ; %信噪比for isnr=1:10X0=A*S;%接收信号X=awgn(X0,snr0(isnr),'measured') ; %添加噪声Rxx=X*X'/n;%计算协方差矩阵[EV,D]=eig(Rxx); %特征值分解EVA=diag(D)';[EVA,I]=sort(EVA); %特征值从小到大排序EVA=fliplr(EVA); %左右翻转,从大到小排序EV=fliplr(EV(:,I));%对应特征矢量排序estimates=(tls_esprit(dd,Rxx,iwave)); %调用子程序doaes(isnr,:)=sort(estimates(1,:));enddisp(doaes);

运行结果如下:

参考文献:

阵列信号处理及MATLAB实现;张小飞,陈华伟,仇小锋(编著)

现代信号处理(第三版);张贤达(编著)