CN103616682A - 一种基于曲面投影的多基线InSAR处理方法 - Google Patents

一种基于曲面投影的多基线InSAR处理方法 Download PDF

Info

Publication number
CN103616682A
CN103616682A CN201310452744.2A CN201310452744A CN103616682A CN 103616682 A CN103616682 A CN 103616682A CN 201310452744 A CN201310452744 A CN 201310452744A CN 103616682 A CN103616682 A CN 103616682A
Authority
CN
China
Prior art keywords
note
antenna
phase
texas tower
done
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.)
Granted
Application number
CN201310452744.2A
Other languages
English (en)
Other versions
CN103616682B (zh
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.)
University of Electronic Science and Technology of China
Original Assignee
University of Electronic Science and Technology of China
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 University of Electronic Science and Technology of China filed Critical University of Electronic Science and Technology of China
Priority to CN201310452744.2A priority Critical patent/CN103616682B/zh
Publication of CN103616682A publication Critical patent/CN103616682A/zh
Application granted granted Critical
Publication of CN103616682B publication Critical patent/CN103616682B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

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
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/9021—SAR image post-processing techniques
    • G01S13/9023—SAR image post-processing techniques combined with interferometric techniques
    • 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

Abstract

本发明公开了一种基于曲面投影的多基线InSAR处理方法，它是通过多基线InSAR系统得到多个通道的测量数据，然后对短基线通道的测量数据进行成像和干涉处理得到低精度的数字高程曲面，使用低精度DEM建立成像空间曲面，对长基线通道的测量数据在低精度高程曲面上进行反向投影成像，经过干涉处理得到反演高程并与低精度数字高程相加，最终得到高精度数字高程。将长基线的测量数据投影到以低精度高程曲面建立的InSAR成像空间中，提高了InSAR成像精度和干涉相位保持精度，同时提取干涉相位以后，长基线干涉图中陡变地形区域的干涉条纹稀疏，后续相位解缠的难度降低，相位解缠的精度提高，使得高程反演的精度提高。

Description

