基于三极化嵌套阵列的部分极化信号参数估计方法
阅读说明：本技术 基于三极化嵌套阵列的部分极化信号参数估计方法 (Partial polarization signal parameter estimation method based on triple polarization nested array ) 是由 潘玉剑 高晓欣 于 2021-06-28 设计创作，主要内容包括：本发明公开了一种基于三极化嵌套阵列的部分极化信号参数估计方法,主要用于解决欠定条件下部分极化信号波达方向及极化参数估计问题。其实现步骤是：根据优化法则布置三极化嵌套阵列；构建三极化阵列的部分极化信号接收模型；利用子协方差矩阵相加构造虚拟均匀线阵；利用子协方差矩阵相加构造虚拟协阵输出；利用Toeplitz方法恢复矩阵的秩；估计信号波达方向；估计噪声方差；估计信号极化参数。该方法优势：采用的三极化阵列及子协方差矩阵相加方法可累积信号所有极化分量的功率；采用了优化的嵌套阵列,可在欠定条件下估计噪声方差,从而实现极化参数的无偏估计；除了能估计信号的波达方向和极化度,还能估计极化指向角和极化椭圆率角。(The invention discloses a partial polarization signal parameter estimation method based on a tri-polarization nested array, which is mainly used for solving the problems of partial polarization signal direction of arrival and polarization parameter estimation under an underdetermined condition. The method comprises the following implementation steps: arranging a tri-polarization nested array according to an optimization rule; constructing a partial polarization signal receiving model of the triple polarization array; constructing a virtual uniform linear array by adding the sub-covariance matrixes; constructing virtual covariance matrix output by adding the sub covariance matrices; recovering the rank of the matrix by using a Toeplitz method; estimating the direction of arrival of the signal; estimating a noise variance; a signal polarization parameter is estimated. The method has the advantages that: the power of all polarization components of the signal can be accumulated by adopting a triple polarization array and sub covariance matrix addition method; by adopting the optimized nested array, the noise variance can be estimated under an underdetermined condition, so that the unbiased estimation of the polarization parameters is realized; in addition to estimating the direction of arrival and the degree of polarization of the signal, the polarization pointing angle and the polarization ellipticity angle can also be estimated.)
技术领域
本发明属于信号处理技术领域，特别涉及对部分极化信号的波达方向和极化参数的估计，具体是一种基于三极化嵌套阵列的部分极化信号参数估计方法。
背景技术
阵列信号处理是利用传感器阵列来进行高维信号处理的技术，广泛应用于雷达、通信、声呐及医学诊断等军用和民用技术领域。嵌套阵列由于具备欠定估计的能力，即估计的信号个数大于阵元数，且其虚拟协阵是均匀线阵，而受到学者的广泛关注。另外，由于常规的完全极化信号只属于部分极化信号的一种特例，研究部分极化信号的参数估计更具普遍意义。因此，本专利关注的是利用嵌套阵列进行部分极化信号的参数估计，涉及的参数包含波达方向和极化参数。
为了能够估计信号的极化参数，阵列必须包含多极化传感器，形成极化敏感阵列。文献：He J,Zhang Z,Shu T,et al.Direction Finding of Multiple PartiallyPolarized Signals With a Nested Cross-Diople Array[J].IEEE Antennas&WirelessPropagation Letters,2017,16:1679-1682提出一种基于正交子空间的方法(简称正交子空间法)来估计多个不同部分极化信号的波达方向。该方法采用的是嵌套阵，每个阵元为双极化阵元。但该方法没有给出估计信号极化参数的方法。文献：Shu T,He J,Han X,etal.Joint DOA and Degree-of-Polarization Estimation of Partially-PolarizedSignals Using Nested Arrays[J].IEEE Communications Letters,2020,24(10):2182-2186基于同上的阵型，提出一种四元数方法(简称四元数法)。该方法根据四元数理论，不仅能估计信号的波达方向还可以估计极化参数中的极化度。但以上方法存在以下问题：第一，采用的都是双极化阵列，无法累积信号所有极化分量的功率；第二，由于欠定条件下无法估计噪声方差，导致四元数法估计的极化度是有偏的；第三，以上两种方法都无法估计信号极化参数中的极化指向角和极化椭圆率角。
发明内容
本发明的目的在于针对上述现有技术的不足，提供一种基于三极化嵌套阵列的部分极化信号参数估计方法。
本发明是通过以下技术方案来实现：
步骤(1)根据优化法则布置三极化嵌套阵列，总阵元个数L＝L1+L2，L1和L2分别为嵌套阵两个子阵列的阵元个数；
步骤(2)构建三极化阵列的部分极化信号接收模型：K个来自θ＝[θ1,θ2,…,θK]方向的窄带不相关信号沿y-z平面入射到三极化嵌套阵列上，θk∈[0,π]为以逆时针方向从y轴正半轴旋转到各个入射信号方向的夹角，k＝1,2,...,K；
阵列在第n个快拍时所有三极化阵元的接收数据用向量形式表示为：n＝1,2,...,N，N表示采样的快拍数，(·)T表示转置操作，表示复数域，ηl[n]为第l个阵元在第n个快拍时接收的数据，l＝1,2,...,L；极化阵列流型表示Kronecker积，第k个信号的导向矢量a(θk)＝[a1(θk),a2(θk),…,aL(θk)]T，al(θk)表示第k个信号在第l个阵元处的响应，极化响应矩阵第n个快拍的信号矢量sk[n]表示第k个信号在第n个快拍时水平和垂直极化分量构成的向量；第n个快拍的噪声矢量υl[n]表示第l个阵元在第n个快拍时的噪声向量；
旋转矩阵极化向量其中，αk为极化指向角，βk为极化椭圆率角，满足-π/2＜αk≤π/2，-π/4≤βk≤π/4，
步骤(3)利用子协方差矩阵相加构造虚拟协阵输出；
步骤(4)利用Toeplitz方法恢复矩阵的秩；
步骤(5)估计信号波达方向：通过子空间类方法估计信号的波达方向，得到其估计值
步骤(6)估计噪声方差：对η[n]的协方差矩阵进行特征值分解，并对其最小的2L-K个特征值取平均得到噪声方差的估计值(·)H表示取共轭转置；
步骤(7)估计信号极化参数：
结合已估计出的和重构出各个信号的协方差矩阵k＝1,2,...,K，即恢复出中的各个参数rk,HH、rk,VV和rk,HV，rk,HH、rk,VV分别表示第k个信号水平和垂直极化分量的功率，rk,HV表示第k个信号的水平和垂直极化分量的相关系数；恢复的rk,HH构成向量rHH，恢复的rk,VV构成向量rVV，恢复的rk,HV构成向量rHV，则重构部分极化信号的协方差矩阵表示为：
对第k个重构的协方差矩阵进行特征分解，设特征值分别为uk,1和uk,2，且uk,1＞uk,2，uk,2对应的特征向量ξk,2，信号的极化度估计为：
解方程[G(αk)W(βk)]Hξk,2＝0，得到信号的极化指向角αk和极化椭圆率角βk的估计值。
进一步，步骤(1)具体是：在三维直角坐标中，沿y轴布置嵌套阵列，满足优化法则L1＜(L2+1)/(L2-2)；各个阵元的位置构成向量μd＝[1,2,…,L1,L1+1,2(L1+1),…,L2(L1+1)]Td，d表示相邻阵元的最小间距，取d＝λ/2，λ为信号波长；每个阵元包含三个极化方向，分别为x方向极化、y方向极化和z方向极化，每个极化方向有单独的输出，即每个阵元有三个输出端口。
更进一步，步骤(2)中，和ηl,z[n]对应于第l个阵元在x、y和z方向极化的输出；
al(θk)＝exp[j(2πμldcosθk)/λ]，μl表示向量μ的第l个元素；
sk[n]＝[sk,H[n],sk,V[n]]T，包含第k个信号的水平和垂直极化分量；
υl[n]＝[υl,x[n],υl,y[n],υl,z[n]]T，包含第l个阵元在x、y和z方向极化的输出包含的噪声，满足其中σ2表示噪声方差，I3表示3阶单位矩阵，表示均值为μ、方差为σ2的复高斯分布；
第k个信号sk[n]由其协方差矩阵描述：
其中，和分别表示第k个信号的随机极化成分的功率和完全极化成分的功率，该信号的极化度表示为ρk∈[0,1]，E{·}表示求期望。
再进一步，步骤(3)具体是：
将各个三极化阵元的x、y和z方向极化输出分别抽取出来，并按顺序构成x方向极化输出向量x[n]、y方向极化输出向量y[n]及z方向极化输出向量z[n]，表示为：其中，η1:3:3L-2[n]表示由η[n]的第1,4,7,…,3L-2个元素构成的向量，其余表示类似；
x[n]、y[n]和z[n]的协方差矩阵表示为： 将三个不同极化方向的协方差矩阵相加，得到新的协方差矩阵IL表示L阶单位矩阵；
定义辅助向量vec(·)表示将矩阵按列向量化；将μv中的元素从小到大排列并去除重复的元素以保证每个元素只出现一次，得到一个连续整数序列，构成长度为L′＝2L2(L1+1)-1的列向量μc；计算转换矩阵满足T的第i行第j列的元素i＝1,2,...,L′，j＝1,2,...,L2；其中，表示实数域，[·]i表示取向量的第i个元素；δi,j为Kronecker Delta函数，只有当i＝j时，δi,j＝1，否则δi,j＝0；虚拟协阵的输出为：
又进一步，其特征在于，步骤(4)中Toeplitz矩阵表示为RT的信号子空间与V＝[v(θ1),v(θ2),…,v(θK)]的列空间相同，其中v(θk)的第l个元素表示为vl(θk)＝exp(j2π(l-1)dcos(θk)/λ)，l＝1,2,…,L2(L1+1)，k＝1,2,...,K，当L2(L1+1)＞K时，满足rank(RT)＞K，rank(·)表示矩阵的秩。
还进一步，步骤(7)中所述恢复出Rsk中的各个参数rk,HH、rk,VV和rk,HV具体是：
恢复参数rk,HH：
首先计算矩阵然后利用公式得到参数rHH＝[r1,HH,r2,HH,…,rK,HH]T；其中，(·)*表示取共轭，为求伪逆；
恢复参数rk,VV：
首先计算矩阵其中 表示矩阵的第k列；然后求其中，则rVV＝[r1,VV,r2,VV,…,rK,VV]T；
恢复参数rk,HV：
通过计算其中，rHV＝[r1,HV,r2,HV,…,rK,HV]T，且有diag(·)表示以向量元素作为对角线元素构成对角矩阵。
本发明与现有技术相比有以下优点：
首先，本发明采用的是三极化阵列，为累积信号所有极化分量的功率和估计噪声方差提供条件；其次，本发明提出了子协方差矩阵相加的方法，一方面可将三极化嵌套阵列模型转换为传统嵌套阵列模型，从而能构造虚拟协阵及提升阵列自由度，另一方面累积了信号所有极化分量的功率；再者，本发明除了能估计信号的波达方向和极化度，还能估计极化指向角和极化椭圆率角；最后，本发明采用的是优化的嵌套阵列，可在欠定条件下估计噪声方差，从而实现极化参数的无偏估计。
附图说明
图1是本发明的方法总体流程框图；
图2是本发明中三极化阵列接收信号示意图；
图3是本发明中部分极化信号的极化椭圆示意图；
图4是本发明方法在欠定条件下多次估计的空间谱示意图；
图5是本发明方法与其他方法在不同信噪比下的波达方向估计性能比较；
图6是本发明方法与其他方法在不同快拍数下的波达方向估计性能比较；
图7是本发明方法与四元数法在不同信噪比下的极化度性能比较。
具体实施方式
以下参照附图，对发明的具体技术方案和效果作进一步地详细说明。
如图1所示，基于三极化嵌套阵列的部分极化信号参数估计方法，具体如下：
步骤(1)根据优化法则布置三极化嵌套阵列：
在三维直角坐标中，沿y轴布置嵌套阵列，满足优化法则L1＜(L2+1)/(L2-2)，L1和L2分别为嵌套阵两个子阵列的阵元个数，总阵元个数L＝L1+L2；各个阵元的位置构成向量μd＝[1,2,…,L1,L1+1,2(L1+1),…,L2(L1+1)]Td，d表示相邻阵元的最小间距，取d＝λ/2，λ为信号波长；每个阵元包含三个极化方向，分别为x方向极化、y方向极化和z方向极化，每个极化方向有单独的输出，即每个阵元有三个输出端口。
该三极化阵列如图2所示。相比与双极化天线，利用三极化天线进行信号接收有潜力完全积累信号所有极化分量的功率，不会造成信号功率的浪费。
步骤(2)构建三极化阵列的部分极化信号接收模型：
假设K个来自θ＝[θ1,θ2,…,θK]方向的窄带不相关信号沿y-z平面入射到三极化嵌套阵列上，波达方向θk∈[0,π]为以逆时针方向从y轴正半轴旋转到各个入射信号方向的夹角，k＝1,2,...,K；第l个阵元在第n个快拍时接收的数据表示为：
其中，表示复数域，ηl,x[n]、ηl,y[n]和ηl,z[n]对应于第l个阵元在x、y和z方向极化的输出，N表示采样的快拍数；表示第k个信号的极化响应矩阵；al(θk)＝exp[j(2πμldcosθk)/λ]表示第k个信号在第l个阵元处的响应，μl表示向量μ的第l个元素；sk[n]＝[sk,H[n],sk,V[n]]T包含第k个信号的水平和垂直极化分量；第l个阵元噪声υl[n]＝[υl,x[n],υl,y[n],υl,z[n]]T包含第l个阵元在x、y和z方向极化的输出包含的噪声，满足其中σ2表示噪声方差，I3表示3阶单位矩阵，表示均值为μ、方差为σ2的复高斯分布，(·)T表示转置操作；
第k个信号sk[n]由其协方差矩阵描述：
其中，rk,HH、rk,VV分别表示第k个信号水平和垂直极化分量的功率，rk,HV表示第k个信号的水平和垂直极化分量的相关系数；和分别表示第k个信号的随机极化成分的功率和完全极化成分的功率，该信号的极化度(DOP)表示为ρk∈[0,1]，(·)H表示取共轭转置，E{·}表示求期望；
旋转矩阵极化向量其中，αk为极化指向角(POA)，βk为极化椭圆率角(PEA)，如图3所示，满足-π/2＜αk≤π/2，-π/4≤βk≤π/4；
阵列在第n个快拍时所有三极化阵元的接收数据用向量形式表示为：
其中，第k个信号的极化阵列流型第k个信号的导向矢量a(θk)＝[a1(θk),a2(θk),…,aL(θk)]T，第n个快拍的信号矢量噪声矢量表示Kronecker积。
步骤(3)利用子协方差矩阵相加构造虚拟协阵输出：
将各个三极化阵元的x、y和z方向极化输出分别抽取出来，并按顺序构成x方向极化输出向量x[n]、y方向极化输出向量y[n]及z方向极化输出向量z[n]，表示为：
其中，η1:3:3L-2[n]表示由η[n]的第1,4,7,…,3L-2个元素构成的向量，其余表示类似；x[n]、y[n]和z[n]的协方差矩阵表示为：
将三个不同极化方向的协方差矩阵相加，得到新的协方差矩阵：
其中，IL表示L阶单位矩阵；当快拍数较大时，式中的约等号视为等号。
定义辅助向量vec(·)表示将矩阵按列向量化；将μv中的元素从小到大排列并去除重复的元素以保证每个元素只出现一次，然后得到一个连续整数序列，构成列向量μc；根据嵌套阵理论，可知μc的长度为2L2(L1+1)-1；计算转换矩阵使其满足：其中，表示实数域，Ti,j表示T的第i行第j列，[·]i表示取向量的第i个元素。δi,j为Kronecker Delta函数，只有当i＝j时，δi,j＝1，否则δi,j＝0；虚拟协阵的输出为：
容易看出，该步骤通过子协方差矩阵相加得到的协方差矩阵的表达式类似于传统的嵌套阵列协方差矩阵模型，且pk为第k个信号的水平和垂直极化分量的功率之和，完整的利用了信号各个极化分量的功率，这一点是双极化阵列无法做到的。
步骤(4)利用Toeplitz方法恢复矩阵的秩：
由于虚拟协阵输出γ为单测量向量，其对应的秩亏问题使得无法同时对多个信号参数进行估计。通过构造Toeplitz矩阵的方法可以恢复矩阵的秩。Toeplitz矩阵表示为：
RT的信号子空间与V＝[v(θ1),v(θ2),…,v(θK)]的列空间相同，其中v(θk)的第l个元素表示为vl(θk)＝exp(j2π(l-1)dcos(θk)/λ)，l＝1,2,…,L2(L1+1)，k＝1,2,...,K，当L2(L1+1)＞K时，满足rank(RT)＞K，rank(·)表示矩阵的秩。
步骤(5)估计信号波达方向：
通过子空间类方法估计信号的波达方向，得到其估计值
以MUSIC算法为例，首先对RT进行特征值分解，其L2(L1+1)-K个最小特征值对应的特征向量构成噪声子空间Un，然后利用V＝[v(θ1),v(θ2),…,v(θK)]与噪声子空间正交的关系得到MUSIC空间谱其中v(θ)的第l个元素vl(θ)＝exp(j2π(l-1)dcos(θ)/λ)。该空间谱的K个最大峰值对应的角度即为信号的波达方向。将信号的估计结果表示为
步骤(6)估计噪声方差：
对总的三极化阵列的输出信号η[n],n＝1,2,...,N的协方差矩阵进行特征值分解，并对其最小的2L-K个特征值取平均得到噪声方差的估计值
需要注意的是，欠定条件下能够实现这一步骤的前提是步骤一中的嵌套阵需满足优化法则L1＜(L2+1)/(L2-2)。这是由于当L1＜(L2+1)/(L2-2)时，由于使用了三极化阵列，即使在欠定条件下也能得到R的噪声子空间维度大于零，从而可用平均对应于噪声子空间的特征值的方法估计噪声方差。若使用的是双极化阵列，则无法使用本步骤进行噪声方差估计。
步骤(7)估计信号极化参数：
结合已估计出的和重构出各个信号的协方差矩阵即恢复出该矩阵中的各个参数rk,HH、rk,VV和rk,HV，k＝1,2,…K；
对rk,HH的估计：首先计算矩阵然后利用公式得到参数rHH＝[r1,HH,r2,HH,…,rK,HH]T；其中，(·)*表示取共轭，为求伪逆；
对rk,VV的估计：首先计算矩阵其中表示矩阵的第k列；然后求其中，则rVV＝[r1,VV,r2,VV,…,rK,VV]T；
对rk,HV的估计，通过计算其中，rHV＝[r1,HV,r2,HV,…,rK,HV]T，且有 diag(·)表示以向量元素作为对角线元素构成对角矩阵；
重构部分极化信号的协方差矩阵表示为：
然后，对每个估计出的部分极化信号的协方差矩阵进行处理。
以第k个重构的协方差矩阵为例，对其进行特征分解，假设特征值分别为uk,1和uk,2，且uk,1＞uk,2，uk,2对应的特征向量ξk,2，于是信号的极化度可估计为：
再解方程[G(αk)W(βk)]Hξk,2＝0，即得到信号的极化指向角αk和极化椭圆率角βk的估计值。
该步骤中，由于在求解rHH和rVV时均使用了估计的噪声方差使得最终得到的是对的无偏估计。若使用的是双极化阵列，由于无法估计噪声方差，最终只能得到对的有偏估计。
下面结合仿真实例对本发明的效果进行验证。
仿真实例1：设置嵌套阵列的两个子阵阵元数分别为L1＝3，L2＝3，满足优化法则，即共L＝L1+L2＝6个三极化阵元，假设10个信号的入射角度均匀分布于[30150]间，极化指向角均匀分布于[-π/3,π/4]之间，极化椭圆率角均匀分布于[-π/5,π/6]之间，信号极化度全部设置为0.9。信噪比设置为20dB，采样快拍数设置为500。将算法运行10次，可得10次叠加的信号的空间功率谱如图4所示，其中虚线部分表示真实信号的角度，峰值位置代表所估计信号的角度。从图中可以看出，本发明能够稳定估计出这10个信号的波达方向，且信号源个数10大于阵元个数6，成功实现了欠定估计。下表给出了各个信号的波达方向、极化度、极化指向角和极化椭圆率角真实值与估计值的比较，估计值对应于10次运行结果的均值。从表中可以看出，本发明能够有效估计出信号的角度和各个极化参数，且估计精度较高。
仿真实例2：为了验证本发明的估计性能，将背景技术中的正交子空间法、四元数法及克拉美罗下限与本发明作对比，其中本发明对应的克拉美罗算法可通过文献：C.Zhou,Y.Gu,X.Fan,Z.Shi,G.Mao,Y.D.Zhang.Direction-of-Arrival Estimation for CoprimeArray via Virtual Array Interpolation[J].IEEE Trans.Signal Process.,vol.66,no.22,pp.5956-5971,2018推得。设置快拍数为300个，信噪比从-10dB扫描至10dB，蒙特卡洛实验次数为500次，信号的极化度统一设置为0.7。其余参数与仿真实例1相同，仿真结果如图5所示。从图中可以看出，本发明方法在各个信噪比的均方根误差达到最小，估计效果最好。正交子空间算法估计效果优于四元数算法。四元数法精度最差主要因为其性能以信号水平极化和垂直极化功率相同为前提。
仿真实例3：设置信噪比为20dB，快拍数从50个扫描至500个，蒙特卡洛实验次数为500次，对比算法和其余参数与仿真实例2相同，仿真结果如图6所示。从图中可以看出，本发明算法在不同快拍数下的均方根误差均达到最小，估计效果最好。
仿真实例4：比较本发明方法和其他方法在极化度估计上的统计性能。由于正交子空间法不能对信号极化度进行估计，于是此例中将本发明与四元数方法作对比。设置信号快拍数为500个，信噪比从-8dB扫描至12dB，蒙特卡洛实验次数为500次，信号极化度统一设置为0.7，其余参数与仿真实例1相同，仿真结果如图7所示，从图中可以看出由于本发明方法由于估计了噪声方差，实现了对极化参数的无偏估计，其均方根误差远小于四元数法。
以上所述仅为本发名的较佳实施范例，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。
- 上一篇：一种医用注射器针头装配设备
- 下一篇：基于耦合张量分解的L型互质阵列波达方向估计方法