CN116202454B - Surface deformation monitoring method and system based on backward sequential least square - Google Patents
Surface deformation monitoring method and system based on backward sequential least squareInfo
- Publication number
- CN116202454B CN116202454B CN202310187871.8A CN202310187871A CN116202454B CN 116202454 B CN116202454 B CN 116202454B CN 202310187871 A CN202310187871 A CN 202310187871A CN 116202454 B CN116202454 B CN 116202454B
- Authority
- CN
- China
- Prior art keywords
- deformation
- time period
- determining
- historical
- matrix
- 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
Links
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01B—MEASURING LENGTH, THICKNESS OR SIMILAR LINEAR DIMENSIONS; MEASURING ANGLES; MEASURING AREAS; MEASURING IRREGULARITIES OF SURFACES OR CONTOURS
- G01B15/00—Measuring arrangements characterised by the use of electromagnetic waves or particle radiation, e.g. by the use of microwaves, X-rays, gamma rays or electrons
- G01B15/06—Measuring arrangements characterised by the use of electromagnetic waves or particle radiation, e.g. by the use of microwaves, X-rays, gamma rays or electrons for measuring the deformation in a solid
-
- 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
Landscapes
- Engineering & Computer Science (AREA)
- Remote Sensing (AREA)
- Radar, Positioning & Navigation (AREA)
- Physics & Mathematics (AREA)
- Electromagnetism (AREA)
- General Physics & Mathematics (AREA)
- Computer Networks & Wireless Communication (AREA)
- Radar Systems Or Details Thereof (AREA)
- Testing And Monitoring For Control Systems (AREA)
Abstract
本发明公开了一种基于后向序贯最小二乘的地表形变监测方法及系统，方法包括：确定监测区域在当前时间段内每一个时间单元的第一形变速率及对应的协因数矩阵；确定监测区域在历史时间段内的SAR影像对应的多幅历史解缠干涉图；根据多个第一形变速率、协因数矩阵和多幅历史解缠干涉图，利用后向序贯最小二乘方法，确定历史时间段内的每一个时间单元的第二形变速率并更新第一形变速率；根据多个更新后的第一形变速率和多个第二形变速率，确定监测区域在当前时间段和历史时间段之间的累积形变。本发明可向后恢复历史时间的形变时间序列，方法简单，且因其计算过程中需要存储的数据较少，故可以极大地节省存储资源。
The present invention discloses a method and system for monitoring surface deformation based on backward sequential least squares, the method comprising: determining the first deformation rate and the corresponding cofactor matrix of each time unit of the monitoring area in the current time period; determining multiple historical unwrapped interference graphs corresponding to the SAR images of the monitoring area in the historical time period; according to the multiple first deformation rates, the cofactor matrix and the multiple historical unwrapped interference graphs, using the backward sequential least squares method, determining the second deformation rate of each time unit in the historical time period and updating the first deformation rate; according to the multiple updated first deformation rates and the multiple second deformation rates, determining the cumulative deformation of the monitoring area between the current time period and the historical time period. The present invention can restore the deformation time series of the historical time backward, the method is simple, and because less data needs to be stored in the calculation process, it can greatly save storage resources.
Description
Technical Field
The invention discloses a ground surface deformation monitoring method and system based on backward sequential least squares, in particular relates to a ground surface deformation monitoring method and system of a backward sequential least squares small baseline set InSAR for recovering historical time sequence deformation, and belongs to the technical field of ground surface deformation monitoring.
Background
The earth surface is influenced by natural factors and human activities, and can generate deformation with different degrees, so that the earth surface building is damaged, the life safety of people is threatened, and economic loss is caused. The method is used for rapidly and accurately monitoring the ground deformation and is important for ground deformation disaster prevention, control and treatment and health and safety assessment of the ground structure.
The synthetic aperture radar interferometry (Interferometric Synthetic Aperture Radar, inSAR) is a ground observation technique taking the phase difference of radar signals as an observation value, has the advantages of high precision, high spatial resolution, wide monitoring range and low cost, and has been widely used for monitoring the surface deformation of urban settlement, glacier movement, mining subsidence in mining areas, volcanic eruption, landslide, debris flow and the like. The differential interferometry short baseline set time sequence analysis technology (Small Baseline Subset InSAR, SBAS-InSAR) in the multi-temporal InSAR technology is an InSAR time sequence method based on multiple primary images, only utilizes the short interference pair of space-time baselines to extract ground surface deformation information, and generates a small baseline set by setting proper time and space baseline thresholds so as to inhibit the influence of incoherence, so that the method has a huge application prospect in the field of ground subsidence monitoring of large-area and long-time sequences.
Many SAR satellites are transmitted in various countries around the world at present, continuous SAR data can be provided, and the large data age of SAR is entered. With the increase of SAR data and the need of dynamic monitoring of important monitoring targets, researchers have conducted research work on InSAR deformation dynamic monitoring. However, the current InSAR deformation dynamic monitoring is to update the deformation time sequence forward. Taking dynamic deformation monitoring by using SBAS-InSAR as an example, when new SAR data is obtained, repeated calculation is needed by combining with historical SAR data, so that a historical accumulated deformation time sequence is obtained.
Disclosure of Invention
The application aims to provide a ground surface deformation monitoring method and system based on backward sequential least square, which are used for solving the technical problem that all unwrapped interferograms in the calculation process need to be stored when a deformation time sequence is calculated by utilizing SBAS-InSAR in the prior art, so that the storage resources of a computer are greatly occupied.
The first aspect of the invention provides a surface deformation monitoring method based on backward sequential least squares, which comprises the following steps:
Determining a first deformation rate of each time unit of the monitoring area in the current time period and a corresponding co-factor matrix;
determining a plurality of historical unwrapping interferograms corresponding to SAR images of a monitoring area in a historical time period;
Determining a second deformation rate of each time unit in the historical time period by using a backward sequential least square method according to a plurality of the first deformation rates, a co-factor matrix and the plurality of historical unwrapping interferograms, and updating the first deformation rates;
And determining the accumulated deformation of the monitoring area between the current time period and the historical time period according to the updated first deformation rate and the updated second deformation rates.
Preferably, determining the first deformation rate and the corresponding co-factor matrix of each time unit in the current time period of the monitoring area specifically includes:
determining a plurality of unwrapped interferograms corresponding to SAR images of a monitoring area in a current time period;
And determining a first deformation rate of each time unit and a corresponding co-factor matrix of the monitoring area in the current time period by using an SBAS-InSAR technology according to the plurality of unwrapping interferograms.
Preferably, the first deformation rate of the monitoring area in each time unit in the current time period is determined according to a first formula:
V1＝A1X1-L1,P1
Wherein V 1 is a first residual vector, A 1 is a design matrix of M 1×N1, the matrix is composed of a plurality of time units in a current time period, M 1 is the number of unwrapped interferograms, N 1 is the number of SAR images in the current time period, L 1 is a vector composed of differential interference phases of the unwrapped interferograms, P 1 is a weight matrix corresponding to L 1, and X 1 is a vector composed of a plurality of first deformation rates ,X1＝(A1 TP1A1)-1A1 TP1L1.
Preferably, the co-factor matrix of the first deformation rate of each time unit of the monitoring area in the current time period is determined according to a second formula:
In the formula, For the co-factor matrix of X 1,Is the transposed matrix of A 1.
Preferably, the second rate of change of each time cell within the historical time period is determined according to a third formula:
Wherein V 2 is a second residual vector, Y is a vector formed by a plurality of second deformation rates, X 1 is a vector formed by a plurality of first deformation rates, X 2 is a vector formed by X 1 and is subjected to adjustment again, that is, an updated first deformation rate, B is a coefficient matrix of Y, a 2 is a coefficient matrix of X 2, L 2 is a vector formed by differential interference phases of the historical unwrapped interferogram, and P 2 is a weight matrix corresponding to L 2.
Preferably, the adjustment criterion corresponding to the third formula is as a fourth formula, where the fourth formula is:
Wherein V 2 is a second residual vector, X 2 is a vector of X 1 after adjustment again, namely an updated first deformation rate, P 2 is a weight matrix corresponding to L 2, Is a co-factor matrix of X 1.
Preferably, the fourth formula corresponds to a normal equation matrix form such as a fifth formula:
In the formula, X 1 is a vector of a plurality of said first deformation rates,As the co-factor matrix of X 1, X 2 is the vector of X 1 after the re-adjustment, a 2 is the coefficient of X 2, L 2 is the vector of the differential interference phases of the historical unwrapped interferogram, P 2 is the weight matrix corresponding to L 2, Y is the vector of a plurality of second deformation rates, B is the coefficient matrix of Y, the solution of the fifth formula is as the sixth formula and the seventh formula, and the sixth formula is:
The seventh formula is:
Where N -1 is the inverse of N, Is a co-factor matrix of [ Y X 2]T ].
Preferably, the N -1 is solved by using a blocking matrix of an eighth formula:
In the formula,
Preferably, determining the cumulative deformation of the monitoring area between the current time period and the historical time period according to the updated first deformation rate and the updated second deformation rates specifically includes:
integrating the updated first deformation rates and the updated second deformation rates in a time domain to obtain the accumulated deformation of the monitoring area between the current time period and the historical time period.
A second aspect of the present invention provides a backward sequential least squares-based surface deformation monitoring system method, comprising:
the first rate determining module is used for determining a first deformation rate of each time unit and a corresponding co-factor matrix of the monitoring area in the current time period;
The interferogram determining module is used for determining a plurality of historical unwrapping interferograms corresponding to SAR images of the monitoring area in a historical time period;
The second rate determining and first rate updating module is used for determining a second deformation rate of each time unit in the historical time period and updating the first deformation rate by using a backward sequential least square method according to a plurality of the first deformation rates, a co-factor matrix and the plurality of historical unwrapping interferograms;
and the deformation determining module is used for determining the accumulated deformation of the monitoring area between the current time period and the historical time period according to the updated first deformation rates and the updated second deformation rates.
Compared with the prior art, the earth surface deformation monitoring method and system based on backward sequential least square have the following beneficial effects:
The method and the system can recover the deformation time sequence of the history time backwards by only combining the first deformation rate of each time unit in the current time period and the co-factor matrix thereof, and the method is simple, and the data required to be stored in the calculation process only has the first deformation rate, the co-factor matrix and the history unwrapping interferogram, so that compared with the SBAS-InSAR technology, all the unwrapping interferograms are required to be stored, and storage resources can be greatly saved.
In addition, when solving the sequential least square equation coefficient matrix, the method utilizes the block matrix inversion to reduce the order of the calculated matrix and improve the calculation efficiency.
Drawings
FIG. 1 is a schematic flow chart of a method for monitoring surface deformation based on backward sequential least squares in an embodiment of the invention;
FIG. 2 is a detailed flow chart of a method for monitoring surface deformation based on backward sequential least squares in an embodiment of the present invention;
FIG. 3 is a schematic diagram of a system for monitoring surface deformation based on backward sequential least squares in an embodiment of the present invention;
FIG. 4 shows the results of the calculation of three deformation models, namely a linear model, a trigonometric function model and a quadratic function model, without adding noise in the embodiment of the invention, wherein (a), (d) and (g) are respectively simulated deformation time sequences, (b), (e) and (h) are respectively the results of the calculation of SBAS-InSAR, and (c), (f) and (i) are respectively the results of the calculation of the backward sequential least square small baseline set InSAR;
FIG. 5 shows the results of the calculation of three deformation models, namely a linear model, a trigonometric function model and a quadratic function model, under the condition of adding noise in the embodiment of the invention, wherein (a), (d) and (g) are respectively simulated deformation time sequences, (b), (e) and (h) are respectively the results of the calculation of SBAS-InSAR, and (c), (f) and (i) are respectively the results of the calculation of the backward sequential least square small baseline set InSAR;
FIG. 6 is a time series of deformations obtained by SBAS-InSAR calculation in an embodiment of the present invention;
FIG. 7 is a deformation time sequence obtained by first solving the backward sequential least squares small baseline set InSAR in the embodiment of the present invention;
FIG. 8 is a deformation time sequence obtained by performing a second solution to the backward sequential least squares small baseline set InSAR in an embodiment of the present invention;
FIG. 9 is a histogram of the difference between the backward sequential least squares small baseline set InSAR results and the SBAS-InSAR results in an embodiment of the present invention;
FIG. 10 is a time series of the deformation at point P1 in FIG. 8, wherein (a) is the result of the calculation of SBAS-InSAR, (b) is the result of the calculation of the backward sequential least squares small baseline set InSAR;
FIG. 11 shows a time series of the deformation at the point P2 in FIG. 8, wherein (a) is the result of the calculation of SBAS-InSAR and (b) is the result of the calculation of the backward sequential least squares small baseline set InSAR.
In the figure, 101 is a first rate determining module, 102 is an interferogram determining module, 103 is a second rate determining and first rate updating module, and 104 is a deformation determining module.
Detailed Description
In the following description, for purposes of explanation and not limitation, specific details are set forth such as the particular system architecture, techniques, etc., in order to provide a thorough understanding of the embodiments of the present invention. It will be apparent, however, to one skilled in the art that the present invention may be practiced in other embodiments that depart from these specific details. In other instances, detailed descriptions of well-known systems, devices, circuits, and methods are omitted so as not to obscure the description of the present invention with unnecessary detail.
The existing InSAR deformation dynamic monitoring is to update the deformation time sequence forward. After a deformation time sequence of a certain period of time is acquired, in order to describe deformation characteristics of the deformation body more accurately, the deformation time sequence needs to be updated backwards, that is, there is a need to recover the historical deformation time sequence before the period of time. Therefore, the invention aims to recover the deformation time sequence of the historical period on the basis of obtaining a certain deformation time sequence, thereby realizing long time sequence monitoring of the key monitoring target.
The first aspect of the present invention provides a method for monitoring surface deformation based on backward sequential least squares, as shown in fig. 1 and fig. 2, including:
Step 1, determining a first deformation rate and a corresponding co-factor matrix of each time unit in a current time period of a monitoring area, wherein the method specifically comprises the following steps:
And 1.1, determining a plurality of unwrapped interferograms corresponding to SAR images of the monitoring area in the current time period.
In the embodiment of the invention, the method is set in the current time periodIn the method, an N 1 +1 scene SLC image (single vision complex SLC image refers to SAR image formed by only using a section of synthetic aperture length) is acquired, then the SLC image is subjected to master-slave image registration, a small baseline set is generated according to a space and time baseline threshold value, then a differential interference pattern is generated by combining external digital elevation model (DigitalElevationModel, DEM) data, then the differential interference pattern is subjected to phase filtering, so that the influence of noise is restrained, and then the interference pattern is subjected to phase unwrapping to obtain M 1 unwrapped interference patterns, wherein:
To improve the accuracy of the unwrapped interferogram, error correction of the unwrapped interferogram is required. In the embodiment of the invention, the track error, the atmospheric delay error and the DEM error are corrected on the unwrapping interferogram, so that the interference of error items is reduced, M 1 unwrapping interferograms after error correction are obtained, and the unwrapping interferograms related later in the invention are all unwrapping interferograms after error correction.
Step 1.2, determining a first deformation rate of each time unit and a corresponding co-factor matrix of the monitoring area in the current time period by utilizing an SBAS-InSAR technology according to the plurality of unwrapping interferograms.
The SBAS-InSAR technology is a commonly used multi-main image MT-InSAR (multi-temporal InSAR) technology, and a differential interference pair is formed by selecting a proper spatial baseline and a time baseline threshold value, so that the influence of incoherence caused by a long-time baseline and a spatial baseline is restrained. The embodiment utilizes a differential interferometry short baseline set timing analysis technology (Small Baseline Subset InSAR, SBAS-InSAR) to determine a first deformation rate of each time unit and a corresponding co-factor matrix of a monitoring area in a current time period, and specifically:
Taking the deformation rate v (delta t i) in the time period t i～ti+1(i＝k,…,k+N1 -1) as an unknown number and taking the differential interference phase delta phi i(i＝1,…,M1) with the DEM error and the atmospheric delay error removed as an observed quantity, the function model of the deformation time sequence is obtained for the first time by the SBAS-InSAR technology:
V1＝A1X1-L1,P1 (2)
Wherein V 1 is a first residual vector, a 1 is a design matrix of M 1×N1, the matrix is composed of a plurality of time units in a current time period, M 1 is a number of unwrapped interferograms, N 1 is a number of SAR images in the current time period, X 1 is a vector composed of a plurality of first deformation rates, L 1 is a vector composed of differential interference phases of the unwrapped interferograms, and P 1 is a weight matrix corresponding to L 1, specifically:
A 1 is a design matrix of M 1×N1, which can be determined by a time base line of the generated small base line set interference pair, X 1 is a parameter to be estimated, L 1 is an observation vector, V 1 is a first residual vector, and P 1 is a weight matrix corresponding to the observation vector L 1. For equation (2), the least squares criterion may be used to obtain the optimal solution:
X1＝(A1 TP1A1)-1A1 TP1L1 (6)
In the formula, The co-factor matrix of X 1 represents the relative accuracy of the deformation rate in different time periods in knowledge,Is the transposed matrix of A 1. And finally, integrating the deformation rates of different time periods in a time domain to obtain a time sequence of accumulated deformation displacement.
And 2, determining a plurality of historical unwrapping interferograms corresponding to SAR images of the monitoring area in a historical time period.
In the embodiment of the invention, a new unwrapping interferogram related to the historical time to be recovered is generated and is recorded as the historical unwrapping interferogram. Registering SAR images in a to-be-recovered historical time period with the main image, generating a new interference pair according to a spatial and temporal baseline threshold, generating a new differential interferogram by combining external DEM data, and then carrying out phase filtering and phase unwrapping on the newly generated interferogram to obtain a newly generated historical unwrapped interferogram.
And then performing error correction on the newly generated historical unwrapped interferogram. The method comprises the steps of carrying out track error, atmosphere delay error and DEM error correction on a newly generated historical unwrapping interferogram, so as to reduce interference of error items.
Step 3, determining a second deformation rate of each time unit in the historical time period by using a backward sequential least square method according to a plurality of first deformation rates, a co-factor matrix and a plurality of historical unwrapping interferograms, and updating the first deformation rates, specifically:
When the deformation time series of the history period needs to be recovered, backward sequential least square may be used, backward solution may be used, the history deformation time series to be recovered may be estimated, and the estimated deformation time series may be updated. Assuming that the deformation timing of N 2 time points needs to be solved backwards on the basis of time t k, namely the acquisition time t i(i＝k-N2,k-N2+1,…,k,k+1,…k+N1+1,k＞N2) relative to the time Only according to M 2 unwrapped interferograms related to a time point before time t k and deformation rates of different time periods calculated before, the deformation rates of all time periods can be recovered by a sequential adjustment method.
In the embodiment of the invention, the function model of the adjustment is as follows:
Wherein V 2 is a second residual vector, Y is a vector composed of a plurality of second deformation rates, X 2 is a vector after the adjustment of X 1, that is, an updated first deformation rate, B is a coefficient matrix of Y, a 2 is a coefficient matrix of X 2, L 2 is a vector composed of differential interference phases of a historical unwrapped interferogram, and P 2 is a weight matrix corresponding to L 2, specifically:
Y is a vector composed of deformation rates of a history period to be recovered, X 2 is a value of X 1 after rescaling, B, A 2 is a coefficient corresponding to a parameter Y, X 2, L 2 is an observation vector composed of differential interference phases related to a time point before a time t k, from which DEM errors and atmospheric delay errors are removed, P 2 is a weight matrix of the observation vector L 2, and V 2 is a second residual vector. The adjustment result X 1 obtained in step 1.2 can be regarded as the prior information of the parameter X 2 to be estimated, so that X 2 and V 2 are both subjected to normal distribution and are independent of each other, namely
Cov(V2,X2)=0 (11)
In the formula,Is the unit weight variance. For equation (9), the adjustment criterion is:
Order the Then
I.e.
Both sides are transposed at the same time to obtain
I.e.
Order the
Then formula (16) can be written as
For the inversion of the matrix N, the calculation efficiency can be improved by using the block matrix inversion, and the calculation efficiency is improved according to a block matrix inversion formula
Wherein, the Thus, the solution of formula (18) is
The unit weight variance of this adjustment is
Where M 1 is the number of interferograms at the first estimate, M 2 is the number of interferograms at the second estimate, and N 1+N2 is the number of total parameters estimated at this time. The variance of the estimated parameters is
And 4, determining the accumulated deformation of the monitoring area between the current time period and the historical time period according to the updated first deformation rates and the updated second deformation rates.
In the embodiment of the invention, the updated first deformation rates and the updated second deformation rates are integrated in a time domain to obtain the accumulated deformation of the monitoring area between the current time period and the historical time period.
A second aspect of the present invention provides a backward sequential least squares based surface deformation monitoring system, as shown in fig. 3, comprising a first rate determination module 101, an interferogram determination module 102, a second rate determination and first rate update module 103, and a deformation determination module 104.
The first rate determining module 101 is configured to determine a first deformation rate of each time unit and a corresponding co-factor matrix in the current time period in the monitoring area;
The interferogram determining module 102 is configured to determine a plurality of historical unwrapped interferograms corresponding to the SAR image of the monitored area in the historical time period;
The second rate determining and first rate updating module 103 is configured to determine a second deformation rate of each time unit in the historical time period and update the first deformation rate according to the plurality of first deformation rates, the corresponding co-factor matrix and the plurality of historical unwrapping interferograms by using a backward sequential least squares method;
The deformation determination module 104 is configured to determine an accumulated deformation of the monitored area between the current time period and the historical time period according to the updated plurality of first deformation rates and the updated plurality of second deformation rates.
The method and system of the present invention will be described in more detail below with reference to specific embodiments.
1. Simulation data experiment
In order to compare the result of the method and the resolving result of SBAS-InSAR, deformation time sequence resolving experiments of three deformation models, namely a linear model, a trigonometric function model and a quadratic function model, are respectively simulated under the conditions of no noise and noise addition. As shown in fig. 4, the interference phases of 84 interference pairs generated from 30-scene SLC data acquired under three deformation models were simulated, and the deformation time series was calculated from these interference phases. The method of the invention comprises three times of resolving, sequentially recovering all deformation time sequences, firstly resolving to obtain deformation time sequences of 21 st to 30 th time points by using an SBAS-InSAR method, and then recovering the deformation time sequences of 11 th to 20 th and 1 st to 10 th historical time by using backward sequential least squares resolving respectively. From fig. 4, it can be seen that the backward sequential least squares small baseline set InSAR solution result and the SBAS-InSAR solution result of the present invention are completely consistent without adding noise. The deformation time sequence resolving results of the three deformation models, namely the linear model, the trigonometric function model and the quadratic function model under the condition of noise addition are shown in figure 5, and after random noise obeying N (0, 1.5 2) is added to the three different deformation models, the result of the backward sequential least squares small baseline set InSAR resolving and the result of the SBAS-InSAR resolving are completely consistent, but have tiny differences with real deformation. In addition, the error in the unit weights of SBAS-InSAR and the backward sequential least squares small baseline set InSAR is also calculated. For linear model deformation, after noise is added, the error in the unit weight calculated by SBAS-InSAR is 1.59mm, and the error in the unit weight calculated by the backward sequential least squares small baseline set InSAR after the first calculation, the second calculation and the third calculation is 1.60mm, 1.48mm and 1.59mm respectively. For the deformation of the trigonometric function model, after noise is added, the error in the unit weight calculated by SBAS-InSAR is 1.79mm, and the error in the unit weight calculated by the backward sequential least squares small baseline set InSAR of the invention after the first calculation, the second calculation and the third calculation is 1.15mm, 1.81mm and 1.79mm respectively. For the deformation of the quadratic function model, after adding noise, the error in the unit weight calculated by SBAS-InSAR is 1.61mm, and the error in the unit weight calculated by the backward sequential least squares small baseline set InSAR of the invention after the first calculation, the second calculation and the third calculation is 1.77mm, 1.57mm and 1.61mm respectively. The error in the final adjustment unit weight of the backward sequential least square small baseline set InSAR is consistent with the error in the unit weight of the SBAS-InSAR.
2. True data experiment
The research area selected in the experiment is a black square platform area of Yongjing county of Hui nationality of Lin Xia in Gansu province, and the used data are 44 interferograms generated by 19-scene derate TerraSAR-X data of the black square platform area 2016 from 24 days of 1 month to 5 days of 11 months of 2016. The deformation time sequence of the black square platform area is obtained by respectively using SBAS-InSAR and the backward sequential least square small baseline set InSAR of the invention, the backward sequential least squares small baseline set InSAR first calculates the deformation time sequence of 10 time points from the year 2016, the month 6, the month 15, the year 2016, the month 11, the month 5, and then recovers the deformation time sequence of the year 2016, the month 1, the month 24, the year 2016, the month 6, the month 15 and the first 9 time points by using the backward sequential least squares. The deformation time sequence result obtained by SBAS-InSAR calculation is shown in figure 6, and the deformation time sequence result obtained by the backward sequential least square small baseline set InSAR calculation of the invention is shown in figures 7 and 8. In order to compare the difference between the deformation time sequence result of the backward sequential least square small baseline set InSAR calculation and the SBAS-InSAR calculation result, the time sequence result is also subjected to difference, and a statistical histogram of the difference is drawn as shown in figure 9. In order to further compare the difference between the results of the solution of the SBAS-InSAR and the result of the solution of the backward sequential least squares small baseline set InSAR, the deformation time sequences of the two unstable points P1 and P2 in the deformation region in FIG. 8 are extracted as shown in FIG. 10 and FIG. 11, and as can be seen from the deformation time sequences of the two unstable points P1 and P2, the result of the solution of the backward sequential least squares small baseline set InSAR is completely consistent with the result of the solution of the SBAS-InSAR, and the backward sequential least squares small baseline set InSAR can be used for recovering the historical deformation time sequences.
According to the earth surface deformation monitoring method and system based on backward sequential least square, the deformation time sequence of the history time can be recovered backward only by combining the newly generated unwrapped interferogram after error correction with the solution calculated last time and the co-factor matrix;
The invention uses the backward sequential least square small baseline set InSAR technology to recover the deformation time sequence of the history time, the data to be stored has the previous unwrapping interference pattern, the co-factor array and the newly generated unwrapping interference pattern, and compared with the SBAS-InSAR technology, the invention needs to store all the unwrapping interference patterns, thereby greatly saving storage resources.
The invention uses the backward sequential least square small baseline set InSAR technology, and when solving the equation coefficient matrix, the method utilizes the block matrix inversion to reduce the order of the calculated matrix and improve the calculation efficiency.
While the application has been described in terms of preferred embodiments, it will be understood by those skilled in the art that various changes and modifications can be made without departing from the scope of the application, and it is intended that the application is not limited to the specific embodiments disclosed.
Claims (7)
1. The earth surface deformation monitoring method based on backward sequential least square is characterized by comprising the following steps of:
Determining a first deformation rate of each time unit of the monitoring area in the current time period and a corresponding co-factor matrix;
determining a plurality of historical unwrapping interferograms corresponding to SAR images of a monitoring area in a historical time period;
Determining a second deformation rate of each time unit in the historical time period by using a backward sequential least square method according to a plurality of the first deformation rates, a co-factor matrix and the plurality of historical unwrapping interferograms, and updating the first deformation rates;
Determining an accumulated deformation of the monitored region between the current time period and the historical time period according to the updated first deformation rates and the updated second deformation rates;
the method for determining the first deformation rate of each time unit and the corresponding co-factor matrix of the monitoring area in the current time period specifically comprises the following steps:
determining a plurality of unwrapped interferograms corresponding to SAR images of a monitoring area in a current time period;
Determining a first deformation rate and a corresponding co-factor matrix of each time unit of a monitoring area in a current time period by using an SBAS-InSAR technology according to the plurality of unwrapping interferograms;
Determining a first deformation rate of each time unit of the monitoring area in the current time period according to a first formula, wherein the first formula is as follows:
V1＝A1X1-L1,P1
Wherein V 1 is a first residual vector, A 1 is a design matrix of M 1×N1, the matrix is composed of a plurality of time units in a current time period, M 1 is the number of unwrapped interferograms, N 1 is the number of SAR images in the current time period, L 1 is a vector composed of differential interference phases of the unwrapped interferograms, P 1 is a weight matrix corresponding to L 1, and X 1 is a vector composed of a plurality of first deformation rates ,X1＝(A1 TP1A1)-1A1 TP1L1;
Determining a co-factor matrix of the first deformation rate of each time unit in the current time period of the monitoring area according to a second formula, wherein the second formula is as follows:
In the formula, For the co-factor matrix of X 1,Is the transposed matrix of A 1.
2. The backward sequential least squares based surface deformation monitoring method of claim 1, wherein the second deformation rate for each time cell within the historical time period is determined according to a third formula:
Wherein V 2 is a second residual vector, Y is a vector formed by a plurality of second deformation rates, X 1 is a vector formed by a plurality of first deformation rates, X 2 is a vector formed by X 1 and is subjected to adjustment again, that is, an updated first deformation rate, B is a coefficient matrix of Y, a 2 is a coefficient matrix of X 2, L 2 is a vector formed by differential interference phases of the historical unwrapped interferogram, and P 2 is a weight matrix corresponding to L 2.
3. The backward sequential least squares-based surface deformation monitoring method according to claim 2, wherein the adjustment criterion corresponding to the third formula is as a fourth formula, and the fourth formula is:
Wherein V 2 is a second residual vector, X 2 is a vector of X 1 after adjustment again, namely an updated first deformation rate, P 2 is a weight matrix corresponding to L 2, Is a co-factor matrix of X 1.
4. The backward sequential least squares-based surface deformation monitoring method of claim 3, wherein the fourth equation corresponds to a normal equation matrix form as a fifth equation:
In the formula, X 1 is a vector formed by a plurality of first deformation rates, Q X1 is a coefficient matrix of X 1, X 2 is a vector formed by X 1 after adjustment, A 2 is a coefficient of X 2, L 2 is a vector formed by differential interference phases of the historical unwrapped interferograms, P 2 is a weight matrix corresponding to L 2, Y is a vector formed by a plurality of second deformation rates, B is a coefficient matrix of Y, solutions of the fifth formula are as a sixth formula and a seventh formula, and the sixth formula is:
The seventh formula is:
Where N -1 is the inverse of N, Is a co-factor matrix of [ Y X 2]T ].
5. The backward sequential least squares-based surface deformation monitoring method of claim 4, wherein the N -1 is solved using a blocking matrix of an eighth formula:
In the formula,
6. The backward sequential least squares-based surface deformation monitoring method of claim 1, wherein determining the cumulative deformation of the monitored region between the current time period and the historical time period based on a plurality of updated first deformation rates and a plurality of second deformation rates, comprises:
integrating the updated first deformation rates and the updated second deformation rates in a time domain to obtain the accumulated deformation of the monitoring area between the current time period and the historical time period.
7. A backward sequential least squares-based surface deformation monitoring system employing the backward sequential least squares-based surface deformation monitoring method of any of claims 1-6, comprising:
the first rate determining module is used for determining a first deformation rate of each time unit and a corresponding co-factor matrix of the monitoring area in the current time period;
The interferogram determining module is used for determining a plurality of historical unwrapping interferograms corresponding to SAR images of the monitoring area in a historical time period;
The second rate determining and first rate updating module is used for determining a second deformation rate of each time unit in the historical time period and updating the first deformation rate by using a backward sequential least square method according to a plurality of the first deformation rates, a co-factor matrix and the plurality of historical unwrapping interferograms;
and the deformation determining module is used for determining the accumulated deformation of the monitoring area between the current time period and the historical time period according to the updated first deformation rates and the updated second deformation rates.
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310187871.8A CN116202454B (en) | 2023-03-02 | 2023-03-02 | Surface deformation monitoring method and system based on backward sequential least square |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| CN202310187871.8A CN116202454B (en) | 2023-03-02 | 2023-03-02 | Surface deformation monitoring method and system based on backward sequential least square |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| CN116202454A CN116202454A (en) | 2023-06-02 |
| CN116202454B true CN116202454B (en) | 2025-07-25 |
Family
ID=86510878
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| CN202310187871.8A Active CN116202454B (en) | 2023-03-02 | 2023-03-02 | Surface deformation monitoring method and system based on backward sequential least square |
Country Status (1)
| Country | Link |
|---|---|
| CN (1) | CN116202454B (en) |
Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106950556A (en) * | 2017-05-03 | 2017-07-14 | 三亚中科遥感研究所 | Heritage area deformation monitoring method based on distributed diffusion body sequential interference SAR technology |
| CN109061641A (en) * | 2018-07-06 | 2018-12-21 | 中南大学 | A kind of InSAR timing earth's surface deformation monitoring method based on sequential adjustment |
Family Cites Families (7)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN111474544B (en) * | 2020-03-04 | 2022-11-18 | 广东明源勘测设计有限公司 | Landslide deformation monitoring and early warning method based on SAR data |
| CN112797886B (en) * | 2021-01-27 | 2022-04-22 | 中南大学 | Winding phase oriented InSAR time sequence three-dimensional deformation monitoring method |
| CN113091596B (en) * | 2021-03-31 | 2022-01-25 | 中国矿业大学 | Surface deformation monitoring method based on multi-polarization time sequence SAR data |
| CN115308743B (en) * | 2022-08-09 | 2025-07-15 | 江苏师范大学 | A sequential BWS-DIE homogeneous point selection method for deformation monitoring in mining areas |
| CN115267779B (en) * | 2022-08-26 | 2025-05-13 | 长安大学 | A method for error estimation of InSAR time series DEM |
| CN115453534B (en) * | 2022-09-19 | 2023-05-16 | 中山大学 | Sequential InSAR time sequence deformation resolving method considering unwrapping error |
| CN115451874B (en) * | 2022-11-11 | 2023-02-28 | 深圳市科荣软件股份有限公司 | Dam deformation monitoring method and device based on InSAR image |
-
2023
- 2023-03-02 CN CN202310187871.8A patent/CN116202454B/en active Active
Patent Citations (2)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| CN106950556A (en) * | 2017-05-03 | 2017-07-14 | 三亚中科遥感研究所 | Heritage area deformation monitoring method based on distributed diffusion body sequential interference SAR technology |
| CN109061641A (en) * | 2018-07-06 | 2018-12-21 | 中南大学 | A kind of InSAR timing earth's surface deformation monitoring method based on sequential adjustment |
Also Published As
| Publication number | Publication date |
|---|---|
| CN116202454A (en) | 2023-06-02 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| Schlögl et al. | Comprehensive time-series analysis of bridge deformation using differential satellite radar interferometry based on Sentinel-1 | |
| CN110412574B (en) | A distributed target InSAR time sequence processing method and device with enhanced spatiotemporal coherence | |
| CN109738892B (en) | Mining area earth surface high-space-time resolution three-dimensional deformation estimation method | |
| CN109061641B (en) | An InSAR time series surface deformation monitoring method based on sequential adjustment | |
| CN111273293B (en) | InSAR residual motion error estimation method and device considering terrain fluctuation | |
| CN111650579B (en) | InSAR mining area three-dimensional deformation estimation method and device for rock migration parameter adaptive acquisition and medium | |
| CN113281742B (en) | SAR landslide early warning method based on landslide deformation information and meteorological data | |
| CN112797886B (en) | Winding phase oriented InSAR time sequence three-dimensional deformation monitoring method | |
| CN113687353A (en) | DS target phase optimization method based on homogeneous pixel time series phase matrix decomposition | |
| CN106066478A (en) | Merge pixel skew to follow the tracks of and the mining area surface deformation calculation method of short baseline set | |
| CN116052009B (en) | GNSS station measurement layout method based on InSAR deformation and application thereof | |
| CN112268517B (en) | Method for monitoring deformation of power transmission tower equipment by using PSInSAR | |
| CN115079172A (en) | A kind of MTInSAR landslide monitoring method, equipment and storage medium | |
| CN108362200A (en) | A kind of method of quick update InSAR Deformation Series results | |
| CN113091600A (en) | Monitoring method for monitoring deformation of soft soil foundation by utilizing time sequence InSAR technology | |
| CN107144213A (en) | The big magnitude three-D sequential deformation method of estimation in mining area and device based on SAR intensity images | |
| CN116202454B (en) | Surface deformation monitoring method and system based on backward sequential least square | |
| Su et al. | An M-estimation method for InSAR nonlinear deformation modeling and inversion | |
| CN119902205B (en) | Surface deformation monitoring method based on intermittent stacking time sequence InSAR | |
| Yan et al. | Sequential SBAS-InSAR backward estimation of deformation time series | |
| CN110146885A (en) | Time-series extraction method of ultra-long deformation in Xincheng land area based on multi-satellite platform MT-InSAR fusion | |
| CN119810740A (en) | MT-InSAR deformation monitoring method combining temporal and spatial characteristics of reservoir bank landslide deformation | |
| CN115792904B (en) | Surface deformation monitoring method, terminal device and computer readable storage medium | |
| CN116502525A (en) | Method for estimating high-time-frequency time sequence deformation of bank slide slope | |
| Li et al. | Fusion of InSAR and GNSS Based on Adaptive Spatio-Temporal Kalman Model for Reconstructing High Spatio-Temporal Resolution Deformation |
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 |