专利名称：一种基于曲面投影的InSAR成像方法
技术领域：
本发明属于雷达技术领域，它特别涉及合成孔径雷达(SAR)成像技术领域。
背景技术：
合成孔径雷达(SAR)是一种高分辨率的微波成像系统。合成孔径雷达利用大时宽带宽信号实现距离向高分辨率，依靠雷达和目标之间的相对运动来合成虚拟阵列获取方位向高分辨率，可以对照射场景进行二维成像。合成孔径雷达干涉测量(InSAR)是一般SAR功能的延伸和扩展，是利用两个或者多个位置不同的天线观测同一个目标场景，根据目标到不同天线的斜距差获得测量数据的干涉相位，再通过平台与地面观测场景的几何关系反演出地面场景的数字高程信息的技术。由于具有全天时、全天候的特点，InSAR已经成为当前提取大面积地表三维图像和地形高程变化信息的一项重要遥感技术，在地形测绘、自然灾害监测和自然资源调查等领域发挥越来越大的作用。高质量的干涉相位是InSAR获取高精度地形数字高程模型(DEM)的基础，随着高程测量精度的不断提高，InSAR数据处理对干涉相位质量的要求也越来越高，一方面需要提高InSAR成像的干涉相位保持精度，另一方面需要稀疏化干涉条纹，以降低后续InSAR数据处理中相位解缠的难度，提高高程反演精度。传统的InSAR成像处理通常采用距离多普勒(RD)算法、变尺度(CS)算法、距离徙动(RM)算法等获取SAR图像，这些成像算法由于基于参考点成像，同时对平台运动轨迹作了近似直线处理，降低了 InSAR成像精度和干涉相位保持精度。在提取干涉相位以后，干涉相位的表现形式是干涉条纹，通常陡变地形区域的干涉条纹密集，使得后续相位解缠的难度增加，同时相位解缠的精度也大大降低，导致最后的高程反演精度也明显降低，因此，如何稀疏化干涉条纹并提高相位解缠精度一直是InSAR技术研究的重点。反向投影算法是一种基于时域相干处理的成像算法，其基本思想是通过计算成像区域内每一散射点到合成孔径长度内雷达平台天线相位中心之间的双程时延，然后将对应的时域回波信号进行相干累加，从而恢复出每个散射点的散射系数信息。曲面投影算法是反向投影算法的改进算法，曲面投影算法通过在低精度的高程曲面上进行反向投影成像，精确计算投影曲面上每个散射点与一个合成孔径长度内雷达平台天线相位中心的斜距信息，实现每个散射点的精确聚焦。将曲面投影算法应用于InSAR数据成像，不仅可以减小由于参考点近似和平台运动轨迹误差对成像的影响，提高InSAR成像精度和干涉相位保持精度，而且可以稀疏化陡变地形区域的干涉条纹，使得提取的干涉相位易于相位解缠，同时提高相位解缠精度，从而提高高程反演精度。
发明内容
为了得到高精度的InSAR干涉相位，同时稀疏化陡变地形区域的干涉条纹，本发明提出了一种基于曲面投影的InSAR成像方法，该方法通过在低精度的高程曲面上进行反向投影成像，将不同雷达天线的测量数据投影到以低精度高程曲面建立的InSAR成像空间中，提高了 InSAR成像精度和干涉相位保持精度，同时提取干涉相位以后，陡变地形区域的干涉条纹稀疏，后续相位解缠的难度降低，相位解缠的精度提高，使得高程反演的精度提闻。为了方便描述本发明的内容，首先作以下术语定义定义1、合成孔径雷达干涉测量(InSAR)合成孔径雷达干涉测量(InSAR)指利用两个或者两个以上的SAR数据中的相位信息进行相干处理，结合雷达参数和雷达几何位置信息反演地表三维及其变化信息的遥感技术，详见文献“合成孔径雷达成像原理”，皮亦鸣等编著，电子科技大学出版社出版。定义2、数字高程模型(DEM)数字高程模型(Digital Elevation Model,DEM)是指利用一组有序数值阵列形式表示地表或地面高程的一种实体地面模型。本发明中DEM表示成一系列地面点的平面坐标X、Y和高程坐标Z组成的数据阵列。对于一个地面区域D，地形DEM表示为DEM=(Di) (xi； Yi, Zi), i e D}其中(Xi, Yi)是第i个地面像素点对应的平面坐标，Zi是对应的高程坐标。定义3、雷达成像空间雷达成像空间是指将场景空间中的散射点投影到距离向一方位向一高度向的三维空间坐标系，该空间由合成孔径雷达成像空间中的三个相互正交的坐标基确定。目前典型的合成孔径雷达成像空间包括距离向一方位向一高度向投影空间。本发明中用以下数学关系表不成像空间M Μ = {Ρ(ν,ζ/)|Ρ(ν,ζ/) = ν· , + · , +H{v,n)-qh,u,v&R}其中孓孓和S表示构成成像空间M的三个相互正交的坐标基，分别表示距离向、方位向和高度向。F(v，W)为成像空间中的采样点向量，U，V分别表示该点的距离向坐标和方位向坐标，H(v，u)表示该点的的高度向坐标，H(v，u)与距离向坐标u和方位向坐标V存在
对应关系，R表不实数。定义4、合成孔径雷达反向投影算法反向投影算法是基于匹配滤波原理的合成孔径雷达成像算法，其主要通过相干累加实现合成孔径雷达数据的聚焦成像。详细内容可参考文献“Research onA novel fastback projection algorithm for strip map bistatic SAR imaging”，Huang Yulin 等。定义5、合成孔径雷达曲面投影算法曲面投影算法是反向投影算法的改进算法，是指利用低精度的场景高程数据，建立雷达成像空间坐标系，并基于此坐标系进行反向投影成像的算法。定义6、合成孔径雷达标准距离压缩方法合成孔径雷达标准距离压缩方法是指利用合成孔径雷达发射参数，采用以下公式生成参考信号，并采用匹配滤波技术对合成孔径雷达的距离向信号进行滤波的过程。
权利要求
1. 一种基于曲面投影的InSAR成像方法，其特征是它包括以下几个步骤步骤1、初始化InSAR成像系统参数InSAR成像空间由InSAR成像空间中的三个相互正交的坐标基确定，定义与雷达平台速度方向平行并在地平面内的单位向量作为InSAR成像空间的第一个坐标基，记做5,，该坐标基方向为方位向；定义在地平面内，并与InSAR成像空间的第一个坐标基I,垂直的单位向量作为InSAR成像空间的第二个坐标基，记做&，该坐标基方向为距离向；定义垂直于地平面向上的单位向量作为InSAR成像空间的第三个坐标基，记做豆,，该坐标基方位为高度向；InSAR雷达平台包含两组天线，即主天线和副天线，两组天线之间的距离为基线长，记做B1，主天线发射脉冲信号，经过Td时间的延迟，主天线和副天线同时接收回波延迟信号； 雷达平台主天线接收的回波数据，记做I，雷达平台副天线接收的回波数据，记做良，其中瓦和反均为二维矩阵，第一维均对应方位向，第二维均对应距离向，即二维矩阵戰,,和瓦的行存储的是方位向数据，二维矩阵1和瓦的列存储的是距离向数据；初始化InSAR成像系统参数包括雷达系统工作的信号波长，记做λ，雷达平台主天线发射信号带宽，记做B，雷达平台主天线发射脉冲时宽，记做 ；，雷达平台接收系统采样频率，记做Fs，雷达系统脉冲重复频率，记做PRF，雷达平台一个合成孔径长度内的慢时刻个数，记做N1,雷达平台速度矢量，记做γ，雷达平台主天线初始位置矢量，记做瓦(O)，雷达平台副天线初始位置矢量，记做艮(O)，场景参考点位置矢量，记做瓦，雷达系统距离向采样点数，记做队，雷达系统方位向采样点数，记做Na，场景距离向散射点间隔，记做4，场景方位向散射点间隔，记做da，场景参考点到雷达平台主天线各慢时刻天线相位中心的最短距离，记做Rm。，场景参考点到雷达平台副天线各慢时刻天线相位中心的最短距离，记做Rs。， 场景参考点在雷达平台主天线回波数据和雷达平台副天线回波数据中的距离门相同，距离门位置记做I。；对于场景散射点P (a，r)，a表示散射点位于场景方位向的第a个位置， a=l,. . . , sa, sa为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r=l,. . . , sr, &为场景距离向总的散射点数，场景散射点P(a，r)的低精度高程，记做 H (a, r);上述参数中，雷达系统工作的信号波长λ，雷达平台主天线发射的信号带宽B，雷达平台主天线发射的脉冲时宽 ；，雷达平台接收系统的采样频率Fs，雷达系统的脉冲重复频率PRF，两组天线之间的基线长B1以及接收系统接收波门相对于发射信号发射波门的延迟Td在InSAR雷达系统设计过程中已经确定；雷达平台一个合成孔径长度内的慢时刻个数N1,雷达平台速度矢量γ ,雷达平台主天线初始位置矢量尺,(0)，雷达平台副天线初始位置矢量场景参考点位置矢量h，雷达系统距离向采样点数队，雷达系统方位向采样点数 Na，场景距离向散射点间隔4，场景方位向散射点间隔da，场景参考点到雷达平台主天线各慢时刻天线相位中心的最短距离Rm。，场景参考点到雷达平台副天线各慢时刻天线相位中心的最短距离Rs。，场景参考点在雷达平台主天线回波数据和雷达平台副天线回波数据中的距离门位置I。以及场景散射点P(a，r)的低精度高程H(a，r)，a表示散射点位于场景方位向的第a个位置，a=l, . . . , sa, sa为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r=l,. . . , sr, sr为场景距离向总的散射点数，在InSAR雷达成像观测方案设计中已经确定；根据InSAR雷达系统方案和InSAR雷达成像观测方案，以上基于曲面投影的 InSAR成像方法需要的初始化成像系统参数均为已知；步骤2 =InSAR原始回波数据进行距离压缩采用传统的合成孔径雷达标准距离压缩方法对雷达平台主天线距离向回波数据瓦进行压缩，得到平台主天线距离压缩后数据，记做;采用传统的合成孔径雷达标准距离压缩方法对雷达平台副天线距离向回波数据瓦进行压缩，得到平台副天线距离压缩后数据，记做瓦;步骤3、计算InSAR成像空间中散射点的距尚史米用公式P,,,(H) = P,,,(O) + V."/PRFtY算得到雷达平台主天线第η个慢时刻的天线相位中心矢量 ,,,(《)，采用公式瓦(H) = Pi(0) + V-h/PRF 计算得到雷达平台副天线第η个慢时刻的天线相位中心矢量天⑷，其中η表示第η个慢时刻，n=l,. . . , Na, Na为步骤I初始化得到的雷达系统方位向采样点数，PRF为步骤I初始化得到的雷达系统脉冲重复频率， 为步骤I初始化得到的雷达平台速度矢量，P,,.10)为步骤I 初始化得到的雷达平台主天线初始位置矢量，F (O』力步骤I初始化得到的雷达平台副天线初始位置矢量；采用公式？(07)=瓦+(『-1)<忑+(^-1)<忑+ //(0^)忑计算得到场景采样点P(a，r)的位置矢量，其中a表示散射点位于场景方位向的第a个位置，a=l, . . . , sa, Sa 为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r=l，…，^为场景距离向总的散射点数，dr为步骤I初始化得到的场景距离向散射点间隔，da为步骤I初始化得到的场景方位向散射点间隔，瓦为步骤I初始化得到的场景参考点位置矢量，H(a，r) 为步骤I初始化得到的场景散射点P (a, r)的低精度高程，I为步骤I中定义的InSAR成像空间的第一个坐标基, ,为步骤I中定义的InSAR成像空间的第二个坐标基，&为步骤I中定义的InSAR成像空间的第三个坐标基；i为雷达平台距离散射点P(a，r)前后半个合成孔径内的慢时刻，i 的取值满足条件i = a-round(队/2)，…，a+round(Nx/2),若 a-round(队/2)〈I, 则 a-round (N1/2) =1, a+round (N1/2) >Na, JjllJ a+round (队/2) =Na,其中 round ( ·)为近似取整函数,N1为步骤I初始化得到的雷达平台一个合成孔径长度内的慢时刻个数，Na为步骤I 初始化得到的雷达系统方位向采样点数，即i只取雷达系统采样范围内的雷达平台距离散射点P (a, r)前后半个合成孔径内的慢时刻；采用公式=p>40- (^)|计算得到雷达平台主天线i时刻处散射点P (a, r)的距离史瓦 (;·; ，Γ)，采用公式兄.(/;fl，r) = ||瓦计算得到雷达平台副天线i时刻处散射点P (a，r)的距离史瓦 ( >，小其中Μ · ||2SL2范数；步骤4、距离压缩后数据插值重采样为雷达平台主天线i时刻处散射点P(a，r)的距离史5 ,dr)对应的距离门，计算公式为= I17 +mund(d(iw)-Rmj/dr)，I ，r)为雷达平台副天线i时刻处散射点 P (a, r)的距离史Rs (/； a,r)对应的距离门，计算公式为 ζ(/； a,r) = Ic + round({Kh (ia,r) — Rx )/,dr),其中iXnmcK ·)为近似取整函数为步骤3中得到的雷达平台主天线i时刻处散射点 P (a, r)的距离史，Λ(/;^为步骤3中得到的雷达平台副天线i时刻处散射点P(a，r)的距离史，Rm。为步骤I初始 化得到的场景参考点到雷达平台主天线各慢时刻天线相位中心的最短距离，Rs。为步骤I初始化得到的场景参考点到雷达平台副天线各慢时刻天线相位中心的最短距离，I。为步骤I初始化得到的场景参考点在雷达平台主天线回波数据和雷达平台副天线回波数据中的距离门位置，4为步骤I初始化得到的场景距离向散射点间隔，i为步骤 3中得到的雷达平台距离采样点P (a，r)前后半个合成孔径内的慢时刻，a表示散射点位于场景方位向的第a个位置，a=l, . . . , sa, sa为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r=l, . . . , sr, sr为场景距离向总的散射点数；从步骤2中得到的平台主天线距离压缩后数据民!的第i行中取第个数据的前Wtl 个数据和后Wtl个数据，采用标准辛格插值方法对这组数据进行插值，得到插值重采样后的数据,从 步骤2中得到的平台副天线距离压缩后数据I的第i行中取第 , 个数据的前Wtl个数据和后Wtl个数据，采用标准辛格插值方法对这组数据进行插值，得到插值重采样后的数据，其中Wtl为标准辛格插值的半个窗长；步骤5、插值重采样后数据相干求和根据补偿相位因子的计算公式i^';a，r) = exP^4_;r·民，得到散射点P(a，r) 在慢时刻i应补偿的相位因子，其中j为虚数单位(即-1开根)，瓦为步骤3 中得到的雷达平台主天线i时刻处散射点P (a，r)的距离史，i为步骤3中得到的雷达平台距离采样点P(a，r)前后半个合成孔径内的慢时刻，λ为步骤I初始化得到的雷达系统工作的信号波长，a表示散射点位于场景方位向的第a个位置，a=l，...，sa，Sa为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r=l,. . . , sr, sr为场景距离向总的散射点数；将步骤4中得到的插值重采样后的数据与散射点P(a，r)在慢时刻i应补偿的相位因子F(;>，r)相乘，得到相位补偿后的数据将步骤4中得到的插值重采样后的数据与散射点P(a，r)在慢时刻i应补偿的相位因子相乘，得到相位补偿后的数据丈^; 7);把雷达平台主天线距离散射点？(3，1·)前后半个合成孔径内的所有慢时刻的相位补偿后的数据I(Hr)相加，得到雷达平台主天线场景散射点P(a，r)成像后的数据S ,ΒΡ5- ^Γ) = ΣΤ ^α>Γ)，把雷达平台副天线距离散射点P (a，r)前后半个合成孔径内的所有慢时刻的相位补偿后的数据相加，得到雷达平台副天线场景散射点 P (a, r)成像后的数据民队小即民一卜写^犰小步骤6、干涉相位提取将步骤5中得到的雷达平台主天线场景散射点P(a，r)成像后的数据i, ( ，r)和雷达平台副天线场景散射点P(a，r)成像后的数据δ,(α^进行共轭相乘，得到InSAR成像的干涉相位 ，完成基于曲面投影的InSAR成像处理，其中a表示散射点位于场景方位向的第a个位置，a=l,. . . , sa, sa为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r=l,. . . , sr, sr为场景距离向总的散射点数。
全文摘要
本发明公开了一种基于曲面投影的InSAR成像方法，该方法通过在低精度的高程曲面上进行反向投影成像，将不同雷达天线的测量数据投影到以低精度高程曲面建立的InSAR成像空间中，提高了InSAR成像精度和干涉相位保持精度，同时提取干涉相位以后，陡变地形区域的干涉条纹稀疏，后续相位解缠的难度降低，相位解缠的精度提高，使得高程反演的精度提高，为InSAR高质量干涉相位生成提供了一种新方法。
文档编号G01S13/90GK103018740SQ20121047516
公开日2013年4月3日 申请日期2012年11月21日 优先权日2012年7月19日
发明者张晓玲, 付涛 申请人:电子科技大学