CN110146883B - Simultaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method based on minimum acceleration - Google Patents

Simultaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method based on minimum acceleration Download PDF

Info

Publication number
CN110146883B
CN110146883B CN201910411920.5A CN201910411920A CN110146883B CN 110146883 B CN110146883 B CN 110146883B CN 201910411920 A CN201910411920 A CN 201910411920A CN 110146883 B CN110146883 B CN 110146883B
Authority
CN
China
Prior art keywords
deformation
date
satellite
data set
time
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201910411920.5A
Other languages
Chinese (zh)
Other versions
CN110146883A (en
Inventor
赵卿
王强
丁敬钊
刘敏
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
East China Normal University
Original Assignee
East China Normal University
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by East China Normal University filed Critical East China Normal University
Priority to CN201910411920.5A priority Critical patent/CN110146883B/en
Publication of CN110146883A publication Critical patent/CN110146883A/en
Application granted granted Critical
Publication of CN110146883B publication Critical patent/CN110146883B/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01S—RADIO DIRECTION-FINDING; RADIO NAVIGATION; DETERMINING DISTANCE OR VELOCITY BY USE OF RADIO WAVES; LOCATING OR PRESENCE-DETECTING BY USE OF THE REFLECTION OR RERADIATION OF RADIO WAVES; ANALOGOUS ARRANGEMENTS USING OTHER WAVES
    • G01S13/00—Systems using the reflection or reradiation of radio waves, e.g. radar systems; Analogous systems using reflection or reradiation of waves whose nature or wavelength is irrelevant or unspecified
    • G01S13/88—Radar or analogous systems specially adapted for specific applications
    • G01S13/89—Radar or analogous systems specially adapted for specific applications for mapping or imaging
    • G01S13/90—Radar or analogous systems specially adapted for specific applications for mapping or imaging using synthetic aperture techniques, e.g. synthetic aperture radar [SAR] techniques
    • G01S13/904—SAR modes
    • G01S13/9058—Bistatic or multistatic SAR
    • 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
    • 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/904—SAR modes
    • G01S13/9064—Inverse SAR [ISAR]

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

The invention provides a minimum acceleration-based contemporaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method, which can decompose the deformation time sequence and the deformation rate of the MT-InSAR to east-west, south-north and vertical directions from a satellite visual line. The method comprises the following steps: step 1, respectively inverting a sight line deformation time sequence of an orbit ascending and descending satellite image set by using MT-InSAR; step 2, extracting a sight direction public deformation time sequence obtained by inverting the orbit ascending and descending satellite platform; step 3, spatially registering and extracting homonymous high coherence points of the image sets of the orbit ascending and descending satellite platforms; and 4, aiming at the characteristics of the mixed occurrence of the data in the same period and the data in the different periods, carrying out mixed arrangement and de-duplication treatment on the ascending and descending time vectors, establishing a deformation decomposition equation set based on minimum acceleration, and solving to obtain a three-dimensional deformation sequence. The invention develops a multi-platform MT-InSAR data three-dimensional decomposition method aiming at the new characteristic of the mixed synchronous and asynchronous data of a new synthetic aperture radar satellite constellation, and obtains a reliable three-dimensional earth surface deformation time sequence.

Description

Simultaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method based on minimum acceleration
Technical Field
The invention belongs to the technical field of remote sensing and geodetic surveying, and relates to a three-dimensional decomposition method for a time sequence of deformation of an MT-InSAR sight line of a synchronous and asynchronous orbit-rising platform to the earth surface, which can decompose the deformation time sequence and the deformation rate of the MT-InSAR sight line of a satellite sight line direction to east, west, south and north and vertical directions. The method can realize the three-dimensional decomposition of the multi-platform MT-InSAR data aiming at the new characteristics of the new synthetic aperture radar satellite constellation in the same period and different periods of data mixture, and obtain the reliable three-dimensional surface deformation time sequence.
Background
The time-series synthetic aperture radar interferometry method (MT-InSAR) can obtain reliable high-precision satellite line-of-sight (LOS) deformation time series and deformation rate, and is widely applied to the fields of landslide, earthquake, bridges, urban deformation monitoring and the like. However, due to the imaging geometry of a single satellite platform, the deformation quantity inverted by the single platform MT-InSAR can only reflect the deformation of the ground object in the satellite sight line direction, and in practice, the deformation quantity in the satellite sight line direction is often required to be converted into the deformation quantities in the east-west direction, the south-north direction and the vertical direction, and for example, the health diagnosis of urban buildings and infrastructure needs to accurately obtain the deformation results in the vertical direction and the horizontal direction.
The current multi-platform MT-InSAR deformation time sequence three-dimensional decomposition method mostly decomposes from the angle of mathematical transformation, and lacks of introducing necessary physical condition constraints, the obtained mathematical solution may not have physical significance, and the robustness degree of the mathematical solution is also required to be improved. On the other hand, along with the in-orbit operation of a new generation of radar satellites, particularly a satellite constellation in which a plurality of satellites carrying the same synthetic aperture radar operate in the same orbit, the same area can be imaged for multiple times in the same day, and when the same area is observed for a long time, the situations of orbit rising and orbit falling satellite data, synchronous or asynchronous data and mixed synchronous and asynchronous data can occur. The method is a new characteristic of acquisition of new generation of synthetic aperture radar satellite constellation data, and the existing three-dimensional deformation decomposition method is not considered.
Disclosure of Invention
The invention aims to develop a multi-platform MT-InSAR data three-dimensional decomposition method aiming at the new characteristic of the mixed occurrence of the contemporaneous and the anormal data of a new synthetic aperture radar satellite constellation, and obtain a reliable three-dimensional earth surface deformation time sequence.
The specific technical scheme for realizing the purpose of the invention is as follows:
a three-dimensional deformation decomposition method (MT-InSAR) based on minimum acceleration is a synchronous and asynchronous multi-satellite platform time sequence radar interferometry method, aiming at new characteristics of new synthetic aperture radar satellite constellation synchronous and asynchronous data mixture, a reliable three-dimensional earth surface deformation time sequence can be obtained; the method specifically comprises the following steps:
step 1: method for respectively inverting sight line direction deformation time sequence of orbit ascending satellite platform and orbit descending satellite platform by using MT-InSAR
Respectively inverting the sight line deformation time sequence and the deformation rate of different orbit ascending and descending synthetic aperture radar satellite platforms by using MT-InSAR;
step 2: extracting sight direction public deformation time sequence obtained by inversion of orbit ascending and descending satellite platform
Respectively extracting common deformation time sequences of MT-InSAR sight line direction deformation time sequences obtained by inversion of the orbit ascending satellite platform and the orbit descending satellite platform in a common time range, further extracting common space ranges of the MT-InSAR sight line direction deformation time sequences and the velocity fields of the orbit ascending satellite platform and the orbit descending satellite platform, and extracting the common range according to the common space ranges;
and step 3: space registration and extraction of homonymous high-coherence points of image sets of orbit ascending and descending satellite platforms
Carrying out spatial registration on high-coherence points extracted by the MT-InSAR of the ascending and descending satellite platforms respectively, and resampling to the same resolution, so that the high-coherence points extracted by the MT-InSAR of the ascending and descending satellite platforms respectively correspond to one another; counting and recording the number and the positions of homonymous high coherent points of the orbit ascending and descending satellite platforms;
and 4, step 4: aiming at the characteristics of the mixed occurrence of the data in the same period and the data in the different periods, the ascending and descending time vectors are arranged in a mixed mode and subjected to de-duplication processing, a deformation decomposition equation set based on minimum acceleration is established, and a three-dimensional deformation sequence is obtained through solving; the method specifically comprises the following steps:
a. after the processing of the step 1-3, the N different synthetic aperture radar satellite data sets have different visual angles and comprise orbit rising data and orbit falling data; each data set has Qi(i ═ 1, 2.., N) time series, the deformation rate for each dataset line-of-sight direction is described by:
Figure BDA0002063077210000021
each row in the above formula represents a time sequence of deformation of a synthetic aperture radar satellite platform data set at a common point, the initial time sequence period of each platform is different,
Figure BDA0002063077210000022
is the gaze direction deformation time series vector from the first data set to the gaze direction deformation time series vector of the nth data set,
Figure BDA0002063077210000023
is the accumulated deformation from the first phase to the last phase in the first data set,
Figure BDA0002063077210000024
is the accumulated amount of deformation for each phase from the first phase to the last phase in the second data set,
Figure BDA0002063077210000025
is the accumulated deformation quantity, Q, of each phase from the first phase to the last phase in the Nth data set1Is the time series number, Q, of the first data set2，…，QNIs the number of time series from the second data set to the Nth data set;
likewise, the date of acquisition for each time series of data sets is represented as:
Figure BDA0002063077210000026
in the above formula, Ti(i 1, 2.., N) is the data acquisition date vector for the first data set to the data acquisition date vector for the nth data set,
Figure BDA0002063077210000027
is the first to the greatest period in the first data setThe date of data acquisition of the latter period,
Figure BDA0002063077210000031
is the date of data acquisition from the first date to the last date in the second data set,
Figure BDA0002063077210000032
is the date of data acquisition from the first to the last phase in the Nth data set, Q1Is the time series number, Q, of the first data set2，…，QNIs the number of time series from the second data set to the Nth data set;
aiming at a synthetic aperture radar satellite constellation, different synthetic aperture radar satellites can image the same area for multiple times in the same day, namely, synchronous data; in the above formula, the acquisition dates corresponding to different single platform time sequences may be the same, i.e. in-phase, or different, i.e. out-of-phase; for the different-period multi-platform MT-InSAR deformation time sequence, the acquisition time is sequenced according to the sequence; for a synchronous multi-platform MT-InSAR deformation time sequence, date deduplication processing is carried out on the time sequence, and the finally obtained total date sequence is unique;
b. the date of the same period and the date of the different period are mixed and checked for duplication before all the acquisition dates are sorted, and the process is described by the following formula:
T＝Sort(Unique(T1，T2，...，TN))＝[t0，t1，t2，...，tQ-1] (3)
in the formula, T is a finally generated non-repeated time vector, Sort is a sorting function and is used for sorting all times in order according to the sequence, Unique is a deduplication function and is used for reserving only one date which appears repeatedly, and T1，T2，...，TNIs a data acquisition date vector, t, from the first data set to the Nth data set0，t1，t2，...，tQ-1All data set date sequences reserved after sorting and de-duplication, wherein Q is the total number of unrepeated dates in N data sets;
after a non-repeating date sequence T of N data set data acquisition dates is obtained, the date lists are subjected to staggered subtraction to obtain Q-1 difference results, and the difference results are expressed by a mathematical expression as follows:
ΔT＝[Δt1，Δt2，...，ΔtQ-1] (4)
where Δ T is the difference date interval obtained by subtracting the date sequence in T one after the other, and Δ T1，Δt2，...，ΔtQ-1Is the first to Q-1 difference date interval, delta T, obtained by the staggered subtraction of the date sequences in T1Is T in T1-t0Result of (1), Δ t2Is T in T2-t1Result of (1), Δ tQ-1Is T in TQ-1-tQ-2As a result, Q is the total number of non-repeating dates in the N data sets;
because the ground surface deformation result obtained by the MT-InSAR method is the deformation of the satellite sight line direction, the deformation of the satellite sight line direction is decomposed into the east-west direction, the south-north direction and the vertical direction according to the satellite imaging geometric relation, and the decomposition process is described by the following formula:
Figure BDA0002063077210000033
in the above equation theta is the satellite view angle,
Figure BDA0002063077210000034
is the satellite flight direction course angle, dLOSRepresenting the amount of deformation of the satellite's line of sight, dEast-WestThe amount of deformation in the east-west direction, dNorth-southIs the north-south direction deformation, dUp-DownIs the deformation amount in the vertical direction; in addition, the flight direction of the satellite is in the near south-north direction, the inversion deformation result is insensitive to the south-north direction, and if the deformation component in the south-north direction is not considered, the inversion deformation result is simplified into a two-dimensional deformation decomposition method expressed as follows: angle of course
Figure BDA0002063077210000035
The distortion of the satellite sight direction is only decomposed into east-west and vertical directions in the above formula;
c. there are Q-1 different date intervals, and the average deformation speeds in the east-west direction, the south-north direction and the vertical direction in each date interval are represented by a matrix, obtaining the following formula:
Figure BDA0002063077210000041
in the above formula VEast-West，VNorth-south，VUp-DownAnd represents the average deformation velocity vector v in the east-west direction, the south-north direction and the vertical directionE1，υE2，...，υEQ-1Is the average deformation velocity component, upsilon, in each time interval in the east-west direction delta TN1，υN2，...，υNQ-1Average deformation velocity component, upsilon, in each time interval in north-south direction delta TU1，υU2，...，υUQ-1Is the average deformation velocity component in each period of time interval in the vertical direction delta T, and Q is the total number of unrepeated dates in N data sets;
with the expression of the above equation, the average deformation rate of each date interval can be combined into a total deformation amount sequence by combining the equation (5), and the rate-deformation conversion equation is expressed as the following equation:
Figure BDA0002063077210000042
Figure BDA0002063077210000043
in the above formula, B is a rate transformation variable reconstruction comprehensive matrix, Bi(i ═ 1, 2.., N) is a rate distortion quantity reconstruction matrix, which is specifically expressed in the following equation (10), VEast-West，VNorth-south，VUp-DownThe average deformation velocity vectors in the east-west direction, the south-north direction and the vertical direction are respectively,
Figure BDA0002063077210000044
is the time-series vector of the line-of-sight deformation from the time-series vector of the line-of-sight deformation of the first data set to the time-series vector of the line-of-sight deformation of the Nth data set, theta is the satellite view,
Figure BDA00020630772100000410
is the satellite flight direction course angle;
to construct Bi(i 1, 2.. times.n), the position of the acquisition date corresponding to the time series of each data set in the sorted time interval series Δ T is set as
Figure BDA0002063077210000045
Expressed as:
Figure BDA0002063077210000046
in the above formula
Figure BDA0002063077210000047
Is the date vector T of the first data set1The position index in the time vector T of all data sets after sorting and deduplication,
Figure BDA0002063077210000048
to
Figure BDA0002063077210000049
Is the date vector T of the second data set2Position sequence numbers in time vector T of all sorted and deduplicated data sets to date vector T of Nth data setNPosition order number, Q, in time vector T of all data sets after sorting and de-duplicationi(i 1, 2.., N) represents the number of time series of 1 to N per data set;
with the help of the above formula, BiCan be constructed using the following formula:
Figure BDA0002063077210000051
in the above formula BiA matrix is reconstructed for the rate-distortion variables,
Figure BDA0002063077210000052
is the position serial number in the formula (9), and Q is the total number of unrepeated dates in the N data sets;
d. introducing a minimum acceleration assumption into the solution of the formula (7), and limiting deformation components in all directions to have minimum acceleration in time sequence, namely the speed difference between continuous time intervals is minimum, so that an equation for solving deformation decomposition has physical significance and a stable solution is obtained; the above description may be specifically expressed as:
Figure BDA0002063077210000053
where C is the minimum acceleration defining equation set, and α is a regularization factor that balances the relative weights of the minimum acceleration to limit the residual norm minima in the equation,
Figure BDA0002063077210000054
the average deformation speed difference results of adjacent date intervals in the east-west direction, the south-north direction and the vertical direction are respectively obtained, and Q is the total number of unrepeated dates in the N data sets;
and (5) integrating the equations (7) to (11) to obtain an equation (12) for solving the average deformation speed in the east-west direction, the south-north direction and the vertical direction.
Figure BDA0002063077210000055
In the above formula, B is shown in formula (8), C is shown in formula (11), and VEast-West，VNorth-south，VUp-DownThe average deformation velocity vectors in the east-west direction, the south-north direction and the vertical direction are respectively,
Figure BDA0002063077210000056
is a line-of-sight shape from a first data setChanging the time sequence vector to the sight line direction of the Nth data set;
e. and solving the linear equation set by using a crack Berg-Marquardt method to obtain the deformation time sequences of each high coherence point in the vertical direction, the south direction, the north direction and the east-west direction.
f. And applying the steps a-e to all homonymous high coherence points to obtain a three-dimensional surface deformation field and a deformation rate based on the MT-InSAR in the lifting and lowering rail platform in a public area.
The invention discloses a minimum acceleration-based contemporaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method which can decompose the deformation time sequence and the deformation rate of the MT-InSAR in the satellite sight direction to east, west, south and north and vertical directions. The method can realize the multi-platform MT-InSAR deformation time sequence three-dimensional decomposition aiming at the new characteristics of the new synthetic aperture radar satellite constellation in the same period and different periods of data mixture, and obtain the reliable three-dimensional surface deformation time sequence.
The invention has the advantages that: (1) the method can acquire a reliable three-dimensional earth surface deformation time sequence and a three-dimensional earth surface deformation field aiming at the mixed situations of synchronous, asynchronous and synchronous phases of new synthetic aperture radar satellite constellation ascending and descending orbit satellite data. (2) The algorithm is robust and efficient. The algorithm fully considers the physical constraint of minimum acceleration of deformation of adjacent time sequences, ensures positive definite of the equation set to be solved, directly solves the deformation time sequences by matrix operation, and is efficient and stable.
The method is suitable for the condition that the data acquisition directions of the multi-source radar are complementary, namely, at least one pair of lifting and lowering rail data sets is needed to ensure that the three-dimensional deformation decomposition result is stable and reliable.
Drawings
FIG. 1 is a flow chart of the present invention;
FIG. 2 is a time distribution graph of the S1A rise, fall, in-orbit, out-of-phase data acquisition;
FIG. 3 is a schematic diagram of an imaging geometry of a raised and lowered synthetic aperture radar satellite platform;
FIG. 4 is a graph of the east-west and vertical average deformation rates after decomposition of the three-dimensional deformation of an embodiment of the present invention;
FIG. 5 is a time-series diagram of the east-west and vertical deformation obtained by decomposing four feature points according to the embodiment of the present invention.
Detailed Description
Examples
The following takes a certain research area as an example to illustrate the specific implementation steps of the present invention in the examples:
step 1: method for respectively inverting sight line direction deformation time sequence of orbit ascending satellite platform and orbit descending satellite platform by using MT-InSAR
Searching the synthetic aperture radar data covering the ascending orbit and the descending orbit of a certain research area, and inverting the sight line direction deformation time sequence of different satellite platforms (including the ascending orbit and the descending orbit) by using MT-InSAR. The present invention will be described in detail herein with reference to the obtained rising and falling orbit data of a certain research area Sentinel-1A as an example.
The information of the lifting track data of the research area is shown in the following table, and the data acquisition time distribution is shown in the attached figure 2:
TABLE 1 Sentinel-1A track data description
Figure BDA0002063077210000061
Step 2: extracting sight direction public deformation time sequence obtained by inversion of orbit ascending and descending satellite platform
Respectively extracting S1A common deformation time sequences of MT-InSAR sight line deformation time sequences obtained by inversion of the orbit ascending satellite platform and the orbit descending satellite platform in a common time range, further extracting common space ranges of the MT-InSAR sight line deformation time sequences and the velocity fields of the orbit ascending satellite platform and the orbit descending satellite platform, and extracting the common range according to the common space ranges; the imaging geometry of the lifting and lowering synthetic aperture radar satellite platform is shown in figure 3.
And step 3: space registration and extraction of homonymous high-coherence points of image sets of orbit ascending and descending satellite platforms
And carrying out spatial registration on the high-coherence points extracted by the MT-InSAR of the orbit-ascending satellite platform and the MT-InSAR of the orbit-descending satellite platform of S1A, and resampling to the same resolution, so that the high-coherence points extracted by the MT-InSAR of the orbit-ascending satellite platform and the MT-InSAR of the orbit-descending satellite platform are in one-to-one correspondence. Counting and recording the number and the positions of homonymous high coherent points of the orbit ascending and descending satellite platforms;
and 4, step 4: establishing a three-dimensional deformation decomposition equation set based on minimum acceleration to generate a three-dimensional deformation sequence
a. After the preprocessing of the first three steps, data sets (up-track and down-track) of 2 different viewing angles with common time and space range are obtained, the two data sets have 23 deformation time sequences and 22 deformation time sequences respectively, and the deformation rate of the line-of-sight directions of the two data sets can be described by the following formula:
Figure BDA0002063077210000071
each row in the above equation represents the deformation time series of the up-track and down-track datasets at some common point,
Figure BDA0002063077210000072
representing the gaze direction deformation time series vector from the gaze direction deformation time series vector of the first data set to the gaze direction deformation time series vector of the second data set,
Figure BDA0002063077210000073
representing the amount of single phase accumulated deformation from the first phase to the last phase in the first data set,
Figure BDA0002063077210000074
a single phase cumulative deformation quantity representing a first phase to a last phase in the second data set;
likewise, the acquisition date corresponding to the lifting rail time series can be expressed as:
Figure BDA0002063077210000075
in the above formula, Ti(i ═ 1, 2) represents the data acquisition date vector for the first data set and the data acquisition date vector for the second data set,
Figure BDA0002063077210000076
representing the first data set from the firstDate of data acquisition from the date of the last date,
Figure BDA0002063077210000077
representing a date of data acquisition from the first date to the last date in the second data set;
the acquisition dates corresponding to different single platform time series in the above formula may be the same (contemporaneous) or different (episodic). For the different-period multi-platform MT-InSAR deformation time sequence, the acquisition time is only required to be sequenced according to the sequence; for the deformation time sequence of the synchronous multi-platform MT-InSAR, date deduplication processing needs to be carried out on the deformation time sequence, and the finally obtained total date sequence is unique.
b. Considering the influence of the synchronization and the non-synchronization, the invention performs a duplicate check before sorting all the acquisition dates, and the process can be described by the following formula:
T＝Sort(Unique(T1，T2))＝[t0，t1，t2，...，t22] (3)
in the formula, T represents a finally generated non-repeated time vector, Sort represents a sorting function for sorting all times in order according to the sequence, Unique represents a deduplication function for reserving only one repeated date, and T1，T2Data acquisition date vector, t, representing a first data set and a second data set0，t1，t2，...，t22Representing all data set date sequences retained after sorting and deduplication;
after the non-repeating date sequence T is obtained, the date list offsets are subtracted to obtain 22 difference results, expressed by the mathematical expression:
ΔT＝[Δt1，Δt2，...，Δt22] (4)
where Δ T represents the difference date interval obtained by subtracting the date sequence in T one after the other, Δ T1，Δt2，...，Δt22Representing the first to 22 nd differential date intervals obtained by subtracting the date sequences in T from offset, e.g. Δ T1Represents T in T1-t0Result of (1), Δ t2Represents T in T2-t1Result of (1), Δ t22Represents T in T22-t21The result of (1);
because the earth surface deformation result obtained by the MT-InSAR method only reflects the component of the satellite line-of-sight (LOS), according to the figure 1 of the satellite imaging geometric relationship, the deformation in the LOS direction is decomposed to the east-west direction and the vertical direction, the south-north direction component is ignored, a simplified two-dimensional deformation decomposition mode is adopted, and the decomposition process can be described by the following formula:
dLOS＝sinθ*dEast-West+cosθ*dUp-Down (5)
in the above formula, θ is the satellite view angle, dLOSRepresenting the amount of deformation of the satellite's line of sight, dEast-WestThe amount of deformation in the east-west direction, dUp-DownIs the amount of deformation in the vertical direction.
c. Considering 22 different date intervals in total, the following equation can be obtained by representing the average deformation speed in the east-west direction and the vertical direction in each date interval by a matrix:
Figure BDA0002063077210000081
in the above formula VEast-West，VNorth-south，VUp-DownAnd represents the average deformation velocity vector v in the east-west direction, the south-north direction and the vertical directionE1，υE2，...，υE22，υU1，υU2，...，υU22Representing the average deformation velocity component in each period of time interval in the east-west direction and the vertical direction delta T;
with the expression of the above equation, the average deformation rate of each date interval can be combined to form the total deformation amount sequence in the formula (5), and the rate-deformation amount conversion equation can be expressed as the following formula:
Figure BDA0002063077210000082
Figure BDA0002063077210000083
in the above formula, B is a rate transformation variable reconstruction comprehensive matrix, Bi(i ═ 1, 2) is a rate distortion quantity reconstruction matrix, which is specifically expressed in the following equation (10), VEast-West，VUp-DownRepresenting the average deformation velocity vector in the east-west direction and the vertical direction,
Figure BDA0002063077210000091
representing the gaze direction deformation time series vector from the gaze direction deformation time series vector of the first data set to the gaze direction deformation time series vector of the second data set, theta being the satellite view,
Figure BDA00020630772100000911
is the satellite flight direction course angle.
To construct BiIt is necessary to find the acquisition dates T corresponding to the time series of the lifting rail data sets respectivelyi(i-1, 2) positions in the sorted time interval sequence Δ T are set as
Figure BDA0002063077210000092
Can be expressed as:
Figure BDA0002063077210000093
in the above formula
Figure BDA0002063077210000094
A date vector T representing the first data set1The position index in the time vector T of all data sets after sorting and deduplication,
Figure BDA0002063077210000095
a date vector T representing the second data set2Position sequence numbers in all sorted and deduplicated data set time vectors T;
with the help of the above formula, BiCan be constructed using the following formula:
Figure BDA0002063077210000096
in the above formula BiA matrix is reconstructed for the rate-distortion variables,
Figure BDA0002063077210000097
represents the position number in formula (9);
d. the invention introduces the minimum acceleration assumption into the solution of the matrix, and by limiting deformation components in all directions to have the minimum acceleration in time sequence, namely the speed difference between continuous time intervals is minimum, the equation for solving the deformation decomposition has physical significance, and a stable solution can be obtained. The above description may be specifically expressed as:
Figure BDA0002063077210000098
where C is the minimum acceleration defining equation set, and α is a regularization factor that balances the relative weights of the minimum acceleration to limit the residual norm minima in the equation,
Figure BDA0002063077210000099
representing the average deformation speed difference result between adjacent date intervals in the east-west direction and the vertical direction;
in conclusion, the complete synchronous and asynchronous lifting and lowering rail multi-platform MT-InSAR three-dimensional deformation decomposition technology can be expressed as follows:
Figure BDA00020630772100000910
in the above formula, B is shown as formula (8), and C is shown as formula (11)Show, VEast-West，VUp-DownRepresenting the average deformation velocity vector in the east-west direction and the vertical direction,
Figure BDA0002063077210000101
representing a gaze direction deformation time series vector from a gaze direction deformation time series vector of a first data set to a gaze direction deformation time series vector of a second data set;
e. and solving the linear equation set by using a crack Berg-Marquardt method to obtain the deformation time sequences of each high coherence point in the vertical direction, the south direction, the north direction and the east-west direction.
f. And applying the steps a-e to all homonymous high coherence points to obtain a three-dimensional surface deformation field and a deformation rate based on the MT-InSAR in the lifting and lowering rail platform in a public area. In this embodiment, the east-west and vertical direction average deformation rate maps after the three-dimensional deformation decomposition are shown in fig. 4, and four feature points are selected to plot the east-west and vertical direction deformation time series obtained by the decomposition, as shown in fig. 5.
Details not described in this specification are within the skill of the art that are well known to those skilled in the art.

Claims (1)

1. A three-dimensional deformation decomposition method for a contemporaneous and an ectopic multi-satellite platform MT-InSAR based on minimum acceleration is characterized by comprising the following specific steps:
step 1: method for respectively inverting sight line direction deformation time sequence of orbit ascending satellite platform and orbit descending satellite platform by using MT-InSAR
Respectively inverting the sight line deformation time sequence and the deformation rate of different orbit ascending and descending synthetic aperture radar satellite platforms by using MT-InSAR;
step 2: extracting sight direction public deformation time sequence obtained by inversion of orbit ascending and descending satellite platform
Respectively extracting common deformation time sequences of MT-InSAR sight line direction deformation time sequences obtained by inversion of the orbit ascending satellite platform and the orbit descending satellite platform in a common time range, further extracting common space ranges of the MT-InSAR sight line direction deformation time sequences and the velocity fields of the orbit ascending satellite platform and the orbit descending satellite platform, and extracting the common range according to the common space ranges;
and step 3: space registration and extraction of homonymous high-coherence points of image sets of orbit ascending and descending satellite platforms
Carrying out spatial registration on high-coherence points extracted by the MT-InSAR of the ascending and descending satellite platforms respectively, and resampling to the same resolution, so that the high-coherence points extracted by the MT-InSAR of the ascending and descending satellite platforms respectively correspond to one another; counting and recording the number and the positions of homonymous high coherent points of the orbit ascending and descending satellite platforms;
and 4, step 4: aiming at the characteristics of the mixed occurrence of the data in the same period and the data in the different periods, the ascending and descending time vectors are arranged in a mixed mode and subjected to de-duplication processing, a deformation decomposition equation set based on minimum acceleration is established, and a three-dimensional deformation sequence is obtained through solving; the method specifically comprises the following steps:
a. after the processing of the step 1-3, the N different synthetic aperture radar satellite data sets have different visual angles and comprise orbit rising data and orbit falling data; each data set has Qi(i ═ 1, 2.., N) time series, the deformation rate for each dataset line-of-sight direction is described by:
Figure FDA0002063077200000011
each row in the above formula represents a time sequence of deformation of a synthetic aperture radar satellite platform data set at a common point, the initial time sequence period of each platform is different,
Figure FDA0002063077200000012
is the gaze direction deformation time series vector from the first data set to the gaze direction deformation time series vector of the nth data set,
Figure FDA0002063077200000013
is the accumulated deformation from the first phase to the last phase in the first data set,
Figure FDA0002063077200000014
from the first phase to the last phase in the second data setThe amount of deformation is accumulated at each stage,
Figure FDA0002063077200000015
is the accumulated deformation quantity, Q, of each phase from the first phase to the last phase in the Nth data set1Is the time series number, Q, of the first data set2，…，QNIs the number of time series from the second data set to the Nth data set;
likewise, the date of acquisition for each time series of data sets is represented as:
Figure FDA0002063077200000021
in the above formula, Ti(i 1, 2.., N) is the data acquisition date vector for the first data set to the data acquisition date vector for the nth data set,
Figure FDA0002063077200000022
is the date of data acquisition from the first date to the last date in the first data set,
Figure FDA0002063077200000023
is the date of data acquisition from the first date to the last date in the second data set,
Figure FDA0002063077200000024
is the date of data acquisition from the first to the last phase in the Nth data set, Q1Is the time series number, Q, of the first data set2，…，QNIs the number of time series from the second data set to the Nth data set;
aiming at a synthetic aperture radar satellite constellation, different synthetic aperture radar satellites can image the same area for multiple times in the same day, namely, synchronous data; in the above formula, the acquisition dates corresponding to different single platform time sequences may be the same, i.e. in-phase, or different, i.e. out-of-phase; for the different-period multi-platform MT-InSAR deformation time sequence, the acquisition time is sequenced according to the sequence; for a synchronous multi-platform MT-InSAR deformation time sequence, date deduplication processing is carried out on the time sequence, and the finally obtained total date sequence is unique;
b. the date of the same period and the date of the different period are mixed and checked for duplication before all the acquisition dates are sorted, and the process is described by the following formula:
T＝Sort(Unique(T1，T2，...，TN))＝[t0，t1，t2，...，tQ-1] (3)
in the formula, T is a finally generated non-repeated time vector, Sort is a sorting function and is used for sorting all times in order according to the sequence, Unique is a deduplication function and is used for reserving only one date which appears repeatedly, and T1，T2，...，TNIs a data acquisition date vector, t, from the first data set to the Nth data set0，t1，t2，...，tQ-1All data set date sequences reserved after sorting and de-duplication, wherein Q is the total number of unrepeated dates in N data sets;
after a non-repeating date sequence T of N data set data acquisition dates is obtained, the date lists are subjected to staggered subtraction to obtain Q-1 difference results, and the difference results are expressed by a mathematical expression as follows:
ΔT＝[Δt1，Δt2，...，ΔtQ-1] (4)
where Δ T is the difference date interval obtained by subtracting the date sequence in T one after the other, and Δ T1，Δt2，...，ΔtQ-1Is the first to Q-1 difference date interval, delta T, obtained by the staggered subtraction of the date sequences in T1Is T in T1-t0Result of (1), Δ t2Is T in T2-t1Result of (1), Δ tQ-1Is T in TQ-1-tQ-2As a result, Q is the total number of non-repeating dates in the N data sets;
because the ground surface deformation result obtained by the MT-InSAR method is the deformation of the satellite sight line direction, the deformation of the satellite sight line direction is decomposed into the east-west direction, the south-north direction and the vertical direction according to the satellite imaging geometric relation, and the decomposition process is described by the following formula:
Figure FDA0002063077200000025
in the above equation theta is the satellite view angle,
Figure FDA0002063077200000026
is the satellite flight direction course angle, dLOsRepresenting the amount of deformation of the satellite's line of sight, dEast-WestThe amount of deformation in the east-west direction, dNorth-SouthIs the north-south direction deformation, dUp-DownIs the deformation amount in the vertical direction; in addition, the flight direction of the satellite is in the near south-north direction, the inversion deformation result is insensitive to the south-north direction, and if the deformation component in the south-north direction is not considered, the inversion deformation result is simplified into a two-dimensional deformation decomposition method expressed as follows: angle of course
Figure FDA0002063077200000038
The distortion of the satellite sight direction is only decomposed into east-west and vertical directions in the above formula;
c. there are Q-1 different date intervals, and the average deformation speeds in the east-west direction, the south-north direction and the vertical direction in each date interval are represented by a matrix, obtaining the following formula:
Figure FDA0002063077200000031
in the above formula VEast-West，
Figure FDA0002063077200000032
VUp-DownRepresents the average deformation velocity vector in the east-west direction, the south-north direction and the vertical direction, vE1，vE2，...，vEQ-1Is the average deformation velocity component, v, in the east-west direction Δ T in each period time intervalN1，vN2，...，vNQ-1Time intervals of respective periods in the north-south direction Δ TInner average deformation velocity component, vU1，vU2，...，vUQ-1Is the average deformation velocity component in each period of time interval in the vertical direction delta T, and Q is the total number of unrepeated dates in N data sets;
with the expression of the above formula, the average deformation rate of each date interval can be combined into a total deformation amount sequence by combining the formula (5), and the rate-deformation conversion equation is expressed as the following formula:
Figure FDA0002063077200000033
Figure FDA0002063077200000034
in the above formula, B is a rate transformation variable reconstruction comprehensive matrix, Bi(i ═ 1, 2.., N) is a rate distortion quantity reconstruction matrix, which is specifically expressed in the following equation (10), VEast-West，
Figure FDA0002063077200000035
VUp-DownThe average deformation velocity vectors in the east-west direction, the south-north direction and the vertical direction are respectively,
Figure FDA0002063077200000036
is the time-series vector of the line-of-sight deformation from the time-series vector of the line-of-sight deformation of the first data set to the time-series vector of the line-of-sight deformation of the Nth data set, theta is the satellite view,
Figure FDA0002063077200000039
is the satellite flight direction course angle;
to construct Bi(i 1, 2.. times.n), the position of the acquisition date corresponding to the time series of each data set in the sorted time interval series Δ T is set as
Figure FDA0002063077200000037
Expressed as:
Figure FDA0002063077200000041
in the above formula
Figure FDA0002063077200000042
Is the date vector T of the first data set1The position index in the time vector T of all data sets after sorting and deduplication,
Figure FDA0002063077200000043
to
Figure FDA0002063077200000044
Is the date vector T of the second data set2Position sequence numbers in time vector T of all sorted and deduplicated data sets to date vector T of Nth data setNPosition order number, Q, in time vector T of all data sets after sorting and de-duplicationi(i 1, 2.., N) represents the number of time series of 1 to N per data set;
with the help of the above formula, BiConstructed using the following formula:
Figure FDA0002063077200000045
in the above formula BiA matrix is reconstructed for the rate-distortion variables,
Figure FDA0002063077200000046
is the position serial number in the formula (9), and Q is the total number of unrepeated dates in the N data sets;
d. introducing a minimum acceleration assumption into the solution of the formula (7), and limiting deformation components in all directions to have minimum acceleration in time sequence, namely the speed difference between continuous time intervals is minimum, so that an equation for solving deformation decomposition has physical significance and a stable solution is obtained; the above description may be specifically expressed as:
Figure FDA0002063077200000047
where C is the minimum acceleration defining equation set, and α is a regularization factor that balances the relative weights of the minimum acceleration to limit the residual norm minima in the equation,
Figure FDA0002063077200000048
the average deformation speed difference results of adjacent date intervals in the east-west direction, the south-north direction and the vertical direction are respectively obtained, and Q is the total number of unrepeated dates in the N data sets;
synthesizing equations (7) - (11) to obtain equation (12) for solving the average deformation speed in the east-west direction, the south-north direction and the vertical direction;
Figure FDA0002063077200000049
in the above formula, B is shown in formula (8), C is shown in formula (11), and VEast-West，
Figure FDA00020630772000000410
VUp-DownThe average deformation velocity vectors in the east-west direction, the south-north direction and the vertical direction are respectively,
Figure FDA0002063077200000051
is a line-of-sight direction deformation time series vector from the line-of-sight direction deformation time series vector of the first data set to the line-of-sight direction deformation time series vector of the Nth data set;
e. solving the linear equation set by using a crack Berge-Marquardt method to obtain deformation time sequences of each high coherence point in the vertical direction, the south direction, the north direction and the east-west direction;
f. and applying the steps a-e to all homonymous high coherence points to obtain a three-dimensional surface deformation field and a deformation rate based on the MT-InSAR in the lifting and lowering rail platform in a public area.
CN201910411920.5A 2019-05-17 2019-05-17 Simultaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method based on minimum acceleration Active CN110146883B (en)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201910411920.5A CN110146883B (en) 2019-05-17 2019-05-17 Simultaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method based on minimum acceleration

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201910411920.5A CN110146883B (en) 2019-05-17 2019-05-17 Simultaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method based on minimum acceleration

Publications (2)

Publication Number Publication Date
CN110146883A CN110146883A (en) 2019-08-20
CN110146883B true CN110146883B (en) 2022-04-05

Family

ID=67595465

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201910411920.5A Active CN110146883B (en) 2019-05-17 2019-05-17 Simultaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method based on minimum acceleration

Country Status (1)

Country Link
CN (1) CN110146883B (en)

Families Citing this family (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111458709B (en) * 2020-06-08 2023-12-22 河南大学 Method and device for monitoring wide-area earth surface two-dimensional deformation field of spaceborne radar
CN113091598B (en) * 2021-04-06 2022-02-08 中国矿业大学 Method for defining stability grade range of goaf building site by InSAR
CN113281748B (en) * 2021-05-24 2022-02-11 西南石油大学 Surface deformation monitoring method
EP4361678A1 (en) * 2022-10-27 2024-05-01 Planetek Italia S.r.l. System for determining the components of the displacement of the persistent scatterers and method thereof

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108983232A (en) * 2018-06-07 2018-12-11 中南大学 A kind of InSAR two dimension earth's surface deformation monitoring method based on adjacent rail data
CN109540159A (en) * 2018-10-11 2019-03-29 同济大学 A kind of quick complete automatic Pilot method for planning track

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US9417323B2 (en) * 2012-11-07 2016-08-16 Neva Ridge Technologies SAR point cloud generation system

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN108983232A (en) * 2018-06-07 2018-12-11 中南大学 A kind of InSAR two dimension earth's surface deformation monitoring method based on adjacent rail data
CN109540159A (en) * 2018-10-11 2019-03-29 同济大学 A kind of quick complete automatic Pilot method for planning track

Non-Patent Citations (5)

* Cited by examiner, † Cited by third party
Title
Investigation of the ground displacement in Saint Petersburg, Russia, using multiple-track differential synthetic aperture radar interferometry;Qiang Wang et al.;《Int J Appl Earth Obs Geoinformation》;20200111;正文全文 *
On the Effects of InSAR Temporal Decorrelation and Its Implications for Land Cover Classification:The Case of the Ocean-Reclaimed Lands of the Shanghai Megacity;Guanyu Ma et al.;《Sensors》;20180904;正文全文 *
基于DInSAR和MAI技术揭示地震三维形变场;杨红磊等;《地球物理学进展》;20141215(第06期);正文全文 *
基于多卫星平台永久散射体雷达干涉提取三维地表形变速度场;刘国祥等;《地球物理学报》;20120815(第08期);正文全文 *
基于多级联合策略的多平台 MT-InSAR 形变时间序列联合分析—以上海市新成陆区为例;禹雷;《中国知网》;20190115;第二章 *

Also Published As

Publication number Publication date
CN110146883A (en) 2019-08-20

Similar Documents

Publication Publication Date Title
CN110146883B (en) Simultaneous and asynchronous multi-satellite platform MT-InSAR three-dimensional deformation decomposition method based on minimum acceleration
Joughin et al. A complete map of Greenland ice velocity derived from satellite data collected over 20 years
CN110058236B (en) InSAR and GNSS weighting method oriented to three-dimensional surface deformation estimation
CN113624122B (en) Bridge deformation monitoring method fusing GNSS data and InSAR technology
CN109738892B (en) Mining area earth surface high-space-time resolution three-dimensional deformation estimation method
Doin et al. Presentation of the Small Baselin NSBAS Processing Chain on a Case Example: The Etan Deformation Monitoring from 2003 to 2010 Using Envisat Data
CN111458709B (en) Method and device for monitoring wide-area earth surface two-dimensional deformation field of spaceborne radar
KR101111689B1 (en) The method for three-dimensional deformation measurement and the apparatus thereof
CN108983232B (en) InSAR two-dimensional surface deformation monitoring method based on adjacent rail data
CN103091676A (en) Mining area surface subsidence synthetic aperture radar interferometry monitoring and calculating method
CN108663678B (en) Multi-baseline InSAR phase unwrapping algorithm based on mixed integer optimization model
CN106842199A (en) It is a kind of to merge the method that different resolution SAR data monitors Ground Deformation
CN103714546A (en) Data processing method of imaging spectrometer
Duffy et al. Surface subsidence in urbanized coastal areas: PSI methods based on Sentinel-1 for Ho Chi Minh City
CN108171732A (en) A kind of detector lunar surface landing absolute fix method based on multi-source image fusion
CN111623749B (en) Novel railway goaf boundary extraction method based on D-InSAR
CN116026283A (en) Method for monitoring and evaluating earth surface subsidence reducing effect of coal mine filling mining
CN115201825A (en) Atmospheric delay correction method in InSAR (interferometric synthetic aperture radar) inter-seismic deformation monitoring
CN112526515A (en) Surface deformation detection method based on synthetic aperture radar interferometry
Wu et al. Detecting the deformation anomalies induced by underground construction using multiplatform MT-InSAR: A case study in to Kwa Wan Station, Hong Kong
CN115096257A (en) Mining area mining subsidence monitoring method and device
CN112835043B (en) Method for monitoring two-dimensional deformation in any direction
CN111505635B (en) GEO SAR two-star formation configuration design method for coherent chromatography
CN112505696A (en) Strip mine slope deformation dynamic monitoring method based on satellite-borne SAR
Nikolakopoulos et al. Preliminary results of using Sentinel-1 SAR data for DSM generation

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