完整word版esprit算法研究.docx

资源描述

完整word版esprit算法研究.docx

《完整word版esprit算法研究.docx》由会员分享，可在线阅读，更多相关《完整word版esprit算法研究.docx（21页珍藏版）》请在冰豆网上搜索。

完整word版esprit算法研究.docx

完整word版esprit算法研究

课程设计报告

实验名称：

ESPRIT算法研究

实验日期:

名:

口

号:

哈尔滨工业大学（威海）

实现空间谱估计算法，并考察算法性能。

方案设计

1）

2）

3）

4）

由均匀线阵形式，确定阵列的导向矢量；

由阵列导向矢量，对接收信号进行建模仿真；

由ESPRIT算法实现信号DOA古计；考察算法性能与信噪比，采样率，观测时间等参数的关系。

3.1空间谱估计数学模型

空间谱估计就是利用空间阵列实现空间信号的参数估计的一项专门技术。

整个空间谱估计系

统应该由三部分组成：

空间信号入射、空间阵列接收及参数估计。

相应地可分为三个空间，即目标空间、观察空间及估计空间，也就是说空间谱估计系统由这三个空间组成，其框图见

图1。

对于上述的系统结构，作以下几点说明。

（1）目标空间是一个由信号源的参数与复杂环境参数张成的空间。

对于空间谱估计系统，

就是利用特定的一些方法从这个复杂的目标空间中估计出信号的未知参数。

（2）观察空间是利用空间按一定方式排列的阵元，来接收目标空间的辐射信号。

由于环

境的复杂性，所以接收数据中包括信号特征（方位、距离、极化等）和空间环境特征（噪声、

杂波、干扰等）。

另外由于空间阵元的影响，接收数据中同样也含有空间阵列的某些特征（互

耦、通道不一致、频带不一致等）。

这里的观察空间是一个多维空间，即系统的接收数据是由多个通道组成，而传统的时域处理方法通常只有一个通道。

特别需要指出的是：

通道与阵

元并不是一一对应，通道是由空间的一个、几个或所有阵元合成的（可用加权或不加权），当

然空间某个特定的阵元可包含在不同的通道内。

（3）估计空间是利用空间谱估计技术（包括阵列信号处理中的一些技术，如阵列校正、空域滤波等技术）从复杂的观察数据中提取信号的特征参数。

从系统框图中可以清晰的看出，估计空间相当于是对目标空间的一个重构过程，这个重构的精度由众多因素决定，如环境的复杂性、空间阵元间的互耦、通道不一致、频带不一

致等。

3.2阵列信号处理

首先，考虑N个远场的窄带信号入射到空间某阵列上，阵列天线由M个阵元组成，这里

假设阵元数等于通道数，即各阵元接收到信号后经过各自的传输信道送到处理器，也就是说

处理器接收来自M个通道的数据。

si（t）=ui（t）ej3細）

根据式（3.2-1）和式（3.2-2），显然有下式成立:

则可以得到第L个阵元接收信号为

「X1（t）■

「e^

ejt0i12…

e血t5n■

「S1（t厂

「n1（t）■

X2（t）

■

』辭1

ej…

ejCCb^N

S2（t）

■

n2（t）

■

LXm（t）-

