CN115390007A - 基于张量链分解的稀疏立方阵列波达方向估计方法 - Google Patents

基于张量链分解的稀疏立方阵列波达方向估计方法 Download PDF

Info

Publication number
CN115390007A
CN115390007A CN202211017720.XA CN202211017720A CN115390007A CN 115390007 A CN115390007 A CN 115390007A CN 202211017720 A CN202211017720 A CN 202211017720A CN 115390007 A CN115390007 A CN 115390007A
Authority
CN
China
Prior art keywords
tensor
matrix
array
dimensional
cubic
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Pending
Application number
CN202211017720.XA
Other languages
English (en)
Inventor
徐子骏
郑航
史治国
周成伟
陈积明
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Zhejiang University ZJU
Original Assignee
Zhejiang University ZJU
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Zhejiang University ZJU filed Critical Zhejiang University ZJU
Priority to CN202211017720.XA priority Critical patent/CN115390007A/zh
Publication of CN115390007A publication Critical patent/CN115390007A/zh
Pending legal-status Critical Current

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S3/00—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received
    • G01S3/02—Direction-finders for determining the direction from which infrasonic, sonic, ultrasonic, or electromagnetic waves, or particle emission, not having a directional significance, are being received using radio waves
    • G01S3/14—Systems for determining direction or deviation from predetermined direction
    • G01S3/16—Systems for determining direction or deviation from predetermined direction using amplitude comparison of signals derived sequentially from receiving antennas or antenna systems having differently-oriented directivity characteristics or from an antenna system having periodically-varied orientation of directivity characteristic
    • G01S3/22—Systems for determining direction or deviation from predetermined direction using amplitude comparison of signals derived sequentially from receiving antennas or antenna systems having differently-oriented directivity characteristics or from an antenna system having periodically-varied orientation of directivity characteristic derived from different combinations of signals from separate antennas, e.g. comparing sum with difference

Landscapes

  • Physics & Mathematics (AREA)
  • Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Remote Sensing (AREA)
  • Measurement Of Velocity Or Position Using Acoustic Or Ultrasonic Waves (AREA)

Abstract

本发明公开了一种基于张量链分解的稀疏立方阵列波达方向估计方法，主要解决现有方法中通过canonical polyadic(CP)分解直接处理高维虚拟域协方差张量时计算复杂度高的问题，其实现步骤是：构建三维嵌套互质立方阵列；嵌套互质立方阵列接收信号的张量建模；通过互相关张量变换构造对应于虚拟均匀立方阵列的虚拟域张量；推导虚拟域协方差张量；对虚拟域协方差张量进行张量链分解；从张量链首矩阵和核张量中提取CP因子矩阵并获取波达方向估计结果。本发明通过张量链分解对高维虚拟域协方差张量进行降维表征，进而从低维度矩阵和核张量中高效提取角度信息，实现了低复杂度、高精度的波达方向估计，可用于目标定位。

Description