一种基于曲面投影的多基线InSAR处理方法
技术领域
本发明属于雷达技术领域，它特别涉及合成孔径雷达（SAR）成像技术领域。
背景技术
合成孔径雷达（SAR）是一种高分辨率的微波成像系统。合成孔径雷达利用大时宽带宽信号实现距离向高分辨率，依靠雷达和目标之间的相对运动来合成虚拟阵列获取方位向高分辨率，可以对照射场景进行二维成像。合成孔径雷达干涉测量（InSAR）是一般SAR功能的延伸和扩展，是利用两个或者多个位置不同的天线观测同一个目标场景，根据目标到不同天线的斜距差获得测量数据的干涉相位，再通过平台与地面观测场景的几何关系反演出地面场景的数字高程信息的技术。由于具有全天时、全天候的特点，InSAR已经成为当前提取大面积地表三维图像和地形高程变化信息的一项重要遥感技术，在地形测绘、自然灾害监测和自然资源调查等领域发挥越来越大的作用。
高质量的干涉相位是InSAR获取高精度地形数字高程模型（DEM）的基础，随着高程测量精度的不断提高，InSAR数据处理对干涉相位质量的要求也越来越高，一方面需要提高InSAR成像的干涉相位保持精度，另一方面需要稀疏化干涉条纹，以降低后续InSAR数据处理中相位解缠的难度，提高高程反演精度。传统的单基线InSAR成像系统仅有两个天线，在提取干涉相位后得到一幅干涉图，干涉相位的表现形式是干涉条纹，通常陡变地形区域的干涉条纹密集，使得后续相位解缠的难度增加，同时相位解缠的精度也大大降低，导致最后的高程反演精度也明显降低（详见“谢先明，InSAR及多基线InSAR关键技术研究，电子科技大学博士学位论文，2011”）。多基线InSAR成像系统利用多个天线得到回波数据进行成像处理，在提取干涉相位后可以得到多幅干涉图。多组回波数据和干涉相位可以有效的提高干涉相位的精度，同时利用长短基线对应干涉相位的差异可以降低解缠难度，提高解缠精度，得到更高精度的反演高程。如何合理的利用多基线InSAR系统提供的多组测量数据提高干涉处理的精度是一个值得研究的问题。
反向投影算法是一种基于时域相干处理的成像算法，其基本思想是通过计算成像区域内每一散射点到合成孔径长度内雷达平台天线相位中心之间的双程时延，然后将对应的时域回波信号进行相干累加，从而恢复出每个散射点的散射系数信息。曲面投影算法是反向投影算法的改进算法，曲面投影算法通过在低精度的高程曲面上进行反向投影成像，精确计算投影曲面上每个散射点与一个合成孔径长度内雷达平台天线相位中心的斜距信息，实现每个散射点的精确聚焦。将曲面投影算法应用于InSAR数据成像，不仅可以减小由于参考点近似和平台运动轨迹误差对成像的影响，提高InSAR成像精度和干涉相位保持精度，而且可以稀疏化陡变地形区域的干涉条纹，使得提取的干涉相位易于相位解缠，同时提高相位解缠精度，从而提高高程反演精度。多基线InSAR系统可以得到多组测量数据和干涉相位，短基线干涉相位的干涉条纹稀疏，易于相位解缠但反演得到的数字高程精度低。长基线干涉相位的干涉条纹密集，相位解缠困难但是反演得到的数字高程精度高。利用短基线对应的干涉相位反演可以得到低精度的高程，在低精度高程曲面上利用曲面投影算法对长基线对应测量数据进行反向投影成像，最终可以得到高精度反演高程。
发明内容
为了得到高精度的InSAR干涉相位，同时稀疏化陡变地形区域的干涉条纹，本发明提出了一种基于曲面投影的多基线InSAR处理方法，本发明方法与传统的方法相比，提高了InSAR成像精度和干涉相位保持精度，同时提取干涉相位以后，长基线干涉图中陡变地形区域的干涉条纹稀疏，后续相位解缠的难度降低，相位解缠的精度提高，使得高程反演的精度提高。
为了方便描述本发明的内容，首先作以下术语定义：
定义1、合成孔径雷达干涉测量（InSAR）
合成孔径雷达干涉测量（InSAR）指利用两个或者两个以上的SAR数据中的相位信息进行相干处理，结合雷达参数和雷达几何位置信息反演地表三维及其变化信息的遥感技术，详见文献“合成孔径雷达成像原理”，皮亦鸣等编著，电子科技大学出版社出版。
定义2、多基线InSAR技术
多基线InSAR技术指利用两个以上的SAR测量数据中的相位信息进行相干处理，结合雷达参数和雷达几何位置信息反演地表三维及其变化信息的遥感技术
定义3、数字高程模型（DEM）
数字高程模型（Digital Elevation Model，DEM）是指利用一组有序数值阵列形式表示地表或地面高程的一种实体地面模型。本发明中DEM表示成一系列地面点的平面坐标X、Y和高程坐标Z组成的数据阵列。对于一个地面区域D，地形DEM表示为
DEM={Di|(xi,yi,zi),i∈D}
其中(xi，yi)是第i个地面像素点对应的平面坐标，zi是对应的高程坐标。
定义4、雷达成像空间
雷达成像空间是指将场景空间中的散射点投影到距离向－方位向－高度向的三维空间坐标系，该空间由合成孔径雷达成像空间中的三个相互正交的坐标基确定。目前典型的合成孔径雷达成像空间包括距离向－方位向－高度向投影空间。本发明中用以下数学关系表示成像空间M
M = { P - ( v , u ) | P - ( v , u ) = v · ς - v + u · ς - u + H ( v , u ) · ς - h , u , v ∈ R }
其中
Figure BDA0000388960630000032
和
Figure BDA0000388960630000033
表示构成成像空间M的三个相互正交的坐标基，分别表示距离向、方位向和高度向。
Figure BDA0000388960630000034
为成像空间中的采样点向量，u，v分别表示该点的距离向坐标和方位向坐标，H(v，u)表示该点的的高度向坐标，H(v，u)与距离向坐标u和方位向坐标v存在一一对应关系，R表示实数。
定义5、合成孔径雷达反向投影算法（BP）
反向投影算法是基于匹配滤波原理的合成孔径雷达成像算法，其主要通过相干累加实现合成孔径雷达数据的聚焦成像。详细内容可参考文献：“Research on A novel fast back projectionalgorithm for strip map bistatic SAR imaging”，Huang Yulin等。
定义6、合成孔径雷达曲面投影算法
曲面投影算法是反向投影算法的改进算法，是指利用低精度的场景高程数据，建立雷达成像空间坐标系，并基于此坐标系进行反向投影成像的算法。
定义7、合成孔径与慢时刻
合成孔径是指对于测绘场景中的一个散射点从进入雷达波束照射范围至离开雷达波束照射范围的这段时间内，雷达波束中心所走过的长度。
雷达平台飞过一个合成孔径所需要的时间称为慢时间，雷达系统以一定的重复周期Tr发射接收脉冲，因此慢时间可以表示为一个以重复周期Tr为步长的离散化时间变量，其中每一个离散时间变量值为一个慢时刻。
详见文献“合成孔径雷达成像原理”，皮亦鸣等编著，电子科技大学出版社出版。
定义8、距离史与距离门
距离史是指不同天线相位中心到场景中散射点的距离组成的序列。
距离门是指对应距离史的回波数据在整个回波数据中的位置，详见文献“合成孔径雷达成像原理”，皮亦鸣等编著，电子科技大学出版社出版。
定义9、合成孔径雷达成像场景参考点
合成孔径雷达成像场景参考点是指合成孔径雷达成像空间中的某个散射点，作为分析和处理场景中其他散射点的参照。
定义10、图像配准
图像配准是指将不同时间、不同传感器(成像设备)或不同条件下(天候、照度、摄像位置和角度等)对同一个地区获取的两幅或多幅图像进行地理坐标的匹配，使不同图像中的相同目标或者特征点位于图像中同一个位置的过程。
本发明中的图像配准特指将不同天线SAR图像空间中对应的相同目标进行匹配的过程，详见文献“合成孔径雷达成像原理”，皮亦鸣等编著，电子科技大学出版社出版。
定义11、天线相位中心
天线相位中心是指雷达天线向外辐射信号的中心，本发明中天线相位中心指雷达平台天线的轨迹位置。
定义12、合成孔径雷达原始回波仿真方法
合成孔径雷达原始回波仿真方法是指基于合成孔径雷达成像原理仿真出一定系统参数条件下具有SAR回波信号特性的原始信号的方法，详细内容可参考文献：“InSAR回波信号与系统仿真研究”，张剑琦，哈尔滨工业大学硕士论文。
定义13、干涉相位提取
干涉相位提取是指将经过配准后的两幅复图像对应的点共轭相乘得到干涉图的过程。详见文献“星载合成孔径雷达干涉成像”，王超等编著，科学出版社出版。
定义14、相位解缠
相位解缠是指将干涉相位由主值或相位差值恢复到真实值的过程，详见文献“星载合成孔径雷达干涉成像”，王超等编著，科学出版社出版。
定义15、相位滤波
相位滤波是指将干涉相位包含的相位噪声进行去除，得到真实干涉相位的过程。
定义16、高程反演
高程反演是指利用解缠相位、成像几何及轨道参数重建得到数字高程模型的过程，详见文献“星载合成孔径雷达干涉成像”，王超等编著，科学出版社出版。
本发明提供的一种基于曲面投影的多基线InSAR处理方法，它包括以下几个步骤：
步骤1、初始化多基线InSAR成像系统参数
InSAR成像空间由InSAR成像空间中的三个相互正交的坐标基确定，定义与雷达平台速度方向平行并在地平面内的单位向量作为InSAR成像空间的第一个坐标基，记做
Figure BDA0000388960630000051
该坐标基方向为方位向；定义在地平面内，并与InSAR成像空间的第一个坐标基
Figure BDA0000388960630000052
垂直的单位向量作为InSAR成像空间的第二个坐标基，记做
Figure BDA0000388960630000053
该坐标基方向为距离向；定义垂直于地平面向上的单位向量作为InSAR成像空间的第三个坐标基，记做
Figure BDA0000388960630000054
该坐标基方位为高度向。
InSAR雷达平台包含三组天线，即主天线、副天线1和副天线2，主天线和副天线1之间的距离为短基线长，记做B1，主天线和副天线2之间的距离为长基线长，记做B2，接收系统接收波门相对于发射信号发射波门的延迟，记为Td，主天线发射脉冲信号，经过Td时间的延迟，主天线、副天线1和副天线2同时接收回波延迟信号。雷达平台主天线接收的回波数据，记做雷达平台副天线1接收的回波数据，记做雷达平台副天线2接收的回波数据，记做
Figure BDA0000388960630000057
其中和均为二维矩阵，第一维均对应方位向，第二维均对应距离向，即二维矩阵
Figure BDA00003889606300000510
和
Figure BDA00003889606300000511
的行存储的是方位向数据，二维矩阵和
Figure BDA00003889606300000513
的列存储的是距离向数据。
初始化多基线InSAR成像系统参数包括：雷达系统工作的信号波长，记做λ，雷达平台主天线发射信号带宽，记做B，雷达平台主天线发射脉冲时宽，记做Tr，雷达平台接收系统采样频率，记做Fs,雷达系统脉冲重复频率，记做PEF，雷达平台一个合成孔径长度内的慢时刻个数，记做Nt，雷达平台速度矢量，记做
Figure BDA00003889606300000514
雷达平台主天线初始位置矢量，记做雷达平台副天线1初始位置矢量，记做雷达平台副天线2初始位置矢量，记做
Figure BDA00003889606300000517
场景参考点位置矢量，记做雷达系统距离向采样点数，记做Nr，雷达系统方位向采样点数，记做Na，场景距离向散射点间隔，记做dr，场景方位向散射点间隔，记做da，场景参考点到雷达平台主天线各慢时刻天线相位中心的最短距离，记做Rmc，场景参考点到雷达平台副天线1各慢时刻天线相位中心的最短距离，记做Rsc1，场景参考点到雷达平台副天线2各慢时刻天线相位中心的最短距离，记做Rsc2，场景参考点在雷达平台主天线回波数据和雷达平台副天线回波数据中的距离门相同，距离门位置记做Ic。对于场景散射点P(a，r)，a表示散射点位于场景方位向的第a个位置，a＝1，...，sa，sa为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r＝1，...，sr，sr为场景距离向总的散射点数。上述参数中，雷达系统工作的信号波长λ，雷达平台主天线发射的信号带宽B，雷达平台主天线发射的脉冲时宽Tr，雷达平台接收系统的采样频率Fs，雷达系统的脉冲重复频率PRF，主天线和副天线1之间的基线长B1，主天线和副天线2之间的基线长B2以及接收系统接收波门相对于发射信号发射波门的延迟Td在InSAR雷达系统设计过程中已经确定；雷达平台一个合成孔径长度内的慢时刻个数Nl，雷达平台速度矢量雷达平台主天线初始位置矢量
Figure BDA0000388960630000062
雷达平台副天线1初始位置矢量雷达平台副天线2初始位置矢量
Figure BDA0000388960630000064
场景参考点位置矢量
Figure BDA0000388960630000065
雷达系统距离向采样点数Nr，雷达系统方位向采样点数Na，场景距离向散射点间隔dr，场景方位向散射点间隔da，场景参考点到雷达平台主天线各慢时刻天线相位中心的最短距离Rmc，场景参考点到雷达平台副天线1各慢时刻天线相位中心的最短距离Rsc1，场景参考点到雷达平台副天线2各慢时刻天线相位中心的最短距离Rsc2，场景参考点在雷达平台主天线回波数据和雷达平台副天线回波数据中的距离门位置Ic，a表示散射点位于场景方位向的第a个位置，a＝1，...，sa，sa为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r1，...，sr，sr为场景距离向总的散射点数，在InSAR雷达成像观测方案设计中已经确定。根据多基线InSAR雷达系统方案和多基线InSAR雷达成像观测方案，以上基于反向投影的多基线InSAR成像方法需要的初始化成像系统参数均为已知。
步骤2：InSAR原始回波数据进行BP成像处理
采用传统的合成孔径雷达BP成像方法对步骤1中雷达平台主天线距离向回波数据
Figure BDA0000388960630000071
进行成像处理，得到平台主天线成像处理后数据，记做
采用传统的合成孔径雷达BP成像方法对步骤1中雷达平台副天线1距离向回波数据进行成像处理，得到平台副天线1成像处理后数据，记做
Figure BDA0000388960630000074
步骤3：相位滤波处理
将步骤2中平台主天线成像处理后数据
Figure BDA0000388960630000075
和平台副天线1成像处理后数据
Figure BDA0000388960630000076
进行共轭相乘，并采用传统的相位提取方法，得到短基线干涉相位，记做
采用传统的相位滤波方法对短基线干涉相位
Figure BDA0000388960630000078
进行降噪处理，得到短基线滤波后的干涉相位，记做
Figure BDA0000388960630000079
步骤4：相位解缠处理
采用传统的相位解缠方法对步骤3中短基线滤波后干涉相位
Figure BDA00003889606300000710
进行相位解缠处理，得到短基线相位解缠后相位，记做
步骤5：高程反演处理
采用传统的高程反演方法对步骤4中短基线相位解缠后相位
Figure BDA00003889606300000712
进行高程反演处理，得到短基线反演高程，记做
Figure BDA00003889606300000713
步骤6：建立曲面BP成像的坐标系
以步骤1中InSAR成像空间的第一个坐标基
Figure BDA00003889606300000714
作为曲面BP成像空间的第一个坐标基，
记做
Figure BDA0000388960630000081
该坐标基方向为方位向；以步骤1中InSAR成像空间的第二个坐标基
Figure BDA0000388960630000082
作为曲面BP成像空间的第二个坐标基，记做
Figure BDA0000388960630000083
该坐标基方向为距离向；将步骤5中短基线反演高程
Figure BDA0000388960630000084
作为曲面BP成像空间的第三个坐标基，记做
Figure BDA0000388960630000085
该坐标基方位为高度向。
步骤7：进行曲面BP成像处理
在步骤5中短基线反演高程
Figure BDA0000388960630000086
形成的曲面上对步骤1中雷达平台主天线距离向回波数据
Figure BDA0000388960630000087
进行反向投影成像即BP成像，得到平台主天线曲面投影成像处理后数据，记做
Figure BDA0000388960630000088
在步骤5中短基线反演高程
Figure BDA0000388960630000089
形成的曲面上对步骤1中雷达平台副天线2距离向回波数据
Figure BDA00003889606300000810
进行反向投影成像即BP成像，得到平台副天线2曲面投影成像处理后数据，记做
Figure BDA00003889606300000811
步骤8：相位滤波处理
将步骤7中平台主天线曲面投影成像处理后数据
Figure BDA00003889606300000812
和平台副天线2曲面投影成像处理后数据
Figure BDA00003889606300000813
进行共轭相乘，并采用传统的相位提取方法，得到长基线干涉相位，记做
Figure BDA00003889606300000814
采用传统的相位滤波方法对长基线干涉相位
Figure BDA00003889606300000815
进行降噪处理，得到长基线滤波后的干涉相位，记做
Figure BDA00003889606300000816
步骤9：相位解缠处理
采用传统的相位解缠方法对步骤8中长基线滤波后干涉相位
Figure BDA00003889606300000817
进行相位解缠处理，得到长基线相位解缠后相位，记做步骤10：高程反演处理
采用传统的高程反演方法对步骤9中长基线相位解缠后相位
Figure BDA0000388960630000091
进行高程反演处理，得到长基线反演高程，记做
Figure BDA0000388960630000092
用步骤5中短基线的反演高程
Figure BDA0000388960630000093
和长基线得到的反演高程相加得到最终反演高程，记做
Figure BDA0000388960630000095
本发明的创新点在于提出了一种基于曲面投影的多基线InSAR处理方法，该方法首先通过多基线InSAR系统得到多个通道的测量数据，然后对短基线通道的测量数据进行成像和干涉处理得到低精度的数字高程曲面，使用低精度DEM建立成像空间曲面，对长基线通道的测量数据在低精度高程曲面上进行反向投影成像，经过干涉处理得到反演高程并与低精度数字高程相加，最终得到高精度数字高程。将长基线的测量数据投影到以低精度高程曲面建立的InSAR成像空间中，提高了InSAR成像精度和干涉相位保持精度，同时提取干涉相位以后，长基线干涉图中陡变地形区域的干涉条纹稀疏，后续相位解缠的难度降低，相位解缠的精度提高，使得高程反演的精度提高。
本发明的优点在于针对多基线InSAR成像系统可以得到多个通道回波数据和干涉相位的优势，利用短基线对应干涉相位易于相位解缠和长基线干涉相位利于高程反演的特点，将长基线对应测量数据在短基线生成的低精度高程曲面上进行曲面投影处理，不仅提高了干涉相位的精度，同时稀疏了长基线对应干涉相位的干涉条纹，提高了相位解缠和高程反演的精度。
附图说明
图1为发明所提供方法的流程示意框图；
图2为InSAR成像空间中长基线实际测绘高程与短基线测绘低精度高程在距离向－高度向的二维坐标系中的关系图；
图3为本发明具体实施方式采用的多基线干涉合成孔径雷达平台飞行几何关系及仿真坡面场景图；
其中，横坐标（X轴）为方位向，纵坐标（Y轴）为距离向，垂直坐标（Z轴）为高度向，B1为雷达平台主天线和副天线1之间的基线，基线的两个端点分别为雷达平台主天线和副天线1，B2为雷达平台主天线和副天线2之间的基线，基线的两个端点分别为雷达平台主天线和副天线2，P(a，r)为场景散射点，a表示散射点位于场景方位向的第a个位置，a＝1，...，sa，sa为仿真场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r＝1，...，sr，sr为仿真场景距离向总的散射点数，da为仿真场景方位向散射点间隔，dr为仿真场景距离向散射点间隔，
Figure BDA0000388960630000101
为雷达平台速度矢量，
Figure BDA0000388960630000102
为雷达平台主天线第i个慢时刻的天线相位中心矢量，
Figure BDA0000388960630000103
为雷达平台副天线1第i个慢时刻的天线相位中心矢量，为雷达平台副天线2第i个慢时刻的天线相位中心矢量，为雷达平台主天线i时刻处散射点P(a，r)的距离史，
Figure BDA0000388960630000106
为雷达平台副天线1在i时刻处散射点P(a，r)的距离史，
Figure BDA0000388960630000107
为雷达平台副天线2在i时刻处散射点P(a，r)的距离史，i为雷达平台距离采样点P(a，r)前后半个合成孔径内的慢时刻，平面ABC′D′为InSAR成像空间中短基线对应低精度高程组成的坡面，平面ABCD为长基线对应高程曲面。
具体实施方式
本发明主要采用仿真实验的方法进行验证，所有步骤、结论都在MATLAB7.0上验证正确。具体实施步骤如下：
步骤1、初始化InSAR成像系统参数
选择与雷达平台速度方向平行并在地平面内的单位向量作为InSAR成像空间的第一个坐标基
Figure BDA0000388960630000108
该坐标基方向为方位向；选择InSAR成像空间的第二个坐标基
Figure BDA0000388960630000109
该坐标基方向为距离向；选择InSAR成像空间的第三个坐标基
Figure BDA00003889606300001010
该坐标基方向为高度向。
初始化多基线InSAR成像系统参数包括：雷达系统工作的信号波长λ＝0.03m，雷达平台主天线发射信号带宽B＝150MHz，雷达平台主天线发射脉冲时宽Tr＝1μs，雷达平台接收系统采样频率Fs＝300MHz,雷达系统脉冲重复频率PRF＝500Hz，雷达平台一个合成孔径长度内的慢时刻个数Nl＝500，雷达平台速度矢量
Figure BDA0000388960630000111
速度的单位为m/s，雷达平台主天线初始位置矢量雷达平台副天线1初始位置矢量
Figure BDA0000388960630000113
雷达平台副天线2初始位置矢量
Figure BDA0000388960630000114
场景参考点位置矢量
Figure BDA0000388960630000115
位置的单位为m，雷达系统距离向采样点数Nr＝1024，雷达系统方位向采样点数Na＝2048，场景距离向散射点间隔dr＝0.5m，场景方位向散射点间隔da＝0.3m，场景参考点到雷达平台主天线各慢时刻天线相位中心的最短距离Rmc＝1000m，场景参考点到雷达平台副天线1各慢时刻天线相位中心的最短距离Rsc1＝9992m，场景参考点到雷达平台副天线2各慢时刻天线相位中心的最短距离Rsc2＝9988m，场景参考点在雷达平台主天线回波数据和雷达平台副天线回波数据中的距离门Ic＝1。
仿真场景为一平面斜坡，仿真场景方位向总的散射点数sa＝512，仿真场景距离向总的散射点数sr＝512，斜坡高度由0m到50m，对于场景散射点P(a，r)，a表示散射点位于场景方位向的第a个位置，a＝1，...，sa，r表示散射点位于场景距离向的第r个位置，r＝1，...，sr，仿真坡面场景如图3所示。
采用传统的合成孔径雷达原始回波仿真方法，生成雷达平台主天线仿真回波数据
Figure BDA0000388960630000116
采用传统的合成孔径雷达原始回波仿真方法，生成雷达平台副天线1仿真回波数据
Figure BDA0000388960630000117
采用传统的合成孔径雷达原始回波仿真方法，生成雷达平台副天线2仿真回波数据
步骤2：InSAR原始回波数据进行BP成像处理
采用传统的合成孔径雷达BP成像方法对步骤1中雷达平台主天线距离向回波数据
Figure BDA0000388960630000119
进行成像处理，得到平台主天线成像处理后数据，记做
采用传统的合成孔径雷达BP成像方法对步骤1中雷达平台副天线1距离向回波数据
Figure BDA0000388960630000122
进行成像处理，得到平台副天线1成像处理后数据，记做
Figure BDA0000388960630000123
步骤3：相位滤波处理
将步骤2中平台主天线成像处理后数据
Figure BDA0000388960630000124
和平台副天线1成像处理后数据
Figure BDA0000388960630000125
进行共轭相乘，采用传统的相位提取方法，得到短基线干涉相位，记做
Figure BDA0000388960630000126
采用传统的相位滤波方法对短基线干涉相位
Figure BDA0000388960630000127
进行降噪处理，得到短基线滤波后的干涉相位，记做
Figure BDA0000388960630000128
步骤4：相位解缠处理
采用传统的相位解缠方法对步骤3中短基线滤波后干涉相位
Figure BDA0000388960630000129
进行相位解缠处理，得到短基线相位解缠后相位，记做
步骤5：高程反演处理
采用传统的高程反演方法对步骤4中短基线相位解缠后相位
Figure BDA00003889606300001211
进行高程反演处理，得到短基线反演高程，记做
步骤6：建立曲面BP成像的坐标系
以步骤1中InSAR成像空间的第一个坐标基作为曲面BP成像空间的第一个坐标基，记做
Figure BDA00003889606300001214
该坐标基方向为方位向；以步骤1中InSAR成像空间的第二个坐标基作为曲面BP成像空间的第二个坐标基，记做
Figure BDA0000388960630000131
该坐标基方向为距离向；将步骤5中短基线反演高程
Figure BDA0000388960630000132
作为曲面BP成像空间的第三个坐标基，记做该坐标基方位为高度向。
步骤7：进行曲面BP成像处理
在步骤5中高程
Figure BDA0000388960630000134
形成的曲面上对步骤1中雷达平台主天线距离向回波数据
Figure BDA0000388960630000135
进行反向投影成像，得到平台主天线曲面投影成像处理后数据，记做
Figure BDA0000388960630000136
在步骤5中高程
Figure BDA0000388960630000137
形成的曲面上对步骤1中雷达平台副天线2距离向回波数据
Figure BDA0000388960630000138
进行反向投影成像，得到平台副天线2曲面投影成像处理后数据，记做
Figure BDA0000388960630000139
步骤8：相位滤波处理
将步骤7中平台主天线曲面投影成像处理后数据
Figure BDA00003889606300001310
和平台副天线2曲面投影成像处理后数据
Figure BDA00003889606300001311
进行共轭相乘，采用传统的相位提取方法，得到长基线干涉相位，记做
Figure BDA00003889606300001312
采用传统的相位滤波方法对长基线干涉相位
Figure BDA00003889606300001313
进行降噪处理，得到长基线滤波后的干涉相位，记做
Figure BDA00003889606300001314
步骤9：相位解缠处理
采用传统的相位解缠方法对步骤8中长基线滤波后干涉相位进行相位解缠处理，得到长基线相位解缠后相位，记做
Figure BDA00003889606300001316
步骤10：高程反演处理
采用传统的高程反演方法对步骤9中长基线相位解缠后相位
Figure BDA00003889606300001317
进行高程反演处理，得到长基线反演高程，记做
Figure BDA0000388960630000141
用步骤5中短基线得到的反演高程
Figure BDA0000388960630000142
和长基线得到的反演高程
Figure BDA0000388960630000143
相加得到最终反演高程
Figure BDA0000388960630000144
通过本发明具体实施方式的仿真及测试，本发明所提供的一种基于曲面投影的多基线InSAR处理方法，对比传统多基线InSAR处理方法，不仅具有理想的干涉相位保持精度，还可以稀疏化长基线中陡变地形区域的干涉条纹，使得后续相位解缠的难度降低，同时提高相位解缠精度和高程反演精度，为多基线InSAR高精度成像提供了一种新方法。

