WO2023142205A1 - 一种InSAR时序相位的优化方法及装置 - Google Patents

一种InSAR时序相位的优化方法及装置 Download PDF

Info

Publication number
WO2023142205A1
WO2023142205A1 PCT/CN2022/077183 CN2022077183W WO2023142205A1 WO 2023142205 A1 WO2023142205 A1 WO 2023142205A1 CN 2022077183 W CN2022077183 W CN 2022077183W WO 2023142205 A1 WO2023142205 A1 WO 2023142205A1
Authority
WO
WIPO (PCT)
Prior art keywords
data set
sar
coordinate system
interference data
image coordinate
Prior art date
Application number
PCT/CN2022/077183
Other languages
English (en)
French (fr)
Inventor
蒋弥
程晓
梁琦
Original Assignee
中山大学
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 中山大学 filed Critical 中山大学
Priority to US17/962,545 priority Critical patent/US20230236311A1/en
Publication of WO2023142205A1 publication Critical patent/WO2023142205A1/zh

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

Definitions

  • the invention relates to the technical field of synthetic aperture radar interferometry, in particular to an optimization method and device for an InSAR time series phase.
  • Time-series phase optimization is one of the most critical steps in InSAR data processing. It can resist radar signal de-correlation and enhance the interference signal-to-noise ratio, thereby achieving high-quality process control and reducing unwrapping error propagation.
  • the mainstream methods of timing phase optimization are divided into two categories. One is to only consider spatial information, and use Goldstein and other noise reduction methods to filter the interferogram one by one. The representative one is the small baseline technology, but this method does not take into account the image details and loses space. Resolution is a method of increasing the deviation in exchange for improving the accuracy. Therefore, the technology represented by SBAS is often efficient in large-scale time-series InSAR deformation monitoring, but it cannot achieve fine monitoring and has low reliability.
  • the representative one is the SqueeSAR technology. This method can take into account the image details, but it needs to search for the single-view complex image sequence pixel by pixel. Homogeneous objectives in order to satisfy statistical model assumptions. Therefore, the calculation efficiency of the technology represented by SqueeSAR is very low, and it is difficult to apply to large-scale deformation monitoring applications.
  • the present invention provides a method and device for optimizing InSAR timing phase to solve the technical problem in the prior art that there is no high-efficiency and high-precision timing phase optimization method.
  • an embodiment of the present invention provides a method for optimizing InSAR timing phase, including:
  • the model and the preset land cover map are registered and geocoded to obtain the digital elevation under the SAR image coordinate system and the land cover map under the SAR image coordinate system respectively; according to the digital elevation under the SAR image coordinate system, the The interference data set is differentially operated to obtain a differential interference data set; according to the land cover map in the SAR image coordinate system and the differential interference data set, the covariance matrix of the differential interference data set is estimated, and according to the preset
  • the timing phase maximum likelihood estimation formula estimates and obtains the optimized timing phase.
  • the acquisition of the time-series SAR data set is performed, and the time-series SAR data set is subjected to L visual processing to obtain the L apparent intensity data set and the interference data set, specifically:
  • n is a positive integer greater than 1.
  • registering the time-series SAR data set and performing L-view processing can improve the SAR intensity sequence and interference intensity obtained after L-view processing.
  • the accuracy of the data set also avoids the need for subsequent processing of time-series SAR single-view images.
  • the multi-view data set obtained after L-view processing can greatly reduce the volume of the data set to be processed.
  • the L apparent intensity data set is used as a reference to obtain and register and geocode the preset digital elevation model and the preset land cover map respectively, and obtain the digital elevation and the digital elevation in the SAR image coordinate system respectively.
  • the land cover map in the SAR image coordinate system specifically:
  • a differential operation is performed on the interference data set to obtain a differential interference data set, specifically:
  • the n(n-1)/2 interference data sets are differentially operated to obtain n(n-1)/2 differential interference data sets; wherein, the The differential interference data set is in one-to-one correspondence with the pixel positions of the land cover map in the SAR image coordinate system.
  • the differential operation of the interferometric data set according to the digital elevation in the SAR image coordinate system can eliminate the contribution of the terrain phase and effectively improve the accuracy of the interferometric phase.
  • the pixel positions of the land cover map below are in one-to-one correspondence, so that in the subsequent step of estimating the covariance matrix of the differential interference data set, the pixels with the same category label as the central pixel in the window can be accurately counted, and at the same time Compared with the prior art, the overall computing efficiency is also improved.
  • the covariance matrix of the differential interference data set is estimated, specifically:
  • Sliding windows are established for the n(n-1)/2 differential interference data sets in turn, and according to the land cover map under the SAR image coordinate system, the sliding windows corresponding to each differential interference data set are marked in turn with preset Referring to the pixels with the same attribute as the pixels, the pixels with the same class label as the central pixel in the sliding window are counted sequentially, and the elements of the covariance matrix are calculated, thereby estimating the covariance matrix of the differential interference data set.
  • the land cover map under the SAR image coordinate system is used Mark the pixels with the same attributes as the reference pixels in the current sliding window, and directly use the pixels marked with the same category, that is, homogeneous pixels, to estimate the covariance matrix, which greatly improves the operation efficiency and has obvious advantages in efficient operations.
  • the optimized timing phase is estimated and obtained according to the preset timing phase maximum likelihood estimation formula, specifically:
  • ⁇ is the positive value of the timing phase, represents the real part operator, and ⁇ is the covariance matrix of the differential interference data set.
  • the optimized timing phase can be accurately estimated through the preset timing phase maximum likelihood estimation formula and the covariance matrix of the differential interference data set, and the accuracy of SAR timing phase optimization can be further improved.
  • the present invention also provides a device for optimizing InSAR timing phase, including: an L-view processing module, a registration encoding module, a differential operation module, and a timing optimization module;
  • the L-view processing module is used to obtain a time-series SAR data set, and perform L-view processing on the time-series SAR data set to obtain a L-view intensity data set and an interference data set respectively;
  • the registration and coding module is used to use the L apparent intensity data set as a reference to obtain and perform registration and geocoding on the preset digital elevation model and the preset land cover map, respectively, to obtain the SAR image coordinate system Land cover map in digital elevation and SAR image coordinate system;
  • the differential operation module is configured to perform a differential operation on the interference data set according to the digital elevation in the SAR image coordinate system to obtain a differential interference data set;
  • the timing optimization module is used for estimating the covariance matrix of each pixel according to the land cover map in the SAR image coordinate system and the differential interference data set, and according to the preset timing phase maximum likelihood estimation formula , to estimate and obtain the optimized timing phase.
  • the L-view processing module is used to obtain a time-series SAR data set, and perform L-view processing on the time-series SAR data set to obtain an L-view intensity data set and an interference data set, specifically:
  • n time-series SAR data sets and register the n time-series SAR data sets, perform L-view processing on the registered n time-series SAR data sets, and obtain n SAR intensity sequences and n(n- 1)/2 interference data sets; wherein, n is a positive integer greater than 1.
  • the registration and encoding module is used to obtain and perform registration and geocoding on the preset digital elevation model and the preset land cover map respectively with reference to the L visual intensity data set to obtain SAR images respectively
  • the digital elevation in the coordinate system and the land cover map in the SAR image coordinate system specifically:
  • the differential operation module is configured to perform a differential operation on the interference data set according to the digital elevation in the SAR image coordinate system to obtain a differential interference data set, specifically:
  • the n(n-1)/2 interference data sets are differentially operated to obtain n(n-1)/2 differential interference data sets; wherein, the The differential interference data set is in one-to-one correspondence with the pixel positions of the land cover map in the SAR image coordinate system.
  • the timing optimization module is used to estimate the covariance matrix of each pixel according to the land cover map in the SAR image coordinate system and the differential interference data set, specifically:
  • Sliding windows are established for the n(n-1)/2 differential interference data sets in turn, and according to the land cover map under the SAR image coordinate system, the sliding windows corresponding to each differential interference data set are marked in turn with preset Refer to the pixels with the same attributes as the pixels, so as to count the pixels in the sliding window with the same class mark as the central pixel in turn, and calculate the elements of the covariance matrix, so as to estimate the covariance matrix of each pixel.
  • the timing optimization module is configured to estimate and obtain the optimized timing phase according to a preset timing phase maximum likelihood estimation formula, specifically:
  • ⁇ is the positive value of the timing phase, represents the real part operator, and ⁇ is the covariance matrix of the differential interference data set.
  • Fig. 1 a flow chart of steps for an optimization method of an InSAR timing phase provided by an embodiment of the present invention
  • Figure 2 a schematic flow diagram of an InSAR timing phase optimization method provided by an embodiment of the present invention
  • Fig. 3 A schematic flow diagram of covariance matrix estimation in an InSAR timing phase optimization method provided by an embodiment of the present invention
  • Fig. 4 A comparison diagram of an InSAR timing phase optimization method provided by the embodiment of the present invention with the original differential interferogram and the optimization result of the existing method;
  • Fig. 5 is a schematic diagram of an InSAR timing phase optimization device provided by an embodiment of the present invention.
  • an InSAR timing phase optimization method provided for an embodiment of the present invention includes the following steps S101-S104:
  • S101 Obtain a time-series SAR data set, and perform L-view processing on the time-series SAR data set to obtain a L-view intensity data set and an interference data set respectively.
  • a time-series SAR dataset with n time-series SAR single-view images is obtained, and the n time-series SAR single-view images are registered, and L-view processing is performed on the registered n time-series SAR single-view images, Obtain n SAR intensity sequences and n(n-1)/2 interference data sets respectively; where n is a positive integer greater than 1.
  • the time-series SAR data set of n time-series SAR single-view images acquired is the SAR data of Sentinel 1.
  • the selection of the number of views L is determined according to actual needs, preferably, in this embodiment, the number of views L is determined by the resolution of the preset land cover map in step S102
  • the satellite is an Earth observation satellite in the European Space Agency's Copernicus program (GMES), consisting of two satellites, carrying a C-band synthetic aperture radar, which can provide continuous images, including but not limited to day, night and various Weather image data information.
  • GMES European Space Agency's Copernicus program
  • registering the time-series SAR data set and performing L-view processing can improve the SAR intensity sequence and interference intensity obtained after L-view processing.
  • the accuracy of the data set also avoids the need for subsequent processing of time-series SAR single-view images.
  • the multi-view data set obtained after L-view processing can greatly reduce the volume of the data set to be processed.
  • S102 Using the L apparent intensity data set as a reference, acquire and perform registration and geocoding on the preset digital elevation model and the preset land cover map respectively, and respectively obtain the digital elevation in the SAR image coordinate system and the SAR image coordinate system The land cover map below.
  • a preset cross-correlation maximization algorithm is used to obtain and register and geocode the preset digital elevation model and the L apparent intensity data set to obtain a SAR image Digital elevation in the coordinate system; Acquire and register and geocode the preset satellite orbit data information with the preset land cover map to obtain the land cover map in the SAR image coordinate system.
  • the preset cross-correlation maximization algorithm has the characteristics of insensitivity to noise and accurate matching positions, while the digital elevation model realizes the digital simulation of ground terrain through limited terrain elevation data, that is, the topography surface shape Digital expression, after registering and geocoding the preset digital elevation model and the L apparent intensity data set, the digital elevation under the SAR image coordinate system can be obtained; preferably, in this embodiment, the preset The satellite orbit data information is provided by the Sentinel-1 satellite, and the preset land cover map reflects the condition of the earth surface. By registering and geocoding the preset satellite orbit data information and the preset land cover map, the SAR can be obtained Land cover map in image coordinate system.
  • S103 Perform a differential operation on the interference data set according to the digital elevation in the SAR image coordinate system to obtain a differential interference data set.
  • the n(n-1)/2 interference data sets are differentially operated to obtain n(n-1)/2 differential interference data sets; wherein , the differential interference data set corresponds to the pixel position of the land cover map in the SAR image coordinate system in a one-to-one correspondence.
  • the differential operation of the interferometric data set according to the digital elevation in the SAR image coordinate system can eliminate the contribution of the terrain phase and effectively improve the accuracy of the interferometric phase.
  • the pixel positions of the land cover map below are in one-to-one correspondence, so that in the subsequent step of estimating the covariance matrix of the differential interference data set, the pixels with the same category label as the central pixel in the window can be accurately counted, and at the same time Compared with the prior art, the overall computing efficiency is also improved.
  • S104 Estimate the covariance matrix of the differential interference data set according to the land cover map in the SAR image coordinate system and the differential interference data set, and estimate and Obtain the optimized timing phase.
  • a sliding window is sequentially established for the n(n-1)/2 differential interference data sets, and the corresponding sliding window of each differential interference data set is marked sequentially according to the land cover map in the SAR image coordinate system. Pixels with the same attributes as the preset reference pixels, so as to sequentially count the pixels with the same category mark as the central pixel in the sliding window, and calculate the elements of the covariance matrix, thereby estimating the covariance matrix of the differential interference data set.
  • ⁇ e j ⁇ ; so as to estimate and obtain the optimized timing phase; where, is the optimized timing phase, ⁇ is the positive value of the timing phase, represents the real part operator, and ⁇ is the covariance matrix of the differential interference data set.
  • the definition of the L-parallel differential interferogram is s i and s k represent the i-th and k-th single-view complex images respectively, * represents the conjugate operator, Represents the terrain phase compensation item, which can be obtained from Sentinel-1 satellite set parameters and a preset digital elevation model in this embodiment.
  • the L-view differential interferogram is simplified as where ⁇ i,k represent the original phase.
  • any element of the covariance matrix of the single-view complex image set can be obtained through the L-view differential interferogram, and the use of M samples for the single-view complex image set is equivalent to the use of Q samples for the L-view data, so a large The memory space of the input data is reduced, and the covariance matrix ⁇ 0 of the original data set can be restored from the covariance matrix ⁇ of the L-view data set through the transformation of the above formula.
  • Figure 3 is a schematic diagram of the covariance matrix estimation of the L visual difference interferogram assisted by the land cover map in the SAR image coordinate system, by obtaining the i and kth L visual average windows in the original data, and by The land cover map in the SAR image coordinate system is used to select homogeneous pixels for covariance estimation.
  • the land cover map under the SAR image coordinate system is used Mark the pixels with the same attributes as the reference pixels in the current sliding window, and directly use the pixels marked with the same category, that is, homogeneous pixels, to estimate the covariance matrix, which greatly improves the operation efficiency and has obvious advantages in efficient operations.
  • the maximum likelihood estimation is a parameter estimation method, through several experiments, observing the results, and using the results to estimate the approximate value of the parameters.
  • the optimized timing phase is estimated by using the positive value ⁇ of the timing phase and the covariance matrix of the differential interference data set
  • the optimized timing phase can be accurately estimated through the preset timing phase maximum likelihood estimation formula and the covariance matrix of the differential interference data set, and the accuracy of SAR timing phase optimization can be further improved.
  • FIG. 2 is a schematic flow chart of an InSAR timing phase optimization method in an embodiment of the present invention.
  • L-view SAR intensity data set and L-view interference data set are obtained after L-view processing, and the SAR coordinate system is obtained by registering and geocoding the digital elevation model and land cover map respectively
  • the digital elevation and land cover map under the SAR coordinate system; the digital elevation under the SAR coordinate system and the L-view interference data set are differentially operated to obtain the L-view differential interference data set, and with the assistance of the land cover map under the SAR coordinate system, Based on the selection of similar targets and the estimation of the covariance matrix for L, and according to the maximum likelihood timing phase estimation, the timing optimization phase is obtained.
  • Figure 4 is a comparison diagram of another embodiment of the present invention, the original difference interferogram and the optimization result of the existing method.
  • Both the single-view complex image difference in the existing method and the embodiment of the present invention can effectively suppress noise.
  • the KS test uses the KS test to select homogeneous samples on the data set under single-view, and then use the same likelihood estimator to estimate the time series phase.
  • the KS test selects a window of 11*21 size.
  • the result under the test of MATLAB2021a version As shown in Table 1, the input and parameter output of multi-view data (that is, the data processed by L-view) occupy about 17% of the workspace space (169.98M) under single-view (978.73M), which proves that multi-view data can significantly
  • the memory space is reduced, and the sample selection of the land cover map under the SAR image coordinates only takes one minute, and the efficiency is 272 times that of the traditional method, which proves the high-efficiency computing capability of this embodiment and the potential of massive data processing.
  • the embodiment of the present invention performs L-view processing on the time-series SAR data set after acquiring the time-series SAR data set.
  • the existing single-view complex data set is converted into an L-view data set.
  • the occupied workspace space of the input and parameter output is much lower than that of single-view, which can greatly reduce the volume of the data set, and estimate the covariance matrix of the differential interference data set through the obtained land cover map in the SAR image coordinate system, greatly The calculation efficiency is improved, and estimation through the homogeneous sample statistical selection method in the prior art is avoided.
  • the present invention also provides an InSAR timing phase optimization device, including: an L-view processing module 201 , a registration encoding module 202 , a differential operation module 203 and a timing optimization module 204 .
  • the L-view processing module 201 is configured to obtain a time-series SAR data set, and perform L-view processing on the time-series SAR data set to obtain a L-view intensity data set and an interference data set respectively.
  • the L-view processing module 201 is configured to acquire n time-series SAR data sets, and perform registration on the n time-series SAR data sets, and perform L-view processing on the registered n time-series SAR data sets , to obtain n SAR intensity sequences and n(n-1)/2 interference data sets respectively; where, n is a positive integer greater than 1.
  • the registration and encoding module 202 is used to obtain and perform registration and geocoding on the preset digital elevation model and the preset land cover map, respectively, using the L apparent intensity data set as a reference, to obtain the SAR image coordinate system respectively.
  • the registration encoding module 202 is configured to use the L visual intensity data set as a reference and adopt a preset cross-correlation maximization algorithm to obtain and compare the preset digital elevation model with the L visual intensity data set. Perform registration and geocoding to obtain the digital elevation in the SAR image coordinate system; obtain and perform registration and geocoding on the preset satellite orbit data information and the preset land cover map to obtain the land cover map in the SAR image coordinate system .
  • the differential operation module 203 is configured to perform a differential operation on the interference data set according to the digital elevation in the SAR image coordinate system to obtain a differential interference data set.
  • the difference operation module 203 is configured to perform a difference operation on the n(n-1)/2 interference data sets according to the digital elevation in the SAR image coordinate system to obtain n(n-1) /2 differential interference data sets; wherein, the differential interference data sets are in one-to-one correspondence with the pixel positions of the land cover map in the SAR image coordinate system.
  • the timing optimization module 204 is used to estimate the covariance matrix of each pixel according to the land cover map in the SAR image coordinate system and the differential interference data set, and estimate the maximum likelihood of each pixel according to the preset timing phase Formula, estimate and obtain the optimized timing phase.
  • the timing optimization module 204 is configured to sequentially establish a sliding window for the n(n-1)/2 differential interference data sets, and mark each frame in turn according to the land cover map in the SAR image coordinate system Pixels with the same attributes as the preset reference pixels in the sliding window corresponding to the differential interference data set, so as to sequentially count the pixels with the same category mark as the central pixel in the sliding window, and calculate the elements of the covariance matrix, thereby estimating each pixel The covariance matrix of .
  • the embodiment of the present invention performs L-view processing on the time-series SAR data set after acquiring the time-series SAR data set.
  • the existing single-view complex data set is converted into an L-view data set.
  • the occupied workspace space of input and parameter output is much lower than that of single-view data, so the data set volume can be greatly reduced, and at the same time, the covariance matrix of the differential interference data set can be estimated through the obtained land cover map in the SAR image coordinate system , which greatly improves the computing efficiency, avoids estimation through the statistical selection method of homogeneous samples in the prior art, and has the ability of high-efficiency computing and the potential of massive data processing.

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)
  • Radar Systems Or Details Thereof (AREA)