基于张量链分解的稀疏立方阵列波达方向估计方法
技术领域
本发明属于阵列信号处理技术领域，尤其涉及基于稀疏阵列虚拟域张量的统计信号处理技术，具体是一种基于张量链分解的稀疏立方阵列波达方向估计方法，可用于目标定位。
背景技术
相比于传统遵循奈奎斯特采样速率的均匀阵列，稀疏阵列具有大孔径、高分辨率的优势，能够突破传统均匀阵列波达方向估计在估计性能与成本开销上的性能瓶颈。具有系统化结构的稀疏阵列基于接收信号的二阶统计量推导增广虚拟阵列，从而基于虚拟域二阶等价信号实现奈奎斯特匹配的波达方向估计。嵌套互质阵列作为一种典型的稀疏阵列架构，能够推导出连续的虚拟阵列，因此得到广泛应用。然而，传统的虚拟域信号处理方法将接收信号表示成矢量，并通过矢量化接收信号的协方差矩阵推导虚拟域二阶等价信号；随着实际应用中稀疏阵列维度的不断扩张，多维稀疏阵列的接收信号涵盖多维度时空信息，而这种矢量化信号的处理方法将破坏多维接收信号的原始结构，造成严重的性能损失。
为了保留多维接收信号的结构化信息，张量作为一种多维的数据类型，开始被应用于阵列信号处理领域，用于表征涵盖多维时空信息的接收信号。现有面向多维稀疏阵列的张量信号处理方法通过推导二阶虚拟域张量，将对应的高维虚拟域协方差张量表示为canonical polyadic(CP)模型，对其进行直接canonical polyadic分解以实现高精度、高分辨的波达方向估计。然而，canonical polyadic分解在处理高维张量的时候，存在收敛速度慢、计算复杂度高的问题。而张量链分解作为另外一种张量分解形式，能够将高维张量高效分解为若干低维度核张量和矩阵的链式连接形式，进而针对低维核张量做后续处理，能够在保障计算效率的前提下有效提取信号特征，从而实现信道估计、谐波检索、高光谱图像去噪等具体功能。由此可见，张量链分解在提升稀疏阵列波达方向估计的处理效率方面具有极大潜力。针对高维虚拟域协方差张量，如何运用张量链分解设计一种低计算复杂度的波达方向估计方法仍是一个亟待解决的问题。
发明内容
本发明的目的在于针对现有技术的不足，提供一种基于张量链分解的稀疏立方阵列波达方向估计方法。
本发明的目的是通过以下技术方案实现的：一种基于张量链分解的稀疏立方阵列波达方向估计方法，包含以下步骤：
(1)接收端使用MxMyMz+NxNyNz-1个物理天线阵元，按照三维嵌套互质立方阵列的结构进行架构；该嵌套互质立方阵列分解为一个均匀立方子阵列
Figure BDA0003811903610000021
和一个稀疏立方子阵列
Figure BDA0003811903610000022
其中
Figure BDA0003811903610000023
包含Mx×My×Mz个天线阵元，阵元间距d为入射窄带信号波长λ的一半，即d＝λ/2，
Figure BDA0003811903610000024
包含Nx×Ny×Nz个天线阵元，在x轴方向、y轴方向和z轴方向上的间距分别为Mxd、Myd和Mzd；
(2)假设有K个来自
Figure BDA0003811903610000025
方向的远场窄带非相关信号源，θk和
Figure BDA0003811903610000026
分别为第k个入射信号源的方位角和俯仰角，k＝1，2，…，K，将嵌套互质立方阵列中均匀立方子阵列
Figure BDA0003811903610000027
的T个采样快拍信号在第四个维度进行叠加后，得到一个四维张量信号
Figure BDA0003811903610000028
建模为：
Figure BDA0003811903610000029
其中，sk＝[sk，1，sk，2，…，sk，T]T为对应第k个入射信号源的多快拍采样信号波形向量，[·]T表示转置操作，
Figure BDA00038119036100000216
表示矢量外积，
Figure BDA00038119036100000210
为与各信号源相互独立的噪声张量，
Figure BDA00038119036100000211
和
Figure BDA00038119036100000212
分别为
Figure BDA00038119036100000213
在x轴、y轴和z轴方向上的导引矢量，表示为：
Figure BDA00038119036100000214
Figure BDA00038119036100000215
Figure BDA0003811903610000031
其中，
Figure BDA0003811903610000032
稀疏立方子阵列
Figure BDA0003811903610000033
的接收信号用四维张量信号
Figure BDA0003811903610000034
表示为：
Figure BDA0003811903610000035
其中，
Figure BDA0003811903610000036
为与各信号源相互独立的噪声张量，
Figure BDA0003811903610000037
和
Figure BDA0003811903610000038
分别为
Figure BDA0003811903610000039
在x轴、y轴和z轴方向上的导引矢量，表示为：
Figure BDA00038119036100000310
Figure BDA00038119036100000311
Figure BDA00038119036100000312
通过求四维张量信号
Figure BDA00038119036100000313
和
Figure BDA00038119036100000314
的互相关统计量，得到二阶互相关张量
Figure BDA00038119036100000315
Figure BDA00038119036100000316
其中，
Figure BDA00038119036100000317
表示第k个入射信号源的功率，
Figure BDA00038119036100000318
表示六维互相关噪声张量，<·，·>r表示两个张量沿着第r维度的张量缩并操作，E[·]表示取数学期望操作，(·)*表示共轭操作；这里，六维互相关噪声张量
Figure BDA00038119036100000319
的第(1，1，1，1，1，1)个元素取值为噪声功率
Figure BDA00038119036100000320
而其他元素取值均为0；
(3)定义维度集合
Figure BDA00038119036100000321
通过对互相关张量
Figure BDA00038119036100000322
进行维度合并的张量变换，得到一个三维张量
Figure BDA00038119036100000323
Figure BDA00038119036100000324
其中，
Figure BDA00038119036100000325
Figure BDA00038119036100000326
通过在指数项上形成差集数组，构造出一个大小为MxNx×MyNy×MzNz的虚拟均匀立方阵列
Figure BDA0003811903610000041
Figure BDA0003811903610000042
表示克罗内克积；噪声张量
Figure BDA0003811903610000043
的第(1，1，1)个元素取值为
Figure BDA0003811903610000044
而其他元素取值均为0；将三维张量
Figure BDA0003811903610000045
中的元素进行排列，以对应
Figure BDA0003811903610000046
中虚拟阵元的位置，得到三维虚拟域张量
Figure BDA0003811903610000047
Figure BDA0003811903610000048
其中，
Figure BDA0003811903610000049
Figure BDA00038119036100000410
和
Figure BDA00038119036100000411
分别为虚拟均匀立方阵列
Figure BDA00038119036100000412
沿x轴、y轴和z轴上的导引矢量，噪声张量
Figure BDA00038119036100000413
的第(MxNx-Mx+1，MyNy-My+1，MzNz-Mz+1)个元素取值为
Figure BDA00038119036100000414
而其他元素取值均为0；
(4)计算对应于虚拟域张量
Figure BDA00038119036100000415
的虚拟域协方差张量
Figure BDA00038119036100000416
Figure BDA00038119036100000417
表示为如下canonical polyadic模型：
Figure BDA00038119036100000418
其中，对角张量
Figure BDA00038119036100000419
的对角线元素为
Figure BDA00038119036100000420
Bx＝
Figure BDA00038119036100000421
Figure BDA00038119036100000422
为canonical polyadic因子矩阵，
Figure BDA00038119036100000423
为噪声张量，×i表示沿着第i维度的张量-矩阵内积；
(5)通过张量链分解，将虚拟域协方差张量
Figure BDA00038119036100000424
表示为矩阵与三维核张量的链式连接形式，即：
Figure BDA00038119036100000425
其中，
Figure BDA00038119036100000426
为首矩阵，
Figure BDA00038119036100000427
Figure BDA00038119036100000428
为四个三维核张量，
Figure BDA00038119036100000429
为尾矩阵，
Figure BDA00038119036100000430
表示沿着张量
Figure BDA00038119036100000431
第r1维度与张量
Figure BDA00038119036100000432
第r2维度的张量缩并操作；
(6)对核张量
Figure BDA00038119036100000433
进行canonical polyadic分解，表示为：
Figure BDA00038119036100000434
其中，对角张量
Figure BDA00038119036100000435
的对角线元素为幅度因子[λ1，λ2，…，λK]，
Figure BDA0003811903610000051
为对应虚拟域协方差张量
Figure BDA0003811903610000052
第一和第二维度的基变换矩阵，
Figure BDA0003811903610000053
为canonical polyadic因子矩阵By的估计值；由此，
Figure BDA0003811903610000054
的第一个canonicalpolyadic因子矩阵Bx通过张量链首矩阵G1的基变换估计得到，即：
Figure BDA0003811903610000055
其中(·)-1表示矩阵求逆操作；从估计的canonical polyadic因子矩阵
Figure BDA0003811903610000056
和
Figure BDA0003811903610000057
中提取角度信息，得到二维波达方向估计结果
Figure BDA0003811903610000058
进一步地，步骤(1)所述的三维嵌套互质立方阵列结构具体描述为：均匀立方子阵列
Figure BDA0003811903610000059
中天线阵元在三维坐标系上的位置坐标为{(mxd，myd，mzd)|mx＝0，1，…，Mx-1，my＝0，1，…，My-1，mz＝0，1，…，Mz-1}；稀疏立方子阵列
Figure BDA00038119036100000510
中天线阵元在三维坐标系上的位置坐标为{(nxMxd，nyMyd，nzMzd)|nx＝0，1，…，Nx-1，ny＝0，1，…，Ny-1，nz＝0，1，…，Nz-1}；{Mx，Nx}、{My，Ny}、{Mz，Nz}分别为一对互质整数；将
Figure BDA00038119036100000525
和
Figure BDA00038119036100000526
按照坐标系原点位置阵元重叠的方式进行子阵列组合，由于
Figure BDA00038119036100000527
和
Figure BDA00038119036100000528
的阵元排布满足互质数条件，因此除坐标系原点位置以外的其余位置阵元均不重叠，获得实际包含MxMyMz+NxNyNz-1个天线阵元的嵌套互质立方阵列。
进一步地，步骤(2)所述的二阶互相关张量
Figure BDA00038119036100000511
推导，在实际中，
Figure BDA00038119036100000512
通过计算四维张量信号
Figure BDA00038119036100000513
和
Figure BDA00038119036100000514
的采样互相关统计量近似得到，即采样互相关张量
Figure BDA00038119036100000515
Figure BDA00038119036100000516
进一步地，步骤(5)中，虚拟域协方差张量的张量链分解按如下步骤进行计算：沿着第一维度对
Figure BDA00038119036100000517
进行展开，得到
Figure BDA00038119036100000518
对
Figure BDA00038119036100000519
进行截断奇异值分解：
Figure BDA00038119036100000520
其中，
Figure BDA00038119036100000521
为奇异值矩阵，
Figure BDA00038119036100000522
为右奇异值矩阵；将
Figure BDA00038119036100000523
重构为一个辅助矩阵
Figure BDA00038119036100000524
其中reshape{·}表示将一个矩阵按指定维度大小重构为另一矩阵或张量的操作；通过对C1使用截断奇异值分解得到：
C1＝U2Λ2V2，
其中，
Figure BDA0003811903610000061
为左奇异值矩阵，
Figure BDA0003811903610000062
为奇异值矩阵，
Figure BDA0003811903610000063
Figure BDA0003811903610000064
为右奇异值矩阵，通过重构U2得到核张量
Figure BDA0003811903610000065
即
Figure BDA0003811903610000066
重复上述步骤，核张量
Figure BDA0003811903610000067
和尾矩阵G6利用左奇异值矩阵{Ur，r＝3，4，5，6}依次对应生成。
进一步地，步骤(6)中，从估计的canonical polyadic因子矩阵
Figure BDA0003811903610000068
和
Figure BDA0003811903610000069
中提取参数μk和vk，表示为：
Figure BDA00038119036100000610
Figure BDA00038119036100000611
其中，
Figure BDA00038119036100000612
与
Figure BDA00038119036100000613
分别表示
Figure BDA00038119036100000614
的第m行元素与
Figure BDA00038119036100000615
的第n行元素，∠(·)表示对一个复数取幅角操作；基于(μk，vk)与二维波达方向
Figure BDA00038119036100000616
的对应关系，得到二维波达方向估计
Figure BDA00038119036100000617
的闭式解为：
Figure BDA00038119036100000618
Figure BDA00038119036100000619
本发明的有益效果是：本发明主要解决现有方法中通过canonical polyadic(CP)分解直接处理高维虚拟域协方差张量时计算复杂度高的问题，其实现步骤是：构建三维嵌套互质立方阵列；嵌套互质立方阵列接收信号的张量建模；通过互相关张量变换构造对应于虚拟均匀立方阵列的虚拟域张量；推导虚拟域协方差张量；对虚拟域协方差张量进行张量链分解；从张量链首矩阵和核张量中提取CP因子矩阵并获取波达方向估计结果。本发明通过张量链分解对高维虚拟域协方差张量进行降维表征，进而从低维度矩阵和核张量中高效提取角度信息，实现了低复杂度、高精度的波达方向估计，可用于目标定位。
附图说明
图1是本发明的总体流程框图。
图2是本发明所构建的嵌套互质立方阵列结构示意图。
图3是本发明所提方法在不同信噪比条件下的波达方向估计精度性能比较图。
图4是本发明所提方法在不同采样快拍数条件下的波达方向估计精度性能比较图。
图5是本发明所提方法在不同信噪比条件下的运行时间比较图。
图6是本发明所提方法在不同采样快拍数条件下的运行时间比较图。
具体实施方式
以下参照附图，对本发明的技术方案作进一步的详细说明。
为了解决现有方法通过canonical polyadic分解直接处理高维虚拟域协方差张量时计算复杂度高的问题，本发明提出了一种基于张量链分解的稀疏立方阵列波达方向估计方法，通过对虚拟域协方差张量进行张量链分解，从降维张量链首矩阵和核张量中提取canonical polyadic因子矩阵，实现高效的波达方向估计。参照图1，本发明的实现步骤如下：
步骤1：构建三维嵌套互质立方阵列。接收端使用MxMyMz+NxNyNz-1个物理天线阵元，按照三维嵌套互质立方阵列的结构进行架构，如图2所示，该嵌套互质立方阵列分解为一个均匀立方子阵列
Figure BDA0003811903610000075
和一个稀疏立方子阵列
Figure BDA0003811903610000071
其中
Figure BDA0003811903610000072
包含Mx×My×Mz个天线阵元，阵元间距d为入射窄带信号波长λ的一半，即d＝λ/2，
Figure BDA0003811903610000073
中天线阵元在三维坐标系上的位置坐标为{(mxd，myd，mzd)|mx＝0，1，…，Mx-1，my＝0，1，…，My-1，mz＝0，1，…，Mz-1}；
Figure BDA0003811903610000074
包含Nx×Ny×Nz个天线阵元，这些天线阵元在x轴方向、y轴方向和z轴方向上的间距分别为Mxd、Myd和Mzd，在三维坐标系上的位置坐标为{(nxMxd，nyMyd，nzMzd)|nx＝0，1，…，Nx-1，ny＝0，1，…，Ny-1，nz＝0，1，…，Nz-1}；{Mx，Nx}、{My，Ny}、{Mz，Nz}分别为一对互质整数；将
Figure BDA0003811903610000081
和
Figure BDA0003811903610000082
按照坐标系原点位置阵元重叠的方式进行子阵列组合，由于
Figure BDA0003811903610000083
和
Figure BDA0003811903610000084
的阵元排布满足互质数条件，因此除坐标系原点位置以外的其余位置阵元均不重叠，获得实际包含MxMyMz+NxNyNz-1个天线阵元的嵌套互质立方阵列；
步骤2：嵌套互质立方阵列接收信号的张量建模。假设有K个来自
Figure BDA0003811903610000085
方向的远场窄带非相关信号源，θk和
Figure BDA0003811903610000086
分别为第k个入射信号源的方位角和俯仰角，k＝1，2，…，K，将嵌套互质立方阵列中均匀立方子阵列
Figure BDA0003811903610000087
的T个采样快拍信号在第四个维度进行叠加后，可以得到一个四维张量信号
Figure BDA0003811903610000088
建模为：
Figure BDA0003811903610000089
其中，sk＝[sk，1，sk，2，…，sk，T]T为对应第k个入射信号源的多快拍采样信号波形向量，[·]T表示转置操作，
Figure BDA00038119036100000810
表示矢量外积，
Figure BDA00038119036100000811
为与各信号源相互独立的噪声张量，
Figure BDA00038119036100000812
和
Figure BDA00038119036100000813
分别为
Figure BDA00038119036100000814
在x轴、y轴和z轴方向上的导引矢量，表示为：
Figure BDA00038119036100000815
Figure BDA00038119036100000816
Figure BDA00038119036100000817
其中，
Figure BDA00038119036100000818
类似地，稀疏立方子阵列
Figure BDA00038119036100000819
的接收信号可以用四维张量信号
Figure BDA00038119036100000820
Figure BDA00038119036100000821
表示为：
Figure BDA00038119036100000822
其中，
Figure BDA00038119036100000823
为与各信号源相互独立的噪声张量，
Figure BDA00038119036100000824
和
Figure BDA00038119036100000825
分别为
Figure BDA00038119036100000826
在x轴、y轴和z轴方向上的导引矢量，表示为：
Figure BDA00038119036100000827
Figure BDA0003811903610000091
Figure BDA0003811903610000092
通过求四维张量信号
Figure BDA0003811903610000093
和
Figure BDA0003811903610000094
的互相关统计量，得到二阶互相关张量
Figure BDA0003811903610000095
Figure BDA0003811903610000096
其中，
Figure BDA0003811903610000097
表示第k个入射信号源的功率，
Figure BDA0003811903610000098
表示六维互相关噪声张量，<·，·>r表示两个张量沿着第r维度的张量缩并操作，E[·]表示取数学期望操作，(·)*表示共轭操作。这里，六维互相关噪声张量
Figure BDA0003811903610000099
的第(1，1，1，1，1，1)个元素取值为噪声功率
Figure BDA00038119036100000910
而其他元素取值均为0。在实际中，
Figure BDA00038119036100000911
通过计算张量信号
Figure BDA00038119036100000912
和
Figure BDA00038119036100000913
的采样互相关统计量近似得到，即采样互相关张量
Figure BDA00038119036100000914
Figure BDA00038119036100000915
步骤3：通过互相关张量变换构造对应于虚拟均匀立方阵列的虚拟域张量。由于互相关张量
Figure BDA00038119036100000916
中包含了均匀立方子阵列
Figure BDA00038119036100000917
和稀疏立方子阵列
Figure BDA00038119036100000918
的空间信息，通过合并
Figure BDA00038119036100000919
中表征同一方向空间信息的维度，可以使两个子阵列对应相同方向的导引矢量在指数项上形成差集数组，从而构造三维的虚拟均匀立方阵列。具体地，互相关张量
Figure BDA00038119036100000920
的第1、4维度表征x轴方向的空间信息，第2、5维度表征y轴方向的空间信息，第3、6维度表征z轴方向的空间信息；为此，定义维度集合
Figure BDA00038119036100000921
通过对互相关张量
Figure BDA00038119036100000922
进行维度合并的张量变换，得到一个三维张量
Figure BDA00038119036100000923
Figure BDA00038119036100000924
其中，
Figure BDA00038119036100000925
Figure BDA0003811903610000101
通过在指数项上形成差集数组，构造出一个大小为MxNx×MyNy×MzNz的虚拟均匀立方阵列
Figure BDA0003811903610000102
Figure BDA0003811903610000103
表示Kronecker积。噪声张量
Figure BDA0003811903610000104
的第(1，1，1)个元素取值为
Figure BDA0003811903610000105
而其他元素取值均为0。将三维张量
Figure BDA0003811903610000106
中的元素进行排列，以对应
Figure BDA0003811903610000107
中虚拟阵元的位置，即可得到三维虚拟域张量
Figure BDA0003811903610000108
Figure BDA0003811903610000109
Figure BDA00038119036100001010
其中，
Figure BDA00038119036100001011
Figure BDA00038119036100001012
和
Figure BDA00038119036100001013
分别为虚拟均匀立方阵列
Figure BDA00038119036100001014
沿x轴、y轴和z轴上的导引矢量，噪声张量
Figure BDA00038119036100001015
的第(MxNx-Mx+1，MyNy-My+1，MzNz-Mz+1)个元素取值为
Figure BDA00038119036100001016
而其他元素取值均为0；
步骤4：推导虚拟域协方差张量。计算对应于虚拟域张量
Figure BDA00038119036100001017
的虚拟域协方差张量
Figure BDA00038119036100001018
表示为如下canonicalpolyadic模型：
Figure BDA00038119036100001019
其中，对角张量
Figure BDA00038119036100001020
的对角线元素为
Figure BDA00038119036100001021
Figure BDA00038119036100001022
Figure BDA00038119036100001023
Figure BDA00038119036100001024
为canonical polyadic因子矩阵，
Figure BDA00038119036100001025
为噪声张量，×i表示沿着第i维度的张量-矩阵内积；
步骤5：对虚拟域协方差张量进行张量链分解。为了确保对虚拟域协方差张量
Figure BDA00038119036100001026
中的canonical polyadic因子进行高效提取，不同于传统方法对其直接进行canonicalpolyadic分解的操作，通过张量链分解对其进行降维处理。具体而言，将虚拟域协方差张量
Figure BDA00038119036100001027
表示为矩阵与三维核张量的链式连接形式，即：
Figure BDA00038119036100001028
其中，
Figure BDA00038119036100001029
为首矩阵，
Figure BDA00038119036100001030
Figure BDA0003811903610000111
为四个三维核张量，
Figure BDA0003811903610000112
为尾矩阵，
Figure BDA0003811903610000113
表示沿着张量
Figure BDA0003811903610000114
第r1维度与张量
Figure BDA0003811903610000115
第r2维度的张量缩并操作。上述虚拟域协方差张量的张量链分解可按如下步骤进行计算：沿着第一维度对
Figure BDA0003811903610000116
进行展开，得到
Figure BDA0003811903610000117
进而对
Figure BDA0003811903610000118
进行截断奇异值分解(Truncatedsingular value decomposition)：
Figure BDA0003811903610000119
其中，
Figure BDA00038119036100001110
为奇异值矩阵，
Figure BDA00038119036100001111
为右奇异值矩阵。将
Figure BDA00038119036100001112
重构为一个辅助矩阵
Figure BDA00038119036100001113
其中reshape{·}表示将一个矩阵按指定维度大小重构为另一矩阵或张量的操作。通过对C1使用截断奇异值分解可得到：
C1＝U2Λ2V2，
其中，
Figure BDA00038119036100001114
为左奇异值矩阵，
Figure BDA00038119036100001115
为奇异值矩阵，
Figure BDA00038119036100001116
Figure BDA00038119036100001117
为右奇异值矩阵，通过重构U2可得到核张量
Figure BDA00038119036100001118
即
Figure BDA00038119036100001119
重复上述步骤，核张量
Figure BDA00038119036100001120
和尾矩阵G6可利用左奇异值矩阵{Ur，r＝3，4，5，6}依次对应生成；
步骤6：从张量链首矩阵和核张量中提取CP因子矩阵并获取波达方向估计结果。由于canonical polyadic因子矩阵Bx与By可用于估计信号源的二维波达方向，因此将虚拟域协方差张量
Figure BDA00038119036100001121
通过张量链分解得到的首矩阵G1和核张量
Figure BDA00038119036100001122
转换为Bx与By。具体过程为：对核张量
Figure BDA00038119036100001123
进行canonical polyadic分解，表示为：
Figure BDA00038119036100001124
其中，对角张量
Figure BDA00038119036100001125
的对角线元素为幅度因子[λ1，λ2，…，λK]，
Figure BDA00038119036100001126
为对应虚拟域协方差张量
Figure BDA00038119036100001127
第一和第二维度的基变换矩阵，
Figure BDA00038119036100001128
为canonical polyadic因子矩阵By的估计值。由此，
Figure BDA00038119036100001129
的第一个canonical polyadic因子矩阵Bx可以通过张量链首矩阵G1的基变换得到，即：
Figure BDA00038119036100001130
其中(·)-1表示矩阵求逆操作。进而，从
Figure BDA0003811903610000121
和
Figure BDA0003811903610000122
中估计μk和νk，表示为：
Figure BDA0003811903610000123
Figure BDA0003811903610000124
其中，
Figure BDA0003811903610000125
与
Figure BDA0003811903610000126
分别表示
Figure BDA0003811903610000127
的第m行元素与
Figure BDA0003811903610000128
的第n行元素，∠(·)表示对一个复数取幅角操作。基于(μk，vk)与二维波达方向
Figure BDA0003811903610000129
的对应关系，得到二维波达方向估计
Figure BDA00038119036100001210
的闭式解为：
Figure BDA00038119036100001211
Figure BDA00038119036100001212
下面结合仿真实例对本发明的效果做进一步的描述。
仿真实例：采用嵌套互质立方阵列接收入射信号，其参数选取为Mx＝My＝Mz＝3，Nx＝Ny＝Nz＝4，即架构的嵌套互质立方阵列共包含MxMyMz+NxNyNz-1＝90个物理阵元。所构造的虚拟均匀立方阵列共有12×12×12个虚拟阵元。将本发明所提基于张量链分解的稀疏立方阵列波达方向估计方法与传统基于canonical polyadic分解的方法和基于旋转不变子空间(estimating signal parameter via rotational invariance techniques，ESPRIT)的方法进行对比。假定有2个入射信号，每次试验的信号方位角和俯仰角均在[15°，60°]范围内随机生成。在采样快拍数T＝150的条件下，绘制均方根误差(root-mean-square error，RMSE)随信噪比(signal-to-noise ratio，SNR)变化的性能对比曲线，如图3所示；在信噪比SNR＝-5dB的条件下，绘制RMSE随采样快拍数T变化的性能对比曲线，如图4所示。进一步地，对比上述方法的计算效率，绘制计算机运行时间随SNR变化的对比曲线，如图5所示；并绘制计算机运行时间随采样快拍数T的对比曲线，如图6所示。从对比结果看出，相比于传统基于ESPRIT的方法，本发明所提方法通过构建虚拟域张量，充分利用所对应高维虚拟域协方差张量的结构化信息，具备更优的波达方向估计性能，且基于张量链分解对虚拟域协方差张量进行高效特征提取，具备更短的运行时间；相比基于canonical polyadic分解的方法，本发明所提方法在虚拟域协方差张量的张量链分解过程中最大程度降低信息损失，不仅表现出更高的估计精度，还使得计算复杂度大幅度降低，显著提升了计算效率。
综上所述，本发明通过对高维虚拟域协方差张量进行张量链分解，并设计了相应canonical polyadic因子的高效提取策略，实现了低复杂度、高精度的稀疏立方阵列波达方向估计。
以上所述仅是本发明的优选实施方式，虽然本发明已以较佳实施例披露如上，然而并非用以限定本发明。任何熟悉本领域的技术人员，在不脱离本发明技术方案范围情况下，都可利用上述揭示的方法和技术内容对本发明技术方案做出许多可能的变动和修饰，或修改为等同变化的等效实施例。因此，凡是未脱离本发明技术方案的内容，依据本发明的技术实质对以上实施例所做的任何的简单修改、等同变化及修饰，均仍属于本发明技术方案保护的范围内。

