[go: up one dir, main page]

CN111239736B - 基于单基线的地表高程校正方法、装置、设备及存储介质 - Google Patents

基于单基线的地表高程校正方法、装置、设备及存储介质 Download PDF

Info

Publication number
CN111239736B
CN111239736B CN202010195430.9A CN202010195430A CN111239736B CN 111239736 B CN111239736 B CN 111239736B CN 202010195430 A CN202010195430 A CN 202010195430A CN 111239736 B CN111239736 B CN 111239736B
Authority
CN
China
Prior art keywords
sub
phase
aperture
full
view
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
Application number
CN202010195430.9A
Other languages
English (en)
Other versions
CN111239736A (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.)
Central South University
Original Assignee
Central South University
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 Central South University filed Critical Central South University
Priority to CN202010195430.9A priority Critical patent/CN111239736B/zh
Publication of CN111239736A publication Critical patent/CN111239736A/zh
Application granted granted Critical
Publication of CN111239736B publication Critical patent/CN111239736B/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
    • 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/882—Radar or analogous systems specially adapted for specific applications for altimeters
    • 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/885—Radar or analogous systems specially adapted for specific applications for ground probing
    • 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/9004—SAR image acquisition techniques
    • G01S13/9011—SAR image acquisition techniques with frequency domain processing of the SAR signals in azimuth
    • 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/40—Means for monitoring or calibrating

Landscapes

  • Engineering & Computer Science (AREA)
  • Remote Sensing (AREA)
  • Radar, Positioning & Navigation (AREA)
  • Physics & Mathematics (AREA)
  • Electromagnetism (AREA)
  • Computer Networks & Wireless Communication (AREA)
  • General Physics & Mathematics (AREA)
  • Signal Processing (AREA)
  • Radar Systems Or Details Thereof (AREA)

Abstract

本发明公开了一种基于单基线的地表高程校正方法、装置、设备及存储介质，其中方法为：获取研究区域的两个全孔径影像，采用子孔径分解技术获取每个全孔径影像分解为两个子视影像；对全孔径影像和子视影像预处理得到对应的干涉图；对两个子视干涉图均使用DEM数据进行去地形处理，以及去平地和轨道误差相位，得到包括相同对流层延迟误差的两个子视差分干涉图；采用小波分解技术，对两个子视差分干涉图的高频小波系数进行相关性分析，提取对流层延迟误差相位；将得到的对流层延迟误差相位从全孔径干涉图中去除，再根据卫星轨道参数计算得到地表高程数据。本发明仅需使用两景全孔径影像即可提取对流层延迟误差相位，并且得到高精度的地表高程数据。

Description