Abstract

一种InSAR时序相位的优化方法及装置，包括：获取时序SAR数据集，并对时序SAR数据集进行L视处理，分别获得L视强度数据集和干涉数据集(S101)；以L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图(S102)；根据SAR图像坐标系下的数字高程，对干涉数据集进行差分操作，得到差分干涉数据集(S103)；根据SAR图像坐标系下的土地覆盖图和差分干涉数据集，估算差分干涉数据集的协方差矩阵，并根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位(S104)。InSAR时序相位的优化方法及装置能够应用于大范围的高效时序相位优化估计的过程。

Description

一种InSAR时序相位的优化方法及装置 技术领域
本发明涉及合成孔径雷达干涉测量技术领域，尤其涉及一种InSAR时序相位的优化方法及装置。
背景技术
随着当代SAR卫星重返频率逐步提高，海量观测数据不仅为更加客观地了解形变的发展动态和规律提供了契机，同时也给快速、高效的时序InSAR带来了新的挑战。尤其是在欧空局哥白尼计划中卫星重返周期12天、幅宽高达200km*200km的TOPS哨兵1数据，传统数据处理方法往往难以满足监测时效性需求。
大范围时序相位优化是InSAR数据处理中最为关键的步骤之一，能够抵抗雷达信号失相关，增强干涉信噪比，从而实现高质量的流程控制，减小解缠误差传播。目前时序相位优化的主流方法分为两类，一是仅考虑空间信息，使用Goldstein等降噪方法对干涉图逐个滤波，具有代表性的是小基线技术，但是该方法不顾及影像细节，损失空间分辨率，是一类以偏差增大换取精度提升的方法。因此，以SBAS为代表的技术在大范围时序InSAR形变监测时往往高效，但是不能实现精细监测，且可靠性低；二是兼顾时空信息，在同类地物具有相似后向散射属性的假设基础上，使用KS(Kolmogorov-Smirnov)检验和BFGS(Broyden Fletcher Goldfarb Shanno)相位三角算法优化时序相位，具有代表性的是SqueeSAR技术，该方法能够顾及影像细节，但是需要对单视复数影像序列逐像素寻找同类目标，以便满足统计模型假设前提。因此，以SqueeSAR为代表的技术计算效率很低，难以适用于大范围形变监测应用。
因此，目前仍需一种如何实现高效、高精度的InSAR时序相位的优化方法， 以解决现有技术中缺乏同时具备高效、高精度的InSAR时序相位优化方法。
发明内容
本发明提供了一种InSAR时序相位的优化方法及装置，以解决现有技术中缺乏同时具备高效、高精度的时序相位优化方法的技术问题。
为了解决上述技术问题，本发明实施例提供了一种InSAR时序相位的优化方法，包括：
获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集合和干涉数据集；以所述L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图；根据所述SAR图像坐标系下的数字高程，对所述干涉数据集进行差分操作，得到差分干涉数据集；根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算所述差分干涉数据集的协方差矩阵，并根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位。
可以理解的是，在获取时序SAR数据集后，对时序SAR数据集进行L视处理，与现有技术相比，将现有的单视复数数据集转换成L视的数据集，L视数据输入和参数输出的占用工作区空间远低于单视，能够大幅降低数据集体积，并通过得到的SAR图像坐标系下的土地覆盖图，来估算所述差分干涉数据集的协方差矩阵，大幅提高了运算效率，避免了现有技术中通过同质样本统计选择方法来进行估计。
作为优选方案，所述获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集合和干涉数据集，具体为：
获取具有n幅时序SAR单视影像的时序SAR数据集，并对所述n幅时序SAR单视影像进行配准，对配准后的n幅时序SAR单视影像进行L视处理，分别获得n幅SAR强度序列和n(n-1)/2个干涉数据集；其中，n为大于1的正整数。
可以理解的是，获取具有n幅时序SAR单视影像的时序SAR数据集后，对时序SAR数据集进行配准，并进行L视处理，能够提高L视处理后所得到的SAR强度序列和干涉数据集的准确性，同时也避免后续需要对时序SAR单视影像进行处理，经过L视处理后所得到多视的数据集，能够大幅降低所要处理的数据集体积。
作为优选方案，所述以所述L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图，具体为：
以所述L视强度数据集为参考，采用预设的互相关最大化算法，获取并对预设数字高程模型与所述L视强度数据集进行配准和地理编码，得到SAR图像坐标系下的数字高程；获取并对预设的卫星轨道数据信息与预设土地覆盖图进行配准和地理编码，得到SAR图像坐标系下的土地覆盖图。
可以理解的是，在以所述L视强度数据集合为参考，对预设数字高程模型与所述L视强度数据集中的SAR强度图像进行配准和地理编码，能够提高刻画地形精确程度，以及对预设的卫星轨道数据信息与预设土地覆盖图进行配准和地理编码，使得能够为后续估算所述差分干涉数据集的协方差矩阵的步骤中提供参考像素，避免了现有技术中采用同质样本选择的方法，从而提高了整体的计算效率。
作为优选方案，所述根据所述SAR图像坐标系下的数字高程，对所述干涉数据集进行差分操作，得到差分干涉数据集，具体为：
根据所述SAR图像坐标系下的数字高程，对所述n(n-1)/2个干涉数据集进行差分操作，得到n(n-1)/2幅差分干涉数据集；其中，所述差分干涉数据集与所述SAR图像坐标系下的土地覆盖图的像素位置一一对应。
可以理解的是，根据SAR图像坐标系下的数字高程来对干涉数据集进行差分操作，能够消除地形相位贡献，可以有效提高干涉相位的精度，同时所得到的差分干涉数据集与SAR图像坐标系下的土地覆盖图的像素位置一一对应，以 使后续对估算所述差分干涉数据集的协方差矩阵的步骤中，能够准确对窗口内与中心像素具有相同类别标签的像素进行统计，同时相比于现有技术，也提高了整体的计算效率。
作为优选方案，根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算所述差分干涉数据集的协方差矩阵，具体为：
依次对所述n(n-1)/2幅差分干涉数据集建立滑动窗口，根据所述SAR图像坐标系下的土地覆盖图，依次标记每幅差分干涉数据集对应的滑动窗口内与预设参考像素相同属性的像素，从而依次统计所述滑动窗口内与中心像素具有相同类别标记的像素，计算协方差矩阵的元素，从而估算出所述差分干涉数据集的协方差矩阵。
可以理解的是，在估算所述差分干涉数据集的协方差矩阵的步骤中，代替传统分布式目标数据处理步骤中同质样本统计选择的方法，而采用通过SAR图像坐标系下的土地覆盖图标记当前滑动窗口内领域像素与参考像素相同属性的像素，并直接用相同类别标记的像素，即同质像素，来进行协方差矩阵估计，大幅提高了运算效率，在高效运算的方面具有明显的优势。
作为优选方案，所述根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位，具体为：
根据预设的时序相位极大似然估计公式：
Figure PCTCN2022077183-appb-000001
ξ＝e jθ；从而估计并获得优化后的时序相位；
其中，
Figure PCTCN2022077183-appb-000002
为优化后的时序相位，θ为时序相位的正值，
Figure PCTCN2022077183-appb-000003
表示取实部算子，∑为所述差分干涉数据集的协方差矩阵。
可以理解的是，通过预设的时序相位极大似然估计公式，以及差分干涉数据集的协方差矩阵，能够准确估计出优化后的时序相位，进一步地提高对SAR时序相位优化的准确性。
相应地，本发明还提供一种InSAR时序相位的优化装置，包括：L视处理模块、配准编码模块、差分操作模块和时序优化模块；
所述L视处理模块，用于获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集合和干涉数据集；
所述配准编码模块，用于以所述L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图；
所述差分操作模块，用于根据所述SAR图像坐标系下的数字高程，对所述干涉数据集进行差分操作，得到差分干涉数据集；
所述时序优化模块，用于根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算每一个像素的协方差矩阵，并根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位。
作为优选方案，所述L视处理模块，用于获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集合和干涉数据集，具体为：
获取n幅时序SAR数据集，并对所述n幅时序SAR数据集进行配准，对配准后的n幅时序SAR数据集进行L视处理，分别获得n幅SAR强度序列和n(n-1)/2个干涉数据集；其中，n为大于1的正整数。
作为优选方案，所述配准编码模块，用于以所述L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图，具体为：
以所述L视强度数据集为参考，采用预设的互相关最大化算法，获取并对预设数字高程模型与所述L视强度数据集进行配准和地理编码，得到SAR图像坐标系下的数字高程；获取并对预设的卫星轨道数据信息与预设土地覆盖图进行配准和地理编码，得到SAR图像坐标系下的土地覆盖图。
作为优选方案，所述差分操作模块，用于根据所述SAR图像坐标系下的数字高程，对所述干涉数据集进行差分操作，得到差分干涉数据集，具体为：
根据所述SAR图像坐标系下的数字高程，对所述n(n-1)/2个干涉数据集进 行差分操作，得到n(n-1)/2幅差分干涉数据集；其中，所述差分干涉数据集与所述SAR图像坐标系下的土地覆盖图的像素位置一一对应。
作为优选方案，所述时序优化模块，用于根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算每一个像素的协方差矩阵，具体为：
依次对所述n(n-1)/2幅差分干涉数据集建立滑动窗口，根据所述SAR图像坐标系下的土地覆盖图，依次标记每幅差分干涉数据集对应的滑动窗口内与预设参考像素相同属性的像素，从而依次统计所述滑动窗口内与中心像素具有相同类别标记的像素，计算协方差矩阵的元素，从而估算出每一个像素的协方差矩阵。
作为优选方案，所述时序优化模块，用于根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位，具体为：
根据预设的时序相位极大似然估计公式：
Figure PCTCN2022077183-appb-000004
ξ＝e jθ；从而估计并获得优化后的时序相位；
其中，
Figure PCTCN2022077183-appb-000005
为优化后的时序相位，θ为时序相位的正值，
Figure PCTCN2022077183-appb-000006
表示取实部算子，∑为所述差分干涉数据集的协方差矩阵。
附图说明
图1：为本发明实施例所提供的一种InSAR时序相位的优化方法的步骤流程图；
图2：为本发明实施例所提供的一种InSAR时序相位的优化方法的流程示意图；
图3：为本发明实施例所提供的一种InSAR时序相位的优化方法的中协方差矩阵估计的流程示意图；
图4：为本发明实施例所提供的一种InSAR时序相位的优化方法与原始差分干涉图和现有方法优化结果的效果对比图；
图5：为本发明实施例所提供的一种InSAR时序相位的优化装置示意图。
具体实施方式
下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。
实施例一
请参照图1和图2，为本发明实施例提供的一种InSAR时序相位的优化方法，包括以下步骤S101-S104：
S101：获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集合和干涉数据集。
具体地，获取具有n幅时序SAR单视影像的时序SAR数据集，并对所述n幅时序SAR单视影像进行配准，对配准后的n幅时序SAR单视影像进行L视处理，分别获得n幅SAR强度序列和n(n-1)/2个干涉数据集；其中，n为大于1的正整数。
需要说明的是，优选地，获取的n幅时序SAR单视影像的时序SAR数据集为哨兵1号的SAR数据，进行配准后，对其进行L视处理，分别获得n幅SAR强度序列和n(n-1)/2个干涉数据集；视数L的选择根据实际的需求来确定，优选地，在本实施例中，视数L由步骤S102中的预设土地覆盖图的分辨率来确定，哨兵1号的SAR数据空间分辨率约为20m*5m(方位向*距离向)，则对于20米的预设土地覆盖图，L＝1*5视；其中，哨兵1号(Sentinel-1)卫星是欧洲航天局哥白尼计划(GMES)中的地球观测卫星，由两颗卫星组成，载有C波段合成孔径雷达，可提供连续图像，包括但不限于白天、夜晚和各种天气的图像数据信息。
可以理解的是，获取具有n幅时序SAR单视影像的时序SAR数据集后，对时序SAR数据集进行配准，并进行L视处理，能够提高L视处理后所得到的SAR 强度序列和干涉数据集的准确性，同时也避免后续需要对时序SAR单视影像进行处理，经过L视处理后所得到多视的数据集，能够大幅降低所要处理的数据集体积。
S102：以所述L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图。
具体地，以所述L视强度数据集为参考，采用预设的互相关最大化算法，获取并对预设数字高程模型与所述L视强度数据集进行配准和地理编码，得到SAR图像坐标系下的数字高程；获取并对预设的卫星轨道数据信息与预设土地覆盖图进行配准和地理编码，得到SAR图像坐标系下的土地覆盖图。
需要说明的是，预设的互相关最大化算法具有对噪声不敏感、匹配位置准确等特点，而数字高程模型则是通过有限的地形高程数据实现对地面地形的数字化模拟，即地形表面形态的数字化表达，在对预设数字高程模型与所述L视强度数据集进行配准和地理编码后，即可得到SAR图像坐标系下的数字高程；优选地，在本实施例中，预设的卫星轨道数据信息则通过哨兵一号卫星提供，预设土地覆盖图则是反映地表的状况，通过对预设的卫星轨道数据信息与预设土地覆盖图进行配准和地理编码，即可得到SAR图像坐标系下的土地覆盖图。
可以理解的是，在以所述L视强度数据集合为参考，对预设数字高程模型与所述L视强度数据集中的SAR强度图像进行配准和地理编码，能够提高刻画地形精确程度，以及对预设的卫星轨道数据信息与预设土地覆盖图进行配准和地理编码，使得能够为后续估算所述差分干涉数据集的协方差矩阵的步骤中提供参考像素，避免了现有技术中采用同质样本选择的方法，从而提高了整体的计算效率。
S103：根据所述SAR图像坐标系下的数字高程，对所述干涉数据集进行差分操作，得到差分干涉数据集。
具体地，根据所述SAR图像坐标系下的数字高程，对所述n(n-1)/2个干涉 数据集进行差分操作，得到n(n-1)/2幅差分干涉数据集；其中，所述差分干涉数据集与所述SAR图像坐标系下的土地覆盖图的像素位置一一对应。
需要说明的是，对n(n-1)/2个干涉数据集进行差分操作，以消除地形相位的贡献，即对干涉数据集进行地形去除以及相位的解缠，有助于提升InSAR的形变监测能力，从而得到n(n-1)/2幅差分干涉数据集。
可以理解的是，根据SAR图像坐标系下的数字高程来对干涉数据集进行差分操作，能够消除地形相位贡献，可以有效提高干涉相位的精度，同时所得到的差分干涉数据集与SAR图像坐标系下的土地覆盖图的像素位置一一对应，以使后续对估算所述差分干涉数据集的协方差矩阵的步骤中，能够准确对窗口内与中心像素具有相同类别标签的像素进行统计，同时相比于现有技术，也提高了整体的计算效率。
S104：根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算所述差分干涉数据集的协方差矩阵，并根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位。
具体地，依次对所述n(n-1)/2幅差分干涉数据集建立滑动窗口，根据所述SAR图像坐标系下的土地覆盖图，依次标记每幅差分干涉数据集对应的滑动窗口内与预设参考像素相同属性的像素，从而依次统计所述滑动窗口内与中心像素具有相同类别标记的像素，计算协方差矩阵的元素，从而估算出所述差分干涉数据集的协方差矩阵。
具体地，根据预设的时序相位极大似然估计公式：
Figure PCTCN2022077183-appb-000007
Figure PCTCN2022077183-appb-000008
ξ＝e jθ；从而估计并获得优化后的时序相位；其中，
Figure PCTCN2022077183-appb-000009
为优化后的时序相位，θ为时序相位的正值，
Figure PCTCN2022077183-appb-000010
表示取实部算子，∑为所述差分干涉数据集的协方差矩阵。
需要说明的是，L视差分干涉图的定义为
Figure PCTCN2022077183-appb-000011
s i和s k分别表示第i和k幅单视复数影像，*表示共轭算子，
Figure PCTCN2022077183-appb-000012
表示地形相位补偿项，在本实施例中，可有哨兵1号卫星集合参数和预设数字高程模型得到。 在本实施例中，若不考虑SAR强度的情况下，L视差分干涉图简化为
Figure PCTCN2022077183-appb-000013
Figure PCTCN2022077183-appb-000014
其中φ i,k表示原始相位。在本实施例中，定义单视复数影像集的协方差矩阵为：
Figure PCTCN2022077183-appb-000015
其中
Figure PCTCN2022077183-appb-000016
将η i,k变换为
Figure PCTCN2022077183-appb-000017
Figure PCTCN2022077183-appb-000018
其中M＝L·Q。从上式可知，单视复数影像集的协方差矩阵的任一元素可通过L视差分干涉图获得，对单视复数影像集使用M个样本等价于L视数据使用Q个样本，因此大幅降低了输入数据的内存空间，同时通过上述的公式变换，可以从L视数据集的协方差矩阵∑恢复原始数据集的协方差矩阵∑ 0。
请参阅图3，其为在SAR图像坐标系下的土地覆盖图辅助下的L视差分干涉图的协方差矩阵估计的示意图，通过获取原始数据中第i，k个L视平均窗口，并通过SAR图像坐标系下的土地覆盖图来选取同质像素进行协方差的估计。
可以理解的是，在估算所述差分干涉数据集的协方差矩阵的步骤中，代替传统分布式目标数据处理步骤中同质样本统计选择的方法，而采用通过SAR图像坐标系下的土地覆盖图标记当前滑动窗口内领域像素与参考像素相同属性的像素，并直接用相同类别标记的像素，即同质像素，来进行协方差矩阵估计，大幅提高了运算效率，在高效运算的方面具有明显的优势。
需要说明的是，极大似然估计是一种参数估计方法，通过若干次试验，观察结果，利用结果估算推出参数的大概值。在本实施例中，利用时序相位的正值θ，以及差分干涉数据集的协方差矩阵，来估算优化后的时序相位
Figure PCTCN2022077183-appb-000019
可以理解的是，通过预设的时序相位极大似然估计公式，以及差分干涉数据集的协方差矩阵，能够准确估计出优化后的时序相位，进一步地提高对SAR时序相位优化的准确性。
请参阅图2，其为本发明实施例中的一种InSAR时序相位的优化方法的流 程示意图。通过获取的时序SAR数据集，进行L视处理后得到L视SAR强度数据集和L视干涉数据集，并通过分别对数字高程模型和土地覆盖图的配准和地理编码，从而得到SAR坐标系下的数字高程和土地覆盖图；将SAR坐标系下的数字高程与L视干涉数据集进行差分操作，从而得到L视差分干涉数据集，并在SAR坐标系下的土地覆盖图的辅助下，对L视同类目标选取以及协方差矩阵的估计，并根据极大似然时序相位估计，得到时序优化相位。
请参阅图4，其为本发明另一实施例与原始差分干涉图和现有方法优化结果的效果对比图，现有方法中的单视复数图像差分和本发明实施例均能有效抑制噪声，获得高信噪比相位，但是由于本发明实施例采用了L＝5视数据，因此图像尺寸较传统方法获得的结果小五倍。在该实施例中，选取2017年Sent inel-1TOPS模式21景单视复数数据集，SRTM30米分辨率DEM以及Esri2020年全球10米分辨率土地覆盖图，将所有数据集配准至20米格网，采用三中的流程估计时间序列相位。同时，采用传统方法，在单视下对数据集使用KS检验选取同质样本，后采用相同的似然估计量估计时序相位，在协方差矩阵估计过程中，KS检验选取了11*21的窗口尺寸。本发明方法使用了斜距向L＝5视对数据预处理后，选取了11*11窗口尺寸估计协方差矩阵,等效到单视下11*55的估计窗口；在MATLAB2021a版本测试下的结果如表1所示，多视数据(即L视处理后的数据)输入和参数输出占用工作区空间(169.98M)约为单视下(978.73M)的17％，证明用多视数据能够大幅降低内存空间，使用SAR图像坐标下的土地覆盖图的样本选取仅需要一分钟，效率是传统方法的272倍，证明该实施例高效运算的能力和海量数据处理方面的潜力。
表1不同方法MATLAB工作区内存消耗与统计选择耗时比较
Figure PCTCN2022077183-appb-000020
Figure PCTCN2022077183-appb-000021
实施本发明实施例，具有如下效果：
本发明实施例通过在获取时序SAR数据集后，对时序SAR数据集进行L视处理，与现有技术相比，将现有的单视复数数据集转换成L视的数据集，L视数据输入和参数输出的占用工作区空间远低于单视，能够大幅降低数据集体积，并通过得到的SAR图像坐标系下的土地覆盖图，来估算所述差分干涉数据集的协方差矩阵，大幅提高了运算效率，避免了现有技术中通过同质样本统计选择方法来进行估计。
实施例二
相应地，本发明还提供一种InSAR时序相位的优化装置，包括：L视处理模块201、配准编码模块202、差分操作模块203和时序优化模块204。
所述L视处理模块201，用于获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集合和干涉数据集。
具体地，所述L视处理模块201，用于获取n幅时序SAR数据集，并对所述n幅时序SAR数据集进行配准，对配准后的n幅时序SAR数据集进行L视处理，分别获得n幅SAR强度序列和n(n-1)/2个干涉数据集；其中，n为大于1的正整数。
所述配准编码模块202，用于以所述L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图。
具体地，所述配准编码模块202，用于以所述L视强度数据集为参考，采用预设的互相关最大化算法，获取并对预设数字高程模型与所述L视强度数据集进行配准和地理编码，得到SAR图像坐标系下的数字高程；获取并对预设的卫星轨道数据信息与预设土地覆盖图进行配准和地理编码，得到SAR图像坐标系下的土地覆盖图。
所述差分操作模块203，用于根据所述SAR图像坐标系下的数字高程，对 所述干涉数据集进行差分操作，得到差分干涉数据集。
具体地，所述差分操作模块203，用于根据所述SAR图像坐标系下的数字高程，对所述n(n-1)/2个干涉数据集进行差分操作，得到n(n-1)/2幅差分干涉数据集；其中，所述差分干涉数据集与所述SAR图像坐标系下的土地覆盖图的像素位置一一对应。
所述时序优化模块204，用于根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算每一个像素的协方差矩阵，并根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位。
具体地，所述时序优化模块204，用于依次对所述n(n-1)/2幅差分干涉数据集建立滑动窗口，根据所述SAR图像坐标系下的土地覆盖图，依次标记每幅差分干涉数据集对应的滑动窗口内与预设参考像素相同属性的像素，从而依次统计所述滑动窗口内与中心像素具有相同类别标记的像素，计算协方差矩阵的元素，从而估算出每一个像素的协方差矩阵。
具体地，所述时序优化模块204，还用于根据预设的时序相位极大似然估计公式：
Figure PCTCN2022077183-appb-000022
ξ＝e jθ；从而估计并获得优化后的时序相位；其中，
Figure PCTCN2022077183-appb-000023
为优化后的时序相位，θ为时序相位的正值，
Figure PCTCN2022077183-appb-000024
表示取实部算子，∑为所述差分干涉数据集的协方差矩阵。
所属领域得技术人员可以清楚地了解到，为描述得方便和简洁，上述描述的具体工作过程，可以参考前述方法实施例中的对应过程，在此不再赘述。
实施以上实施例，具有如下效果：
本发明实施例通过在获取时序SAR数据集后，对时序SAR数据集进行L视处理，与现有技术相比，将现有的单视复数数据集转换成L视的数据集，L视数据输入和参数输出的占用工作区空间远低于单视数据，因此能够大幅降低数据集体积，同时通过得到的SAR图像坐标系下的土地覆盖图，来估算所述差分干涉数据集的协方差矩阵，大幅提高了运算效率，避免了现有技术中通过同质样本统计选择方法来进行估计，具备高效运算的能力和海量数据处理方面的潜 力。
以上所述的具体实施例，对本发明的目的、技术方案和有益效果进行了进一步的详细说明，应当理解，以上所述仅为本发明的具体实施例而已，并不用于限定本发明的保护范围。特别指出，对于本领域技术人员来说，凡在本发明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (10)

  1. 一种InSAR时序相位的优化方法，其特征在于，包括：
    获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集和干涉数据集；
    以所述L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图；
    根据所述SAR图像坐标系下的数字高程，对所述干涉数据集进行差分操作，得到差分干涉数据集；
    根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算所述差分干涉数据集的协方差矩阵，并根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位。
  2. 如权利要求1所述的一种InSAR时序相位的优化方法，其特征在于，所述获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集和干涉数据集，具体为：
    获取具有n幅时序SAR单视影像的时序SAR数据集，并对所述n幅时序SAR单视影像进行配准，对配准后的n幅时序SAR单视影像进行L视处理，分别获得n幅SAR强度序列和n(n-1)/2个干涉数据集；其中，n为大于1的正整数。
  3. 如权利要求1所述的一种InSAR时序相位的优化方法，其特征在于，所述以所述L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图，具体为：
    以所述L视强度数据集为参考，采用预设的互相关最大化算法，获取并对预设数字高程模型与所述L视强度数据集进行配准和地理编码，得到SAR图像 坐标系下的数字高程；获取并对预设的卫星轨道数据信息与预设土地覆盖图进行配准和地理编码，得到SAR图像坐标系下的土地覆盖图。
  4. 如权利要求2所述的一种InSAR时序相位的优化方法，其特征在于，所述根据所述SAR图像坐标系下的数字高程，对所述干涉数据集进行差分操作，得到差分干涉数据集，具体为：
    根据所述SAR图像坐标系下的数字高程，对所述n(n-1)/2个干涉数据集进行差分操作，得到n(n-1)/2幅差分干涉数据集；其中，所述差分干涉数据集与所述SAR图像坐标系下的土地覆盖图的像素位置一一对应。
  5. 如权利要求4所述的一种InSAR时序相位的优化方法，其特征在于，所述根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算所述差分干涉数据集的协方差矩阵，具体为：
    依次对所述n(n-1)/2幅差分干涉数据集建立滑动窗口，根据所述SAR图像坐标系下的土地覆盖图，依次标记每幅差分干涉数据集对应的滑动窗口内与预设参考像素相同属性的像素，从而依次统计所述滑动窗口内与中心像素具有相同类别标记的像素，计算协方差矩阵的元素，从而估算出所述差分干涉数据集的协方差矩阵。
  6. 如权利要求1所述的一种InSAR时序相位的优化方法，其特征在于，所述根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位，具体为：
    根据预设的时序相位极大似然估计公式：
    Figure PCTCN2022077183-appb-100001
    ξ＝e jθ；从而估计并获得优化后的时序相位；
    其中，
    Figure PCTCN2022077183-appb-100002
    为优化后的时序相位，θ为时序相位的正值，
    Figure PCTCN2022077183-appb-100003
    表示取实部算子，∑为所述差分干涉数据集的协方差矩阵。
  7. 一种InSAR时序相位的优化装置，其特征在于，包括：L视处理模块、配准编码模块、差分操作模块和时序优化模块；
    所述L视处理模块，用于获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集合和干涉数据集；
    所述配准编码模块，用于以所述L视强度数据集为参考，获取并分别对预设数字高程模型和预设土地覆盖图进行配准与地理编码，分别得到SAR图像坐标系下的数字高程和SAR图像坐标系下的土地覆盖图；
    所述差分操作模块，用于根据所述SAR图像坐标系下的数字高程，对所述干涉数据集进行差分操作，得到差分干涉数据集；
    所述时序优化模块，用于根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算所述差分干涉数据集的协方差矩阵，并根据预设的时序相位极大似然估计公式，估计并获得优化后的时序相位。
  8. 如权利要求7所述的一种InSAR时序相位的优化装置，其特征在于，所述L视处理模块，用于获取时序SAR数据集，并对所述时序SAR数据集进行L视处理，分别获得L视强度数据集合和干涉数据集，具体为：
    获取具有n幅时序SAR单视影像的时序SAR数据集，并对所述n幅时序SAR单视影像进行配准，对配准后的n幅时序SAR单视影像进行L视处理，分别获得n幅SAR强度序列和n(n-1)/2个干涉数据集；其中，n为大于1的正整数。
  9. 如权利要求8所述的一种InSAR时序相位的优化装置，其特征在于，所述差分操作模块，用于根据所述SAR图像坐标系下的数字高程，对所述干涉数据集进行差分操作，得到差分干涉数据集，具体为：
    根据所述SAR图像坐标系下的数字高程，对所述n(n-1)/2个干涉数据集进行差分操作，得到n(n-1)/2幅差分干涉数据集；其中，所述差分干涉数据集与 所述SAR图像坐标系下的土地覆盖图的像素位置一一对应。
  10. 如权利要求9所述的一种InSAR时序相位的优化装置，其特征在于，所述时序优化模块，用于根据所述SAR图像坐标系下的土地覆盖图和所述差分干涉数据集，估算所述差分干涉数据集的协方差矩阵，具体为：
    依次对所述n(n-1)/2幅差分干涉数据集建立滑动窗口，根据所述SAR图像坐标系下的土地覆盖图，依次标记每幅差分干涉数据集对应的滑动窗口内与预设参考像素相同属性的像素，从而依次统计所述滑动窗口内与中心像素具有相同类别标记的像素，计算协方差矩阵的元素，从而估算出所述差分干涉数据集的协方差矩阵。
