《基于L阵列的目标波达方向估计算法改进(共25页).doc》由会员分享，可在线阅读，更多相关《基于L阵列的目标波达方向估计算法改进(共25页).doc（25页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、精选优质文档-倾情为你奉上基于L型阵列的目标波达方向估计算法改进摘 要：从目标波达方向估计的工程应用角度分析，普通的阵列布局已不能满足现代人们的需求。本文主要阐述了一种采用次最小冗余阵的L型阵列的波达方向估计，与传统L型阵列相比，采用次最小冗余阵能够减少阵元的使用同时提高角度分辨率。本文采用9个阵元的L型阵列，通过采用次最小冗余阵能够实现只使用9个阵元就能达到19个阵元的效果。实验表明，结合传统的MUSIC算法，在进行DOA判断时比普通的L型阵列的分辨率具有明显的提高。关键词：目标波达方向；次最小冗余阵；阵元；MUSIC专心-专注-专业Research on directionofarriva
2、l estimation algorithm based on L - type arrayAbstract：From the perspective of engineering application of reach direction estimation, the general array layout can no longer meet the needs of modern people. In this paper, a method for estimating the direction of arrival of wave in an l-type array wit
3、h subminimum redundancy matrix is presented. In this paper, the l-type array with 9 elements is adopted, and the effect of 19 elements can be achieved by using only 9 elements with sub-minimal redundant matrix. The experimental results show that the resolution of DOA with traditional MUSIC algorithm
4、 is significantly higher than that of ordinary l-type arrayKey words: directionofarrival; Subminimum redundant matrix;Pixelarray; MUSIC目 录1 引言1.1 研究背景及意义技术不断革新的当代生活，人们对信息化以及智能化的设备需求也在不断提升，表现形式也愈加丰富多彩，波达方向（Direction of arrival, DOA）估计主要用于判断信源方向的一种手段。随着近年来关于阵列信号理论和研究的不断发展与进步，加之越来越多的地方应用到智能化的处理，DOA技术又被
5、研究学者们推上研究的重心。阵列信号也就是将一或多组的阵元（传感器）按照一定的规律进行排布组合，构成一个阵列，形成多个阵元接收信号的系统。通过系统不断的在空间内进行对信号的搜寻以及发射接收等动作，将有用的信号保留并放大，对无用的信号或者噪声等信号进行抑制甚至消除，将目标信号中的特征信息也就是有用信息进行提取1。对于传统技术而言，DOA估计算法就必须是符合空间理论的满阵，也就是平常所说的阵元之间的距离必须小于等于入射信号的半波长，以此用来避免产生测量时产生的模糊误差。但是在现代所有产品都趋向微型化和智能化的社会，且在现实生活中，各种干扰因素过多，信号复杂且交叠，参数范围变化幅度愈加不明显，这对传统
6、的满阵DOA方法充满了挑战，导致传统满阵的DOA系统显然已经逐渐被淘汰，所以如何实现空间方位识别率更准，鲁棒性更高，实时性更快的DOA估计方法成为全球相关科研人员及研究机构研究的重点内容2,3。通过查阅资料了解到，阵列的孔径直接对阵列对空间的空间分别率起着重要的作用，且成正比关系，也就是孔径越大则系统的空间分辨率越低。同时，阵元增加不仅增加了硬件成本，对多通道信号的的处理也提出了更高的要求，不适合工程应用4。稀疏阵列指的是阵元之间的距离不小于入射信源半波长的阵列排布。和传统的满阵相比，稀疏阵列孔径变得更大，自由度提高，在方向测量精准度、空间分辨率和多信源监测时拥有更好的表现。本系统采用最小冗余
7、阵和MUSIC算法结合的方法来提高系统分辨率。1.2 研究历史和现状在初期的研究中，最常用的估计方向的方法现在被称为波束形成算法，由于多种因素的原因，该算法在很多应用环境下受限严重，虽然现在技术已经很发达了，但是还是很难在该算法上做出新的突破，影响最严重的是因为该算法受到瑞利的制约，然而根据现代科研学者的研究，发现通过增大天线的孔径可以实现对精度的提高，但是如何合理的增大天线孔径又成了一个难题5。知道20世纪60年代，关于如何提高分辨率的问题再次出现在人们的视野，其中以最大熵法和最小方差法在这方面应用最为突出，但是这两种算法的局限性也很强，使用这两种算法在对线性函数进行预测时准确率很高6，但针
8、对非线性预测时却效果很不理想。而且由于在实际进行DOA估计时，条件非常苛刻，理想环境下的情况几乎不存在，两种线性算法在这时使用的效果很不理想，所以在DOA估计方面很大程度限制了这两种线性算法的发展7。后续出现了更多的优秀算例，例如有些研究学者将时频分布和小波分析加入进了超分辨算法当中。综上所述，窄带信号不是该领域的局限，而当代研究学者的研究方向更倾向于再实际情况下使用，因此对宽带信号得研究逐渐被挖掘出来9,10。信源的数量估计有着不可替代的关键地位，在目前已知的算法中其大部分需要事先知道入射信号源的数量，如果对信源数量的估计不正确的话，不仅是虚警和漏报的问题，同时也会影响子空间之间的正交特性。
9、信源数量估计的错误还会对子空间和信号维度的判断造成影响，致使对方位进行判断时产生的角度差过大。在这类算法的初期阶段，通过增加惩罚函数的方法避免了主观门限的弊端。上述准则的局限性也非常大，仅仅适应在白噪声和非相干信号的前提下，而且这种情况对快拍数的要求也非常的高，然而后面通过Wu等人针对于盖氏圆盘定理的信源数目估算方法，该算法在对信源数量进行估计时不需用对协方差矩阵的特征值，而是利用噪声和信号圆盘半径差距的特性。通过上述分析可以知道，超分辨测向技术存在这非常多的缺点，最知名的缺点就是该技术对环境要求非常高，而且对阵列矢量的精度要求非常高，必须在没有任何误差的情况下，才能最大限度的保持其精准度，这
10、样就限制了其在实际生活中的应用，因为在正常的实际使用中，实际的阵列形状会受到多种外界环境的影响，会对测量的结果形成偏差和扰动，从而对方向的判定效果和性能造成极大影响。因此，在实际的生活应用中，阵列误差必须被考虑在内，并且必须要通过手段来对其进行调整，确保其始终保持在一定的范围内，如何保证让误差始终保持在一定的范围内，这个问题一直到上世纪90年代末才被研究的一种新方法解决，该方法通过将阵列校正通过一种参数估计的方法，使得其波动范围始终在一定的范围内，这样有效的解决了上述的问题，而且，使用该方法硬件成本低，设备兼容性高，同时还能够针对不同的环境进行自适应，从而将该项技术又向前推进了一步。2 DOA
11、基本概念2.1 DOA估计阵列数学模型在进行DOA空间谱估的数学建模时，可在理想状态下记性假设，其假设如下： 每个待检测的信源极化相同且互不相关，通常我们认为信源为窄带信号，且每个信源有着一样的中心频率，假设待检测的信源数量为D。 天线阵列看成是由M（MD）个阵元构成的同距线阵，阵元特性相同、各项同性，设d为阵元之间的距离，根据阵元间距规则阵元之间的距离不能大于最高频率波长信号的一半。 在每个远场的信源中必须要放置天线，因为对方向的估计就是通过信源中天线的波动来实现断定的。 对于除了主接收线路外，其他各支路在理论和实际应用中的特征完全一样。图1为DOA估计的一般模型图。图 21 DOA估计的一
12、般模型假使使用第个信源用于对天线阵列的辐射，其中辐射信号为，此时我们令作为这次辐射信号中的窄带信号，这样我们可以将表示为 通过一系列的变化可以得到如式2.2所示： 式中，为光速；为信号波长。在我们常规的测量中使用来表示经过天线阵列时电磁波所需用的时间，同时使用窄带假设理论对其进行计算，得到如下的结果： 因此延迟后的波前信号如下所示： 因此，假设选择的参考点为第一个阵元，那么在t时刻同距线阵中的第个阵元能够对第k个信号源进行感应，感应的信号为 通过上式，我们可以知道是第k个信源受第m个阵元影响产生的。上面已经对阵元方向进行了假设，即所有阵元均无方向，所以当时第k个信源的方位角我们使用来进行表示。
13、式2.6是我们通过计算得到的第m个阵元信号的输出方程，想要得到该方程就必须要对所有的声源进行监测，通过要去除所采集的声源中噪音和无用的部分。 式中，为测量噪声。标号即代表着第多少个阵元。设 通过矩阵运算，我们可以变换出更加简单明了的表达式： 式中，在一般的信号处理中，一段真实的语音信号可以分解成若干信号的叠加，各种各样的声音通过叠加形成了我们采集到的声音，在这些叠加的声音中，有我们需要的特征声音，当然也包含了许多杂乱的无用信号，如何去除这些无用信号也是语音信号处理的关键要解决的问题，在DOA估计上，我们可以通过对X进行阵列输出的相关处理，从而得到 式中，表示矩阵共轭转置。 2.2 DOA估计方
14、法波束形成和电子引导作为DOA估计的传统方法的理论概念，这两种方法并非是单一的对信号进行接收，然后通过信号的特征进行识别然后建模而实现的估计，在这里我们需要用这两种方法对信号进行采集，与此同时对阵列进行一定程度上的引导，通过原始的DOA估计方法结合现有的技术电子引导，这样可以使我们在进行DOA估计时任意的调节我们需要的方向，并且能够针对峰值信息进行寻找。下面我们将对传统的DOA估计方法进行描述和讨论。2.2.1 延迟相加法延时-相加法其采用的是最经典的波束形成和傅里叶变换的方法。图2给出了经典窄带波束形成器的结构。图 22 窄带波束形成器构造输出信号y(t)可以看做是传感器阵元的输出的线性加权
15、和，即 式中，中针对对个方面的因素进行了包含和容纳。原始的波束形成器的总输出功率可以用下式进行表示 式中， 是关于的导引向量；是阵列输入端的噪声向量；和分别是信号和噪声的功率。可以看出，当时，输出功率最大。2.2.2 Capon最小方差法Capon最小方差法的初衷是为了解决延时法针对分辨率过低这一缺点而提出的一种新的方法，这种方法既然被提取出来使因为其算法的优越性主要体现在仅使用一个波束就可对多方面的问题进行监测和估计，使用这种方法可以将输出功率降低，达到增益观测方向保持为零数的前提下降低非期望干扰的贡献，增益观测方向通常为1。 其约束条件为。利用拉格朗日乘子，权向量解为 利用此方法，可以通过
16、Capon空间谱得到输出功率关于DOA的函数，即 2.2.3 MUSIC算法MUSIC算法和其他算法相比起来具有很强的适应性，可以对一维和二维的阵列进行DOA估计，在知道阵列分布的情况下即。可得 (2.19)根据上式，得到矩阵为满秩矩，非奇异，所以有逆存在。所以上式通过变化可转换为，这种变化足以说明中噪声子空间和各个列矢量之间存在正交关系，所以有根据信号我们通过对上式的变换，采用寻找波峰的方法用来估计到达角，上式还可进行变换为 (2.20)式中：，为的投影矩阵。通过上述理论知识的讲解，我们可以将MUSIC算法归结为以下几个步骤： 协方差矩阵和估计值我们可以通过N个接收信号矢量计算得到，即 (2
17、.21)我们通过特征分解的方法对式2.21进行分解能够得到 通过对所得到的特征值的大小进行对比，将特征值和K个信号与对应的矢量信号看成一个信号的子空间然后把剩下的全部看做噪声子空间，有 (2.22)3 次最小冗余线阵及其数学模型3.1 非均匀线阵阵列模型先以线阵为例，然后根据线阵对L阵采用同样方法进行处理。我们先假设一个由M个阵元组成的非均匀的线阵，将左边第一个阵元作为系统的参考阵元，我们可以使用D=（d1,d2,dM）来表示不同阵元的位置，分别从发射P个波长为的窄带平面波，其第i个阵元接收到的信号为 (3.1)式中，作为复包络表示第m个信号的。看做加性噪声是第i个阵元接收到的。将式（1）以向
18、量的形式进行表示： (3.2)采用数学的方法结合自己的DOA只是对式3.2进行分析，阵列信号的输出使用表示，阵列信号的输入采用表示接收的噪声向量用表示。定义为方向向量构成的阵列流形，其中。记，则信号在方向上的方向向量可以表示成： (3.3)3.2 次最小冗余线阵通过3.1节的分析与讨论，我们分别采用了理论知识与数学推导结合的方法，对可能存在的线阵阵元中的输出共轭循环函数进行了分析，如式3.4。 (3.4)其中，为信号的循环频率。设第个信源共轭循环自相关函数为，将式（1）代入式（4），由的循环独立性，有 (3.5)从上市可以得到，共轭循环相关函数对于非均匀线阵来说至于和 有关。 (3.6)上式对
19、的所有可能进行的列出，当A与B等价时，即是一组存在一样孔径的阵列，通过式3.6也可以看出，在均匀的线阵模式下，由于线阵不同组合存在不同的冗余度，多以很难获得具有相同共轭循环相关的函数值。通过上述的众多理论分析，我们开始对最小冗余线阵进行配置，最小冗余线阵也就是采用最少的阵元来实现相同的效果，最小冗余线阵就是对传统的线阵进行了一系列的优化，其配置是不固定的，我们可以通过计算机穷尽搜索的办法得到，表1给出了其中一种配置，其中M、N分别代表了阵元数和孔径大小，集合中的数代表每个阵元距离第一个阵元的距离数。表3-1 最小冗余线阵归一化配置MN6120,1,5,7,97160,1,3,6,9,13820
20、0,1,4,11,16,19,219260,1,4,11,17,18,25,3710320,1,4,12,19,28,31,30,4011410,1,3,9,15,29,37,41,45,4612590,1,3,6,15,20,28,36,38,44,49上面表格只是我们通过穷举法穷举出来的其中一种阵列排布方式，虽然我们通过穷举法可以穷举出无穷多种可能，但是阵元的排布方式并不是随便的更改的，我们通过表1能够清晰的看出，最小冗余阵阵元间的距离并不是对称的，由于非对称的线阵并不适应于MUSIC算法，所以我们采用同等变换的方法，构建出新的等价阵列，这种阵列被称为次最小冗余阵。4 仿真与结果分析4.1
21、 3个信号源L线阵模型对比实验中设置信号是频率为f的调幅信号，采样率为2f，快拍数为1024，其他条件相同的情况下，对比可以看出普通L阵的阵列在采用9个阵元和采用最小冗余阵列采用9个阵元在角度分辨率上明显不足。 图4-1 3个信号源情况下L线阵阵元分布图在相同的实验条件及设置下，可以明显看出，使用次最小冗余阵结合MUSIC算法的系统比传统MUSIC算法从空间谱的角度分析，分辨率明显提高，精确度也显著提高。图4-2左为处理前右为处理后4.2 5个信号源L线阵模型对比实验中设置信号是频率为f的调幅信号，采样率为2f，快拍数为1024，其他条件相同的情况下，对比可以看出普通L阵的阵列在采用9个阵元和
22、采用最小冗余阵列采用9个阵元在角度分辨率上明显不足。图4-1 5个信号源情况下L线阵阵元分布图在相同的实验条件及设置下，可以明显看出，使用次最小冗余阵结合MUSIC算法的系统比传统MUSIC算法从空间谱的角度分析，分辨率明显提高，精确度也显著提高。图4-2 左为普通music算法右为最小冗余处理后4.3 最小均方误差利用本文所提出的改进型算法除了可以通过对角度判定的分辨率以外，我们还可以通过最小均方误差来看图4-3 最小均方误差计算综上所述，无论是从分辨率还是从最小均方误差，本文所提出的算法都比传统的MUSIC算法要优越。完全符合设计要求。5 结论本文以结合工程实际应用为出发点，提出了一种采用
23、次最小冗余线阵的DOA估计方法，并进行的matlab仿真，通过实验结果我们可以看出，MUSIC与之结合实现的效果非常理想，并且采用次最小冗余线阵，其一是为了使用相同的阵元实现更高的分辨率，提高天线口径；其二由于次最小冗余矩阵拥有很强的灵活性，在进行选择时能够提供多种布阵方式，我们可以根据不同的实际应用来进行选择，通过选择的方法在一定程度上对误差较大的通道进行了筛选和避免，并且这种方法还可以在一定程度上保证了相同的孔径的分辨率不变得理想环境下，通过对阵元数量得降低来对算法和多通道信号处理得复杂度进行一定程度的降低。该系统使用的次最小冗余矩阵虽然带来了很多优点，但是正是由于其阵列的变化多端，所以在
24、阵列选择的优先顺序上，仍有很长的路要走。参考文献1 张小飞，汪飞,陈伟华等阵列信号处理的理论与应用M.北京国防工业出版，2013:1-6.2 Radar Emitter Signals Recognition and Classification with Feedforward NetworksJ.Nedyalko Petrov,Ivan J ordanov, Jon Roe.Procedia Computer Sci ence.2013.3 王猛.基于高阶累积量的声矢量阵列信号DOA估计算法研究D.吉林大学，2017.4 王姝.宽带信号DOA估计和DBF技术算法研究D.成都:电子科技大学硕
25、士学位论文.20185 Cands E. Compressive samplingC. Proceedings of the International Congress of Mathematicians, Madrid, Spain, 2016, 3: 1433- 1452.6 Cands E,Romberg J,Tao T.Stable signal recovery from incomplete and inaccurate measurementsJ.Comm.on Pure and Applied Math, 2014, 59(8): 1207-1223. DonohoDL. C
26、ompressed sensingJ.IEEE Trans.on Information Theory,2017,52(4):1289-13067 Cands E, Wakin M B. An introduction to compressive samplingJ. IEEE Signal Processing Magazine, 2018, 25(2): 21-30.8 Sofotasios PC, Mohjazi L, Muhaidat S, Al-Qutayri M, Karagiannidis GK. Energy detection of unknown signals over
27、 cascaded fading channels. IEEE Anten Wirel Propag Lett.2016;15:1359 Wu Y, Black RA, Jeffs BD, Warnick KF. Cancelling non-linear processing products due to strong out-of-band interference in radio astronomical arrays. In: 2015 IEEE Signal Processing and Signal Processing Education Workshop (SP/SPE).
28、 Salt Lake City, Utah, USA; 2015: 272-277.10 Hayashi H, Ohtsuki T. DOA estimation for wideband signals based on weighted squared TOPS. EURASIP J Wireless Commun Networking.2016;2016(1):243.致谢本论文在陈浩老师的悉心指导下终于成功完成。回首这三个月左右的时间，从刚刚拿到这篇论文时的迷茫无力感，到自己真正成功写出一篇自己的毕业论文的成就感，可以称得上是为大学四年的落幕添了一道彩虹。在此我要特别感谢陈浩老师，在论
29、文编写过程中认真负责，严格而又不失亲切，多次对论文提出建设性的意见，为我们组逐一批注我们每个人在格式和内容上的缺点。也是陈浩老师的努力才使得我们的论文进度进行的如此顺利。在此还要感谢我组的每一位组员。我们一直在不断交流不断分享我们之间点点滴滴从这次论文编写中获得的经验，为我们每个人都的成长都有不可或缺的作用。时光匆匆如流水转眼已是大学毕业时节，犹如白驹过隙，离校时间已近在眼前，在此再次由衷的感谢我遇到的老师同学，是你们构筑了我前进道路上的每一片砖瓦。附录1clc % 初始化matlabclearclose all% 参数设置% 常数设置c = 3.0e8; % 光速% 工作参数设置f = 28
30、e9; % 频率lambda = c/f; % 波长SN = 1024; % 快拍数k = 2*pi/lambda; % 波数snr = 20; % 信噪比% 设置参数M = 5;N = 5;dx = 0.5 * lambda;dy = 0.5 * lambda;% 得到阵列单元位置,阵列完整的数据保存在cy cx内coorsx = (0:N-1)*dx, zeros(1,M-1) ;coorsy = zeros(1,N), (1:M-1)*dy ; % 转化为一维数据cy = coorsy(:);cx = coorsx(:);% 画图figureplot( cx/lambda, cy/lam
31、bda, r., markersize,20 )xlabel(x/)ylabel(y/)% 设置来波方向doas_theta = 10, 20, 30, 40, 50 *pi/180; % 范围是090doas_phi = 10, 40, 70, 100, 150 *pi/180; % 范围是090% 提取信源数量source_number = length( doas_theta );% 模拟接收数据% 得到A 6.3.3中S前面的矩阵temp = cx*( sin(doas_theta).*cos(doas_phi) ) + cy*( sin(doas_theta).*sin(doas_p
32、hi) );A = exp(1i.*k.*temp);% 生成随机信源st = rand(source_number,SN) + 1i * rand(source_number,SN);% 生成接收信号x,并添加噪声x = A * st;x = awgn( x, snr, measured );% MUSIC算法doa估计% 计算协方差矩阵RR = x*x/SN;% 特征值分解V,D = eig(R);D = diag( D );% 对D进行排序D,index = sort( D, descend );V = V(:,index);% 根据D判断信源数source_number_music 如
33、果与source_number不符 则报错temp = D(1:end-1)./D(2:end);source_number_music = find( temp 2, 1,last);% 判断if source_number_music = source_number fprintf( music结果信源数为%dn,source_number);else% error(信源数估计错误n);end% 求空间谱theta = linspace(0,pi/2,1001);phi = linspace(0,2*pi,1001);phph,thth = meshgrid( phi, theta );%
34、 初始化空间谱PP = zeros( size(phph) );% 提取噪声子空间EN = V(:,source_number+1:end);% 开始求空间谱for i1 = 1:size(P,1) for i2 = 1:size(P,2) % 提取当前的theta和phi theta = thth(i1,i2); phi = phph(i1,i2); % 计算方向矢量 temp = cx*( sin(theta).*cos(phi) ) + cy*( sin(theta).*sin(phi) ); An = exp(1i.*k.*temp); % 计算P(i1,i2) temp = An*E
35、N; temp = temp * temp; P(i1,i2) = real( 1./temp ); endend% 画空间谱PP = P./max(max(P);P = 10.*log10(P);figuremesh(thth*180./pi,phph*180./pi,P);xlabel(theta(); %xlim(0,90)ylabel(phi(); %ylim(0,90)title(空间谱P(dB)% 求波达角temp = imregionalmax( P );aa = find( temp = 1 );% 提取所有角度 去除重复的all = thth(aa)*180./pi, php
36、h(aa)*180./pi;m = 1;tol = 3;while m= size(all,1)-1 % 计算距离 index = m+1:size(all,1); dist = pdist2( all(index,:), all(m,:) ); loc = find( dist 2, 1,last);% 判断if source_number_music = source_number fprintf( music结果信源数为%dn,source_number);else error(信源数估计错误n);end% 求空间谱theta = linspace(0,pi/2,1001);phi =
37、linspace(0,2*pi,1001);phph,thth = meshgrid( phi, theta );% 初始化空间谱PP = zeros( size(phph) );% 提取噪声子空间EN = V(:,source_number+1:end);% 开始求空间谱for i1 = 1:size(P,1) for i2 = 1:size(P,2) % 提取当前的theta和phi theta = thth(i1,i2); phi = phph(i1,i2); % 计算方向矢量 temp = cx*( sin(theta).*cos(phi) ) + cy*( sin(theta).*s
38、in(phi) ); An = exp(1i.*k.*temp); % 计算P(i1,i2) temp = An*EN; temp = temp * temp; P(i1,i2) = real( 1./temp ); endend% 画空间谱PP = P./max(max(P);P = 10.*log10(P);figuremesh(thth*180./pi,phph*180./pi,P);xlabel(theta(); %xlim(0,90)ylabel(phi(); %ylim(0,90)title(空间谱P(dB)% 求波达角temp = imregionalmax( P );aa =
39、find( temp = 1 );% 提取所有角度 去除重复的all = thth(aa)*180./pi, phph(aa)*180./pi;m = 1;tol = 3;while m= size(all,1)-1 % 计算距离 index = m+1:size(all,1); dist = pdist2( all(index,:), all(m,:) ); loc = find( dist tol ); loc = index(loc); all(loc,:) = ; aa(loc,:) = ; % 更新m m = m + 1;end% 找到loc中P值最大的几个PV = P(aa);PV,index = sort(PV,descend);aa = aa(index);aa = aa(1:source_number_music);% 得到波达角度doas_theta_music = thth(aa)*180./pi;doas_phi_music = phph(aa)*180./pi;% 排序doas_theta_musi