Claims (1)

1.一种基于曲面投影的多基线InSAR处理方法，其特征是它包括以下几个步骤：
步骤1、初始化多基线InSAR成像系统参数
InSAR成像空间由InSAR成像空间中的三个相互正交的坐标基确定，定义与雷达平台速度方向平行并在地平面内的单位向量作为InSAR成像空间的第一个坐标基，记做
Figure FDA0000388960620000011
该坐标基方向为方位向；定义在地平面内，并与InSAR成像空间的第一个坐标基
Figure FDA0000388960620000012
垂直的单位向量作为InSAR成像空间的第二个坐标基，记做
Figure FDA0000388960620000013
该坐标基方向为距离向；定义垂直于地平面向上的单位向量作为InSAR成像空间的第三个坐标基，记做
Figure FDA0000388960620000014
该坐标基方位为高度向；
InSAR雷达平台包含三组天线，即主天线、副天线1和副天线2，主天线和副天线1之间的距离为短基线长，记做Ｂ1，主天线和副天线2之间的距离为长基线长，记做Ｂ2，接收系统接收波门相对于发射信号发射波门的延迟，记为Td，主天线发射脉冲信号，经过Ｔd时间的延迟，主天线、副天线1和副天线2同时接收回波延迟信号；雷达平台主天线接收的回波数据，记做雷达平台副天线1接收的回波数据，记做
Figure FDA0000388960620000016
雷达平台副天线2接收的回波数据，记做
Figure FDA0000388960620000017
其中
Figure FDA0000388960620000018
Figure FDA0000388960620000019
和
Figure FDA00003889606200000110
均为二维矩阵，第一维均对应方位向，第二维均对应距离向，即二维矩阵
Figure FDA00003889606200000111
Figure FDA00003889606200000112
和
Figure FDA00003889606200000113
的行存储的是方位向数据，二维矩阵
Figure FDA00003889606200000114
和
Figure FDA00003889606200000116
的列存储的是距离向数据；
初始化多基线InSAR成像系统参数包括：雷达系统工作的信号波长，记做λ，雷达平台主天线发射信号带宽，记做B，雷达平台主天线发射脉冲时宽，记做Tr，雷达平台接收系统采样频率，记做Fs,雷达系统脉冲重复频率，记做PRF，雷达平台一个合成孔径长度内的慢时刻个数，记做Nl，雷达平台速度矢量，记做雷达平台主天线初始位置矢量，记做
Figure FDA00003889606200000118
雷达平台副天线1初始位置矢量，记做雷达平台副天线2初始位置矢量，记做
Figure FDA00003889606200000120
场景参考点位置矢量，记做
Figure FDA0000388960620000021
雷达系统距离向采样点数，记做Nr，雷达系统方位向采样点数，记做Na，场景距离向散射点间隔，记做dr，场景方位向散射点间隔，记做da，场景参考点到雷达平台主天线各慢时刻天线相位中心的最短距离，记做Rmc，场景参考点到雷达平台副天线1各慢时刻天线相位中心的最短距离，记做Rsc1，场景参考点到雷达平台副天线2各慢时刻天线相位中心的最短距离，记做Rsc2，场景参考点在雷达平台主天线回波数据和雷达平台副天线回波数据中的距离门相同，距离门位置记做Ic；对于场景散射点P(a，r)，a表示散射点位于场景方位向的第a个位置，a＝1，...，sa，sa为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r＝1，...，sr，sr为场景距离向总的散射点数；上述参数中，雷达系统工作的信号波长λ，雷达平台主天线发射的信号带宽B，雷达平台主天线发射的脉冲时宽Tr，雷达平台接收系统的采样频率Fs，雷达系统的脉冲重复频率PRF，主天线和副天线1之间的基线长B1，主天线和副天线2之间的基线长B2以及接收系统接收波门相对于发射信号发射波门的延迟Td在InSAR雷达系统设计过程中已经确定；雷达平台一个合成孔径长度内的慢时刻个数Nl，雷达平台速度矢量
Figure FDA0000388960620000022
雷达平台主天线初始位置矢量
Figure FDA0000388960620000023
雷达平台副天线1初始位置矢量
Figure FDA0000388960620000024
雷达平台副天线2初始位置矢量
Figure FDA0000388960620000025
场景参考点位置矢量
Figure FDA0000388960620000026
雷达系统距离向采样点数Nr，雷达系统方位向采样点数Na，场景距离向散射点间隔dr，场景方位向散射点间隔da，场景参考点到雷达平台主天线各慢时刻天线相位中心的最短距离Rmc，场景参考点到雷达平台副天线1各慢时刻天线相位中心的最短距离Rsc1，场景参考点到雷达平台副天线2各慢时刻天线相位中心的最短距离Rsc2，场景参考点在雷达平台主天线回波数据和雷达平台副天线回波数据中的距离门位置Ic，a表示散射点位于场景方位向的第a个位置，a＝1，...，sa，sa为场景方位向总的散射点数，r表示散射点位于场景距离向的第r个位置，r＝1，...，sr，sr为场景距离向总的散射点数，在InSAR雷达成像观测方案设计中已经确定；根据多基线InSAR雷达系统方案和多基线InSAR雷达成像观测方案，以上基于反向投影的多基线InSAR成像方法需要的初始化成像系统参数均为已知；
步骤2：InSAR原始回波数据进行BP成像处理
采用传统的合成孔径雷达BP成像方法对步骤1中雷达平台主天线距离向回波数据
Figure FDA0000388960620000031
进行成像处理，得到平台主天线成像处理后数据，记做
Figure FDA0000388960620000032
采用传统的合成孔径雷达BP成像方法对步骤1中雷达平台副天线1距离向回波数据
Figure FDA0000388960620000033
进行成像处理，得到平台副天线1成像处理后数据，记做
Figure FDA0000388960620000034
步骤3：相位滤波处理
将步骤2中平台主天线成像处理后数据
Figure FDA0000388960620000035
和平台副天线1成像处理后数据
Figure FDA0000388960620000036
进行共轭相乘，并采用传统的相位提取方法，得到短基线干涉相位，记做
Figure FDA0000388960620000037
采用传统的相位滤波方法对短基线干涉相位
Figure FDA0000388960620000038
进行降噪处理，得到短基线滤波后的干涉相位，记做
Figure FDA0000388960620000039
步骤4：相位解缠处理
采用传统的相位解缠方法对步骤3中短基线滤波后干涉相位
Figure FDA00003889606200000310
得到短基线相位解缠后相位，记做
Figure FDA00003889606200000311
步骤5：高程反演处理
采用传统的高程反演方法对步骤4中短基线相位解缠后相位
Figure FDA00003889606200000312
进行高程反演处理，得到短基线反演高程，记做
Figure FDA00003889606200000313
步骤6：建立曲面BP成像的坐标系
以步骤1中InSAR成像空间的第一个坐标基
Figure FDA00003889606200000417
作为曲面BP成像空间的第一个坐标基，记做
Figure FDA0000388960620000041
该坐标基方向为方位向；以步骤1中InSAR成像空间的第二个坐标基
Figure FDA0000388960620000042
作为曲面BP成像空间的第二个坐标基，记做
Figure FDA0000388960620000043
该坐标基方向为距离向；将步骤5中短基线反演高程
Figure FDA00003889606200000419
Figure FDA0000388960620000045
该坐标基方位为高度向；
步骤7：进行曲面BP成像处理
在步骤5中短基线反演高程
Figure FDA0000388960620000046
形成的曲面上对步骤1中雷达平台主天线距离向回波数据
Figure FDA0000388960620000047
进行反向投影成像即BP成像，得到平台主天线曲面投影成像处理后数据，记做
Figure FDA0000388960620000048
在步骤5中短基线反演高程
Figure FDA0000388960620000049
形成的曲面上对步骤1中雷达平台副天线2距离向回波数据进行反向投影成像即BP成像，得到平台副天线2曲面投影成像处理后数据，记做
Figure FDA00003889606200000418
步骤8：相位滤波处理
将步骤7中平台主天线曲面投影成像处理后数据和平台副天线2曲面投影成像处理后数据
Figure FDA00003889606200000412
进行共轭相乘，并采用传统的相位提取方法，得到长基线干涉相位，记做
采用传统的相位滤波方法对长基线干涉相位
Figure FDA00003889606200000414
进行降噪处理，得到长基线滤波后的干涉相位，记做
Figure FDA00003889606200000415
步骤9：相位解缠处理
采用传统的相位解缠方法对步骤8中长基线滤波后干涉相位进行相位解缠处理，得到长基线相位解缠后相位，记做
Figure FDA0000388960620000051
步骤10：高程反演处理
采用传统的高程反演方法对步骤9中长基线相位解缠后相位
Figure FDA0000388960620000052
进行高程反演处理，得到长基线反演高程，记做
Figure FDA0000388960620000053
用步骤5中短基线的反演高程和长基线得到的反演高程
Figure FDA0000388960620000055
相加得到最终反演高程，记做
Figure FDA0000388960620000056
CN201310452744.2A 2013-09-27 2013-09-27 一种基于曲面投影的多基线InSAR处理方法 Active CN103616682B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201310452744.2A CN103616682B (zh) 2013-09-27 2013-09-27 一种基于曲面投影的多基线InSAR处理方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201310452744.2A CN103616682B (zh) 2013-09-27 2013-09-27 一种基于曲面投影的多基线InSAR处理方法