基于单基线的地表高程校正方法、装置、设备及存储介质
技术领域
本发明涉及合成孔径雷达干涉测量(InSAR)技术在大范围地形测量中的应用，属于星载合成孔径雷达提取地表高程技术领域，具体是指一种基于单基线的地表高程校正方法、装置、设备及存储介质。
背景技术
数字高程模型(digital elevation model,DEM)是国家经济建设、社会发展和国家安全必不可少的基础图件，是基础设施规划、设计和施工及资源调查、开发所必须的基础资料。合成孔径雷达(InSAR)技术，具备全天候、全天时、观测周期短、大范围的优势，成为了一种获取全球DEM的强有力工具。然而，在进行InSAR对地观测时，雷达信号会受到多种因素的影响，这也误差严重影响了采用InSAR技术进行地形测量的精度。对于目前常用的重轨星载SAR系统而言，对流层延迟误差是影响其测高精度最主要的误差源之一。因此，为了获取大范围、高精度的InSAR数字高程模型，必须对干涉相位中包含的对流层延迟误差进行有效的估计和去除。
目前，InSAR测量中常用的对流层延迟误差改正方法主要包括以下三类：1)相位-高程的相关性分析分析法；2)外部水汽数据辅助法；3)时序InSAR分析法。其中，第一类方法通常需要依赖于一定的先验信息，即对流层延迟误差相位与地面高度呈现一定的相关性，且该方法只能去除与地形相关的对流层延迟误差相位，当干涉图中与地形不相关的湍流误差占优时，该方法不能得到较好的改正结果；相较于第一类方法，第二类方法要求SAR数据与外部水汽数据之间具有较好的时空一致性，这就限制了该方法的推广使用；与前两类方法不同的是，为了使用第三类方法得到较好的改正结果，需要覆盖同一地区的大量SAR数据量。然而，在实际应用中，由于外部数据的缺乏或SAR数据量较少，上述方法很难取得较好的对流层延迟误差改正效果。
因此，有必要设计一种无需依赖外部水汽辅助数据或过量SAR数据的对流层延迟误差校正方法。
发明内容
本发明所要解决的技术问题在于，提供一种基于单基线的地表高程校正方法、装置、设备及存储介质，通过提取对流层延迟误差相位以对地表相位进行校准，从而得到更高精度的地表高程数据。
为实现上述技术目的，本发明采用如下技术方案：
一种基于单基线的地表高程校正方法，包括以下步骤：
步骤1，获取研究区域的两个全孔径SAR影像，采用子孔径分解技术将每个全孔径SAR影像分解为两个与方位向入射角对应的子视影像；
步骤2，对两个全孔径SAR影像进行预处理得到全孔径干涉图；对同一方位向入射角的两个子视影像进行预处理，得到与方位向入射角对应的子视干涉图；
步骤3，对两个子视干涉图，均使用DEM数据并采用差分干涉技术进行去地形处理，然后再去除平地相位和轨道误差相位处理，得到两个如以下公式所示的子视差分干涉图：
Figure BDA0002417446130000021
式中，α1和α2代表两个不同的方位向入射角，φdiff为子视差分干涉图中的差分干涉相位，△φtopo为地形误差相位，φatmo为对流层延迟误差相位，φnoise为随机噪声误差相位；
步骤4，对两个子视差分干涉图进行小波分解，对两个子视差分干涉图在每个分解尺度下的高频小波系数进行相关性分析，得到对应分解尺度下两个高频小波系数之间的相关值，根据相关值和任一高频小波系数计算得到对流层延迟误差相位相关的小波系数；然后对所述与对流层延迟误差相位相关的小波系数进行小波逆变换，即可得到对流层延迟误差相位；
步骤5，对全孔径干涉图进行解缠处理，将步骤4得到的对流层延迟误差相位从全孔径干涉图的解缠相位中去除；最后，根据卫星轨道参数以及去除对流层延迟误差相位的全孔径干涉图解缠相位，计算得到地表高程数据。
在更优的技术方案中，对两个子视干涉图，使用两种由不同技术获得的DEM数据进行去地形处理。
在更优的技术方案中，对子视差分干涉图进行小波分解表示为：
Figure BDA0002417446130000022
式中，m和n分别表示干涉图的行数和列数，x,y是像素点在雷达坐标系的坐标，ix和iy分别表示干涉图的第ix行和第iy列，Φ和Ψ分别表示小波分解的平移函数和母小波函数，J表示最优分解尺度，j代表不同的分解尺度，v,w分别表示低频小波系数和高频小波系数，ε＝1,2,3分别对应水平、垂直和对角方向；
步骤4中的相关性分析方法为，取两个子视差分干涉图在每个分解尺度下的高频小波系数
Figure BDA0002417446130000031
和
Figure BDA0002417446130000032
按以下公式进行相关值计算：
Figure BDA0002417446130000033
在小波分解得到ε＝1,2,3时的高频小波系数
Figure BDA0002417446130000034
和
Figure BDA0002417446130000035
时，代入上述相关性计算公式中，即可通过方程组计算得到两个子视差分干涉图在分解尺度j时的相关值fj和整体偏移系数cj；
再使用两个子视差分干涉图在分解尺度j时的相关性fj和整体偏移系数cj，按以下公式计算对流层延迟误差相位对应的高频小波系数：
Figure BDA0002417446130000036
式中，
Figure BDA0002417446130000037
代表对流层延迟误差相位所对应的高频小波系数；
Figure BDA0002417446130000038
表示方位向入射角为αi子视差分干涉图在分解尺度j时高频小波系数。
在更优的技术方案中，去除平地相位时，平地相位的计算公式为：
Figure BDA0002417446130000039
式中，φflat为平地相位，λ为雷达信号波长，B//是平行基线长度。
在更优的技术方案中，轨道误差相位的去除方法为：使用融入地表高程h的多项式模型拟合轨道误差相位，通过最小二乘法解算多项式模型中的未知参数，然后将轨道误差相位从子视差分干涉图中对应像素点的解缠相位中去除，得到新的子视差分干涉图。
最后，将轨道误差相位φorbit从子视差分干涉图中对应像素点的解缠相位中去除。
在更优的技术方案中，步骤1采用子孔径分解技术将全孔径SAR影像分解的过程为：
首先，对全孔径SAR影像进行一维傅利叶变换，得到SAR影像在方位向上的频谱；
然后，对SAR影像在方位向上的频谱按不同的方位向入射角进行分割，得到与方位向入射角对应的频谱；
最后，对所有与方位向入射角对应的频谱均进行傅里叶逆变换，得到与方位向入射角对应的影像，即为与方位向入射角对应的子视影像。
在更优的技术方案中，对全孔径SAR影像进行预处理包括：对2个全孔径SAR影像进行共配准，然后逐像素共轭相乘生成全孔径干涉图；
对主、辅影像每个方位向入射角对应的2个子视影像进行预处理包括：将当前方位向入射角对应的2个子视影像逐像素共轭相乘得到当前方位向入射角对应的子视干涉图。
本发明还提供一种基于单基线的地表高程校正装置，包括：
子孔径分解模块，用于：获取两个全孔径SAR影像，采用子孔径分解技术将每个全孔径SAR影像分解为两个与方位向入射角对应的子视影像；
数据预处理模块，用于：对两个全孔径SAR影像进行预处理得到全孔径干涉图；对同一方位向入射角的两个子视影像进行预处理，得到与方位向入射角对应的子视干涉图；
差分干涉相位建模模块，用于：对两个子视干涉图，均使用DEM数据并采用差分干涉技术进行去地形处理，然后再去除平地相位和轨道误差相位处理，得到两个如以下公式所示的子视差分干涉图：
Figure BDA0002417446130000041
式中，α1和α2分别代表两个不同的方位向入射角，φdiff为子视差分干涉图中的差分干涉相位，△φtopo为地形误差相位，φatmo为对流层延迟误差相位，φnoise为随机噪声误差相位；
对流层延迟误差相位提取模块，用于：对两个子视差分干涉图进行小波分解，对两个子视差分干涉图在每个分解尺度下的高频小波系数进行相关性分析，得到对应分解尺度下两个高频小波系数之间的相关值，根据相关值和任一高频小波系数计算得到对流层延迟误差相位相关的小波系数；然后对所述与对流层延迟误差相位相关的小波系数进行小波逆变换，即可得到对流层延迟误差相位；
地表高程数据校正模块，用于：对全孔径干涉图进行解缠处理，将得到的对流层延迟误差相位从全孔径干涉图的解缠相位中去除；最后，根据卫星轨道参数以及去除对流层延迟误差相位的全孔径干涉图解缠相位，计算得到地表高程数据。
本发明还提供一种设备，包括处理器和存储器；其中：所述存储器用于存储计算机指令；所述处理器用于执行所述存储器存储的计算机指令，具体执行上述任一所述的方法。
本发明还提供一种计算机存储介质，用于存储程序，所述程序被执行时，用于实现上述任一所述的方法。
有益效果
本发明一种基于子孔径相关性分析的对流层延迟误差改正方法，其优点是：
(1)本发明基于子孔径分解技术获取不同方位向入射角下的子视干涉图，通过相关性分析提取相同的对流层延迟误差相位并予以剔除，削弱对流层延迟误差对InSAR测高的影响，从而提高了InSAR技术的测高精度；
(2)本发明无需过多SAR数据量，只需两景全孔径SAR影像形成单基线即可，且无需外部水汽辅助数据的情况下，通过对流层延迟误差相位在与方位向入射角的独立性，即可获取对流层延迟相位的分布特征；
(3)本发明利用相关性分析，从两个子视差分干涉图中提取相同的对流层延迟误差相位成分，即可提取消除干涉图中混叠的对流层延迟误差相位，相较于传统方法，本发明可以同时去除与地形相关的垂直分层大气以及与地面高度不相关的湍流大气延迟成分，从而极大地拓宽了InSAR测高的应用范围。
附图说明
图1为本发明实施例所述方法的流程图。
图2为本发明实施例所述子孔径相关性分析对流层延迟误差改正算法详细流程；
图3为本发明实施例所述SAR影像子孔径分解的示意图。
图4为本发明实施例所述方法得到的DEM结果。
具体实施方式
为了更好的说明本发明的方法与步骤，利用位于洛杉矶地区两景间隔92天的ALOS1PALSAR雷达数据，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施仅用以解释本发明，并不用于限定本发明。
本发明实施例提供一种基于单基线的地表高程校正方法，如图1、2所示，包括以下步骤：
步骤1，获取研究区域的两个全孔径SAR影像，采用子孔径分解技术将每个全孔径SAR影像分解为两个与方位向入射角对应的子视影像，基本原理如图3所示；每个全孔径SAR影像对应得到两个子视影像，称为前视影像和后视影像。
其中的两个全孔径SAR影像：1个全孔径SAR影像为2010年7月11日获取的ALOS1PALSAR雷达数据，本实施例中以该影像为主影像；另1个全孔径SAR影像为2010年10月11日获取的ALOS1 PALSAR雷达数据，本实施例以该影像为辅影像。
采用子孔径分解技术将全孔径SAR影像分解的过程具体为：
首先，对全孔径SAR影像进行一维傅利叶变换，得到SAR影像在方位向上的频谱；
然后，对SAR影像在方位向上的频谱按不同的方位向入射角进行分割，得到与方位向入射角对应的频谱；
最后，对所有与方位向入射角对应的频谱均进行傅里叶逆变换，得到与方位向入射角对应的影像，即为与方位向入射角对应的子视影像。
在具体实施中，采用子孔径分解技术对每个全孔径SAR影像均分解为两个不同方位向入射角的子视影像，具体通过采用两个不同方位向入射角的带通滤波器，生成不同方位向观测角度的子视影像。
步骤2，对两个全孔径SAR影像进行预处理得到全孔径干涉图；对同一方位向入射角的两个子视影像进行预处理，得到与方位向入射角对应的子视干涉图；
对全孔径SAR影像进行预处理包括：对2个全孔径SAR影像进行共配准，然后逐像素共轭相乘生成全孔径干涉图；
对每个方位向入射角的两个子视影像进行预处理包括：将当前方位向入射角对应的2个子视影像进行干涉处理，得到当前方位向入射角对应的子视干涉图。
步骤3，对两个子视干涉图，均使用DEM数据并采用差分干涉技术进行去地形处理，然后再去除平地相位和轨道误差相位处理，得到两个如以下公式所示的子视差分干涉图：
Figure BDA0002417446130000063
式中，α1和α2代表两个不同的方位向入射角，φdiff为子视差分干涉图中的差分干涉相位，△φtopo为地形误差相位，φatmo为对流层延迟误差相位，φnoise为随机噪声误差相位。
具体地，星载重复轨道InSAR在进行对地观测时，会受到大气、地表形变等多种因素的影响，因此，InSAR干涉相位表现为多种相位成分的叠加。
为了保证干涉图的相位质量，本实施例选择时间基线较短的两个SAR影像，地表形变几乎为零，因此可以不考虑地表形变的影响；此外，相较于对流层延迟而言，电离层延迟的量级一般都比较小，所以忽略电离层延迟的影响。因此，本发明实施例针对InSAR地形测量，在不考虑地表形变和电离层延迟的影响下，将不同观测角度下干涉图的相位φint可以表示为：
Figure BDA0002417446130000061
式中：
α1和α2代表两个不同的方位向入射角；
φtopo是地形测量相关的地形相位，其与裸地表高程h和地表散射目标的高度△h有关；
φflat为平地相位，可以表示为：
Figure BDA0002417446130000062
可以根据卫星的轨道参数将其进行有效的去除；λ是雷达波长，B//是平行基线；
φorbit为轨道误差相位，是由不精确的轨道参数导致，在常规InSAR处理中，通常采用多项式模型拟合方法去除；
φnoise为随机噪声误差相位，可以采用干涉图滤波处理将其削弱；
φatmo为对流层延迟误差相位，可以表示为：
Figure BDA0002417446130000071
主要产生原因为两次成像过程中微波信号受到不同水汽环境的干扰，干涉后以相位的形式体现在干涉图中。
根据以上对干涉图中相位成分的分析，干涉图中只有对流层延迟误差很难通过固定的模型进行有效的去除，其成为了星载重轨SAR地形测量的主要误差来源。
目前已有研究表明，方位向入射角的变化会导致后向散射特性的变化，因此不同子视干涉图具有不同的干涉相位中心。另外，由于不同方位向入射角下的子视影像共用了近似相同的轨道参数，由上述平地相位和轨道误差相位的分析可得，各子视干涉相位中包含了相同的平地相位和轨道误差相位。由于SAR传感器合成孔径时间一般都很短(一般为几秒)，并且考虑到星载SAR传感器较小的波束宽度角，可以认为不同子孔径雷达波在传输过程中穿过的大气成分近似相同，进而导致不同子孔径影像含有相同的对流层延迟误差，这些相同的对流层延迟误差经过干涉处理，就转换为相同的对流层延迟误差相位。尽管各孔径影像只包含全孔径SAR影像频谱的一部分，但是其包含的对流层延迟误差相位是相同的。该特征为相同对流层延迟误差的提取和去除提供了可能。
在本步骤3中，为了提取相同的对流层延迟误差相位，需要先将各子视干涉图中的裸地高程对应的地形相位、平地相位和轨道误差相位成分进行去除。
首先采用差分干涉技术进行去地形处理，即采用外部DEM数据模拟地形相位，将地形相位从干涉相位中去除：
Figure BDA0002417446130000072
Figure BDA0002417446130000073
式中，
Figure BDA0002417446130000074
为采用外部DEM数据h模拟得到的地形相位，△φtopo为地形误差相位，由外部DEM数据的不准确所导致；B⊥代表垂直基线的长度，λ为雷达信号波长，R1和R2分别表示主、辅影像的天线中心到地面目标点的距离。
尽管不同方位向观测角度下的干涉图对应着不同的相位中心，在差分干涉进行去地形处理过程中，如果采用同一个外部DEM数据来去除不同极化干涉图中的地形相位，会导致各子视影像对应的子视差分干涉图中含有相似的地形误差相位，这将严重干扰后续相同对流层延迟误差相位的提取。因此，本实施例分别采用不同技术手段获取的外部DEM数据进行去地形处理，分别去除相应子视干涉图中的地形相位，可以增强最终得到的两个子视差分干涉图的区分度，从而极大提高对流层延迟误差的估计精度。
然后，使用平地相位计算公式
Figure BDA0002417446130000081
根据卫星的轨道参数λ、R2、R1计算平地相位，进而将平地相位从子视干涉相位中去除：
Figure BDA0002417446130000082
再后，使用融入地表高程h的多项式模型拟合轨道误差相位，通过最小二乘法解算多项式中的未知参数，然后将轨道误差相位从子视差分干涉图中对应像素点的解缠相位中去除，得到新的子视差分干涉图。
具体在本实施例中，使用多项式a0+a1x+a2y+a3h拟合轨道误差相位φorbit，得到轨道误差相位模型为：
Figure BDA0002417446130000083
其中x,y是像素点在雷达坐标系的坐标，hi为第i个像素对应的地表高程，a1,a2,a3为待求参数；
为了准确估计未知参数a1,a2,a3的值，可以采用上述轨道误差相位模型对解缠后的子视差分干涉相位进行拟合，然后采用最小二乘方法解算未知参数a1,a2,a3；随后即可以根据出每一像素对应的坐标x,y与地表高程hi，通过上述轨道误差相位模型计算其轨道误差相位φorbit。
另外需要注意的是：对于非线性的轨道误差相位，在实际应用中可以通过选择高次多项式来替代本实施例中的多项a0+a1x+a2y+a3h进行计算。
将轨道误差相位φorbit从子视差分干涉图中对应像素点的解缠相位中去除，即可得到新的子视差分干涉图：
Figure BDA0002417446130000084
其中，φdiff为子视差分干涉图中的差分干涉相位。
步骤4，对两个子视差分干涉图进行小波分解，对两个子视差分干涉图在每个分解尺度下的高频小波系数进行相关性分析，得到对应分解尺度下两个高频小波系数之间的相关值，根据相关值和任一高频小波系数计算得到对流层延迟误差相位相关的小波系数；然后对所述与对流层延迟误差相位相关的小波系数进行小波逆变换，即可得到对流层延迟误差相位。
对于对流层延迟误差在空间上分布的复杂性，使得我们难以通过简单的空间滤波将其进行有效的分离。然而，差分干涉图中不同的相位成分在频率域表现为了不同的波长成分，因此小波变换可以将差分干涉图分解为不同的波长成分，允许我们在不同尺度下对其进行研究。因此，本发明实施例采用小波变换将差分干涉图分解到不同的尺度，并在每一分解尺度下估计对流层延迟误差相位的大小。
对子视差分干涉图进行小波分解表示为：
Figure BDA0002417446130000091
式中，m和n分别表示干涉图的行数和列数，x,y是像素点在雷达坐标系的坐标，ix和iy分别表示干涉图的第ix行和第iy列，Φ和Ψ分别表示小波分解的平移函数和母小波函数，J表示最优分解尺度，j代表不同的分解尺度，v,w分别表示低频小波系数和高频小波系数，ε＝1,2,3分别对应水平、垂直和对角方向；
在两个子视差分干涉图中，如前述分析对流层延迟误差相位相同，且对流层延迟误差在频率域通常表现为长波长的高频信息，因此它们的高频小波系数间具有一定的相似性，然后对不同子视干涉图每一分解尺度下的高频系数(分别表示为
Figure BDA0002417446130000092
)分别进行如下的整体相关性分析计算：
Figure BDA0002417446130000093
在小波分解得到ε＝1,2,3时的高频小波系数
Figure BDA0002417446130000094
和
Figure BDA0002417446130000095
时，代入上述相关性计算公式中，即可通过方程组计算得到两个子视差分干涉图在分解尺度j时的相关值fj和整体偏移系数cj；
再使用两个子视差分干涉图在分解尺度j时的相关性fj和整体偏移系数cj，按以下公式计算对流层延迟误差相位对应的高频小波系数：
Figure BDA0002417446130000096
式中，
Figure BDA0002417446130000097
代表对流层延迟误差相位所对应的高频小波系数；
Figure BDA0002417446130000098
表示方位向入射角为αi子视差分干涉图的在分解尺度j时高频小波系数。对流层延迟误差相位对应的高频小波系数时，具体可取两个方位向入射角中任意一个对应的子视差分干涉图得到的高频小波系数
Figure BDA0002417446130000101
进行计算。
最后，为了提取子视差分干涉图中包含的对流层延迟误差相位，对上述得到的与对流层延迟误差相位相关的小波系数进行小波逆变换，即可得到全孔径干涉图相位中包含的对流层延迟误差相位。
步骤5，对全孔径干涉图进行解缠处理，将步骤4得到的对流层延迟误差相位从全孔径干涉图的解缠相位中去除；最后，根据卫星轨道参数以及去除对流层延迟误差相位的全孔径干涉图解缠相位，按下以公式计算得到地表高程数据h：
Figure BDA0002417446130000102
式中，θ为主影像各像素点的天线入射角，α为基线倾角，φ'topo表示去除对流层延迟误差相位的全孔径干涉图解缠相位。
在更优的实施例中，还包括对全孔径干涉图解缠相位φ'topo进行噪声滤除处理，即去除随机噪声误差相位，再将随机噪声误差相位和对流层延迟误差相位从全孔径干涉图的解缠相位中去除，将此时的地形相位记为φ”topo，然后使用公式计算得到更高精度的地表高程数据：
Figure BDA0002417446130000103
如图4所示，其中图4(a)-(c)分别表示直接使用全孔径SAR影像处理(未进行对流层延迟误差改正)得到的地表高程数据InSAR DEM1、使用本发明实施例去除对流层延迟误差相位后再计算得到的地表高程数据InSAR DEM2以及USGS提供的地表高程数据SRTM DEM；然后将SRTM DEM作为参考数据分别对InSAR DEM1和InSAR DEM2进行了定性分析，得到InSAR DEM1和InSAR DEM2的高程误差分别如图4(d)和4(e)分别所示，可以看到本发明实施方法可以有效的去除InSAR地形测量中的对流层延迟误差。
本发明还提供一种基于单基线的地表高程校正装置，与上述方法实施例相对应，包括：
子孔径分解模块，用于：获取两个全孔径SAR影像，采用子孔径分解技术将每个全孔径SAR影像分解为两个与方位向入射角对应的子视影像；
数据预处理模块，用于：对两个全孔径SAR影像进行预处理得到全孔径干涉图；对同一方位向入射角的两个子视影像进行预处理，得到与方位向入射角对应的子视干涉图；
差分干涉相位建模模块，用于：对两个子视干涉图，均使用DEM数据并采用差分干涉技术进行去地形处理，然后再去除平地相位和轨道误差相位处理，得到两个如以下公式所示的子视差分干涉图：
Figure BDA0002417446130000111
式中，α1和α2分别代表两个不同的方位向入射角，φdiff为子视差分干涉图中的差分干涉相位，△φtopo为地形误差相位，φatmo为对流层延迟误差相位，φnoise为随机噪声误差相位；
对流层延迟误差相位提取模块，用于：对两个子视差分干涉图进行小波分解，对两个子视差分干涉图在每个分解尺度下的高频小波系数进行相关性分析，取相关性大于预设阈值的高频小波系数，均作为对流层延迟误差相位相关的小波系数；然后对所述与对流层延迟误差相位相关的小波系数进行小波逆变换，即可得到对流层延迟误差相位；
地表高程数据校正模块，用于：将得到的对流层延迟误差相位从全孔径干涉图中去除；最后，根据卫星轨道参数以及去除对流层延迟误差相位的全孔径干涉图，计算得到地表高程数据。
本发明还提供一种设备，包括处理器和存储器；其中：所述存储器用于存储计算机指令；所述处理器用于执行所述存储器存储的计算机指令，具体执行上述方法实施例所述的步骤。
本发明还提供一种计算机存储介质，用于存储程序，所述程序被执行时，用于实现上述方法实施例所述的方法。
以上实施例为本申请的优选实施例，本领域的普通技术人员还可以在此基础上进行各种变换或改进，在不脱离本申请总的构思的前提下，这些变换或改进都应当属于本申请要求保护的范围之内。

