一种双基地前视合成孔径雷达的运动误差估计方法与流程

本发明属于雷达
技术领域：
，特别涉及合成孔径雷达(SAR)成像技术中的双基地前视合成孔径雷达运动补偿方法。
背景技术：
：合成孔径雷达是一种具有高分辨率的成像雷达，与光学传感器相比，SAR具有全天时全天候工作能力的独特优点。随着SAR技术的发展和提高，其分辨率越来越高，目前已接近或超光学成像的分辨率，因而被广泛的应用于地球遥感、海洋研究、资源勘探、灾情预报和军事侦察等领域。双基地前视SAR的发射机和接收机分别工作于侧视和前视，不同的收发装置因而具有作用距离远、抗干扰性强、能对前视场景成像等特点。对于双基地前视SAR的高分辨成像，自聚焦是关键步骤之一。传统自聚焦方法均假设场景内各点目标具有相同的相位误差，忽略了运动误差的空变性。随着双基前视SAR分辨率要求的越来越高，这种假设不能成立。针对这一问题，在传统单基地正侧视SAR系统中，产生了一系列运动误差估计方法，先估计载机平台的运动误差，再使用运动补偿方法将运动误差精确补偿，以此得到良好聚焦的SAR图像。文献H.M.CantalloubeandC.E.Nahum,“Multiscalelocalmap-drift-drivenmultilaterationSARautofocususingfastpolarformatimagesynthesis,”IEEETrans.Geosci.RemoteSens.,vol.49,no.10,pp.3730–3736,2011提出了运动误差的参数化模型，并利用多边法来估计机载平台运动误差参数。文献Y.Li,C.Liu,Y.WangandQ.Wang,“Arobustmotionerrorestimationmethodbasedonrawdata,”IEEETrans.Geosci.RemoteSens.,vol.50,no.7,pp.2780–2790,2012利用PGA来估计不同距离单元的强散射点的相位误差，并使用最小加权平方法根据强点目标相位误差反解运动误差。然而以上文献都只考虑了运动误差的距离向空变，相较于单基地正侧视SAR，双基地前视SAR接收机前视使得运动误差的方位空变性无法忽略，因此以上述方法运用到双基地前视SAR模式中，运动误差估计精度较低。技术实现要素：本发明的目的是针对现有的运动误差算法中存在的缺陷，本发明提出了一种双基地前视合成孔径雷达的运动误差估计方法。本发明的技术方案为：一种双基地前视合成孔径雷达的运动误差估计方法，具体包括如下步骤：S0：系统参数初始化以场景参考目标点P0＝[000]T为原点建立坐标系，场景中任意点目标记为PA＝[xyz]T，系统发射机工作于侧视，接收机工作于前视，理想位置分别记于与其中，η为方位向慢时间变量；xT,yR,zT,zR分别对应着接收机或发射机的轴坐标位置；下标T和R分别对应着发射机和接收机；发射机和接收机都沿着y轴以理想速度V运动，当平台存在轨迹误差时其位置分别记为其中，eT＝[ΔxT(η)ΔyT(η)ΔzT(η)]T，eR＝[ΔxR(η)ΔyR(η)ΔzR(η)]T，Δx(η),Δy(η),Δz(η)分别对应着载机于方位向慢时间为η时的x、y、z轴的位置偏差，下标T和R分别对应着发射机和接收机，点目标的双基距离和r＝rT+rR，其中，rT＝||PT-PA||,rR＝||PR-PA||，理想双基距离和其中距离差S1：距离向脉冲压缩获取二维回波数据，并进行距离向脉冲压缩，脉冲压缩后点目标回波数据记为s0(τ,η)：s0(τ,η)=sinc(τ-rc)waz(η)exp(-j2πλr)---(1)]]>其中，τ为距离向快时间变量，waz(η)为方位向时域包络，λ为波长，c为电磁波传播速度；S2：距离走动校正对距离压缩后的双基前视SAR回波数据进行距离走动校正，将s0(τ,η)傅里叶变换到距离频域得到S0(fτ,η)并乘以走动校正相位：s1(τ,η)＝IFFT{S0(fτ,η)exp[-j2πλfτ(VTcosθT+VRcosθR)η/c]}(2)其中，VT、VR为发射、接收平台运动速度，θT、θR为发射、接收平台斜视角，fτ表示快时间频率，走动校正后结果记为s1(τ,η)；S3：方位时间慢解耦引入新的时间轴t，时间轴t垂直于数据τ-η平面，对信号s1(τ,η)进行数据扩展，将个M×N维的零矩阵沿时间轴t置于数据s1(τ,η)上方，将个M×N维的零矩阵沿时间轴t置于数据s1(τ,η)下方，其中，M,N和K分别是τ，η，t轴的采样点数，将新产生的三维数据记为s(τ,η,t)；对信号s(τ,η,t)进行二维傅里叶变换，表示成s(τ,fη,ft)，fη表示方位频率，ft表示慢时间频率，对信号s(τ,fη,ft)乘以二阶距离压缩因子和距离徙动校正因子的共轭，如式(3)所示，结果记为信号s1(τ,fη,ft):s1(τ,fη,ft)＝s(τ,fη,ft)×exp[-jφrcm(fη,ft)-jφsrc(fη,ft)](3)其中φsrc(fη,ft)≈2πc4k2f0[(fηf0)2-(fηf0)3]ft2+k38k23[3k1(cf0)((fηf0)2-(fηf0)3)ft2+(cf0)2(3(fηf0)2-4(fηf0)3)ft3]+9k32-4k2k464k25×6k12(cf0)((fηf0)2-(fηf0)3)ft2+4k1(cf0)2(3(fηf0)2-4(fηf0)3)ft3+(cf0)3(6(fηf0)2-10(fηf0)3)ft4---(4)]]>φrcm(fη,ft)≈2πfη-Rcenc+14k2[k12c-cf02ft2]+k38k23[k13c-3k1cf02ft2-2c2f03ft3]+9k32-4k2k464k25k14c-6k12cf02ft2-8k1c2f03ft3-3c3f04ft4---(5)]]>其中，f0表示发射脉冲的中心频率，系数k1、k2、k3、k4以及孔径中心距离Rcen由当方位时间为η时的距离方程扩展得到，为本领域的现有技术，不再详细说明。对信号s1(τ,fη,ft)沿t轴进行傅里叶反变换，得到信号s1(τ,fη,t)；方位慢时间解耦合操作后的信号表示为s2(τ,fη,t)s2(τ,fη,t)＝s1(τ,fη,t)exp(jφ(fη,t))(6)其中，φ(fη,t)=-2πfηcy2R+z2RVyRt---(7)]]>t=VyRcy2R+z2Rη---(8)]]>把s2(τ,fη,t)在方位域进行傅里叶反变换，得到s3(τ,η,t)。S4：运动误差估计记离散化的数据信号为sm,n,k，m,n,k分别与τ,η,t相对应；发射机和接收机的运动误差矩阵分别表示为eT和eR，慢时间为k时的发射机和接收机的运动误差(运动误差矩阵第k列)记为(eT)k＝[(ΔxT)k(ΔyT)k(ΔzT)k]T，(eR)k＝[(ΔxR)k(ΔyR)k(ΔzR)k]T，其中，(ΔxT)k、(ΔyT)k、(ΔzT)k、(ΔxR)k、(ΔyR)k、(ΔzR)k分别表示发射机和接收机在慢时间为k时的位置偏差；记为慢时间为k时的第m，n个像素点的双基距离和，其中，(rT)m,n,k、(rR)m,n,k分别表示慢时间为k时的第m，n个像素点到发射机、接收机的距离；和分别表示发射机和接收机在慢时间k时的位置矢量，估计步骤如下：步骤S41：令i＝0，(eT)0＝0，(eR)0＝0；步骤S42：若已得到(eT)i,(eR)i，则步骤S43：步骤S44：若停止；否则令i＝i+1，并返回步骤S42。其中，(eT)0和(eR)0分别表示发射机和接收机的初始估计误差；i表示迭代次数；(eT)i、(eR)i分别表示当迭代次数为i时的发射机和接收机的误差；阈值β为预先设定的常数。在S42中，第i次迭代时的目标函数梯度表示为：▿f(eT,(eR)i)=∂f∂(ΔxT)1∂f∂(ΔyT)1∂f∂(ΔzT)1∂f∂(ΔxT)2∂f∂(ΔyT)2∂f∂(ΔzT)2.........∂f∂(ΔxT)k∂f∂(ΔyT)k∂f∂(ΔzT)kT---(9)]]>其中，∂f∂(ΔθT)k=ΣmΣn∂zm,nz*m,n∂(ΔθT)k---(10)]]>θ＝x,y,z(11)zm,n=Σksm,n,kexp(j2πλrk;m,n)---(12)]]>∂zm,nz*m,n∂(ΔθT)k=2Re{z*m,n∂zm,n∂(ΔθT)k}---(13)]]>其中，z*m,n表示zm,n的共轭；把式(12)、(13)带入式(10)得：∂f∂(ΔθT)k=ΣmΣn2Re{z*m,n∂zm,n∂(ΔθT)k}---(14)]]>其中，∂zm,n∂(ΔθT)k=j2πλsm,n,kexp(j2πλrk;m,n)∂(rT)m,n,k∂(ΔθT)k---(15)]]>∂(rT)m,n,k∂(ΔθT)k=(θ~T)k+(ΔθT)k-θm,n||(P~T)k+(eT)k-Pm,n||---(16)]]>θ如式(11)所示，θm,n表示第m,n个像素点的θ坐标值，Pm,n表示第m,n个像素点的坐标位置；把式(15)、(16)带入式(14)得：▿f(eT,(eR)i)=v1v2....vkT---(17)]]>vk=∂f∂(ΔxT)k∂f∂(ΔyT)k∂f∂(ΔzT)kT=ΣmΣnγm,n,kvm,n,k---(18)]]>其中，γm,n,k=2Im{2πλz*m,nsm,n,kexp(j2πλrk;m,n)}||(P~T)k+(eT)k-Pm,n||---(19)]]>vm,n,k=(x~T)k+(ΔxT)k-xm,n(y~T)k+(ΔyT)k-ym,n(z~T)k+(ΔzT)k-zm,nT---(20)]]>其中，为慢时间为k时的理想情况下发射机的坐标值，同理可以得到：▿f((eT)i,eR)=v1v2....vkT---(21)]]>vk=∂f∂(ΔxR)k∂f∂(ΔyR)k∂f∂(ΔzR)kT=ΣmΣnγm,n,kvm,n,k---(22)]]>γm,n,k=2Im{2πλz*m,nsm,n,kexp(j2πλrk;m,n)}||(P~R)k+(eR)k-Pm,n||---(23)]]>vm,n,k=(x~R)k+(ΔxR)k-xm,n(y~R)k+(ΔyR)k-ym,n(z~R)k+(ΔzR)k-zm,nT---(24)]]>其中，为慢时间为k时的理想情况下接收机的坐标值；通过上述步骤，可以精确估计双基地前视SAR运动误差。本发明的有益效果：本发明提供的方法首先对距离向处理后的数据进行方位-慢时间解耦的预处理，有效解决了运动误差的空变性；其次基于最大图像强度准则建立优化模型；最后利用块坐标下降方法(BCD)来求解该优化模型，与现有技术相比，本发明的方法能够更加精确地估计运动误差和得到良好聚焦的双基地前视SAR图像，解决了现有运动误差估计方法未能考虑运动方位空变的问题，从而实现双基地前视合成孔径雷达精确的运动误差估计和良好聚焦。本发明的方法克服了现有的合成孔径雷达运动误差估计方法忽略方位空变，导致误差估计精度恶化的问题。附图说明图1是本发明提供方法的流程图。图2是本发明实施例采用的双基前视雷达系统结构图。图3是本发明实施例采样的双基前视雷达仿真参数表图。图4是本发明实施例的场景点目标分布图。图5是本发明实施例的发射站x方向误差估计示意图。图6是本发明实施例的发射站z方向误差估计示意图。图7是本发明实施例的接收站x方向误差估计示意图。图8是本发明实施例的接收站z方向误差估计示意图。图9是本发明实施例的成像结果示意图。具体实施方式本发明主要采用仿真实验的方法来进行验证，所有仿真都是基于matlab平台实现的。现在就具体的实施方式进行说明。S0：初始化成像系统参数；本实施例采用的双基地前视SAR几何结构如图2所示。其参数如图3参数表所示。本实施例采用的目标场景如图4所示，图中的黑色圆点为布置于地面上的9个点目标，这9个点分布于x轴和y轴上，沿x方向(切航迹)间隔100m，沿y方向(沿航迹)间隔100m，平台沿y轴运动。波束中心位于场景坐标原点处时记为零时刻，平台速度为v发射站的位置坐标为(-800,-500,1000)m，接收站零时刻位置坐标为(-600,-300,800)m，场景中心坐标为(0,0,0)；波速中心位于场景坐标原点处时记为零时刻，平台速度为50m/s，场景中任一点目标的位置坐标为P[xy0]T；用MATLAB仿真出回波数据。对成像区域内的任意目标计算其距离历程r，作为BFSAR仿真回波数据，记为s'(τ,η)系统仿真参数如图3所示。在此步骤中对回波数据加入误差eT,eR：eT＝[ΔxT(η)ΔyT(η)ΔzT(η)]T，eR＝[ΔxR(η)ΔyR(η)ΔzR(η)]T，其中，ΔxT(η)＝(1.2sin(0.05πη)-0.7sin(0.07πη)+Δran)mΔzT(η)＝(-1.2sin(0.05πη)+0.7sin(0.07πη)+Δran)mΔyT(η)＝0ΔxR(η)＝(1.2sin(0.05πη)-0.7cos(0.07πη)+Δran)mΔzR(η)＝(-1.2sin(0.05πη)+0.7cos(0.07πη)+Δran)mΔyR(η)＝0其中，Δran表示(-0.5,0.5)之间随机序列的均值。S1：对回波数据s'(τ,η)进行距离压缩。s0(τ,η)＝IFFT(FFT(s'(τ,η))×FFT(fτ(τ)))；其中，Kr为距离向调频斜率，取Kr＝2.0e+13。S2：对距离压缩之后的数据进行距离徙动校正。把s0(τ,η)代入(2)式，计算出s1(τ,η)。其中VT＝160,VR＝300,θT＝45°,θR＝39.8°。S3：将s1(τ,η)按照S3的方式进行数据扩展得到s(τ,η,t)。对s(τ,η,t)进行二维傅里叶变换之后代入(3)式得到s1(τ,fη,ft)；对s1(τ,fη,ft)进行傅里叶反变换得到s1(τ,fη,t)；将s1(τ,fη,t)代入(6)式得到s2(τ,fη,t)，对其进行傅里叶反变换得到s3(τ,η,t)。S4：设置初始误差值(eT)0＝0、(eT)0＝0，将eT、eR代入式(19)、(20)计算出(17)；代入(23)、(24)计算(21)；使用Armijo线搜索法来获得步长因子βT、βR；按照S4-2、S4-3计算下一次迭代的eT和eR；设定阈值β＝10-3，若满足S4-4的条件则返回上述步骤并重复进行，否则停止迭代，求得估计误差eT、eR。运动误差估计结果如图5-图8所示。通过以上步骤对运动误差进行精确估计后，利用BP算法来对场景进行成像，以验证所提出的办法是否能得到图像的精确聚焦。验证成像结果如图9所示。当前第1页1 2 3