PCT/CN2022/077183 2022-01-26 2022-02-22 一种InSAR时序相位的优化方法及装置 WO2023142205A1 (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
US17/962,545 US20230236311A1 (en) 2022-01-26 2022-10-10 OPTIMIZATION METHOD AND APPARATUS FOR INTERFEROMETRIC SYNTHETIC APERTURE RADAR (InSAR) TIME-SERIES PHASE

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
CN202210097356.6 2022-01-26
CN202210097356.6A CN114488150A (zh) 2022-01-26 2022-01-26 一种InSAR时序相位的优化方法及装置

Related Child Applications (1)

Application Number Title Priority Date Filing Date
US17/962,545 Continuation-In-Part US20230236311A1 (en) 2022-01-26 2022-10-10 OPTIMIZATION METHOD AND APPARATUS FOR INTERFEROMETRIC SYNTHETIC APERTURE RADAR (InSAR) TIME-SERIES PHASE

Publications (1)

Publication Number Publication Date
WO2023142205A1 true WO2023142205A1 (zh) 2023-08-03

Family

ID=81476685

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2022/077183 WO2023142205A1 (zh) 2022-01-26 2022-02-22 一种InSAR时序相位的优化方法及装置

Country Status (2)

Country Link
CN (1) CN114488150A (zh)
WO (1) WO2023142205A1 (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN115267781B (zh) * 2022-09-28 2022-12-16 中山大学 一种基于多视SAR数据集的InSAR相干性估计方法

Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB2446615A (en) * 2007-02-15 2008-08-20 Selex Sensors & Airborne Sys Interleaved bi-static interferometric synthetic aperture radar technique for determining height information for an imaged area
CN108051810A (zh) * 2017-12-01 2018-05-18 南京市测绘勘察研究院股份有限公司 一种InSAR分布式散射体相位优化方法
CN110412574A (zh) * 2019-09-05 2019-11-05 河海大学 一种时空相干性增强的分布式目标InSAR时序处理方法和装置
CN111598929A (zh) * 2020-04-26 2020-08-28 云南电网有限责任公司电力科学研究院 基于时序差分干涉合成孔径雷达数据的二维解缠方法
CN112711021A (zh) * 2020-12-08 2021-04-27 中国自然资源航空物探遥感中心 一种多分辨率InSAR交互干涉时序分析方法
CN112797886A (zh) * 2021-01-27 2021-05-14 中南大学 面向缠绕相位的InSAR时序三维形变监测方法
CN113687353A (zh) * 2021-08-09 2021-11-23 中国矿业大学 基于同质像元时序相位矩阵分解的ds目标相位优化方法

Patent Citations (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
GB2446615A (en) * 2007-02-15 2008-08-20 Selex Sensors & Airborne Sys Interleaved bi-static interferometric synthetic aperture radar technique for determining height information for an imaged area
CN108051810A (zh) * 2017-12-01 2018-05-18 南京市测绘勘察研究院股份有限公司 一种InSAR分布式散射体相位优化方法
CN110412574A (zh) * 2019-09-05 2019-11-05 河海大学 一种时空相干性增强的分布式目标InSAR时序处理方法和装置
CN111598929A (zh) * 2020-04-26 2020-08-28 云南电网有限责任公司电力科学研究院 基于时序差分干涉合成孔径雷达数据的二维解缠方法
CN112711021A (zh) * 2020-12-08 2021-04-27 中国自然资源航空物探遥感中心 一种多分辨率InSAR交互干涉时序分析方法
CN112797886A (zh) * 2021-01-27 2021-05-14 中南大学 面向缠绕相位的InSAR时序三维形变监测方法
CN113687353A (zh) * 2021-08-09 2021-11-23 中国矿业大学 基于同质像元时序相位矩阵分解的ds目标相位优化方法

Non-Patent Citations (1)

* Cited by examiner, † Cited by third party
Title
MA, ZHANGFENG; YUE, DONGJIE; JIANG, MI; LIU, LIAN: "Co-registration of Image Stacks for Sentinel-1A TOPS Mode Based on Dijkstra's Algorithm", GEOMATICS AND INFORMATION SCIENCE OF WUHAN UNIVERSITY, vol. 45, no. 6, 30 June 2020 (2020-06-30), CN , pages 904 - 913, XP009548111, ISSN: 1671-8860, DOI: 10.13203/j.whugis20180412 *

Also Published As

Publication number Publication date
CN114488150A (zh) 2022-05-13

Similar Documents

Publication Publication Date Title
CN110058237B (zh) 面向高分辨率SAR影像的InSAR点云融合及三维形变监测方法
US11238602B2 (en) Method for estimating high-quality depth maps based on depth prediction and enhancement subnetworks
WO2021227423A1 (zh) 基于动态基线的InSAR数字高程模型构建方法及系统
CN107817465B (zh) 超高斯噪声背景下的基于无网格压缩感知的doa估计方法
US20230236311A1 (en) OPTIMIZATION METHOD AND APPARATUS FOR INTERFEROMETRIC SYNTHETIC APERTURE RADAR (InSAR) TIME-SERIES PHASE
CN109212522B (zh) 一种基于双基星载sar的高精度dem反演方法及装置
CN111650579B (zh) 一种岩移参数自适应获取的InSAR矿区三维形变估计方法、装置及介质
CN113065467B (zh) 一种基于深度学习的卫星图像低相干区域识别方法及装置
CN108051810A (zh) 一种InSAR分布式散射体相位优化方法
Wang et al. A new likelihood function for consistent phase series estimation in distributed scatterer interferometry
CN111964667B (zh) 一种基于粒子滤波算法的地磁_ins组合导航方法
CN104657962A (zh) 基于级联线性回归的图像超分辨重建方法
WO2023142205A1 (zh) 一种InSAR时序相位的优化方法及装置
CN107742133A (zh) 一种用于极化sar图像的分类方法
CN109507704A (zh) 一种基于互模糊函数的双星定位频差估计方法
CN103454636A (zh) 基于多像素协方差矩阵的差分干涉相位估计方法
CN111583330B (zh) 一种多尺度时空马尔可夫遥感影像亚像元定位方法及系统
CN117788296B (zh) 基于异构组合深度网络的红外遥感图像超分辨率重建方法
CN111242999A (zh) 基于上采样及精确重匹配的视差估计优化方法
CN108876829B (zh) 基于非线性尺度空间及径向基函数的sar高精度配准方法
CN112184785B (zh) 基于mcd度量和vtm的多模态遥感图像配准方法
CN111062972B (zh) 基于图像频率域转换的图像跟踪方法
CN111948652B (zh) 一种基于深度学习的sar智能参数化超分辨成像方法
CN111832399A (zh) 一种融合注意力机制的跨域道路航标配准算法
CN115546264A (zh) 一种星载InSAR图像精细配准与立体测量方法

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 22923026

Country of ref document: EP

Kind code of ref document: A1