Publications (2)

Publication Number Publication Date
CN103616682A true CN103616682A (zh) 2014-03-05
CN103616682B CN103616682B (zh) 2015-11-18

Family

ID=50167386

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201310452744.2A Active CN103616682B (zh) 2013-09-27 2013-09-27 一种基于曲面投影的多基线InSAR处理方法

Country Status (1)

Country Link
CN (1) CN103616682B (zh)

Cited By (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104698457A (zh) * 2014-09-02 2015-06-10 电子科技大学 一种迭代曲面预测InSAR成像及高度估计方法
CN105182337A (zh) * 2015-09-09 2015-12-23 北京航空航天大学 一种基于曲面后向投影算法的形变反演方法
CN105842695A (zh) * 2016-03-25 2016-08-10 北京理工大学 一种基于混合基线的斜视InSAR动静目标分离方法
CN106249236A (zh) * 2016-07-12 2016-12-21 北京航空航天大学 一种星载InSAR长短基线图像联合配准方法
CN106707281A (zh) * 2017-01-05 2017-05-24 北京航空航天大学 一种基于多频数据处理的机载D‑InSAR形变检测方法
CN113640797A (zh) * 2021-08-09 2021-11-12 北京航空航天大学 一种用于参考条带模式InSAR的前斜视测高方法
CN113835090A (zh) * 2021-08-31 2021-12-24 电子科技大学 一种基于多通道sar系统的高精度干涉相位获取方法

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101866002A (zh) * 2010-06-01 2010-10-20 中国人民解放军信息工程大学 基于中国余数定理的多基线、多波段InSAR相位解缠方法
CN103018741A (zh) * 2012-12-11 2013-04-03 电子科技大学 一种基于后向投影的InSAR成像去平地一体化方法
CN103018740A (zh) * 2012-07-19 2013-04-03 电子科技大学 一种基于曲面投影的InSAR成像方法

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101866002A (zh) * 2010-06-01 2010-10-20 中国人民解放军信息工程大学 基于中国余数定理的多基线、多波段InSAR相位解缠方法
CN103018740A (zh) * 2012-07-19 2013-04-03 电子科技大学 一种基于曲面投影的InSAR成像方法
CN103018741A (zh) * 2012-12-11 2013-04-03 电子科技大学 一种基于后向投影的InSAR成像去平地一体化方法

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
JIE LI ET AL.: "A Multi-baseline InSAR DEM Reconstruction Approach without Ground Control Points", 《GEOSCIENCE AND REMOTE SENSING SYMPOSIUM, 2007. IGARSS 2007. IEEE INTERNATIONAL》 *
ZHANG HONGMIN ET AL.: "Simulation of Interferograms for Multi-Baseline InSAR Study", 《SYNTHETIC APERTURE RADAR, 2009. APSAR 2009. 2ND ASIAN-PACIFIC CONFERENCE ON》 *

Cited By (14)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN104698457A (zh) * 2014-09-02 2015-06-10 电子科技大学 一种迭代曲面预测InSAR成像及高度估计方法
CN104698457B (zh) * 2014-09-02 2017-02-15 电子科技大学 一种迭代曲面预测InSAR成像及高度估计方法
CN105182337B (zh) * 2015-09-09 2018-04-17 北京航空航天大学 一种基于曲面后向投影算法的形变反演方法
CN105182337A (zh) * 2015-09-09 2015-12-23 北京航空航天大学 一种基于曲面后向投影算法的形变反演方法
CN105842695A (zh) * 2016-03-25 2016-08-10 北京理工大学 一种基于混合基线的斜视InSAR动静目标分离方法
CN105842695B (zh) * 2016-03-25 2018-05-15 北京理工大学 一种基于混合基线的斜视InSAR动静目标分离方法
CN106249236A (zh) * 2016-07-12 2016-12-21 北京航空航天大学 一种星载InSAR长短基线图像联合配准方法
CN106249236B (zh) * 2016-07-12 2019-01-22 北京航空航天大学 一种星载InSAR长短基线图像联合配准方法
CN106707281A (zh) * 2017-01-05 2017-05-24 北京航空航天大学 一种基于多频数据处理的机载D‑InSAR形变检测方法
CN106707281B (zh) * 2017-01-05 2019-01-11 北京航空航天大学 一种基于多频数据处理的机载D-InSAR形变检测方法
CN113640797A (zh) * 2021-08-09 2021-11-12 北京航空航天大学 一种用于参考条带模式InSAR的前斜视测高方法
CN113640797B (zh) * 2021-08-09 2022-04-12 北京航空航天大学 一种用于参考条带模式InSAR的前斜视测高方法
CN113835090A (zh) * 2021-08-31 2021-12-24 电子科技大学 一种基于多通道sar系统的高精度干涉相位获取方法
CN113835090B (zh) * 2021-08-31 2024-04-12 电子科技大学 一种基于多通道sar系统的高精度干涉相位获取方法

Also Published As

Publication number Publication date
CN103616682B (zh) 2015-11-18

Similar Documents

Publication Publication Date Title
CN103018740B (zh) 一种基于曲面投影的InSAR成像方法
CN103616682B (zh) 一种基于曲面投影的多基线InSAR处理方法
CN103439693B (zh) 一种线阵sar稀疏重构成像与相位误差校正方法
KR101490981B1 (ko) 위성레이더 간섭도의 이온왜곡 보정방법 및 그 장치
Sansosti et al. Geometrical SAR image registration
CN103018741A (zh) 一种基于后向投影的InSAR成像去平地一体化方法
CN103698763B (zh) 基于硬阈值正交匹配追踪的线阵sar稀疏成像方法
Palm et al. Radargrammetric DEM extraction over urban area using circular SAR imagery
CN103885059B (zh) 一种多基线干涉合成孔径雷达三维重建方法
CN102788979B (zh) 一种基于后向投影InSAR成像配准的GPU实现方法
US20160259046A1 (en) Method and system for rendering a synthetic aperture radar image
CN103576137B (zh) 一种基于成像策略的多传感器多目标定位方法
CN103941243A (zh) 一种基于sar三维成像的自旋式飞行器测高方法
CN104698457A (zh) 一种迭代曲面预测InSAR成像及高度估计方法
CN105677942A (zh) 一种重复轨道星载自然场景sar复图像数据快速仿真方法
CN102662171A (zh) 一种sar层析三维成像方法
CN104007439B (zh) 一种干涉圆迹sar高程估计处理方法
CN103336278A (zh) 多视角观测下前视三维sar成像方法
Allstadt et al. Observations of seasonal and diurnal glacier velocities at Mount Rainier, Washington, using terrestrial radar interferometry
CN103630903B (zh) 基于顺轨干涉sar测量海面流场径向速度的方法
CN105447867A (zh) 基于isar图像的空间目标姿态估计方法
Anghel et al. Compact ground-based interferometric synthetic aperture radar: Short-range structural monitoring
CN103454633A (zh) 一种基于后向投影算法的干涉sar动基线处理方法
Floricioiu et al. Velocities of major outlet glaciers of the Patagonia Icefield observed by TerraSAR-X
CN104849769B (zh) 一种三维测距成像的系统及方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
C10 Entry into substantive examination
SE01 Entry into force of request for substantive examination
C14 Grant of patent or utility model
GR01 Patent grant