Claims (10)

1.一种基于单基线的地表高程校正方法，其特征在于，包括以下步骤：
步骤1，获取研究区域的两个全孔径SAR影像，采用子孔径分解技术将每个全孔径SAR影像分解为两个与方位向入射角对应的子视影像；
步骤2，对两个全孔径SAR影像进行预处理得到全孔径干涉图；对同一方位向入射角的两个子视影像进行预处理，得到与方位向入射角对应的子视干涉图；
步骤3，对两个子视干涉图，均使用DEM数据并采用差分干涉技术进行去地形处理，然后再去除平地相位和轨道误差相位处理，得到两个如以下公式所示的子视差分干涉图：
Figure FDA0002417446120000011
式中，α1和α2代表两个不同的方位向入射角，φdiff为子视差分干涉图中的差分干涉相位，△φtopo为地形误差相位，φatmo为对流层延迟误差相位，φnoise为随机噪声误差相位；
步骤4，对两个子视差分干涉图进行小波分解，对两个子视差分干涉图在每个分解尺度下的高频小波系数进行相关性分析，得到对应分解尺度下两个高频小波系数之间的相关值，根据相关值和任一高频小波系数计算得到对流层延迟误差相位相关的小波系数；然后对所述与对流层延迟误差相位相关的小波系数进行小波逆变换，即可得到对流层延迟误差相位；
步骤5，对全孔径干涉图进行解缠处理，将步骤4得到的对流层延迟误差相位从全孔径干涉图的解缠相位中去除；最后，根据卫星轨道参数以及去除对流层延迟误差相位的全孔径干涉图解缠相位，计算得到地表高程数据。
2.根据权利要求1所述的方法，其特征在于，对两个子视干涉图，使用两种由不同技术获得的DEM数据进行去地形处理。
3.根据权利要求1所述的方法，其特征在于，对子视差分干涉图进行小波分解表示为：
Figure FDA0002417446120000012
式中，m和n分别表示干涉图的行数和列数，x,y是像素点在雷达坐标系的坐标，ix和iy分别表示干涉图的第ix行和第iy列，Φ和Ψ分别表示小波分解的平移函数和母小波函数，J表示最优分解尺度，j代表不同的分解尺度，v,w分别表示低频小波系数和高频小波系数，ε＝1,2,3分别对应水平、垂直和对角方向；
步骤4中的相关性分析方法为，取两个子视差分干涉图在每个分解尺度下的高频小波系数
Figure FDA0002417446120000021
和
Figure FDA0002417446120000022
按以下公式进行相关值计算：
Figure FDA0002417446120000023
在小波分解得到ε＝1,2,3时的高频小波系数
Figure FDA0002417446120000024
和
Figure FDA0002417446120000025
时，代入上述相关性计算公式中，即可通过方程组计算得到两个子视差分干涉图在分解尺度j时的相关值fj和整体偏移系数cj；
再使用两个子视差分干涉图在分解尺度j时的相关性fj和整体偏移系数cj，按以下公式计算对流层延迟误差相位对应的高频小波系数：
Figure FDA0002417446120000026
式中，
Figure FDA0002417446120000027
代表对流层延迟误差相位所对应的高频小波系数；
Figure FDA0002417446120000028
表示方位向入射角为αi子视差分干涉图在分解尺度j时高频小波系数。
4.根据权利要求1所述的方法，其特征在于，去除平地相位时，平地相位的计算公式为：
Figure FDA0002417446120000029
式中，φflat为平地相位，λ为雷达信号波长，B//是平行基线长度。
5.根据权利要求1所述的方法，其特征在于，轨道误差相位的去除方法为：使用融入地表高程h的多项式模型拟合轨道误差相位，通过最小二乘法解算多项式模型中的未知参数，然后将轨道误差相位从子视差分干涉图中对应像素点的解缠相位中去除，得到新的子视差分干涉图。
6.根据权利要求1所述的方法，其特征在于，步骤1采用子孔径分解技术将全孔径SAR影像分解的过程为：
首先，对全孔径SAR影像进行一维傅利叶变换，得到SAR影像在方位向上的频谱；
然后，对SAR影像在方位向上的频谱按不同的方位向入射角进行分割，得到与方位向入射角对应的频谱；
最后，对所有与方位向入射角对应的频谱均进行傅里叶逆变换，得到与方位向入射角对应的影像，即为与方位向入射角对应的子视影像。
7.根据权利要求1所述的方法，其特征在于，对全孔径SAR影像进行预处理包括：对2个全孔径SAR影像进行共配准，然后逐像素共轭相乘生成全孔径干涉图；
对主、辅影像每个方位向入射角对应的2个子视影像进行预处理包括：将当前方位向入射角对应的2个子视影像逐像素共轭相乘得到当前方位向入射角对应的子视干涉图。
8.一种基于单基线的地表高程校正装置，其特征在于，包括：
子孔径分解模块，用于：获取两个全孔径SAR影像，采用子孔径分解技术将每个全孔径SAR影像分解为两个与方位向入射角对应的子视影像；
数据预处理模块，用于：对两个全孔径SAR影像进行预处理得到全孔径干涉图；对同一方位向入射角的两个子视影像进行预处理，得到与方位向入射角对应的子视干涉图；
差分干涉相位建模模块，用于：对两个子视干涉图，均使用DEM数据并采用差分干涉技术进行去地形处理，然后再去除平地相位和轨道误差相位处理，得到两个如以下公式所示的子视差分干涉图：
Figure FDA0002417446120000031
式中，α1和α2分别代表两个不同的方位向入射角，φdiff为子视差分干涉图中的差分干涉相位，△φtopo为地形误差相位，φatmo为对流层延迟误差相位，φnoise为随机噪声误差相位；
对流层延迟误差相位提取模块，用于：对两个子视差分干涉图进行小波分解，对两个子视差分干涉图在每个分解尺度下的高频小波系数进行相关性分析，得到对应分解尺度下两个高频小波系数之间的相关值，根据相关值和任一高频小波系数计算得到对流层延迟误差相位相关的小波系数；然后对所述与对流层延迟误差相位相关的小波系数进行小波逆变换，即可得到对流层延迟误差相位；
地表高程数据校正模块，用于：对全孔径干涉图进行解缠处理，将得到的对流层延迟误差相位从全孔径干涉图的解缠相位中去除；最后，根据卫星轨道参数以及去除对流层延迟误差相位的全孔径干涉图解缠相位，计算得到地表高程数据。
9.一种设备，其特征在于，包括处理器和存储器；其中：所述存储器用于存储计算机指令；所述处理器用于执行所述存储器存储的计算机指令，具体执行如权利要求1-7任一所述的方法。
10.一种计算机存储介质，其特征在于，用于存储程序，所述程序被执行时，用于实现如权利要求1-7任一所述的方法。
CN202010195430.9A 2020-03-19 2020-03-19 基于单基线的地表高程校正方法、装置、设备及存储介质 Active CN111239736B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202010195430.9A CN111239736B (zh) 2020-03-19 2020-03-19 基于单基线的地表高程校正方法、装置、设备及存储介质

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202010195430.9A CN111239736B (zh) 2020-03-19 2020-03-19 基于单基线的地表高程校正方法、装置、设备及存储介质