[

0血爲1

ejCi為2---

e曲也N

■

LSN（t）.

[

『M（t）.

（3.2-6）

假设阵列中各阵元是各向同性的且不存在通道不一致、互耦等因素的中的增益gli可以省略（即归一化

将式（3.2-6）写成矢量形式如下:

S（t）为空间信号的Nx1维矢量，A为空间阵列的MxN维流型矩阵（导向矢量阵），且

A=01（时0）0^（时0）…Fn（⑷0）】（3.2-8）

其中导向矢量

[expjoiMi）.

式中©0=2时=2;ic/几，c为光速，几为波长。

由上述的知识可知，一旦知道阵元间的延迟表达式T,就很容易得出待定空间阵列的

导向矢量或阵列流型。

下面推导一下空间阵元间的延迟表达式。

假设空间任意两个阵元，其

中一个为参考阵元（位于原点），另一个阵元的坐标为（X,y,z），两阵元的几何关系见图，图中“X”表示阵元。

图2空间任意两阵元的几何关系

由几何关系可以推导出两阵元的波程差为

（3.2-10）

T=—（xcos日COS®+ysinScos®+zsinW）c

这里的波程差其实就是位于X轴上两阵元间的延迟、位于y轴上两阵元间的延迟和位于z

轴上两阵元间的延迟之和。

根据式（3.2-10）的结论，下面给出实际环境中常用的几种阵列及阵元间的相互延迟表达式。

（1）平面阵设阵元的位置为（Xk，yk）（k=1,2，…,M），以原点为参考点，另假设信号入射参数为（q,®i）（i=1,2，…，N），分别表示方位角与俯仰角，其中方位角表示与X轴的夹

角。

⑵线阵设阵元的位置为Xk（k=1,2，…,M），以原点为参考点，另假设信号入射参

图3-1均匀线阵的数学模型示意图

XM（t）

数为Ti（i=1,2,…，N），表示方位角，其中方位角表示与则有

y轴的夹角（即与线阵法线的夹角）

Tki=—（XkSin®）

（3）均匀圆阵设以均匀圆阵的圆心为参考点，则有

（3.2-11）

5=1fcos^M"1）]cos%

clIM丿丿

（3.2-12）

其中方位角表示与x轴的夹角，r为圆半径。

3.3旋转不变子空间算法原理

3.3.1信号模型

算法介绍前，首先对信号进行建模。

为了推导分析的方便，将波达方向的

数学模型做如下理想状态的假设:

1）阵列形式为线性均匀阵，阵元间距不大于信号波长的二分之一。

2）存生两个完全相同的子阵，且两个子阵的间距△是己知的。

3）噪声序列为一零均值高斯过程，各阵元间噪声相互独立，噪声与信号也

相互独立。

4）空间信号为零均值平稳随机过程，通常为窄带远场信号。

5）信号源数小于子阵阵列元数，信号取样数大于子阵阵列元数，以确保子

阵阵列流型的各列线性独立。

6）组成阵列的各传感器为各向同性阵元，且无互耦以及通道不一致的干扰。

下图给出了均匀线阵的数学模型示意图

332算法原理

对于均匀线阵，相邻子阵间存在一个固定间距，这个固定间距反映出各相邻

子阵间的一个固定关系，即子阵间的旋转不变性，而ESPRIT算法正是利用了这

个子阵间的旋转不变性实现阵列的DOA估计。

ESPRIT算法最基本的假设是存在两个完全相同的子阵,

且两个子阵的间距

△是已知的。

由于两个子阵的结构完全相同，且子阵的阵元数为m,对于同一个

信号而言，两个子阵的输出只有一个相位差q，i=1,2，…N。

F面假设第一个子阵的接收数据为X1，第二个子阵的接收数据为X2，根据

前面所述的阵列模型可知

X1=[ae）川a（0N）]S+N“AS+N1

（3.1）

X2=[a（3）ej打||a（TN）ej*]S+N2=Aes+N2

（3.2）

式中，子阵1的阵列流型A1=A，子阵2的阵列流型人=A①，且式中

①=diag[ejL.e血]

（3.3）

从上面的数学模型可知，需要求解的是信号的方向，而信号的方向信息包含

在A和①中，由于①是一个对角阵，所以下面只考虑这个矩阵，即

（3.4）

由上可知。

只要得到两个子阵间的旋转不变关系①，就可以方便地得到关于

信号到达角的信息。

下面的任务就是从式（

3.1）和式（3.2）中得到两个子阵间

的关系。

先将两个子阵的模型进行合并，即

显然上式中得到的特征值有如下关打＞-＞ZN＞耳申=•••=◎，Us为大特征

值对应的特征矢量张成的信号子空间，UN为小特征值对应矢里张成的噪声子空

间。

对于实际的快拍数据，式（3.7）应修正如下:

由前面的知识可知，上述的特征分解中大特征矢量张成的信号子空间与阵列流型张成的信号子空间是相等的。

即

此时，存在一个惟一的非奇异矩阵T,使得

很显然，由子阵1的大特征矢量张成的子空间Usi、由子阵2的大特征矢量张成

的子空间Us2与阵列流型A张成的子空间三者相等，即

另外，由两个子阵列在阵列流型上的关系可知

（3.13）

再利用式（3.11）可知两个子阵列的信号子空间的关系如下

两个子阵的阵列接收数据的信号子空间的旋转不变性。

如果阵列流型A是满秩矩阵，则由式（3.14）可以得到

（3.15）

所以上式中空的特征值组成的对角阵一定等于①，而矩阵T的各列就是矩阵屮特

征矢量。

所以一旦得到上述的旋转不变关系矩阵甲，就可以直接利用式（3.4）

得到信号的入射角度。

3.4标准的旋转不变子空间算法

有上节的知识可知，ESPRIT算法的基本原理就是利用式（3.14）的旋转不变性，常规的旋转不变子空间算法就是利用上述的基本原理求解信号的入射角度

信息。

下面就分析解这个等式的两种最经典、应用最广泛方法：

最小二乘（LS）

法和总体最小二乘（TLS）法。

3.4.1最小二乘法

由最小二乘的数学知识，我们知道式（3.14）的最小二乘解的方法等价于

因此最小二乘法的基本思想就是使校正项iUS2尽可能小，而同时保证满足

约束条件。

为了得到LS解，将式（3.14）代入式（3.16）即得

对上式进行展开可得

f（普）=Us1普-Us2

（1）当Us1满秩时，也就是子阵1的信号子空间的维数等于信号源数时，则上

式的解是唯一的，可得上式的最小二乘解

⑵当US1不满秩，即rank（US1）吒N时，也就是信号源间存在相干或相差时,

则甲存在很多解，但我们却无法区别对应于方程的各个不同的解，可以称这些

解是不可辨识的，解的不可辨识性是我们需要解相干的原因所在。

F面给出LS-ESPRIT算法的求解步骤:

1.由两个子阵的接收数据X1,X2，分别得到两个子阵的数据协方差矩阵;

2•对矩阵对｛R,Rn｝进行特征分解，从而得到两个数据矩阵的信号子空间Usi和

US2；

3.按式（3.20）得到矩阵甲LS,然后对其进行特征分解.得到N个特征值，就可得到对应的N个信号的到达角。

当考虑嗓声影响时，上述基于最小二乘算法的估计都是有偏的，这就是为什么需要考虑总体最小二乘ESPRIT算法的原因。

342总体最小二乘法

我们知道，普通最小二乘的基本思想是用一个范数平方为最小的扰动AUS2

去于扰信号子空间Us2，目的是校正Us2中存在的嗓声。

显然这就存在一个问题:

如果同时扰动Usi和Us2，并使扰动范数的平方保持最小，是否可以同时校正

Usi和Us2中存在的嗓声？

答案是肯定的，这就是总休最小二乘（TLS）的思想。

它考虑的是如下矩阵方程的解：

（3.27）

（Us1+也Us1）曾=Us2+也Us2

显然上式可以改写成

所以TLS的解等价于

（3.29）

fmin|8u||

［约束条件：

（AU+U）z=0

定义如下一个矩阵Usi2=［Usi|Us2］，再结合上述分析过程。

我们发现其实

就是寻找一个2Nxm的酉矩阵F,便得矩阵F与Usi2正交，也就说明了由F张

成的空间与Us1或Us2列矢量张成的空间正交。

所以矩阵F可从UH12Us12的特征分

解中得到。

因为

（3.30）

曲丄12=E^Eh

式中的2是由特征值构成的对角矩阵，e是与其相应的特征矢量构成的矩阵。

即

令EN=f12]是由对应特征值为0的特征矢量构成的矩阵.它属子噪声子空

LE220>N

如果令空=-FiF2-，则

（3.34）

屮就可得到

（3.35）

上式说明空的特征值即①是对角线元素。

这说明通过构造一个矩阵

有关信号角度的信息.而这个矩阵的构造可通过式（3.30）得到，即

甲TLS=-E2E；2

F面直接给出TLS-ESPRIT算法的求解步骤:

1.由两个子阵的接收数据X1，X2,由式（3.8）得到数据协方差矩阵R;

2.通过矩阵对于｛r,Rn｝的广义特征分解，得到维数为2M沢N的信号子空间

Us；

3.由Us构造矩阵Us12，并按式防（3.30）进行特征分解得到矩阵E，然后再

按式（3.31）将矩阵分为四个小的矩阵;

4.按式（3.35）得到矩阵甲tls，然后对其进行特征分解,得到N个特征值,

就可得到对应的N个信号的到达角。

通过分析，我们可以得到标准ESPRIT算法的计算过程如下:

（1）通过特征值或奇异值分解（EVD或SVD）分别估计两个存在旋转不变关系的子阵的信号子空；

（2）用上述的LS、TLS等方法求解式（3.14）所示的不变等式;

（3）计算甲k=T弋的特征值，其中①如式（3.3）所示。

然后利用式（3.4）求

解人射信号的角度信息。

就ESPRIT算法而言，TLS算法与LS算法性能基本一致，只是在低信噪比情况下TLS算法性能略好。

四、仿真结果

180

主要分析各个参数对估计误差的影响，误差函数定义如式

（1）：

4.1信噪比SNR对估计误差的影响分析首先对信噪比SNR离散化取值，

方位角相差比较小时，估计的误差也会相应的增大。

另外，若两信号为相干信号,

则此方法将不能对其进行正确的估计。

□>5

*4*.十

0510

SNR

SNR对怙汁误绘的辫

-1'5'-10-5

4.2阵元数L对估计误差的影响分析

与4.1节类似，首先对阵元数L

离散化取值，然后求得不同阵元数下的误

差，从而绘制出误差随阵元数改变的函数曲线如图3所示，图3中阵元数从

K+1取到K+25，间隔为1,运行次数为100次，其余条件如题中所述。

由于

阵元数L需大于信号个数K才能正确估计，故取值中含有信号个数K。

由图3

可知，随着阵元数的增加，估计误差会越来越小，即估计精度会越来越高，但当

阵元数大到一定程度后，对估计精度的影响则会慢慢的减小。

2.5

0.5

*爭哄血*』

510152025

阵元数L

图3阵元對（L对估计误蓋的軫向

3'0|

4.3采样点数N对估计误差的影响分析

与4.1节类似，首先对采样点数N离散化取值，然后求得不同采样点数下

的误差,

从而绘制出误差随采样点数改变的函数曲线如图4所示，图4中采样

点数从

10取到200,间隔为5,运行次数为100次，其余条件如题中所述。

由图4

可知，随着采样点数的增加，估计误差会越来越小，即估计精度会越来

越咼。

圈斗采样点■散N对佶吴差的夥甸

然能得出估计结果，但估计的精度会大受影响。

因此，为了分析两信号之间的不同间隔会对估计精度造成多大的影响，绘制不同GAP下的估计误差曲线如图5所示。

处理方法与4.1节类似，图5中GAP（单位为度）从0.1取到5,间

隔为0.1,独立运行次数为100次，其余条件如题中所述。

由图5可知，GAP

越大估计越准确，但当GAP大到一定程度后则估计精度趋于稳定。

►

522.53354

GAPO

图5两倍号之间的甬度连2AP〉对怙计俣羞的議响

4.5单信号DOA不同分布对估计误差的影响分析

信号波达方向（DOA）的取值区间为-90度到90度，若只考虑只有一个信号的情况，则当信号的DOA不同时，估计误差也会不一样。

因此，为了分析不同的DOA会对估计精度造

成多大的影响，绘制不同DOA下的估计误差曲线如图6所示。

处理方法与4.1节类似,

图6中GAP从-80度取到80度，间隔为5度，独立运行次数为100次，其余条件如题中所述。

由图6可知，DOA越靠近0度估计越准确，越靠近正负90度估计误差越大。

且仿真结果表明，当DOA在正负90附近时，估计误差太大，因此，为了不影响估计结果显示效果，故在图中未绘制正负90度附近的估计误差。

Q5M-

O.B

■

O4£

0.4

0.3*

0.3

V>■■

0：

0.2

0.15k

0.1-

uwL

则■SO^CJO02040soao

DOAC）

團$单信号DOA不同分布对怙计误差的嶷响

2.6减与不减噪声方差（Rn）对估计误差的影响分析

由于有噪声的影响，因此在估计信号自相关矩阵R时，若将无信号时的自相关矩阵Rn减

去，即相当与减去估计出噪声方差，则估计的精度会有所提高。

结合信噪比SNR对估计误

差的影响，绘制减与不减噪声方差两种情况下估计误差随SNR的变化曲线如图7所示，图

7中SNR从-15dB到5dB,间隔为1dB，独立运行次数为100次。

仿真结果表明，若减Rn，主要是在低信噪比时对估计精度的改善较大，当信噪比较大时二者几乎一样。

五、程序清单

%%%%%文件名为drawTLSesprit.m%%%%%

%%%分析基于总体最小二乘的ESPRIT算法（TLS-ESPRIT）的DOA估计的性能

%%%%%

%信号个数：

%阵元数:

%%%%%阵%元数L对估计误差的影响分析%%%%%

%Ln=p+1:

p+25;%阵元数L需大于信号个数p才能正确估计

%forn=1:

length（Ln）

%L=Ln（n）;

%fork=1:

%[estimated,error]=TLSesprit（p,L,K,SNR,DOA）;

%errorm（k）=error;%将每次的估计误差存入变量errorm中，便于求均值

%end

%errorn（n）=sum（errorm）/M;%求多次运行后的估计误差的均值

%绘制曲线，并适当标

阵元数L对估计误差的影

%end

%figure

（2）;plot（Ln,errorn*180/pi,'r:

*','LineWidth',2）;

%xlabel（'阵元数L'）;ylabel（'估计误差（°）'）;title（'

%%结论：

阵元数L越大估计越准确，但当L大到一定程度后则估计精度趋于稳定

%%结论：

快拍数K越大估计越准确，但当K大到一定程度后则估计精度趋于稳定

%%%%信燥比SNR对估计误差的影响分析%%%%%

%SNRn=-15:

15;%对信噪比SNR离散化取值

%forn=1:

length（SNRn）

%SNR=SNRn（n）;

%fork=1:

%[estimated,error]=TLSesprit（p丄,K,SNR,DOA）;

%errorm（k）=error;%将每次的估计误差存入变量errorm中，便于求

均值

%end

%errorn（n）=sum（errorm）/M;%求多次运行后的估计误差的均值

%end

%figure（4）;plot（SNRn,errorn*180/pi,'r:

*','LineWidth',2）;%

注（误差：

角度）

绘制曲线并适当标

%xlabel（'SNR'）;ylabel（'估计误差（°）'）;title（'SNR

%%吉论：

信噪比SNR越大估计越准确，但当信噪比

趋于稳定

对估计误差的影响'）;

SNR大到一定程度后则估计精度

%%%两信号之间的角度差（GAP的大小对估计误差的影响分析%%%%%

GAPn=0.1:

0.1:

5;%对两信号之间的角度差（GAP离散化取值

forn=1:

length（GAPn）

GAP=GAPn（n）;%每次循环只取其中一个值

DOA=[pi*（0/180）pi*（GAP/180）];fork=1:

M%M为独立重复运行次数

[estimated,error]=TLSesprit（p丄,K,SNR,DOA）;

errorm（k）=error;%将每次的估计误差存入变量errorm中，便于求均值

end

errorn（n）=sum（errorm）/M;%求多次运行后的估计误差的均值

end

figure（5）;plot（GAPn,errorn*180/pi,'r:

*','LineWidth',2）;%

（误差：

角度）xlabel（'GAP（°）'）;ylabel（'估计误差（°）'）;

%title（'两信号之间的角度差（GAP对估计误差的影响'）;

%吉论：

GAP越大估计越准确，但当GAP大到一定程度后则估计精度趋于稳定

绘制曲线并适当标注

%%%%单个信号时，信号波达方向分布不同时对估计误差的影响分析对信号波达方向离散化取值（80度到

%%%%%

%DOAn=pi*（-80/180）:

（5/180）:

pi*（80/180）;%

90度时误差太大，因此未取）

%forn=1:

length（DOAn）

%DOA=DOAn（n）;p=1;%每次循环只取其中一个值

为1

，信号个数P设为

%fork=1:

M%M为独立重复运行次数

%[estimated,error]=TLSesprit1（p,L,K,SNR,DOA）;

情况）

%调用TLSespritl（—个信号的

%end

%errorn（n）=sum（errorm）/M;%求多次运行后的估计误差的均值

%end

%forn=1:

length（SNRn）

%SNR=SNRn（n）;

%end

%errorn（n）=sum（errorm）/M;errornRn（n）=sum（errormRn）/M;

%%%本文件名为TLSesprit.m%%%%%

%%%基于总体最小二乘的ESPRIT算法（TLS-ESPRIT）的DOA估计函数％%%%i%iction[estimated,error]=TLSesprit（p,L,K,SNR,DOA）

%调用格式:

[estimated,error]=TLSesprit（p丄,K,SNR,DOA）;

%估计结果（弧度，矢量：

P行1列）：

estimated

%估计误差（弧度，标量：

均方误差）：

error

%信号个数：

%阵元数：

%快拍数：

Rn=Rn+Noise（:

i）*Noise（:

i）：

R=R+Y（:

i）*Y（:

i）';

行P列，

U12=U（p+1:

2*p,1:

p）;%p行p列

TLS=-U11*inv（U12）;%p行p列，U11和U12构成U的噪声子空间

d=eig（TLS）;%特征分解，求TLS的特征值

estimated=（sort（asin（angle（d）/pi）））';%

输出估计（已从小到大排序）结果（弧度）

求出估计误差（弧度）:

均方误差

error=sqrt（sum（（estimated-sort（DOA））.A2）/p）;%

展开阅读全文