（信息与通信工程专业论文）高分辨率频率测量方法研究.pdf

国防科学技术大学研究生院学位论文 摘要 本文的写作是针对在噪声背景下，尤其信噪比比较低时，估计信号的频率参 量。首先回顾了频率测量的历史，总结了经典法和现代法信号频率估计中主要的 频率测量方法，主要包括傅立叶变换( f f t ) 法、最大似然( m l ) 法，修正协方 差法、多信号分类( m u s i c ) 法和子空间旋转不变法( e s p r j t ) 等。针对加噪信 号的两个接近的频率，应用以上提到的各种方法进行频率估计，并对各种方法下 得到的估计结果进行比较，进而得出一种比较理想的频率测量方法- - e s p r i t 法。 这种方法在较低信噪比情况下，仍具有较高的分辨率，在保证精度的基础上。极 大的减少了计算量，且不需作扫描处理。 文章针对对比得到的e s p r i t 方法，就其算法尝试在硬件上进行实现，并对算 法中的关键部分奇异值分解( s v d ) 应用c o r d i c 算法在f p g a 上实现。所有这些 仿真设计都是在m a t l a b 7 0 和x i l i n x 的f p g a 芯片及其i s e 7 1 为基础进行的。 最后，对本文的工作做出了总结。 关键字：频率测量复高斯白噪声信噪比m a t l a b 软件e s p r i tf p g a c o r d i c s y s t o l i c 第i 页 国防科学技术大学研究生院学位论文 a b s t r a c t n 坨w r i t i n go ft h i sp a p e ri su n d e rt h eb a c k g r o u n do fa d d i t i v en o i s e ，e s p e c i a l l yt h e l o ws n rt oe s t i m a t et h ep a r a m e t e ro f 丘铷u 吼c y a tf i r s t , t h eh i s t o r yo ff r e q u e n c y e s t i m a t i o ni sr e v i e w e di nt h i sp a p e r , a n ds o m eo f t h em a i ne s t i m a t i o na r i t h m e t i c s ，w h i c h i n c l u d ef f tm e t h o d ：m a x i m u mi i k e l i h o o dm e t h o d ；m o d i f i e dc o v a r i a n c em e t h o d ； m u l t i p l es i g n a lc l a s s i f i c a t i o nm e t h o d ；e s p r 】tm e t h o da n ds oo n , a r es u m m a r i z e di nt h e s e c o n dc h a p t e r t h e n , t w os i g n a l sw i t l lc l o s ef r e q u e n c i e si nw g na r ep r e s e n t e d , t h e f r e q u e n c i e sa r ee s t i m a t e db yt h em e t h o d sm e n t i o n e da b o v e c o m p u t e rs i m u l a t i o n sh a v e b e e nm a d ei nt h i sp a p e ra n dt h er e s u l t so ft h e s em e t h o d sh a v ec o m p a r e d a tl a s t , a f r e q u e n c ye s t i m a t i o nm e t h o do f b e t t e rc a p a b i l i t y - e s p ti ss e l e c t e d ，w h i c hh a sh i 曲 r e s o l u t i o na b i l i t 3 , i nt h ec o n d i t i o no f l o ws n i 乙g r e a t l yr e d u c e st h ea m o u n to f c a l c u l a t i o n w i t h o u tt h ed e c r e a s eo f p r e c i s i o n , a n dn e e d sn os c a n n l em e t h o d e s p r j ts e l e c t e da b o v ei st r i e dt or e a l i z eo nh a r d w a r e ，a n dt h ec o r e o fi ts v di si m p l e m e n t e db yt h eu s i n go fc o r d i ca l g o r i t h mo nf p g a a l lt h er e s u l t s i nt h i sp a p e ra r es i m u l a t e di nm 魄t l a b7 0a n dd e s i g n e db y 廿l ed e v i c e so fx i l i n xa n d i t sd e v e l o p m e n ts y s t e mi s e7 1 as u m m a r i z a t i o no ft h ew o r ko ft h i sp a p e ri sg i v e n a t t h ee n do f t h i sp a p e r k e yw o r d s ：f r e q u e n c ye s t i m a t i o n ，w g n ，s n r ，m a t l a bs o f t w a r e ，e s p r i t c o r d i c ，s y s t o l i c 国防科学技术大学研究生院学位论文 表目录 表2 1 不同m 值和不同方法下频率的测量值2 3 表2 2 各方法仿真时间2 5 表2 3 调制信号载波频率估计2 6 表4 1 迭代角度3 6 表4 2c o r d i c 算法模式4 2 表4 3c o r d i c 算法功能4 2 表4 4c o r d i c 算法收敛范围4 4 表5 - 1c o r d i c 与m a t l a b 计算结果对比5 3 第1 i l 页 国防科学技术大学研究生院学位论文 图目录 图1 1 侦察接收模型2 图2 1f f t 功率谱图1 7 图2 2 最大似然功率谱图。1 8 图2 3e s p r i t 法功率谱图1 8 图2 , 4m u s i c 法功率谱图1 9 图2 5 修正协方差法功率谱图1 9 图2 6 最大似然频率分辨性能功率谱图2 0 图2 7e s p r i t 频率分辨性能功率谱图2 l 图2 8m u s i c 频率分辨性能功率谱圉2 1 图2 , 9 修正协方差频率估计性能功率谱图2 2 图2 。1 0 最大似然法频率估计精度2 3 图2 1 1e s p r i t 法频率估计精度2 4 图2 1 2m u s t c 法频率估计精度2 4 图2 1 3 修j 下协方差法频率估计精度2 5 图3 1e s p r i t 算法结构| _ 2 7 图3 ，2j a c o b i 并行排列规则3 0 图3 3s y s t o l i c 阵列结构3 l 图3 4 矩阵s v d 的b r e n t - l u k v a nl o a n 阵列结构3 2 图4 1 矢量旋转3 5 图4 2 圆周模式( m = 1 ) 矢量旋转轨迹3 7 图4 j 3 线性模式( m = o ) 矢量旋转轨迹3 7 图4 4 双曲线模式( m - 一1 ) 矢量旋转轨迹3 7 图4 , 5 一次c o r d i c 迭代的基本实现结构( 注：a = 2 17 ) 3 9 图4 , 6 旋转模式。圆周模式( m = 1 ) 矢量旋转轨迹4 0 图4 7 矢量模式，圆周模式( i n _ 1 ) 矢量旋转轨迹4 1 图5 1c o r d i c s v d 阵列结构图4 6 图5 2 数据调度规则4 7 图5 3 6 x6 阵列调度4 7 图5 4 调度过程4 8 图5 5 伪s y s t o l i c 运行过程一4 9 图5 6 加初始数据的串行数据调度4 9 第1 v 页 国防科学技术大学研究生院学位论文 图5 7 旋转角计算单元实现5 0 图5 8 双边旋转实现结构，5 l 图5 9 延时单元5 l 图5 1 02 x 2 矩阵的奇异值分解5 2 图5 ，1 l2 x 2 矩阵s v d 输入5 2 图5 1 22 x 2 矩阵s v d 输出5 3 图5 1 34 x 4 矩阵s v d 结构5 3 图5 ，1 44 x 4 矩阵s v d 输入5 4 图5 1 54 x 4 矩阵s v d 输出5 4 第v 页 独创性声明 本人声明所呈交的学位论文是我本入在导师指导下进行的研究工作及取得的研 究成果尽我所知，除了文中特别加以标注和致谢的地方外，论文中不包含其他人已 经发表和撰写过的研究成果，也不包含为获得国防科学技术大学或其它教育机构的学 位或证书而使用过的材料与我一同工作的同志对本研究所做的任何贡献均已在论文 中作了明确的说明并表示谢意 学位论文题目：矗l 毛趱壅塑窒型量左洼丑峦 学位论文版权使用授权书 本人完全了解国防科学技术大学有关保留，使用学位论文的规定本人授权国 防科学技术大学可以保留并向国家有关部门或机构送交论文的复印件和电子文档，允 许论文被查阅和借阅；可以将学位论文的全部或部分内容编入有关数据库进行检索， 可以采用影印，缩印或扫描等复制手段保存、汇编学位论文 ( 保密学位论文在解密后适用本授权书) 学位论文题目：直l 艺蕴空夔空型量左洼盈究 学位论文作者签名： 作者指导教师签名：互釜监 日甄：o 奄菇年f f 觋f d b a 期：矽弼年f f 蜀i h 国防科学技术大学研究生院学位论文 第一章绪论 本章介绍论文写作的课题背景，对频率测量的意义和研究的历史做简要的概 括，并在此基础上对频率测量技术的关键问题做出一定分析。最后，对论文写作 的主要内容进行阐述。 i i 课题的背景和意义 信号与信息处理是信息科学中近十几年发展最为迅速的学科之一。参数估计 是信号与信息处理学科的重要组成部分，也是近年来非常活跃、发展迅速的一个 研究领域。这是因为参数估计作为一种重要的信号分析手段，有着广泛的应用背 景和重要的应用价值，其应用涉及雷达、声纳、导航、通信、成像、地质勘探， 生物医学工程等诸多军事及民用领域。随着信号处理的快速发展，人们对信号的 检测和参数估计的要求也越来越高，提出精度更高、分辨能力更好、稳健性和容 错性更强，计算更简单的估计方法一直是人们不断努力的方向。 在时间序列中，我们关心的一个重要的信号参数就是信号的频率。人们之所 以投入大量的精力关注到频率估计的领域，是因为这个领域的研究具有非常重要 的现实意义。估计淹没在噪声中的正弦波的频率，是现代信号处理中最有实用价 值的技术之一，是测试所有谱估计性能的基础，也是信号处理技术的基础。该技 术在雷达、通信干扰、声纳、震动测量、地震信号处理及电子监测等领域有着广 泛的应用，若该技术取得一定的突破，必然推动以上领域的飞跃式发展。随着现 代通信及信息处理技术的迅猛发展，对这个领域的研究也必然会提出更高的要求。 本论文设定在一个复杂空间电磁环境下，利用仿真技术对各种侦察方法的性 能进行定量分析。侦察接收机包括天线、接收机前端和信号处理部分，天线和接 收机完成信号截获和参数测量，而信号处理部分完成参数估计、信号分选、调制 识别和解调。而本文所要完成的任务就是a d 采样和参数估计部分，也即主要对 信号频率的测量。其中侦察机的模型如图1 1 。 参数测量是信号截获的目的，也是终端信号处理的基础，被测参数的种类和 测量精度直接影响到终端信号的功能和性能。其中，信号到达角和载波频率是两 个最基本、最重要的参数。它们既是实现侦察所必须的，也是进行信号分选的重 要参数。载波频率与辐射源的特征和属性密切相关，是进行辐射源识别的重要依 据。传统的载波频率估计方法都是假设背景为白噪声。由于实际空间环境噪声的 复杂性，基于信号总是淹没在一定噪声中的实际情况，因此，有必要研究噪声中 的信号载波估计方法。同时，系统面临着多个信号同时到达的问题，因此也要对 第1 页 国防科学技术大学研究生院学位论文 信号源的数目进行估计，这也是利用高分辨率方法实现频率估计的前提条件。这 两个方面也正是本文重要研究的内容。 1 2 1 频率估计技术 接收机前端信号处理模块 图l 1 侦察接收模型 1 2 研究现状和发展趋势 自从2 0 世纪3 0 年代，特别是第二次世界大战以来，由于在战争中雷达和声 纳技术的广泛应用，使人们开始注意到以上两种技术中的关键技术一也就是淹没 于噪声中的正弦信号频率估计问题的研究的重要性。频率估计是参量估计中的一 个重要物理量，有着广泛的应用领域。电子对抗在现代战争中的地位越来越重要， 电子对抗中一种常用手段就是截获敌方情报信号，然后以相同的频率调制方式发 送有利于自己的假情报，干扰敌方部署，或使敌方的通信系统不能正常工作。自 此，频率估计技术开始真正的拓展成为一个相对独立的研究领域。从上个世纪6 0 年代初，为了满足通信技术及电子测量技术的大规模发展应用的要求，人们开始 更广泛的关注这个领域的研究。研究频率估计方法以及各种方法的性能对国民经 济的发展和国防建设都有重要的意义。 1 9 5 8 年b l a c k m a n 和t u r k e y 提出了用取样自相关函数的傅立叶变换作为功率 谱的估计，这种方法称为谱估计的自相关法。应用它进而可以估计正弦波的频率。 1 9 6 5 年由j w c o o l e y 和j w t u r k e y 完善了著名的快速傅立叶变换( t h ef a s t f o u r i e rt r a n s f o r m 简称f f n 算法，把离散傅立叶变换的计算量减少了两个数量 级，从而使傅立叶变换走向工程实用，也使周期图的谱估计方法很快流行起来。 周期图和自相关法以及它们的改进方法称为谱估计的经典方法。 经典谱估计方法都是以d f t 为基础，它隐含着对无限长数据序列进行加窗处 理，而窗函数导致了频谱主瓣能量向旁瓣“泄露”，使主瓣产生很大的失真。为 了克服以上缺点，人们提出了平均、加窗平滑等方法，在一定程度上改善了经典 第2 页 国防科学技术大学研究生院学位论文 谱估计的性能，快速傅立叶变换( f f t ) 的发现，使经典法一度得到广泛的应用。 但是，经典法的局限性在于其分辨率( 其分辨率为1 n ，n 为数据长度) ，所以 要提高分辨率则需要加强数据长度，从而增加了运算时间；反之，要想减少运算 时间，只能采用短数据，这样经典方法始终不能解决频率分辨率和数据长度之间 的矛盾。这就促使了现代谱估计方法研究的开展。 1 9 6 7 年j p b u r g 提出的最大熵谱估计方法和j c a p o n ( 1 9 6 9 年) 的最大似然 谱估计( 最小方差谱估计) 方法是高分辨谱估计方法的开始。随后，高分辨谱估 计方法迅猛发展，出现了大量优秀算法，1 9 7 3 年，由p i s a r e n k o 提出的谐波分解方 法提供了可靠的频率估计方法，而s c h m i d t ( 1 9 7 9 年) 提出的多重信号分类( m u s i c ) 算法开创了子空间信号处理算法的研究先河，可视为信号处理算法研究中的一个 重要里程碑。子空间算法通过对数据矩阵的特征值分解或协方差矩阵的特征分解 来获得信号子空间和噪声予空间。后来，r o y 等人( 1 9 8 6 ) 提出的子空间旋转不变法 ( e s p r i t ) 又为特征结构法的研究翻开了新的一页，该方法利用信号子空间的旋转 不变特性来估计信号参数，避免了m u s i c 算法因需要一维搜索而带来的大量计 算，同时在一定程度上降低了算法对硬件的要求，因而在保证精度的基础上极大 的减少了运算量，而且有较好的鲁棒性和分辨率，其在某些特定的应用领域的应 用可望获得突破。常用的基于特征结构的谱估计方法还有最小范数法( m i n - n o r m ) ， r o o t m u s i c 法，加权子空间拟合法( w s f ) 以及基于多维搜索的最大似然( m l ) 方法 等。 时至今日，可以说关于谱估计的基本理论和基本算法己经相当成熟。但是， 实际的应用系统还十分罕见，其原因是大部分超分辨方法的优良性能都是在理想 的信号及噪声模型下取得的，而在实际应用中，理想条件总是难以满足。现有的 高分辨谱估计方法都假定背景噪声是统计独立的白高斯噪声，然而在实际环境中， 白噪声的假设并不总是成立，当实际噪声模型不满足白噪声时，传统的高分辨方 法则将弓l 起严重的性能恶化甚至失效。将现代谱理论和基于高阶累积量的谱估计 方法有效地结合起来，能更为有效地抑制噪声，为提高频率分辨率开辟了一条新 的途径。但高阶谱估计方法算法复杂，难以在工程上上实现，如何在各种噪声中 提高谱估计精度和谱估计稳定性以及估计速度是频率测量中需要解决的关键问 题。同时应用频率测量方法实现对多频率点、多调制方式未知信号频率参数的测 量也是频率参量估计需要克服的问题。 1 2 2 信号源数目的估计 信号源数目估计问题在信号处理中占有重要的地位，它是基于特征分解等高 分辨率算法得以实现的前提条件。信号处理中的信号源数目估计问题可以追溯到 第3 页 国防科学技术大学研究生院学位论文 上个世纪五、六十年代。当时的信号源数目估计方法主要是采用假设检验的方法。 由于假设检验方法需要认为确定一个检测门限，用以与似然比检验统计量进行比 较，所以基于假设检验的信号源数目估计方法容易受到认为的主观因素的影响， 使得估计的信号数目与实际的信号源数目差别比较大。m w a x 和t k a i l a t h 首先 将信息论准则引入，用来解决信号源数目估计问题。m w a x 和t k a i l a t h 给出了 阵列协方差矩阵特征值形式的a i c 准则和m d l 准则，既方便使用又克服了基于 假设检验的信号源数目估计方法的缺点。在m w a x 和t k a i l a t h 之后，多种改进 的基于信息论准则的信号源数目估计方法相继被提出【2 】。 1 3 本文研究主要内容 根据所总结的关于频率测量方法的知识，以及基于在低信噪比环境中和测频 速度的考虑，对各种频率测量的方法的应用有了整体认识和数学模型的理解。有 针对地对多频率点，噪声环境的干扰、低信噪比、多种调制方式下等；对各种频 率测量方法分别进行m a t l a b 仿真实验。对各种方法性能进行对比和评估，选 择具有良好谱估计精度和谱估计稳定性的方法；进而对该频率测量方法的核心算 法部分在f p g a 上实现。 论文第一章主要介绍了频率测量的背景知识和其重要性，回顾频率测量的历 史并指出当前频率测量方法存在的一些问题。 第二章主要介绍目前应用比较多的几种频率测量方法，包括f f t 、最大似然、 m u s i c 法、e s p r i t 法等，并对其特性做一定分析。针对这几种频率测量方法， 在一定噪声环境中的不同信噪比，不同频率点，多调制方式环境下进行m a t l a b 仿真实验，对比不同方法得到一定结论。 第三章基于e s p r i t 方法中的奇异值分解的计算，讨论了奇异值分解的性质， j a c o b i 计算方法以及相应的s y s t o l i c 结构。 第四章重点讨论c o r d i c 原理及其算法模式和两种操作模式，简单介绍了其 收敛性和精度问题。为下章奠定理论基础。 第五章主要运用伪s y s t o l i c 阵列结构和c o r d i c 处理单元结合的方法来实现 矩阵的s v d 分解。 最后对本课题研究所取得成果进行小结，同时指出存在的问题，以及对汲取 的经验教训进行总结。 第4 页 国防科学技术大学研究生院学位论文 第二章频率测量方法与研究 信号频率估计方法的研究已有多年的历史，一般分为经典法和现代法两大类。 经典法主要以傅立叶为基础的周期图法等。以随机过程的参数模型为基础的现代 谱估计方法，可以由信号的谱峰位置得到，其主要分类两大类，一类是参数模型 法，主要有自回归模型，滑动模型和自回归滑动模型等；另一类现代谱估计方式 是非参数模型法，包括最大似然法，多信号分类法和子空间旋转不变法等【3 1 。下面 对实际中应用比较多且效果较好的频率估计方法进行原理介绍、仿真实现以及结 果性能分析。 2 1 1f f t 法 2 1 原理介绍 基于f f t 主要测频原理为：设单一频率正弦信号表示为： s q ) = a c o s ( 2 n f o r + o o ) ( 2 1 ) 其中口、a 和岛分别为正弦信号的幅度、频率和初相。按等间隔= 丁对j 8 在 o r 区间内进行采样得到长度为n 的序列j ( 功。s ( 以) 的点d 刀记为j ( 七) ，鉴 于实序列的d f t 的对称性，忽略d f t 频谱的负频率成分，只考虑离散频谱的前点 n 2 。有 s ( j i ) ：! ：墅堕生三厶型e j t o o - - 等( t ；- o r ， 、7 2 s i n t c ( k 一五t ) n k = 0 ，1 2 n 2 一l ( 2 2 ) s ( t ) 幅度最大值处的离散频率索引值记做k l ，k t = i n t ( f o t ) ，正弦信号的真实频率 f o = k l ( 恤) ，必然对应于连续功率谱的最大值。故可通过搜索正弦信号的频谱的 最大值处的谱线来确定正弦信号的频率。 2 1 2 最大似然估计法 针对最大似然估计法1 4 ，我们设复高斯自噪声v ( 力中混有复正弦信号为 s ( 坊= a 1e x p ( j w l n + 伊1 ) ( 2 3 ) 设x ( n ) = s ( n ) + v ( ”) 的一次实现的个取样值为： 工( 月) = a ie x p ( j w l n + 伊1 ) + v ( 珂) 作= 0 , 1 2 n 一1 ( 2 4 ) 第5 页 国防科学技术大学研究生院学位论文 式中a ，w ，纯分别是正弦信号的振幅，频率和相位；v ( ”) 为复高斯白噪声。s ( n ) 和 石o ) 用矢量分别表示为： 5 = 删5 ( 一1 ) 】 一 ( 2 5 ) = a le x p u 仍) 1e x p ( j w l ) e x p o ( n 1 ) w 1 ) 】7 和 x = 【x ( o ) ，x ( 1 ) ，x ( n 一1 ) 】丁 ( 2 6 ) 白噪声v ( 胛) 的自相关矩阵为： r 。= i ( 2 7 ) 式中，盯：为v ( n ) 的方差，是阶单位矩阵冷： a c t = a ie x p u 识) ( 2 8 ) 又令e j = 【l ，e x p ( j w t ) ，, e x p u ( n - 0 w 0 7 ( 2 9 ) 则式( 2 5 ) 可表示为： s = a g e ； ( 2 1 0 ) 为求a ，w ，仍的最大似然估计，需求解： l ( a q ，w 1 ) = ( x s ) ”( x s ) ( 2 1 1 ) 如果已知w 。，则可求出使( 2 ，1 1 ) 式最小化解4 c 。，式中m 如果用最大似然估计取代， 则a c 就是正弦波复振幅的最大似然估计( 2 1 1 ) 式的求解即为x ( n ) 的周期图，故 有结论：在复高斯白噪声中混有单个复正弦信号，其频率的最大似然估计，可根据数 据的周期图的最大值所在的频率位置求出来。 下面讨论复高斯白噪声中混有多个复正弦信号的频率估计问题。设s 由m 个 正弦波构成 s = a c l e i ( 2 1 2 ) t = l 为求正弦波参数的最大似然估计，必须使以下的函数最小化 mm l = ( x 一a c l e i ) “( x 一a c l e i ) ( 2 1 3 ) t = lf l l 令e = k p 2 ，a e r ，a e = a c i ，a c 2 ，a a c u 7 ，q = w 1 ，w 2 ，a 】则式( 2 - 1 3 ) 可写成 以下形式 l ( a e ，q ) = ( x e a e ) ”( x e “a e ) ( 2 1 4 ) 假设e 是已知矩阵，这上式的最小化问题实际上是一个典型的线形最小二乘方问题 为求频率的最大似然估计，必须求解上式右端二项的最大化问题，即： l ( a e ，q ) = i 片( x e 7a c ) = x x x 珂e r ( e e ，) 一1 e x ( 2 1 5 ) 这是一个高度非线形的方程，人们提出了用所谓迭代最佳化方法求解，如果各正弦 国防科学技术大学研究生院学位论文 波的频率用周期图能进行分辨，那么式( 2 1 5 ) 的最大解将对应于周期图中对大值 所在的频率本实验中所用仿真数据m = 2 ，用周期图可以分辨这两个频率，故采用这 种方法来估计频率。 2 1 3 互谱e s p r i t 法 2 1 3 1 基本定义 e s p r i t ( e s t i m a t i o no fs i g n a lp a r a m e t e rr o t a t i o ni n v a r i a n c et e c h n i q u e s ) 借助旋 转不变技术估计信号参数方法 6 1 现已成为信源入射方向估计和谐波分析的重要方 法。与传统的m u s i c ( 多信号分类) 方法 7 1 相比较，e s p r i t 方法不需要阵列元响 应特性的具体描述，无需进行阵列校准，而且不必进行频谱搜索，具有比m u s i c 算法更好的稳健性。由于这些特点，e s p r i t 算法得到越来越广泛的应用。算法的 关键在于对数据协方差进行特征分解或对数据进行奇异值分解，从中提取信号子 空间的信息。l i u 等人提出的自适应e s p r i t 算法 s ，以及c h a m p a g n e 提出的自适 应特征分解方法 9 1 ，都大大简化了算法的运算量，使之可以有效的得以工程实现。 我们考虑这样一类复正弦加噪声信号，为如下时间序列： 上 x ( n ) = qe x p j ( c o j n ) + c o , ( n ) ( 2 1 6 ) ，= l j l y ( 胛) = 届e x p j ( c o , n + r , ) + c o y ( n ) ( 2 1 7 ) ，；1 式中口，、( i 1 , 2 ，垡) 为复数谐波信号幅值( 均为实数) 。根据情况不同 口，、屈或是已知的( 如雷达回波等) ，或是q 、层两者之一是已知的，或两者均 为未知的。够( i = 1 , 2 ，q ) 为谐波信号频率。因本文主要讨论互谱的e s p r i t i s 法，所以规定z ( 胛) 、y ( n ) 的所有谐波频率是完全相同的。实际上不存在x ( 栉) y ( 甩) 的 谐波频率完全相同这偶然现象，但由于不同频率的谐波是互不相关的，互不相关 的的谐波信号的互相关函数为零，因此就互谱估计而言，完全可以等效为本文规 定的这种情况。只为信号工0 ) y ( n ) 在各个谐波分量间的相位差。r o a n ) , f o y ( 珂) 为附 加的互相独立的，谱密度未知的，其均值为零的可观测噪声。q ( n ) ，0 ) 可以是 复的，实的；可以是高斯的，也可以是非高斯的。 另外有y ( 聆) 构造一个新过程： 卫 y i ( 厅) = y 伽+ 1 ) = 届e x p j ( c o , ( n + 1 ) ) + o j + c o y ( n + 1 ) ( 2 1 8 ) i = 1 m ( 9 1 ) 的约束条件和物理意义同y ( n ) 。 对于( 2 1 6 ) ( 2 1 7 ) 式所示的信号模型可算出其互相关函数： 第7 页 国防科学技术大学研究生院学位论文 ，0 ( 脚) = 目工( 疗) y ( 厅+ 肌) 】 ：e 杰口，e x p j ( c o , 蝴+ 致o ) 】【妻矗e x p u ( c o , ( ”+ 确+ 只) 】+ 出如+ 搬) 】 lt = l ，= l j ：研芝e - j 哪n 至屈e j o ，( m m ) + o d 】+ 研q ( 行) 羔屈e 柏( 卜日1 】 + e 【a ，口叫嘞”哆( ，l + m ) 】+ e h ( 珂) 哆( 雄+ m ) 】 因为a y x ( n ) ，o j ，( n + 肌) 是零均值，且互不相关的，所以有： ( 研) = 口，f l f e a , p ” ( 2 1 9 ) 将( 2 1 9 ) 代) k p x p 维( p g ) 互相关函数矩阵，其对应的互相关矩阵为： r f = ( o )勺( 一1 )( - p + i ) 勺( 1 )勺( 0 )( p 十2 ) 7 0 ( p 1 ) 勺( p 一2 ) 勺( o ) ：羔q f i r e j mf ；f t ( 2 2 0 ) 1 = 1 f = 【1 ，e 硒，e j 2 。t ，e m 。* r ，其中符号表示共轭转置，丁表示转置。由于向量f j 仅包含了正弦信号的第i 个谐波的频率信息，所以称为f 为信号谐波矢量。显然 r 。也可写成下列矩阵形式： r= f e 9 p f ( 2 2 1 ) x y 其中f = 啊，e 只】为p x q 维复矩阵，e 一= d i a g e 确，p ”，一咚】为复数对角矩 阵。p = d i a g a l 届，口2 屈，口。成】为一实数对角矩阵。由信号( 2 1 6 ) ( 2 1 8 ) 式的 信号模型可算出其互相关函数： n = 【x ( 彬y 1 ( 托+ m = a ，卢，e 嗥e 。嗥“ ( 2 2 2 ) f i i 将( 2 2 2 ) 式代入p x q 维c p g ) 互相关矩阵： r m = ( o )( - 1 )( _ 尸+ 1 ) ( 1 )( o ) r & ( - p + 2 ) j；i； p 一1 ) ( p 一2 ) ( o ) ( 1 )( o ) r ( - e + 0 1 勺( 2 )( 1 )r ，( - e + 2 ) i ii；i r ，, ( p - 2 ) p 一3 ) 勺( 1 ) j ( 2 2 3 a ) 第8 页 国防科学技术大学研究生院学位论文 = f e 一m p f 其中f e j 9 p 同上面含义同。 = d i a g e i * 1 ，e7 4 g q 1 ( 2 2 3 b ) 2 1 3 2 算法基本原理 算法基本思想是：x 经过旋转后得到y ，但这种空间旋转却保持了x 和y 所 对应的信号子空间的不变性。由e s p r i t 方法可知，构造矩阵束和对矩阵束的广义 特征分解是两个关键步骤。同理对互e s p r i t 方法可同样适用。 因此，由r ，和r 矾构造诹。，r 叫j 矩阵束，中为 r ，r 。拍q 旋转因子，使r ，、 r 。张成的子空间是同一空间。 定理：设r 为矩阵束 r 。，r 。 相对应的广义特征值矩阵，若p 是非奇异的， 则r 和q 有如下的关系： r = i 硼 旺m ， 证明：因为e 声、p 均为含有q 个元素对角阵，又由于f 的秩为q ，所以r 。的 秩为q ，类似知r 。的秩也为q ，考虑下式： r ，一r r m = f e 7 9 p f 一，f e7 e ( p p f = f e j a ( i r q ) p f ( 2 2 5 ) 由( 2 2 5 ) 式给出的矩阵束及由秩恒等关系( 一个矩阵，每当它左乘( 右乘) 一个 满列( 满行) 秩的矩阵时，其秩不变) 有下列关系式： ( 1 ) 当r p 的一般情况有： r a n k ( r w r r 。) = r a n k ( i 一，中) = g ( 2 2 6 ) 这是因为矩阵f 显然为非奇异的，q 是由( 2 2 3 b ) 式定义的对角矩阵。因此 ( i 一御) 的任何对角元素都不可能等于零。说明式( 2 2 6 ) 成立。 ( 2 ) 当r = p ”时，对角元素( i 一卿) 的第i 行元素全部为零，则其秩降低l 。因此： r a n k ( r ，一r r 。) = r a n k ( i 一，p ) = q 一1 ( 2 2 7 ) 由( 2 2 5 ) ( 2 2 6 ) 式可知，= p 4 是矩阵对缸，r m j 的广义特征值。又由于 矩阵对的两个子阵列所张成的子空间是同一子空间，所以对应于两个矩阵公共零 空间的广义特征值将为零。即有g 个位于单位圆上并等于旋转符9 的对角元素，其 余尸一口个广义特征值为零。因此( 2 2 4 ) 式成立。这就是基于旋转不变的互谱估 计法。 2 1 3 3 算法基本过程 下面介绍拓广互谱e s p r i t 算法和t l s - - e s p r i t 算法的基本过程： 第9 页 算法一：拓广的互谱e s p r i t 算法： 1 ) 利用已知信息( 2 1 6 ) ( 2 1 7 ) ( 2 1 8 ) 式求互相关函数、，及其互 相关函数矩阵r 。，r 。 2 ) 计算r ，的特征分解式，r ，= f i g p f 及r ，与r 。之间的关系推导a 算 出r m = f e 7 9 q p f 。 3 ) 构造矩阵对僻，r 。 并求矩阵的广义特征分解，位于单位圆上的q 个广 义特征值确定子空问旋转符9 ，由此确定谐波频率略，剩余( p g ) 个广义特征值 为零。 拓广的互e s p r i t 方法可看做是一个p x p 算予的最小二乘估计子。其作用是 将原观测空间约束到一个q 维子空间。文献【1 0 】指出最小二乘算子会导致求广义特 征问题时的某些潜在的数值困难。而另一种e s p r i t 思想是采用总体最小二乘；能 将病态的p x p 维广义特征值问题转化为一个无病态的q x q 维广义特征问题1 。 算法二：t l s 一拓广互谱e s p r i t 算法： 1 ) 利用( 2 1 9 ) ( 2 2 2 ) 两式算出矩阵r 。，、r 。 z ) 应用r ，= u zv - - u 。u ：( 弓o ：业v j 确定r ，的有效秩，并存 储q 个主奇异值及与q 个主奇异值对应的左右主奇异矢量u 。、v l 。 3 ) 用u ，左乘和v l 右乘矩阵束 r ，r 。j 即： u f r ，v l 一声，r 料v l = 一声f r 科v 】。计算矩阵束的广义特征值互。 4 ) 估计谐波频率q = 一 i m 喘。( z i ) 5 ) 利用e ：p = ( f f ) - f r 。f ( f f ) - 1 可求出q 屈及易。 后面仿真实验中，将采用算法二进行计算。 2 1 4i v i u s i c 谱估计技术 m u s i c 是多重信号分类( m u l t i p l es i g n a lc l a s s i f i c a t i o n ) 的英文缩写，近年受到 人们广泛重视的一种优秀的高分辨率算法。这种方法是由r 0 s c h m i d t 于1 9 7 9 年 提出来的【2 8 j ，并公开发表，由于该文章第一次是在一次小规模的谱估计的专题讨 论会上发表的，而该会议的论文集的发行范围非常有限。所以在1 9 8 6 年i e e e t r a n s a c t i o no na n t e n n aa n dp r o p a g a t i o n 的自适应天线系统专刊上又重新发表。 m u s i c 算法1 1 2 1 利用了信号子空问和噪声子空间的正交性，构造空间谱函数，通过 第1 0 页 国防科学技术大学研究生院学位论文 谱峰搜索，检测到信号频谱最大值。它属于特征结构子空间方法，子空间方法建 立在这样一个观察之上。 2 1 4 1 基本定义 我们考虑( 2 1 6 ) ( 2 1 7 ) 中提到的信号模型。其中各符号代表的意思同上。 此模型是比较有代表性的；在雷达、声纳、机械振动分析等许多问题中的信号大 多可以用此信号模型描述。在很多问题中，有时仅能得到单一的正弦信号观测序 列，这时可根据实际情况采用参考信号法( 如锁定放大器) 、双通道隔离放大器 法、时间延迟法【 3 1 均可将原序列等效地用式( 2 1 6 ) ( 2 1 7 ) 来描述。注意到以上 方法比较容易在工程上实现，因此本方法中才用信号仍使用上节模型是有一定普 遍意义的。 由式( 2 1 6 ) ( 2 1 7 ) 定义的信号的互相关函数为( 2 1 9 ) 所示即： 勺( 肌) = e k 0 涉0 + 肌) 】= 口，属p 鸬p 4 4 ( 2 2 8 ) 扭l 将) 代入p p 维q ) 的扩阶互相关函数矩阵： r i y = ( o )勺( 一1 )卜p + 1 ) r y 0 ) ( o )( p + 2 ) 勺( p 一1 ) ( p 一2 ) ( o ) 司得： r q = f e # p f h 式中： f j = 【l ，p 确，e j “，e j ( h 】。 f = 【f l ，f 2 ，l j 为一p q 维复矩阵。 e ，f = 幽曙k m ，e m ，p 填j ，为一复指数对角矩阵。 p = 讲昭k 。届，厦，口i 岛j 为一实对角矩阵。显然r ，的秩为q ， 对称矩阵，因此r 。的奇异值分解式为： r ，= u 陌妒 ( 2 2 9 ) ( 2 3 0 ) 且为非共轭 式中：矩阵u 、v 的列向量分别为矩阵的左右奇异矢量，他们由下式定义： r 掣r p “i = 露2 “ 1 j s p r 吖胃r 廿v i = v 1 | i s p 而= 幽昭b ，正，乞】为矩阵r ，的非零奇异值矩阵。 ( 2 3 1 ) 第1 1 页 国防科学技术大学研究生院学位论文 v 2 是v 中后q q ) 个奇异矢量所构成，即： ii v 1 = l v l ，v z ，j 而v 2 = h 。，v ，j 用同样方法将左奇异矢量矩阵分块成u = 【u 。，u ：】，其中： ( 2 3 2 ) ( 2 3 3 ) u l = 扣l ，u 2 ，u gj ( 2 3 4 ) 而u ：= k 州，。，j ( 2 3 5 ) 将式( 2 3 2 ) ( 2 3 5 ) 代入式( 2 31 ) ，得： r ，= u i ，u 2 样捌 亿s s ， 则有： r ，= u 。v 1 ( 2 3 7 和r 。v 2 = 0 ，r 。v 2 = 0 ( 2 3 8 ) 为以后叙述的方便，称f 为信号矢量空间，称u ，、v ，为噪声特征矢量空间。 定理