Publications (2)

Publication Number Publication Date
CN111239736A CN111239736A (zh) 2020-06-05
CN111239736B true CN111239736B (zh) 2022-02-11

Family

ID=70869555

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202010195430.9A Active CN111239736B (zh) 2020-03-19 2020-03-19 基于单基线的地表高程校正方法、装置、设备及存储介质

Country Status (1)

Country Link
CN (1) CN111239736B (zh)

Families Citing this family (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111833446B (zh) * 2020-06-24 2024-06-25 浙江省测绘科学技术研究院 一种基于特征线提取的架空地物快速纠正方法
CN111812600B (zh) * 2020-06-29 2023-09-08 中南林业科技大学 一种自适应地形相关的srtm-dem校正方法
CN112816983B (zh) * 2021-01-06 2023-09-19 中南大学 基于优化干涉图集的时序InSAR湍流大气延迟校正方法
CN116736306B (zh) * 2023-08-15 2023-10-24 成都理工大学 一种基于高分三号的时序雷达干涉监测方法
CN117665809B (zh) * 2023-12-21 2024-06-21 西南林业大学 反演森林冠层高度方法
CN119828138B (zh) * 2024-11-22 2025-10-17 武汉工程大学 一种基于平地先验信息的多基线InSAR残余相位误差校正方法和系统
CN119805416B (zh) * 2024-12-13 2025-10-28 北京空间机电研究所 一种星载单光子激光雷达海面测距误差修正方法

Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101339245A (zh) * 2008-08-08 2009-01-07 西安电子科技大学 多基线干涉合成孔径雷达干涉相位展开方法
CN101609151A (zh) * 2009-07-17 2009-12-23 重庆大学 一种基于单通道合成孔径雷达(sar)图像序列特征值分解的动目标检测方法
CN104459696A (zh) * 2014-12-24 2015-03-25 中南大学 一种基于平地相位的sar干涉基线精确估计方法
EP3229038A1 (en) * 2014-12-01 2017-10-11 Institute of Electronics, Chinese Academy of Sciences Wavelet domain insar interferometric phase filtering method in combination with local frequency estimation
CN109375222A (zh) * 2018-12-17 2019-02-22 中国国土资源航空物探遥感中心 一种合成孔径雷达干涉测量电离层相位估计和补偿方法
CN110058237A (zh) * 2019-05-22 2019-07-26 中南大学 面向高分辨率SAR影像的InSAR点云融合及三维形变监测方法
CN110487241A (zh) * 2019-08-15 2019-11-22 中国测绘科学研究院 卫星激光测高提取建筑区高程控制点方法
CN110658521A (zh) * 2019-10-16 2020-01-07 中国地质大学（北京） 一种基于缠绕相位的GBInSAR大气校正方法及系统

Family Cites Families (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
AU3797297A (en) * 1996-07-11 1998-02-09 Science Applications International Corporation Terrain elevation measurement by interferometric synthetic aperture radar (ifsar)
AU5721900A (en) * 1999-03-08 2000-09-28 Lockheed Martin Corporation Single-pass interferometric synthetic aperture radar
ITRM20070399A1 (it) * 2007-07-19 2009-01-20 Consiglio Nazionale Ricerche Metodo di elaborazione di dati rilevati mediante radar ad apertura sintetica (synthetic aperture radar - sar) e relativo sistema di telerilevamento.
IT1394733B1 (it) * 2009-07-08 2012-07-13 Milano Politecnico Procedimento per il filtraggio di interferogrammi generati da immagini sar acquisite sulla stessa area.
CN101706577B (zh) * 2009-12-01 2012-01-18 中南大学 InSAR监测高速公路路面沉降方法
CN103323830B (zh) * 2013-05-20 2016-03-09 中国科学院电子学研究所 基于极化干涉合成孔径雷达的三元素分解方法及装置
CN103675790B (zh) * 2013-12-23 2016-01-20 中国国土资源航空物探遥感中心 一种基于高精度DEM提高InSAR技术监测地表形变精度的方法
CN106842200A (zh) * 2017-01-11 2017-06-13 中国科学院电子学研究所 一种双基合成孔径雷达成像方法和装置
CN106772342B (zh) * 2017-01-11 2020-02-07 西南石油大学 一种适用于大梯度地表沉降监测的时序差分雷达干涉方法
CN107991676B (zh) * 2017-12-01 2019-09-13 中国人民解放军国防科技大学 星载单航过InSAR系统对流层误差校正方法
CN108362200B (zh) * 2018-02-24 2020-09-22 深圳市北斗智星勘测科技有限公司 一种快速更新InSAR变形序列结果的方法
CN108761397B (zh) * 2018-05-30 2022-05-27 中南大学 基于电磁散射模拟的极化sar模型分解评价方法
CN109254270A (zh) * 2018-11-01 2019-01-22 西南交通大学 一种星载x波段合成孔径雷达干涉定标方法

Patent Citations (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101339245A (zh) * 2008-08-08 2009-01-07 西安电子科技大学 多基线干涉合成孔径雷达干涉相位展开方法
CN101609151A (zh) * 2009-07-17 2009-12-23 重庆大学 一种基于单通道合成孔径雷达(sar)图像序列特征值分解的动目标检测方法
EP3229038A1 (en) * 2014-12-01 2017-10-11 Institute of Electronics, Chinese Academy of Sciences Wavelet domain insar interferometric phase filtering method in combination with local frequency estimation
CN104459696A (zh) * 2014-12-24 2015-03-25 中南大学 一种基于平地相位的sar干涉基线精确估计方法
CN109375222A (zh) * 2018-12-17 2019-02-22 中国国土资源航空物探遥感中心 一种合成孔径雷达干涉测量电离层相位估计和补偿方法
CN110058237A (zh) * 2019-05-22 2019-07-26 中南大学 面向高分辨率SAR影像的InSAR点云融合及三维形变监测方法
CN110487241A (zh) * 2019-08-15 2019-11-22 中国测绘科学研究院 卫星激光测高提取建筑区高程控制点方法
CN110658521A (zh) * 2019-10-16 2020-01-07 中国地质大学（北京） 一种基于缠绕相位的GBInSAR大气校正方法及系统

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
Topography correlated atmospheric delay correction in radar interferometry using wavelet transforms;M. Shirzaei;《GEOPHYSICAL RESEARCH LETTERS》;20120106;第39卷(第1期);第1305页 *
Underlying Topography Estimation Over Forest Areas Using Single-Baseline InSAR Data;H. Q. Fu, J. J. Zhu, C. C. Wang, R. Zhao and Q. H. Xie;《IEEE Transactions on Geoscience and Remote Sensing》;20181116;第57卷;第2876-2888页 *
多分辨率分析的干涉技术大气改正算法;赵蓉;胡俊;章浙涛;占文俊;付海强;《测绘科学》;20140930;第40卷(第3期);第57-62页 *

Also Published As

Publication number Publication date
CN111239736A (zh) 2020-06-05

Similar Documents

Publication Publication Date Title
CN111239736B (zh) 基于单基线的地表高程校正方法、装置、设备及存储介质
EP3896482B1 (en) Method for the computer-implemented generation of a synthetic data set for training a convolutional neural network for an interferometric sar
Lanari et al. Generation of digital elevation models by using SIR-C/X-SAR multifrequency two-pass interferometry: The Etna case study
Shirzaei A wavelet-based multitemporal DInSAR algorithm for monitoring ground surface motion
CN102144174B (zh) Sar图像序列中的永久散射体的识别和分析
CN111273293A (zh) 一种顾及地形起伏的InSAR残余运动误差估计方法及装置
Liao et al. Improved topographic mapping through high-resolution SAR interferometry with atmospheric effect removal
CN110018476B (zh) 时间差分基线集时序干涉sar处理方法
Pinheiro et al. Generation of highly accurate DEMs over flat areas by means of dual-frequency and dual-baseline airborne SAR interferometry
CN115060208A (zh) 基于多源卫星融合的输变电线路地质灾害监测方法及系统
CN115062260B (zh) 一种适用于异质森林的森林生物量PolInSAR估算方法和系统、存储介质
CN119780923B (zh) 基于干涉雷达影像的林下地形提取方法、系统、设备及介质
CN119199850B (zh) 煤矿采空区沉陷盆地边界的测量方法及装置
CN118191841A (zh) 一种基于相关性分析的地表沉降形变测量校正方法、装置、设备及介质
CN112363165B (zh) 一种林下地形反演方法、装置、设备及介质
CN118244268A (zh) 基于ICESat-2修正轨道误差的DEM生成方法、装置、设备及介质
CN110297243B (zh) 合成孔径雷达层析三维成像中的相位误差补偿方法及装置
CN112711021A (zh) 一种多分辨率InSAR交互干涉时序分析方法
CN113341410B (zh) 一种大范围林下地形估计方法、装置、设备及介质
Yamashita et al. Mitigation of ionospheric noise in azimuth offset based on the split-spectrum method
CN114966681B (zh) 一种基于大气校正C波段InSAR数据的土壤湿度估算方法
CN112649807B (zh) 一种基于小波多尺度相关性分析的机载InSAR轨道误差去除方法
CN119415833A (zh) 一种多尺度自适应时序InSAR地表形变提取方法
CN110297242A (zh) 基于压缩感知的合成孔径雷达层析三维成像方法及装置
Belhadj-Aissa et al. Contextual filtering methods based on the subbands and subspaces decomposition of complex SAR interferograms

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