CN116540203B - 基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法 - Google Patents
基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法 Download PDFInfo
- Publication number
- CN116540203B CN116540203B CN202310811676.8A CN202310811676A CN116540203B CN 116540203 B CN116540203 B CN 116540203B CN 202310811676 A CN202310811676 A CN 202310811676A CN 116540203 B CN116540203 B CN 116540203B
- Authority
- CN
- China
- Prior art keywords
- hrrp
- matrix
- representing
- posterior
- precision
- 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.)
- Active
Links
- 238000000034 method Methods 0.000 title claims abstract description 46
- 238000009825 accumulation Methods 0.000 title claims abstract description 28
- 230000001427 coherent effect Effects 0.000 title claims abstract description 25
- 239000011159 matrix material Substances 0.000 claims abstract description 208
- 239000013598 vector Substances 0.000 claims abstract description 68
- 238000004422 calculation algorithm Methods 0.000 claims abstract description 32
- 238000000354 decomposition reaction Methods 0.000 claims description 33
- 238000012545 processing Methods 0.000 claims description 7
- 238000006073 displacement reaction Methods 0.000 claims description 6
- 230000001131 transforming effect Effects 0.000 claims description 6
- 238000010276 construction Methods 0.000 claims description 4
- 230000007613 environmental effect Effects 0.000 claims description 3
- 238000006243 chemical reaction Methods 0.000 claims description 2
- 230000010354 integration Effects 0.000 claims 1
- 230000006870 function Effects 0.000 description 46
- 238000004364 calculation method Methods 0.000 description 9
- 238000005259 measurement Methods 0.000 description 5
- 238000005457 optimization Methods 0.000 description 5
- 230000008569 process Effects 0.000 description 5
- 230000008859 change Effects 0.000 description 3
- 238000011156 evaluation Methods 0.000 description 3
- 238000009795 derivation Methods 0.000 description 2
- 230000000694 effects Effects 0.000 description 2
- 238000003384 imaging method Methods 0.000 description 2
- 238000012804 iterative process Methods 0.000 description 2
- 238000005070 sampling Methods 0.000 description 2
- 241000764238 Isis Species 0.000 description 1
- 238000007476 Maximum Likelihood Methods 0.000 description 1
- 230000001133 acceleration Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 238000013459 approach Methods 0.000 description 1
- 238000013398 bayesian method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 230000036461 convulsion Effects 0.000 description 1
- 230000007123 defense Effects 0.000 description 1
- 230000001419 dependent effect Effects 0.000 description 1
- 238000002592 echocardiography Methods 0.000 description 1
- 238000005516 engineering process Methods 0.000 description 1
- 230000014509 gene expression Effects 0.000 description 1
- 230000002068 genetic effect Effects 0.000 description 1
- 238000011478 gradient descent method Methods 0.000 description 1
- 238000012544 monitoring process Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
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
- G01S7/00—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00
- G01S7/02—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00
- G01S7/41—Details of systems according to groups G01S13/00, G01S15/00, G01S17/00 of systems according to group G01S13/00 using analysis of echo signal for target characterisation; Target signature; Target cross-section
- G01S7/415—Identification of targets based on measurements of movement associated with the target
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F18/00—Pattern recognition
- G06F18/20—Analysing
- G06F18/24—Classification techniques
- G06F18/241—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches
- G06F18/2415—Classification techniques relating to the classification model, e.g. parametric or non-parametric approaches based on parametric or probabilistic models, e.g. based on likelihood ratio or false acceptance rate versus a false rejection rate
- G06F18/24155—Bayesian classification
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Abstract
本发明公开了一种基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法,包括:建立雷达回波观测模型;基于雷达回波观测模型得到矩阵形式的雷达回波观测模型;基于矩阵形式的雷达回波观测模型的HRRP得到分层先验联合概率分布;基于分层先验联合概率分布得到待重构HRRP的第二后验均值和第二后验协方差矩阵;基于SBL框架和期望最大化算法更新离散形式的HRRP的精度向量和观测噪声精度;基于SBL框架和期望最大化算法更新高阶相位误差矩阵;得到最终的后验均值、后验协方差矩阵、离散形式的HRRP的精度向量、观测噪声精度和高阶相位误差矩阵,以获得最终的高速目标相参积累结果。本发明能够在速度未知的情况下对高速目标进行快速相参积累。
Description
技术领域
本发明属于雷达信号处理技术领域,具体涉及一种基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法。
背景技术
在新体制宽带雷达系统中,高分辨率距离像(HRRP,High Resolution RangeProfile)是目标距离范围内不同距离单元回波的相干叠加,反映了目标沿雷达视距(LOS,Line of Sight)的一维投影分布。HRRP包含了丰富的目标电磁散射特性、结构分布和形状轮廓等信息。同时,HRRP重构是目标识别、逆合成孔径雷达(ISAR,Inverse SyntheticAperture Radar)和三维成像的第一步。与高维目标成像相比,HRRP具有测量简单、处理效率高、运动依赖性低等优点。这些独特的优点使HRRP重构成为新体制宽带雷达系统中不可缺少的功能,在空间目标监视和弹道目标防御等领域发挥着重要作用。然而,对于弹道导弹、人造卫星和空间碎片等高超声速空间目标,由于径向速度在回波中产生二次相位项,导致目标的HRRP拉伸和分裂。为了不影响后续基于HRRP的应用,有必要对高速目标的速度进行测量或估计,以补偿脉冲内高阶相位项。
在实际应用中,利用窄带测量数据来获得目标速度并补偿目标运动。然而,这种方法在处理空间目标时面临着重大挑战:(1)捕获高超音速目标很困难;(2)窄带测量的速度不能反映宽带观测过程中的目标速度;提高测速精度需要更多的系统资源。此外,传统基于信号处理的速度补偿方法大多以波形熵为评价准则,对速度进行搜索,该类方法计算量大并且对信噪比要求较高。
因此,如何提供一种高速运动目标的快速HRRP重构及自聚焦技术成为了亟待解决的问题。
发明内容
为了解决现有技术中所存在的上述问题,本发明提供了一种基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法。
本发明要解决的技术问题通过以下技术方案实现:
一种基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法,所述相参积累方法包括:
步骤1、建立具有高阶相位误差的雷达回波观测模型;
步骤2、基于所述雷达回波观测模型得到矩阵形式的雷达回波观测模型;
步骤3、基于所述矩阵形式的雷达回波观测模型的HRRP得到待重构HRRP的分层先验联合概率分布;
步骤4、基于所述待重构HRRP的分层先验联合概率分布和雷达回波观测数据,利用基于GS分解和PCG的快速算法得到待重构HRRP的第二后验均值和第二后验协方差矩阵;
步骤5、基于SBL框架和期望最大化算法更新当前迭代的离散形式的HRRP的精度向量和观测噪声精度;
步骤6、基于SBL框架和期望最大化算法更新高阶相位误差矩阵;
步骤7、重复步骤4-步骤6,直至满足收敛条件停止迭代,得到最终的后验均值、后验协方差矩阵、离散形式的HRRP的精度向量、观测噪声精度和高阶相位误差矩阵,以获得最终的高速目标相参积累结果。
可选地,所述步骤1包括:
步骤1.1、获取宽带雷达回波信号,所述宽带雷达回波信号表示为:
其中,sr(t)表示宽带雷达回波信号,K表示目标散射体的总数量,k=1,2,…,K,Ak表示第k个目标散射体后向的散射系数,rect(·)表示窗函数,t表示快时间,τk表示时间延迟,Tp表示脉冲宽度,j表示虚数单位,fc表示工作频率,γ表示调频率;
步骤1.2、对所述宽带雷达回波信号进行解线性调频处理,得到雷达回波观测初始模型,所述雷达回波观测初始模型表示为:
其中,sd0(t)表示雷达回波观测初始模型,Rk表示第k个目标散射体和雷达的初始距离,c表示光速,v表示目标的速度;
步骤1.3、根据所述雷达回波观测初始模型得到所述雷达回波观测模型,所述雷达回波观测模型表示为:
其中,sd(t)表示雷达回波观测模型。
可选地,所述步骤2包括:
步骤2.1、对所述雷达回波观测模型进行相位补偿,得到相位补偿后的雷达回波观测模型;
步骤2.2、对所述相位补偿后的雷达回波观测模型进行关于快时间t的快速傅里叶变换并离散化,得到离散形式的HRRP,所述离散形式的HRRP表示为:
其中,h(m)表示第m+1个离散形式的HRRP,n表示回波数据的时间索引,n=0,1,…,N-1,N表示时域样本个数,sd(n)表示第n+1个时域观测数据,j表示虚数单位,m表示距离单元格的索引,m=0,1,…,M-1,M表示距离单元的总数,γ表示调频率,v表示目标的速度,c表示光速,w(m)表示第m个距离单元的环境噪声和系统噪声;
步骤2.3、基于所述离散形式的HRRP,得到矩阵形式的雷达回波观测模型,所述矩阵形式的雷达回波观测模型表示为:
s=EFh+w;
其中,s表示接收到的观测信号,s=[sd(0),sd(1),...,sd(n),...,sd(N-1)],E表示高阶相位误差矩阵,F表示逆傅里叶矩阵,h表示离散形式的HRRP,h=[h(0),h(1),...,h(m),...,h(M-1)],w表示高斯白噪声。
可选地,所述步骤3包括:
步骤3.1、获取HRRP服从的多元复高斯分布和观测信号的先验模型,所述HRRP服从的多元复高斯分布表示为:
所述观测信号的先验模型表示为:
其中,p(h|Λ-1)和p(h|α)表示HRRP服从的多元复高斯分布,p(s|h,β;E)表示观测信号的先验模型,s表示接收到的观测信号,E表示高阶相位误差矩阵,F表示逆傅里叶矩阵,h表示离散形式的HRRP,Λ=diag(α),diag(·)表示对角矩阵,α表示离散形式的HRRP的精度向量,α=[α0,α1,...,αm,...,αM-1],αm表示第m+1个离散形式的HRRP的精度,m=0,1,…,M-1,M表示距离单元的总数,表示服从复高斯分布,β表示观测噪声精度,I表示单位矩阵;
步骤3.2、对超参数施加Gamma先验以诱导稀疏性,得到第一概率密度函数和第二概率密度函数,所述第一概率密度函数表示为:
所述第二概率密度函数表示为:
其中,p(α|a1,b1)表示第一概率密度函数,p(β|a2,b2)表示第二概率密度函数,a1、b1、a2、b2均表示一个正数,a1和b1为α的分布参数,a2和b2为β的分布参数,Γ(·)表示Gamma函数,∏表示求积运算;
步骤3.3、基于所述观测信号的先验模型、所述HRRP服从的多元复高斯分布、所述第一概率密度函数和所述第二概率密度函数得到待重构HRRP的分层先验联合概率分布,所述待重构HRRP的分层先验联合概率分布表示为:
p(s,h,α,β;E)=p(s|h,β;E)p(h|α)p(α|a1,b1)p(β|a2,b2);
其中,p(s,h,α,β;E)表示待重构HRRP的分层先验联合概率分布。
可选地,所述步骤4包括:
步骤4.1、基于观测信号的先验分布和所述待重构HRRP的分层先验联合概率分布得到所述离散形式的HRRP的后验模型,所述离散形式的HRRP的后验模型表示为:
其中,p(h|s,α,β;E)表示离散形式的HRRP的后验模型,p(s,h,α,β;E)表示待重构HRRP的分层先验联合概率分布,p(s|α,β;E)=∫p(s|h,β;E)p(h|α)dh,p(s|h,β;E)表示观测信号的先验模型,p(h|α)表示HRRP服从的多元复高斯分布,h表示离散形式的HRRP,s表示接收到的观测信号,E表示高阶相位误差矩阵,α表示离散形式的HRRP的精度向量,β表示观测噪声精度;
步骤4.2、获取所述离散形式的HRRP的后验模型的第一后验协方差矩阵和第一后验均值,所述第一后验协方差矩阵表示为:
Σ1=(βFHF+Λ)-1;
所述第一后验均值表示为:
μ1=βΣ1FHEHs;
其中,Σ1表示第一后验协方差矩阵,μ1表示第一后验均值,H表示共轭转置,F表示逆傅里叶矩阵,Λ=diag(α),diag(·)表示对角矩阵;
步骤4.3、根据Woodbury矩阵恒等式将所述第一后验协方差矩阵变换为第二后验协方差矩阵,所述第二后验协方差矩阵表示为:
其中,Σ2表示第二后验协方差矩阵,QN表示参数矩阵,QN=I+βFΛ-1FH,I表示单位矩阵;
步骤4.4、根据所述第二后验协方差矩阵将所述第一后验均值变换为第二后验均值,所述第二后验均值表示为:
其中,μ2表示第二后验均值;
步骤4.5、基于GS分解得到参数矩阵QN的逆以求解所述步骤4.4所得到的第二后验均值。
可选地,所述步骤4.5包括:
步骤4.51、将所述参数矩阵QN转换为块格式,所述块格式的参数矩阵QN表示为:
其中,q0表示QN中的第(1,1)个元素,qN-1表示QN中第1列元素中除q0外的其他元素所组成的向量,QN-1表示QN中的一个(N-1)×(N-1)的子矩阵,T表示转置操作,JN-1表示(N-1)×(N-1)的翻转矩阵,*表示取共轭;
步骤4.52、将矩阵反演引理应用于所述步骤4.51中的所述块格式的参数矩阵QN,得到所述参数矩阵QN的逆所述参数矩阵QN的逆/>表示为:
其中,
步骤4.53、根据下移移位算子和所述步骤4.52得到的参数矩阵QN的逆得到所述参数矩阵QN的逆/>的位移秩,所述参数矩阵QN的逆/>的位移秩/>表示为:
其中,gN表示第一低秩分解因子向量,gn表示gN中的第n+1个元素,pN表示第二低秩分解因子向量,/>pn表示pN中的第n+1个元素,ZN表示下移移位算子,/> IN-1表示N-1阶的单位矩阵;
步骤4.54、根据所述参数矩阵QN的逆的位移秩构造所述参数矩阵QN的逆/>
步骤4.55、根据所述参数矩阵QN的逆的位移秩、所述第一低秩分解因子向量和所述第二低秩分解因子向量得到参数矩阵QN的逆/>的GS分解形式,所述参数矩阵QN的逆的GS分解形式表示为:
其中, 表示托普利兹矩阵生成算子;
步骤4.56、基于所述步骤4.55得到的参数矩阵QN的逆的GS分解形式求解所述第二后验均值,其中,/>μ2=βΛθM。
可选地,所述步骤5包括:
步骤5.1、获取所述离散形式的HRRP的精度向量的目标函数和所述观测噪声精度的目标函数,所述离散形式的HRRP的精度向量的目标函数和所述观测噪声精度的目标函数表示为:
其中,Q(α;αi)表示离散形式的HRRP的精度向量的目标函数,Q(β;βi)表示观测噪声精度的目标函数,α表示离散形式的HRRP的精度向量,αi表示第i次迭代时的离散形式的HRRP的精度向量,αm表示第m+1个离散形式的HRRP的精度,h(m)表示第m+1个离散形式的HRRP,m表示距离单元格的索引,m=0,1,…,M-1,M表示距离单元的总数,βi表示第i次迭代时的观测噪声精度,β表示观测噪声精度,N表示时域样本个数,s表示接收到的观测信号,E表示高阶相位误差矩阵,F表示逆傅里叶矩阵,h表示离散形式的HRRP,||·||2表示向量的欧氏范数,<·>表示求期望运算符,const表示与待估计参数无关的常数项,a1、b1、a2、b2均表示一个正数;
步骤5.2、基于所述离散形式的HRRP的精度向量的目标函数和观测噪声精度的目标函数,通过最大化目标函数更新当前迭代的离散形式的HRRP的精度向量和观测噪声精度。
可选地,第i+1次迭代时的第m+1个离散形式的HRRP的精度的更新规则为:
第i+1次迭代时的观测噪声精度的更新规则为:
其中,表示第i+1次迭代时的第m+1个离散形式的HRRP的精度,βi+1表示第i+1次迭代时的观测噪声精度,γm=1-αmΣm,m,Σm,m表示待重构HRRP的第二后验协方差矩阵对角线的第m+1个元素,μm表示待重构HRRP的第二后验均值的第m+1个元素,μ2表示第二后验均值。
可选地,所述步骤6包括:
利用EM算法对所述高阶相位误差矩阵进行估计,以更新所述高阶相位误差矩阵。
可选地,第i+1次迭代时的高阶相位误差矩阵的更新规则为:
其中,表示第i+1次迭代时的高阶相位误差矩阵对角线的第n+1个元素,sn表示观测向量的第n+1个元素,Fn,·表示逆傅里叶矩阵的第n+1行,trace(·)表示对矩阵求迹,μ2表示第二后验均值,Σ2表示第二后验协方差矩阵,H表示共轭转置。
与现有技术相比,本发明的有益效果:
本发明将高速目标相参积累和高阶相位补偿问题整合到一个稀疏贝叶斯框架中,并提出了一种快速的SBL(Sparse Bayesian Learning,稀疏贝叶斯)计算方法。该方法通过对高阶相位误差的高分辨率距离像的观测模型进行贝叶斯推理得到HRRP的后验期望,利用期望最大化算法估计了迭代过程中分层先验参数和高阶相位误差,通过迭代计算实现了HRRP的重构和自动对焦。
利用Hermitian-Toeplitz矩阵可被Gohberg-Semencul(GS)分解的特性,将SBL算法迭代过程中矩阵求逆操作用快速傅里叶变换和快速逆傅里叶变换实现,提高了基于SBL框架的HRRP重构及高阶相位补偿的效率,并使计算复杂度大幅降低。
以下将结合附图及对本发明做进一步详细说明。
附图说明
图1是本发明实施例提供的一种基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法的流程示意图。
具体实施方式
下面结合具体实施例对本发明做进一步详细的描述,但本发明的实施方式不限于此。
实施例一
近些年来,对空间目标进行相参积累的方法主要有三类:第一类是基于信号参数估计的方法,该方法利用高速目标回波为线性调频信号的特性,对目标速度相关的调频项进行估计,得到目标速度,并对二次相位误差进行补偿,然后采用快速傅里叶变换(FFT,Fast Fourier Transform)对补偿后的高速目标雷达回波进行相参积累,得到目标的HRRP。第二类方法是基于聚焦性能评价的搜索优化方法。该方法以波形熵或对比度作为目标函数,通过优化算法搜索目标的准确速度以补偿高阶相位项,然后采用FFT实现目标的相参积累。第三类是利用HRRP所具有的较强的稀疏特性,通过压缩感知(CS,Compressed Sensing)与聚焦性能评价的搜索优化方法相结合,实现空间目标的速度补偿和相参积累。
但是,基于信号参数估计的方法在对高速目标HRRP重构时具有局限性,例如分数阶傅里叶变换(FrFT,Fractional Fourier Transform)可用于线性调频信号的相参积累,但这种方法需要选择最优的FrFT阶数,通常采用搜索的方法,这需要较高的计算资源,之后基于先进时频分析技术的方法相继被采用,如积分三次相函数法(ICPF,Integral CubicPhase Function)和吕氏分布法(LVD,Lv’sDistribution),这些方法对信号采样率和质量要求较高,增加了雷达系统的硬件成本。基于聚焦性能指标的搜索优化方法已在ISAR应用中用于运动补偿和自聚焦处理。而优化问题通常采用梯度下降法、牛顿法、遗传算法等求解,这些算法收敛速度慢,可能陷入局部最优值,目标相参积累失败。在已报道的CS算法中,贝叶斯方法已被证明可以获得更稀疏的解,并且在参数学习和统计信息获取方面具有优势。然而,它涉及到高维矩阵的反演运算,由此带来了巨大的计算负担。
因此,请参见图1,图1是本发明实施例提供的一种基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法的流程示意图,本发明提供一种基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法,该相参积累方法包括:
步骤1、建立具有高阶相位误差的雷达回波观测模型。
步骤1.1、获取宽带雷达回波信号。
在宽带雷达系统中,高距离分辨率通常是通过发射大带宽线性调频信号来实现的。单个LFM(Linear Frequency Modulation,线性调频)脉冲信号可以表示为:
其中,st(t)表示LFM脉冲信号,Tp表示时宽,t表示快时间,rect(·)表示窗函数,j表示虚数单位,fc表示工作频率,γ表示调频率。
通常情况下,观测时间内目标与雷达之间的瞬时径向距离变化可以看作是时间的多项式函数。对于大多数场景,目标的运动可以用三阶多项式函数来描述。然而,加速度和加加速度对距离变化的贡献在毫秒级脉冲时间可以忽略不计。因此,目标的运动可以表示为:
R(t)=R0+vt (2)
其中,R(t)表示目标散射点的距离变化,R0表示散射点和雷达之间的初始距离,v表示目标的速度。
假设目标以速度v接近雷达,目标由K个散射体组成,则宽带雷达回波信号可以表示为:
其中,sr(t)表示宽带雷达回波信号,K表示目标散射体的总数量,k=1,2,…,K,Ak表示第k个目标散射体后向的散射系数,τk表示时间延迟,τk=2(Rk+vt)/c,Rk表示第k个目标散射体和雷达的初始距离,c表示光速。
步骤1.2、对宽带雷达回波信号进行解线性调频(Dechirping)处理,得到雷达回波观测初始模型,雷达回波观测初始模型表示为:
其中,sd0(t)表示雷达回波观测初始模型,v表示目标的速度。
步骤1.3、根据雷达回波观测初始模型得到雷达回波观测模型。
观察公式(4),雷达回波观测初始模型(也即观测信号)的相位由七项组成。第一项是距离测量项,对应于最终HRRP的时域数据。第二项和第六项是线性项,它们引起HRRP的移动,但不影响其分布。在ISAR应用中,可以通过包络对齐来补偿。第四项和第五项为恒相位项,对HRRP没有影响,可以忽略。第三项和第七项是二次相位项,引起HRRP的展宽和分裂。本发明的目标是在稀疏贝叶斯(SBL,Sparse Bayesian Learning)框架下对二次相位项进行补偿,忽略对目标HRRP没有影响的相位项后,HRRP的重构观测信号(即雷达回波观测模型)可以表示为:
其中,sd(t)表示雷达回波观测模型。
步骤2、基于雷达回波观测模型得到矩阵形式的雷达回波观测模型。
步骤2.1、对雷达回波观测模型进行相位补偿,得到相位补偿后的雷达回波观测模型。
步骤2.2、对相位补偿后的雷达回波观测模型进行关于快时间t的快速傅里叶变换并离散化,得到离散形式的HRRP,离散形式的HRRP表示为:
其中,h(m)表示第m+1个离散形式的HRRP,N表示时域样本个数,n表示回波数据的时间索引,n=0,1,…,N-1,sd(n)表示第n+1个时域观测数据,m表示距离单元格的索引,m=0,1,…,M-1,M表示距离单元的总数,w(m)表示第m个距离单元的环境噪声和系统噪声。
步骤2.3、基于离散形式的HRRP,得到矩阵形式的雷达回波观测模型。
即,对于时域观测数据sd(n),对应于公式(6)的逆问题,具有高阶相位误差的雷达回波观测模型可表示为如下矩阵形式:
s=EFh+w (7)
其中,s表示接收到的观测信号, 表示复数域,s=[sd(0),sd(1),...,sd(n),...,sd(N-1)],E表示高阶相位误差矩阵,/>F表示逆傅里叶矩阵,h表示离散形式的HRRP,h=[h(0),h(1),...,h(m),...,h(M-1)],w表示高斯白噪声。在次奈奎斯特采样中,设置M>>N,以达到高的距离分辨率。
到目前为止,具有高阶相位误差的HRRP观测模型已被表述为易于处理的表达式。在高频段中,高速目标的HRRP由几个强散射点构成,具有很强的稀疏特性。因此重构精确的HRRP,需要找到向量h的最稀疏解,同时对于高速运动目标还需要估计高阶相位误差E。
步骤3、基于矩阵形式的雷达回波观测模型的HRRP得到待重构HRRP的分层先验联合概率分布。
步骤3.1、获取HRRP服从的多元复高斯分布和观测信号的先验模型。
在SBL框架中,对高阶相位补偿的雷达回波观测模型的矩阵形式的随机变量进行分层建模,以推导出h的稀疏性。具有高阶相位误差的高速目标HRRP建模过程如下:
第一层假设稀疏向量h服从多元复高斯分布,则HRRP服从的多元复高斯分布表示为:
其中,p(h|Λ-1)和p(h|α)表示HRRP服从的多元复高斯分布,Λ=diag(α),diag(·)表示对角矩阵,即diag(α)表示对角线元素是α的对角矩阵,α表示离散形式的HRRP的精度向量,α=[α0,α1,...,αm,...,αM-1],αm表示第m+1个离散形式的HRRP的精度(方差的倒数),表示服从复高斯分布。
在雷达信号采集过程中,高斯白噪声w通常被建模为具有相同噪声方差的独立同分布复高斯分布,因此观测信号的先验模型可表示为:
其中,p(s|h,β;E)表示观测信号的先验模型,β表示观测噪声精度,I表示单位矩阵。
这里,αm和β为超参数,其中αm控制h的稀疏度,β用于学习观测中的噪声。
步骤3.2、对超参数施加Gamma先验以诱导稀疏性,得到第一概率密度函数和第二概率密度函数。
也就是说,在层次先验的第二层,对超参数施加Gamma先验来诱导稀疏性,假设α中每个元素均服从独立同分布,则第一概率密度函数表示为:
第二概率密度函数表示为:
其中,p(α|a1,b1)表示第一概率密度函数,p(β|a2,b2)表示第二概率密度函数,a1、b1、a2、b2均表示一个正数,且均为接近于0的正数,具体的,a1=b1=a2=b2=10-4,Γ(·)表示Gamma函数,∏表示求积运算。
步骤3.3、基于观测信号的先验模型、HRRP服从的多元复高斯分布、第一概率密度函数和第二概率密度函数得到待重构HRRP的分层先验联合概率分布。
也就是说,根据公式(8)-(11),具有未知参数E的联合概率分布表示为:
p(s,h,α,β;E)=p(s|h,β;E)p(h|α)p(α|a1,b1)p(β|a2,b2) (12)
其中,p(s,h,α,β;E)表示待重构HRRP的分层先验联合概率分布。
步骤4、基于待重构HRRP的分层先验联合概率分布和雷达回波观测数据,利用基于GS分解和PCG的快速算法得到待重构HRRP的第二后验均值和第二后验协方差矩阵。
步骤4.1、基于观测信号的先验分布和待重构HRRP的分层先验联合概率分布得到离散形式的HRRP的后验模型。
这里,根据贝叶斯定理和分层先验模型,离散形式的HRRP的后验模型表示为:
其中,p(h|s,α,β;E)表示离散形式的HRRP的后验模型,p(s|α,β;E)=∫p(s|h,β;E)p(h|α)dh,p(s|h,β;E)表示观测信号的先验模型,p(h|α)表示HRRP服从的多元复高斯分布。
步骤4.2、获取离散形式的HRRP的后验模型的第一后验协方差矩阵和第一后验均值。
这里,根据第二类最大似然估计,p(h|s,α,β;E)服从多元复高斯分布,则p(h|s,α,β;E)的第一后验协方差矩阵表示为:
Σ1=(βFHF+Λ)-1 (14)
p(h|s,α,β;E)的第一后验均值为:
μ1=βΣ1FHEHs (15)
其中,Σ1表示第一后验协方差矩阵,μ1表示第一后验均值。
步骤4.3、根据Woodbury(伍德伯里)矩阵恒等式将第一后验协方差矩阵变换为第二后验协方差矩阵。
具体而言,公式(14)涉及计算协方差矩阵,这需要对大小为N×N的矩阵求逆,其中N是观测数据的长度。众所周知,矩阵反演的计算复杂度为 表示时间复杂度,即使是中等大小的矩阵,也会消耗大量的计算成本,限制了其工程应用。在HRRP观测模型中,感知矩阵是逆傅里叶变换矩阵,这使得迭代过程中的相关矩阵具有特殊的结构。该结构可用于加速矩阵反演计算和矩阵向量乘法计算。因此,根据Woodbury矩阵恒等式,公式(14)可以改写为:
其中,Σ2表示第二后验协方差矩阵,QN表示参数矩阵。
步骤4.4、根据第二后验协方差矩阵将第一后验均值变换为第二后验均值,即将公式(16)代入公式(15)得到:
QN=I+βFΛ-1FH (18)
其中,μ2表示第二后验均值。
步骤4.5、基于GS分解得到参数矩阵QN的逆以求解步骤4.4所得到的第二后验均值。/>
步骤4.51、将参数矩阵QN转换为块格式。
具体而言,令RN=FΛ-1FH,RN是一个Hermitian-Toeplitz矩阵,即对角元素相等且共轭对称。其中(·)N表示大小为N×N的矩阵或长度为N的向量,根据公式(18)可以看出参数矩阵QN也是一个Hermitian-Toeplitz矩阵。因此,构造矩阵RN只需要第一列元素。第一列元素可以通过逆傅里叶变换得到,即:
其中,rn表示RN中第一列向量中的第n+1个元素。
参考快速逆傅里叶变换(IFFT,Inverse Fast Fourier Transform),可以通过对1/αm进行M点IFFT并取前N个值来获得RN的第一列元素。
与任意矩阵不同,Hermitian-Toeplitz矩阵的逆可以通过GS分解来实现,它将矩阵的逆运算转化为对GS分解因子的移位运算。此外,逆矩阵与矢量的乘法可以通过对GS分解因子和矢量进行FFT和IFFT来实现,这大大降低了计算成本。为了说明这一点,本实施例将QN写成块格式,如下所示:
其中,q0表示QN中的第(1,1)个元素,qN-1表示QN中第1列元素中除q0外的其他元素所组成的向量,QN-1表示参数矩阵QN中的一个(N-1)×(N-1)的子矩阵,即QN-1表示的(N-1)×(N-1)的参数矩阵,T表示转置操作,JN-1表示(N-1)×(N-1)的翻转矩阵,*表示取共轭。
设J为一个翻转矩阵,定义如下:
其中,JN表示N×N的翻转矩阵。
步骤4.52、将矩阵反演引理应用于步骤4.51中的块格式的参数矩阵QN,得到参数矩阵QN的逆参数矩阵QN的逆/>表示为:
其中:
步骤4.53、根据下移移位算子和步骤4.52得到的参数矩阵QN的逆得到参数矩阵QN的逆/>的位移秩/>
这里,根据上述公式可以看出:
其中,ZN表示下移移位算子,IN-1表示阶数为N-1的单位矩阵。
根据公式(22)和公式(25),得到参数矩阵QN的逆的位移秩/>参数矩阵QN的逆/>的位移秩/>表示为:
为方便起见,定义以下两个低秩分解因子向量,即:
其中,gN表示第一低秩分解因子向量,pN表示第二低秩分解因子向量,gn表示gN中的第n+1个元素,pn表示pN中的第n+1个元素。
步骤4.54、根据参数矩阵QN的逆的位移秩构造参数矩阵QN的逆/>
具体而言,根据Hermitian-Toeplitz矩阵的特殊结构,可以通过参数矩阵QN的逆/>的位移秩/>来构造,即:
步骤4.55、根据参数矩阵QN的逆的位移秩、第一低秩分解因子向量和第二低秩分解因子向量得到参数矩阵QN的逆/>的GS分解形式。
具体而言,将公式(26)-公式(28)代入公式(29),可以得到的GS分解形式,/>的GS分解形式可以表示为:
其中,G和P是下三角Hermitian-Toeplitz矩阵,分别以gN和pN为第一列元素。它们的结构由下式给出:
公式(30)即为的GS分解公式,其中gN和pN为GS分解因子,可通过求解公式(23)获得,/>表示托普利兹矩阵生成算子。将公式(23)重新表示为/>很明显,上式是一个Toeplitz(托普利兹)线性系统,/>为该系统的解。本发明采用预处理共轭梯度(PCG,Preconditioned Conjugate Gradient)方法迭代求解公式(23),同时考虑到Toeplitz线性系统的性质,PCG算法中的矩阵乘向量操作可通过FFT和IFFT算法实现,此外通过选择合适的预处理矩阵可以加速算法的收敛速度。该算法的时间复杂度为其中u为PCG算法收敛时的迭代次数。
步骤4.56、基于步骤4.55得到的参数矩阵QN的逆的GS分解形式求解第二后验均值。
也就是说,第二后验均值的计算可以分为以下三个公式去求解,即:
μ2=βΛθM (35)
通过PCG算法得到的GS分解形式并带入公式(33),公式(33)的右侧是托普利兹矩阵与向量的乘法运算,可以通过FFT和IFFT有效实现。公式(34)可以通过FFT进行计算。最后根据公式(35),通过θM、β和Λ对角线元素的点积计算μ2。
步骤5、基于SBL框架和期望最大化算法更新当前迭代的离散形式的HRRP的精度向量和观测噪声精度。
步骤5.1、获取离散形式的HRRP的精度向量的目标函数和观测噪声精度的目标函数。
在步骤4中,需要注意的是,对于公式(14)、公式(15)而言,由于EHE=I,因此在公式中省略了与E相关的因素。为了更新剩余的未知隐变量和参数,本发明采用期望最大化算法(EM,Expectation-Maximization algorithm)进行推理,这是模型中存在隐变量和未知参数时,对其进行估计的标准框架。用θ={α,β,E}表示待估计的参数集。在EM算法中,期望步(E-step)和最大化步(M-step)交替迭代进行。设参数在第i次迭代中的估计值为θi。在E-step中,计算完整数据对数似然相对于h后验分布的期望,即 表示求相对于后验分布p(h|s,θi)的期望,然后,在M-step中,通过最大化目标函数Q(θ;θi)得到参数在第(i+1)次迭代中的估计值。
在E-step中,对于待估计参数集{αi,βi,Ei}而言,αi表示第i次迭代时的离散形式的HRRP的精度向量,βi表示第i次迭代时的观测噪声精度,Ei表示第i次迭代时的高阶相位误差矩阵,根据待估计参数集{αi,βi,Ei}的当前迭代值和观测信号s,计算{α,β,E}的完全对数后验对h的期望,即要最大化目标函数Q,因此:
其中,Q(α;αi)表示离散形式的HRRP的精度向量的目标函数,Q(β;βi)表示观测噪声精度的目标函数,||·||2表示向量的欧氏范数,const表示与待估计参数无关的常数项,<·>表示求期望运算符,其具体形式为:
其中,h(m)表示h中的第m+1个元素,μm表示待重构HRRP的第二后验均值的第m+1个元素,trace(·)表示求矩阵的迹,Σm,m表示待重构HRRP的第二后验协方差矩阵对角线的第m+1个元素。
步骤5.2、基于离散形式的HRRP的精度向量的目标函数和观测噪声精度的目标函数,通过最大化目标函数更新当前迭代的离散形式的HRRP的精度向量和观测噪声精度。
在M-step中,通过最大化Q函数来更新当前迭代的隐变量{α,β}和参数E。下边是隐变量{α,β}的更新规则。
α的更新规则为:将公式(38)代入公式(36),对得到的结果求αm的导数,令所求导数为零,第m+1个离散形式的HRRP的精度的更新规则为:
其中,表示第i+1次迭代时的第m+1个离散形式的HRRP的精度,γm=1-αmΣm,m。
β的更新规则为:将公式(39)代入公式(37),对β求导,令求导的结果为零,则第i+1次迭代时的观测噪声精度的更新规则为:
其中,βi+1表示第i+1次迭代时的观测噪声精度。
α和β两个隐变量的更新依赖于第二后验协方差矩阵的对角元素,本发明利用的GS分解,快速求解α和β。
因此,将公式(16)写为:
Σ2=Λ-1-βΛ-1HMΛ-1 (42)
/>
因为Λ是一个对角矩阵,Σ2的对角元素仅依赖于的对角元素。定义一个向量 表示/>的主对角元素构成的向量。根据公式(43),/>的第m+1个元素/>为:
其中,fN(m)表示逆傅里叶矩阵的第m+1列向量,为/>的第l1条对角线上所有元素的和,l1=0表示主对角线,l1>0在主对角线以上,l1<0在主对角线以下,l1=-(N-1),-(N-2),…0,…,(N-2),(N-1)。由主对角元素和下对角元素之和组成的向量表示为e=[e-N+1,e-N+2,…,e-1,e0]T。
根据Hermitian-Toeplitz矩阵的性质和傅里叶变换的周期性,也可以写成:
其中,表示向量/>中的元素,/>表示向量/>中的元素,l2=0,1,…,M-1,/>和/>通过e进行构造,/>0M-N表示长度为M-N的零向量,0M-N+1表示长度为M-N+1的零向量,e(0:N-2)表示e的第1个到第N-1个元素。
公式(45)等价于:
其中,表示/>的主对角元素构成的向量,FFT(·)M是M点傅里叶变换算子。
若有的GS分解,则e为:
其中, gn表示gN中的第n+1个元素,pn表示PN中的第n+1个元素,n=0,1,…,N-1。
明显,公式(47)可以通过FFT和IFFT有效地计算。在公式(42)中,β为标量,Λ-1为对角矩阵,因此通过公式(46)得到的对角线元素/>最后通过/>β和Λ-1的对角线元素的点积得到Σ2的对角线元素。进一步地,将Σ2的对角线元素代入γm=1-αmΣm,m,并根据公式(40)和公式(41)以计算α和β。
步骤6、基于SBL框架和期望最大化算法更新高阶相位误差矩阵。
具体而言,利用EM算法对高阶相位误差矩阵进行估计,以更新高阶相位误差矩阵。
本发明利用EM算法对未知的高阶相位误差矩阵进行估计。根据步骤5有:
其中,Q(E,Ei)表示高阶相位误差矩阵的目标函数,Ei表示第i次迭代时的高阶相位误差矩阵。
E的更新规则为:将公式(39)代入公式(48),对E的主要对角线元素求导,并使求导结果为零,即可得到E的更新规则。E的更新规则可以表示为:
其中,表示第i+1次迭代时的高阶相位误差矩阵的第n+1个对角线元素,sn(·)表示观测信号的第n+1个元素,Fn,·表示逆傅里叶矩阵的第n+1行,trace(·)表示求矩阵的迹,μ2表示第二后验均值,Σ2表示第二后验协方差矩阵,H表示共轭转置。
步骤7、重复步骤4-步骤6,直至满足收敛条件停止迭代,得到最终的后验均值、后验协方差矩阵、离散形式的HRRP的精度向量、观测噪声精度和高阶相位误差矩阵,以获得最终的高速目标相参积累结果。
本发明将快速算法应用到后验均值、后验协方差矩阵、隐变量和高阶相位误差矩阵的计算中,并进行迭代计算,以获得最终的高速目标相参积累结果。
具体的,重复步骤4-步骤6,直至满足收敛条件停止迭代,完成高速目标的相参积累。本实施例设置收敛精度为δ,根据公式(50)判断每次迭代得到的μ2值是否满足收敛条件,如果满足公式(50)则停止迭代,否则继续迭代,收敛条件表示为:
其中,μ2 i+1表示第i+1次迭代时的第二后验均值,μ2 i表示第i次迭代时的第二后验均值。
在大脉宽宽带雷达中,传统空间高速目标相参积累算法对硬件资源和计算资源要求高,高阶相位补偿容易陷入局部最优值导致相参积累失败,而本发明将高速目标的高阶相位补偿及相参积累统一到稀疏贝叶斯框架中,并利用Hermitian-Toeplitz矩阵的逆可以通过GS分解,将稀疏贝叶斯迭代过程中矩阵求逆转化成FFT和IFFT的操作,提高SBL框架下高速目标的高阶相位补偿和相参积累的效率。本发明能够在速度未知的情况下对高速目标进行快速相参积累。
需要说明的是,术语“第一”、“第二”仅用于描述目的,而不能理解为指示或暗示相对重要性或者隐含指明所指示的技术特征的数量。由此,限定有“第一”、“第二”的特征可以明示或者隐含地包括一个或者更多个特征。在本发明的描述中,“多个”的含义是两个或两个以上,除非另有明确具体的限定。
在本说明书的描述中,参考术语“一个实施例”、“一些实施例”、“示例”、“具体示例”、或“一些示例”等的描述意指结合该实施例或示例描述的具体特征或者特点包含于本发明的至少一个实施例或示例中。在本说明书中,对上述术语的示意性表述不必须针对的是相同的实施例或示例。而且,描述的具体特征或者特点可以在任何的一个或多个实施例或示例中以合适的方式结合。此外,本领域的技术人员可以将本说明书中描述的不同实施例或示例进行接合和组合。
尽管在此结合各实施例对本发明进行了描述,然而,在实施所要求保护的本发明过程中,本领域技术人员通过查看所述附图以及公开内容,可理解并实现所述公开实施例的其他变化。在说明书中,“包括”一词不排除其他组成部分或步骤,“一”或“一个”不排除多个的情况。相互不同的实施例中记载了某些措施,但这并不表示这些措施不能组合起来产生良好的效果。
以上内容是结合具体的优选实施方式对本发明所作的进一步详细说明,不能认定本发明的具体实施只局限于这些说明。对于本发明所属技术领域的普通技术人员来说,在不脱离本发明构思的前提下,还可以做出若干简单推演或替换,都应当视为属于本发明的保护范围。
Claims (10)
1.一种基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法,其特征在于,所述相参积累方法包括:
步骤1、建立具有高阶相位误差的雷达回波观测模型;
步骤2、基于所述雷达回波观测模型得到矩阵形式的雷达回波观测模型;
步骤3、基于所述矩阵形式的雷达回波观测模型的HRRP得到待重构HRRP的分层先验联合概率分布;
步骤4、基于所述待重构HRRP的分层先验联合概率分布和雷达回波观测数据,利用基于GS分解和PCG的快速算法得到待重构HRRP的第二后验均值和第二后验协方差矩阵;
步骤5、基于SBL框架和期望最大化算法更新当前迭代的离散形式的HRRP的精度向量和观测噪声精度;
步骤6、基于SBL框架和期望最大化算法更新高阶相位误差矩阵;
步骤7、重复步骤4-步骤6,直至满足收敛条件停止迭代,得到最终的后验均值、后验协方差矩阵、离散形式的HRRP的精度向量、观测噪声精度和高阶相位误差矩阵,以获得最终的高速目标相参积累结果。
2.根据权利要求1所述的相参积累方法,其特征在于,所述步骤1包括:
步骤1.1、获取宽带雷达回波信号,所述宽带雷达回波信号表示为:
其中,sr(t)表示宽带雷达回波信号,K表示目标散射体的总数量,k=1,2,…,K,Ak表示第k个目标散射体后向的散射系数,rect(·)表示窗函数,t表示快时间,表示时间延迟,Tp表示脉冲宽度,j表示虚数单位,fc表示工作频率,γ表示调频率;
步骤1.2、对所述宽带雷达回波信号进行解线性调频处理,得到雷达回波观测初始模型,所述雷达回波观测初始模型表示为:
其中,sd0(t)表示雷达回波观测初始模型,Rk表示第k个目标散射体和雷达的初始距离,c表示光速,v表示目标的速度;
步骤1.3、根据所述雷达回波观测初始模型得到所述雷达回波观测模型,所述雷达回波观测模型表示为:
其中,sd(t)表示雷达回波观测模型。
3.根据权利要求1所述的相参积累方法,其特征在于,所述步骤2包括:
步骤2.1、对所述雷达回波观测模型进行相位补偿,得到相位补偿后的雷达回波观测模型;
步骤2.2、对所述相位补偿后的雷达回波观测模型进行关于快时间t的快速傅里叶变换并离散化,得到离散形式的HRRP,所述离散形式的HRRP表示为:
其中,h(m)表示第m+1个离散形式的HRRP,n表示回波数据的时间索引,n=0,1,…,N-1,N表示时域样本个数,sd(n)表示第n+1个时域观测数据,j表示虚数单位,m表示距离单元格的索引,m=0,1,…,M-1,M表示距离单元的总数,γ表示调频率,v表示目标的速度,c表示光速,w(m)表示第m个距离单元的环境噪声和系统噪声;
步骤2.3、基于所述离散形式的HRRP,得到矩阵形式的雷达回波观测模型,所述矩阵形式的雷达回波观测模型表示为:
s=EFh+w;
其中,s表示接收到的观测信号,s=[sd(0),sd(1),...,sd(n),...,sd(N-1)],E表示高阶相位误差矩阵,F表示逆傅里叶矩阵,h表示离散形式的HRRP,h=[h(0),h(1),...,h(m),...,h(M-1)],w表示高斯白噪声。
4.根据权利要求1所述的相参积累方法,其特征在于,所述步骤3包括:
步骤3.1、获取HRRP服从的多元复高斯分布和观测信号的先验模型,所述HRRP服从的多元复高斯分布表示为:
所述观测信号的先验模型表示为:
其中,p(h|Λ-1)和p(h|α)表示HRRP服从的多元复高斯分布,p(s|h,β;E)表示观测信号的先验模型,s表示接收到的观测信号,E表示高阶相位误差矩阵,F表示逆傅里叶矩阵,h表示离散形式的HRRP,Λ=diag(α),diag(·)表示对角矩阵,α表示离散形式的HRRP的精度向量,α=[α0,α1,...,αm,...,αM-1],αm表示第m+1个离散形式的HRRP的精度,m=0,1,…,M-1,M表示距离单元的总数,表示服从复高斯分布,β表示观测噪声精度,I表示单位矩阵;
步骤3.2、对超参数施加Gamma先验以诱导稀疏性,得到第一概率密度函数和第二概率密度函数,所述第一概率密度函数表示为:
所述第二概率密度函数表示为:
其中,p(α|a1,b1)表示第一概率密度函数,p(β|a2,b2)表示第二概率密度函数,a1、b1、a2、b2均表示一个正数,a1和b1为α的分布参数,a2和b2为β的分布参数,Γ(·)表示Gamma函数,∏表示求积运算;
步骤3.3、基于所述观测信号的先验模型、所述HRRP服从的多元复高斯分布、所述第一概率密度函数和所述第二概率密度函数得到待重构HRRP的分层先验联合概率分布,所述待重构HRRP的分层先验联合概率分布表示为:
p(s,h,α,β;E)=p(s|h,β;E)p(h|α)p(α|a1,b1)p(β|a2,b2);
其中,p(s,h,α,β;E)表示待重构HRRP的分层先验联合概率分布。
5.根据权利要求1所述的相参积累方法,其特征在于,所述步骤4包括:
步骤4.1、基于观测信号的先验分布和所述待重构HRRP的分层先验联合概率分布得到所述离散形式的HRRP的后验模型,所述离散形式的HRRP的后验模型表示为:
其中,p(h|s,α,β;E)表示离散形式的HRRP的后验模型,p(s,h,α,β;E)表示待重构HRRP的分层先验联合概率分布,p(s|α,β;E)=∫p(s|h,β;E)p(h|α)dh,p(s|h,β;E)表示观测信号的先验模型,p(h|α)表示HRRP服从的多元复高斯分布,h表示离散形式的HRRP,s表示接收到的观测信号,E表示高阶相位误差矩阵,α表示离散形式的HRRP的精度向量,β表示观测噪声精度;
步骤4.2、获取所述离散形式的HRRP的后验模型的第一后验协方差矩阵和第一后验均值,所述第一后验协方差矩阵表示为:
Σ1=(βFHF+Λ)-1;
所述第一后验均值表示为:
μ1=βΣ1FHEHs;
其中,Σ1表示第一后验协方差矩阵,μ1表示第一后验均值,H表示共轭转置,F表示逆傅里叶矩阵,Λ=diag(α),diag(·)表示对角矩阵;
步骤4.3、根据Woodbury矩阵恒等式将所述第一后验协方差矩阵变换为第二后验协方差矩阵,所述第二后验协方差矩阵表示为:
其中,Σ2表示第二后验协方差矩阵,QN表示参数矩阵,QN=I+βFΛ-1FH,I表示单位矩阵;
步骤4.4、根据所述第二后验协方差矩阵将所述第一后验均值变换为第二后验均值,所述第二后验均值表示为:
其中,μ2表示第二后验均值;
步骤4.5、基于GS分解得到参数矩阵QN的逆以求解所述步骤4.4所得到的第二后验均值。
6.根据权利要求5所述的相参积累方法,其特征在于,所述步骤4.5包括:
步骤4.51、将所述参数矩阵QN转换为块格式,所述块格式的参数矩阵QN表示为:
其中,q0表示QN中的第(1,1)个元素,qN-1表示QN中第1列元素中除q0外的其他元素所组成的向量,QN-1表示QN中的一个(N-1)×(N-1)的子矩阵,T表示转置操作,JN-1表示(N-1)×(N-1)的翻转矩阵,*表示取共轭,N表示时域样本个数;
步骤4.52、将矩阵反演引理应用于所述步骤4.51中的所述块格式的参数矩阵QN,得到所述参数矩阵QN的逆所述参数矩阵QN的逆/>表示为:
其中,
步骤4.53、根据下移移位算子和所述步骤4.52得到的参数矩阵QN的逆得到所述参数矩阵QN的逆/>的位移秩,所述参数矩阵QN的逆/>的位移秩/>表示为:
其中,gN表示第一低秩分解因子向量,gn表示gN中的第n+1个元素,pN表示第二低秩分解因子向量,/>pn表示pN中的第n+1个元素,ZN表示下移移位算子,/> IN-1表示N-1阶的单位矩阵;
步骤4.54、根据所述参数矩阵QN的逆的位移秩构造所述参数矩阵QN的逆/>
步骤4.55、根据所述参数矩阵QN的逆的位移秩、所述第一低秩分解因子向量和所述第二低秩分解因子向量得到参数矩阵QN的逆/>的GS分解形式,所述参数矩阵QN的逆/>的GS分解形式表示为:
其中, 表示托普利兹矩阵生成算子;
步骤4.56、基于所述步骤4.55得到的参数矩阵QN的逆的GS分解形式求解所述第二后验均值,其中,/>μ2=βΛθM。
7.根据权利要求1所述的相参积累方法,其特征在于,所述步骤5包括:
步骤5.1、获取所述离散形式的HRRP的精度向量的目标函数和所述观测噪声精度的目标函数,所述离散形式的HRRP的精度向量的目标函数和所述观测噪声精度的目标函数表示为:
其中,Q(α;αi)表示离散形式的HRRP的精度向量的目标函数,Q(β;βi)表示观测噪声精度的目标函数,α表示离散形式的HRRP的精度向量,αi表示第i次迭代时的离散形式的HRRP的精度向量,αm表示第m+1个离散形式的HRRP的精度,h(m)表示第m+1个离散形式的HRRP,m表示距离单元格的索引,m=0,1,…,M-1,M表示距离单元的总数,βi表示第i次迭代时的观测噪声精度,β表示观测噪声精度,N表示时域样本个数,s表示接收到的观测信号,E表示高阶相位误差矩阵,F表示逆傅里叶矩阵,h表示离散形式的HRRP,||·||2表示向量的欧氏范数,<·>表示求期望运算符,const表示与待估计参数无关的常数项,a1、b1、a2、b2均表示一个正数;
步骤5.2、基于所述离散形式的HRRP的精度向量的目标函数和观测噪声精度的目标函数,通过最大化目标函数更新当前迭代的离散形式的HRRP的精度向量和观测噪声精度。
8.根据权利要求7所述的相参积累方法,其特征在于,第i+1次迭代时的第m+1个离散形式的HRRP的精度的更新规则为:
第i+1次迭代时的观测噪声精度的更新规则为:
其中,表示第i+1次迭代时的第m+1个离散形式的HRRP的精度,βi+1表示第i+1次迭代时的观测噪声精度,γm=1-αmΣm,m,Σm,m表示待重构HRRP的第二后验协方差矩阵对角线的第m+1个元素,μm表示待重构HRRP的第二后验均值的第m+1个元素,μ2表示第二后验均值。
9.根据权利要求1所述的相参积累方法,其特征在于,所述步骤6包括:
利用EM算法对所述高阶相位误差矩阵进行估计,以更新所述高阶相位误差矩阵。
10.根据权利要求9所述的相参积累方法,其特征在于,第i+1次迭代时的高阶相位误差矩阵的更新规则为:
其中,表示第i+1次迭代时的高阶相位误差矩阵对角线的第n+1个元素,sn表示观测向量的第n+1个元素,Fn,·表示逆傅里叶矩阵的第n+1行,trace(·)表示对矩阵求迹,μ2表示第二后验均值,Σ2表示第二后验协方差矩阵,H表示共轭转置。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202310811676.8A CN116540203B (zh) | 2023-07-04 | 2023-07-04 | 基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202310811676.8A CN116540203B (zh) | 2023-07-04 | 2023-07-04 | 基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN116540203A CN116540203A (zh) | 2023-08-04 |
CN116540203B true CN116540203B (zh) | 2023-09-22 |
Family
ID=87456304
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202310811676.8A Active CN116540203B (zh) | 2023-07-04 | 2023-07-04 | 基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN116540203B (zh) |
Citations (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103713288A (zh) * | 2013-12-31 | 2014-04-09 | 电子科技大学 | 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法 |
CN110068805A (zh) * | 2019-05-05 | 2019-07-30 | 中国人民解放军国防科技大学 | 基于变分贝叶斯推论的高速目标hrrp重构方法 |
CN113030972A (zh) * | 2021-04-26 | 2021-06-25 | 西安电子科技大学 | 基于快速稀疏贝叶斯学习的机动目标isar成像方法 |
WO2022017801A2 (en) * | 2020-07-21 | 2022-01-27 | Veoneer Sweden Ab | Adaptive ofdm radar operation based on time variable subcarrier assignments |
CN115421115A (zh) * | 2022-05-23 | 2022-12-02 | 中国人民解放军空军预警学院 | 一种用于联合相位校正与isar成像的重赋权交替方向乘子法 |
CN115453528A (zh) * | 2022-08-05 | 2022-12-09 | 西安电子科技大学 | 基于快速sbl算法实现分段观测isar高分辨成像方法及装置 |
CN115453527A (zh) * | 2022-08-02 | 2022-12-09 | 西安电子科技大学 | 一种周期性分段观测isar高分辨成像方法 |
CN115616518A (zh) * | 2022-11-03 | 2023-01-17 | 中国人民解放军国防科技大学 | 宽带雷达微弱目标运动参数估计与高分辨距离像重构方法 |
WO2023015623A1 (zh) * | 2021-08-13 | 2023-02-16 | 复旦大学 | 一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法 |
Family Cites Families (1)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
IT1391472B1 (it) * | 2008-09-26 | 2011-12-23 | Mbda italia spa | Procedimento di elaborazione di un segnale di eco radar, prodotto da un bersaglio, per compensare effetti di degradazione di detto segnale dovuti al moto di detto bersaglio. |
-
2023
- 2023-07-04 CN CN202310811676.8A patent/CN116540203B/zh active Active
Patent Citations (9)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN103713288A (zh) * | 2013-12-31 | 2014-04-09 | 电子科技大学 | 基于迭代最小化稀疏贝叶斯重构线阵sar成像方法 |
CN110068805A (zh) * | 2019-05-05 | 2019-07-30 | 中国人民解放军国防科技大学 | 基于变分贝叶斯推论的高速目标hrrp重构方法 |
WO2022017801A2 (en) * | 2020-07-21 | 2022-01-27 | Veoneer Sweden Ab | Adaptive ofdm radar operation based on time variable subcarrier assignments |
CN113030972A (zh) * | 2021-04-26 | 2021-06-25 | 西安电子科技大学 | 基于快速稀疏贝叶斯学习的机动目标isar成像方法 |
WO2023015623A1 (zh) * | 2021-08-13 | 2023-02-16 | 复旦大学 | 一种多旋翼无人机载合成孔径雷达分段孔径成像及定位方法 |
CN115421115A (zh) * | 2022-05-23 | 2022-12-02 | 中国人民解放军空军预警学院 | 一种用于联合相位校正与isar成像的重赋权交替方向乘子法 |
CN115453527A (zh) * | 2022-08-02 | 2022-12-09 | 西安电子科技大学 | 一种周期性分段观测isar高分辨成像方法 |
CN115453528A (zh) * | 2022-08-05 | 2022-12-09 | 西安电子科技大学 | 基于快速sbl算法实现分段观测isar高分辨成像方法及装置 |
CN115616518A (zh) * | 2022-11-03 | 2023-01-17 | 中国人民解放军国防科技大学 | 宽带雷达微弱目标运动参数估计与高分辨距离像重构方法 |
Non-Patent Citations (4)
Title |
---|
High-Speed Targets Detection with Large Pulse Width Wideband LFM Signal Based on IAA-CICPF Algorithm;Hang Dong et al.;《2022 5th International Conference on Information Communication and Signal Processing (ICICSP)》;80-85 * |
Integrated Tracking and ISAR Imaging Using an Integrated Kalman Filter With Wideband Radar;SHAOPENG WEI et al.;《IEEE Transactions on Aerospace and Electronic Systems》;第58卷(第5期);4639-4655 * |
Scattering Centers Analysis via Sparse Bayesian Learning for non-Coherent Short Pulse Radar;Haibo Wang et al.;《2019 International Conference on Microwave and Millimeter Wave Technology (ICMMT)》;1-3 * |
基于迭代加权L2/L1范数块稀疏信号重构的ISAR成像算法;冯俊杰等;《计算机工程》;第44卷(第11期);234-238 * |
Also Published As
Publication number | Publication date |
---|---|
CN116540203A (zh) | 2023-08-04 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN104749553B (zh) | 基于快速稀疏贝叶斯学习的波达方向角估计方法 | |
CN110208735B (zh) | 一种基于稀疏贝叶斯学习的相干信号doa估计方法 | |
Pan et al. | Multi-task hidden Markov modeling of spectrogram feature from radar high-resolution range profiles | |
CN112163373A (zh) | 基于贝叶斯机器学习的雷达系统性能指标动态评估方法 | |
CN110146886A (zh) | 非均匀旋转目标运动参数的快速估计方法 | |
Shi et al. | Independent component analysis | |
CN110726992A (zh) | 基于结构稀疏和熵联合约束的sa-isar自聚焦法 | |
CN104950297A (zh) | 基于矩阵1范数拟合的阵元误差估计方法 | |
Ding et al. | Super‐resolution 3D imaging in MIMO radar using spectrum estimation theory | |
CN115453528A (zh) | 基于快速sbl算法实现分段观测isar高分辨成像方法及装置 | |
CN108614235B (zh) | 一种多鸽群信息交互的单快拍测向方法 | |
CN109783960A (zh) | 一种基于网格部分细化的波达方向估计方法 | |
CN115236584A (zh) | 基于深度学习的米波雷达低仰角估计方法 | |
CN113534065B (zh) | 一种雷达目标微动特征提取与智能分类方法及系统 | |
CN112904298B (zh) | 一种基于局部网格分裂的网格偏离空时自适应处理方法 | |
CN116540203B (zh) | 基于快速稀疏贝叶斯的宽带雷达高速目标的相参积累方法 | |
Huang et al. | Off-grid DOA estimation in real spherical harmonics domain using sparse Bayesian inference | |
Sit et al. | Characterizing evaporation ducts within the marine atmospheric boundary layer using artificial neural networks | |
CN109541567B (zh) | 基于深度学习的高速机动目标检测方法 | |
Cabanes et al. | Non-supervised high resolution Doppler machine learning for pathological radar clutter | |
CN110136167A (zh) | 面向监视系统的多群目标跟踪方法及跟踪系统 | |
CN116165619A (zh) | 一种高机动目标运动参数估计和相参积累检测方法 | |
CN113109760B (zh) | 一种基于组稀疏的多线谱联合doa估计和聚类方法及系统 | |
CN115754896A (zh) | 基于变分推理稳健稀疏贝叶斯学习的波达方向估计方法 | |
CN114624646A (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 | ||
GR01 | Patent grant | ||
GR01 | Patent grant |