Claims (5)

1.一种基于张量链分解的稀疏立方阵列波达方向估计方法，其特征在于，包含以下步骤：
(1)接收端使用MxMyMz+NxNyNz-1个物理天线阵元，按照三维嵌套互质立方阵列的结构进行架构；该嵌套互质立方阵列分解为一个均匀立方子阵列
Figure FDA0003811903600000011
和一个稀疏立方子阵列
Figure FDA0003811903600000012
其中
Figure FDA0003811903600000013
包含Mx×My×Mz个天线阵元，阵元间距d为入射窄带信号波长λ的一半，即d＝λ/2，
Figure FDA00038119036000000122
包含Nx×Ny×Nz个天线阵元，在x轴方向、y轴方向和z轴方向上的间距分别为Mxd、Myd和Mzd；
(2)假设有K个来自
Figure FDA0003811903600000014
方向的远场窄带非相关信号源，θk和
Figure FDA0003811903600000015
分别为第k个入射信号源的方位角和俯仰角，k＝1,2,…,K，将嵌套互质立方阵列中均匀立方子阵列
Figure FDA0003811903600000016
的T个采样快拍信号在第四个维度进行叠加后，得到一个四维张量信号
Figure FDA0003811903600000017
建模为：
Figure FDA0003811903600000018
其中，sk＝[sk,1,sk,2,…,sk,T]T为对应第k个入射信号源的多快拍采样信号波形向量，
Figure FDA0003811903600000019
表示转置操作，
Figure FDA00038119036000000110
表示矢量外积，
Figure FDA00038119036000000111
为与各信号源相互独立的噪声张量，
Figure FDA00038119036000000112
和
Figure FDA00038119036000000113
分别为
Figure FDA00038119036000000114
在x轴、y轴和z轴方向上的导引矢量，表示为：
Figure FDA00038119036000000115
Figure FDA00038119036000000116
Figure FDA00038119036000000117
其中，
Figure FDA00038119036000000118
稀疏立方子阵列
Figure FDA00038119036000000119
的接收信号用四维张量信号
Figure FDA00038119036000000120
表示为：
Figure FDA00038119036000000121
其中，
Figure FDA0003811903600000021
为与各信号源相互独立的噪声张量，
Figure FDA0003811903600000022
和
Figure FDA0003811903600000023
分别为
Figure FDA0003811903600000024
在x轴、y轴和z轴方向上的导引矢量，表示为：
Figure FDA0003811903600000025
Figure FDA0003811903600000026
Figure FDA0003811903600000027
通过求四维张量信号
Figure FDA0003811903600000028
和
Figure FDA0003811903600000029
的互相关统计量，得到二阶互相关张量
Figure FDA00038119036000000210
Figure FDA00038119036000000211
其中，
Figure FDA00038119036000000212
表示第k个入射信号源的功率，
Figure FDA00038119036000000213
表示六维互相关噪声张量，<·,·>r表示两个张量沿着第r维度的张量缩并操作，E[·]表示取数学期望操作，(·)*表示共轭操作；这里，六维互相关噪声张量
Figure FDA00038119036000000214
的第(1,1,1,1,1,1)个元素取值为噪声功率
Figure FDA00038119036000000215
而其他元素取值均为0；
(3)定义维度集合
Figure FDA00038119036000000216
通过对互相关张量
Figure FDA00038119036000000217
进行维度合并的张量变换，得到一个三维张量
Figure FDA00038119036000000218
Figure FDA00038119036000000219
其中，
Figure FDA00038119036000000220
Figure FDA00038119036000000221
通过在指数项上形成差集数组，构造出一个大小为MxNx×MyNy×MzNz的虚拟均匀立方阵列
Figure FDA00038119036000000222
Figure FDA00038119036000000223
表示克罗内克积；噪声张量
Figure FDA00038119036000000224
的第(1,1,1)个元素取值为
Figure FDA00038119036000000225
而其他元素取值均为0；将三维张量
Figure FDA00038119036000000226
中的元素进行排列，以对应
Figure FDA00038119036000000227
中虚拟阵元的位置，得到三维虚拟域张量
Figure FDA00038119036000000228
Figure FDA0003811903600000031
其中，
Figure FDA0003811903600000032
Figure FDA0003811903600000033
和
Figure FDA0003811903600000034
分别为虚拟均匀立方阵列
Figure FDA0003811903600000035
沿x轴、y轴和z轴上的导引矢量，噪声张量
Figure FDA0003811903600000036
的第(MxNx-Mx+1,MyNy-My+1,MzNz-Mz+1)个元素取值为
Figure FDA0003811903600000037
而其他元素取值均为0；
(4)计算对应于虚拟域张量
Figure FDA0003811903600000038
的虚拟域协方差张量
Figure FDA0003811903600000039
Figure FDA00038119036000000310
表示为如下canonical polyadic模型：
Figure FDA00038119036000000311
其中，对角张量
Figure FDA00038119036000000312
的对角线元素为
Figure FDA00038119036000000313
Figure FDA00038119036000000314
Figure FDA00038119036000000315
Figure FDA00038119036000000316
为canonical polyadic因子矩阵，
Figure FDA00038119036000000317
为噪声张量，×i表示沿着第i维度的张量-矩阵内积；
(5)通过张量链分解，将虚拟域协方差张量
Figure FDA00038119036000000318
表示为矩阵与三维核张量的链式连接形式，即：
Figure FDA00038119036000000319
其中，
Figure FDA00038119036000000320
为首矩阵，
Figure FDA00038119036000000321
Figure FDA00038119036000000322
为四个三维核张量，
Figure FDA00038119036000000323
为尾矩阵，
Figure FDA00038119036000000324
表示沿着张量
Figure FDA00038119036000000325
第r1维度与张量
Figure FDA00038119036000000326
第r2维度的张量缩并操作；
(6)对核张量
Figure FDA00038119036000000327
进行canonical polyadic分解，表示为：
Figure FDA00038119036000000328
其中，对角张量
Figure FDA00038119036000000329
的对角线元素为幅度因子[λ1,λ2,…,λK]，
Figure FDA00038119036000000330
为对应虚拟域协方差张量
Figure FDA00038119036000000331
第一和第二维度的基变换矩阵，
Figure FDA00038119036000000332
为canonical polyadic因子矩阵By的估计值；由此，
Figure FDA00038119036000000333
的第一个canonical polyadic因子矩阵Bx通过张量链首矩阵G1的基变换估计得到，即：
Figure FDA00038119036000000334
其中(·)-1表示矩阵求逆操作；从估计的canonical polyadic因子矩阵
Figure FDA0003811903600000041
和
Figure FDA0003811903600000042
中提取角度信息，得到二维波达方向估计结果
Figure FDA0003811903600000043
2.根据权利要求1所述的基于张量链分解的稀疏立方阵列波达方向估计方法，其特征在于，步骤(1)所述的三维嵌套互质立方阵列结构具体描述为：均匀立方子阵列
Figure FDA0003811903600000044
中天线阵元在三维坐标系上的位置坐标为{(mxd,myd,mzd)|mx＝0,1,…,Mx-1,my＝0,1,…,My-1,mz＝0,1,…,Mz-1}；稀疏立方子阵列
Figure FDA0003811903600000045
中天线阵元在三维坐标系上的位置坐标为{(nxMxd,nyMyd,nzMzd)|nx＝0,1,…,Nx-1,ny＝0,1,…,Ny-1,nz＝0,1,…,Nz-1}；{Mx,Nx}、{My,Ny}、{Mz,Nz}分别为一对互质整数；将
Figure FDA0003811903600000046
和
Figure FDA0003811903600000047
按照坐标系原点位置阵元重叠的方式进行子阵列组合，由于
Figure FDA0003811903600000048
和
Figure FDA0003811903600000049
的阵元排布满足互质数条件，因此除坐标系原点位置以外的其余位置阵元均不重叠，获得实际包含MxMyMz+NxNyNz-1个天线阵元的嵌套互质立方阵列。
3.根据权利要求1所述的基于张量链分解的稀疏立方阵列波达方向估计方法，其特征在于，步骤(2)所述的二阶互相关张量
Figure FDA00038119036000000410
推导，在实际中，
Figure FDA00038119036000000411
通过计算四维张量信号
Figure FDA00038119036000000412
和
Figure FDA00038119036000000413
的采样互相关统计量近似得到，即采样互相关张量
Figure FDA00038119036000000414
Figure FDA00038119036000000415
4.根据权利要求1所述的基于张量链分解的稀疏立方阵列波达方向估计方法，其特征在于，步骤(5)中，虚拟域协方差张量的张量链分解按如下步骤进行计算：沿着第一维度对
Figure FDA00038119036000000416
进行展开，得到
Figure FDA00038119036000000417
对
Figure FDA00038119036000000418
进行截断奇异值分解：
Figure FDA00038119036000000419
其中，
Figure FDA00038119036000000420
为奇异值矩阵，
Figure FDA00038119036000000421
为右奇异值矩阵；将
Figure FDA00038119036000000422
重构为一个辅助矩阵
Figure FDA00038119036000000423
其中reshape{·}表示将一个矩阵按指定维度大小重构为另一矩阵或张量的操作；通过对C1使用截断奇异值分解得到：
C1＝U2Λ2V2,
其中，
Figure FDA0003811903600000051
为左奇异值矩阵，
Figure FDA0003811903600000052
为奇异值矩阵，
Figure FDA0003811903600000053
Figure FDA0003811903600000054
为右奇异值矩阵，通过重构U2得到核张量
Figure FDA0003811903600000055
即
Figure FDA0003811903600000056
重复上述步骤，核张量
Figure FDA0003811903600000057
和尾矩阵G6利用左奇异值矩阵{Ur,r＝3,4,5,6}依次对应生成。
5.根据权利要求1所述的基于张量链分解的稀疏立方阵列波达方向估计方法，其特征在于，步骤(6)中，从估计的canonical polyadic因子矩阵
Figure FDA0003811903600000058
和
Figure FDA0003811903600000059
中提取参数μk和vk，表示为：
Figure FDA00038119036000000510
Figure FDA00038119036000000511
其中，
Figure FDA00038119036000000512
与
Figure FDA00038119036000000513
分别表示
Figure FDA00038119036000000514
的第m行元素与
Figure FDA00038119036000000515
的第n行元素，∠(·)表示对一个复数取幅角操作；基于(μk,vk)与二维波达方向
Figure FDA00038119036000000516
的对应关系，得到二维波达方向估计
Figure FDA00038119036000000517
的闭式解为：
Figure FDA00038119036000000518
Figure FDA00038119036000000519
CN202211017720.XA 2022-08-23 2022-08-23 基于张量链分解的稀疏立方阵列波达方向估计方法 Pending CN115390007A (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211017720.XA CN115390007A (zh) 2022-08-23 2022-08-23 基于张量链分解的稀疏立方阵列波达方向估计方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211017720.XA CN115390007A (zh) 2022-08-23 2022-08-23 基于张量链分解的稀疏立方阵列波达方向估计方法

Publications (1)

Publication Number Publication Date
CN115390007A true CN115390007A (zh) 2022-11-25

Family

ID=84121547

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211017720.XA Pending CN115390007A (zh) 2022-08-23 2022-08-23 基于张量链分解的稀疏立方阵列波达方向估计方法

Country Status (1)

Country Link
CN (1) CN115390007A (zh)

Similar Documents

Publication Publication Date Title
CN108872929B (zh) 基于内插虚拟阵列协方差矩阵子空间旋转不变性的互质阵列波达方向估计方法
CN109655799B (zh) 基于iaa的协方差矩阵向量化的非均匀稀疏阵列测向方法
CN107092004B (zh) 基于信号子空间旋转不变性的互质阵列波达方向估计方法
CN109471082B (zh) 基于信号子空间重构的阵元缺损mimo雷达角度估计方法
CN107037392B (zh) 一种基于压缩感知的自由度增加型互质阵列波达方向估计方法
CN107589399B (zh) 基于多采样虚拟信号奇异值分解的互质阵列波达方向估计方法
CN110109050B (zh) 嵌套阵列下基于稀疏贝叶斯的未知互耦的doa估计方法
CN110161452B (zh) 基于互质式l型电磁矢量传感器阵列的波达方向估计方法
WO2021068496A1 (zh) 基于结构化虚拟域张量信号处理的互质面阵二维波达方向估计方法
CN111239678A (zh) 一种基于l型阵列的二维doa估计方法
CN112731278B (zh) 一种部分极化信号的角度与极化参数欠定联合估计方法
CN109917329B (zh) 一种基于协方差匹配准则的l型阵列波达方向估计方法
CN111624545A (zh) 基于结构化虚拟域张量信号处理的互质面阵二维波达方向估计方法
JP7327840B2 (ja) 相互相関テンソルに基づく三次元の互いに素のキュービックアレイの到来方向推定方法
CN109507636B (zh) 基于虚拟域信号重构的波达方向估计方法
CN110837076A (zh) 一种基于张量分解的矢量水听器阵列方位估计方法
CN112733327A (zh) 一种面向非高斯信号的连续和阵稀疏阵列及其设计方法
CN113093144A (zh) 基于采样数据矩阵重构的mimo雷达doa估计方法
CN111983554A (zh) 非均匀l阵下的高精度二维doa估计
CN114624646B (zh) 一种基于模型驱动复数神经网络的doa估计方法
CN106980105B (zh) 电磁矢量传感器阵列空间旋转解相干测向方法
CN113552532B (zh) 基于耦合张量分解的l型互质阵列波达方向估计方法
CN114200388A (zh) 基于四阶采样协方差张量去噪的子阵分置式l型互质阵列波达方向估计方法
Xiao et al. A weighted forward-backward spatial smoothing DOA estimation algorithm based on TLS-ESPRIT
CN112327292A (zh) 一种二维稀疏阵列doa估计方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination