Earthquake fracture process inversion method combined with InSAR (Interferometric Synthetic Aperture Radar) relative phase gradient

Abstract

The invention discloses an InSAR relative phase gradient-combined seismic fracture process inversion method, which comprises the following steps of: performing gradient operation on an InSAR phase difference to obtain an InSAR relative gradient; the relative gradient serves as observation data, the minimum residual error of the relative gradient serves as a target function, the curvature and radian of the interference fringes are restrained, and therefore seismic co-seismic sliding distribution inversion based on the interference fringes which cannot be unwound on the dive zone near-field island is achieved. By adopting the technical scheme of the invention, the inversion precision of the same-seismic sliding distribution of the subduction zone earthquake under the observation constraint of the unconstrained InSAR is improved.

Description

Technical Field

[0001] This invention belongs to the field of seismic inversion technology, specifically relating to a seismic rupture process inversion method based on joint InSAR relative phase gradient. Background Technology

[0002] Since the 1990s, when Interferometric Synthetic Aperture Radar (InSAR) technology was first used for monitoring coseismic deformation, InSAR observation data has gradually become an important data source for studying coseismic deformation. InSAR coseismic deformation observations, combined with source model inversion algorithms, can effectively infer the distribution of coseismic slip, helping us to better understand the physical processes involved in earthquakes.

[0003] When inverting coseismic slip distribution in subduction zones, incoherence in sea surface InSAR observations can lead to unwrapped InSAR interferometric fringes on islands near the epicenter. The applicant has proposed an invention for inverting coseismic slip distribution based on constraints of the relative magnitude of coseismic deformation at various points on the island. However, this invention cannot constrain the curvature and arc of the interferometric fringes, thus limiting the ability to invert earthquake rupture processes using near-field InSAR observations. Summary of the Invention

[0004] To address the problems existing in the prior art, this invention provides a seismic rupture process inversion method based on the relative phase gradient of InSAR, which improves the inversion accuracy of coseismic slip distribution in subduction zones under unconstrained InSAR observation constraints.

[0005] To achieve the above objectives, the present invention provides the following solution: A seismic rupture process inversion method based on joint InSAR relative phase gradient includes: The InSAR relative gradient is obtained by performing gradient calculation on the InSAR phase difference; By using the relative gradient as the observation data and minimizing the residual of the relative gradient as the objective function, the curvature and arc of the interference fringes are constrained to achieve the inversion of the coseismic slip distribution of near-field islands in the subduction zone.

[0006] As a preferred approach, the mathematical expression of the observation values ​​is optimized by removing the flatland effect and topographic undulation influence of the unreferenced InSAR deformation observation values ​​in the area to be measured, and the InSAR relative gradient between the reference point and any point is obtained.

[0007] As a preferred approach, by combining static GNSS, seismic wave, strong ground motion, and high-frequency GNSS waveform observation data, the Green's function between InSAR observations and sub-fault slip is calculated using dislocation theory. Then, a simulated annealing nonlinear inversion algorithm is employed, with the minimum residual of the relative gradient as the objective function, to search and determine the coseismic slip amount and direction on each sub-fault, thereby achieving the inversion of the coseismic slip distribution of the earthquake.

[0008] Compared with the prior art, the beneficial effects of the present invention are as follows: This invention optimizes the mathematical expression of observations by removing the flat-land effect and topographic undulation influence of unreferenced InSAR deformation observations in the test area, and obtains the relative phase difference between the reference point and any other point. By performing gradient calculations on the relative phase difference, the relative gradient is incorporated as observation data into the inversion, and the minimum residual of the relative gradient is used as the objective function to constrain the curvature and arc of the interferometric fringes. Furthermore, this invention can achieve the inversion of coseismic slip distributions of seismic fringes on near-field islands where direct unwrapping of interferometric fringes is not possible, solving the problems of missing near-field geodetic data, poor accuracy of results, and insufficient spatial constraint capability of the slip model in subduction zones. Attached Figure Description

[0009] To more clearly illustrate the technical solution of the present invention, the drawings used in the embodiments are briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0010] Figure 1 This is a flowchart of the seismic rupture process inversion method based on InSAR relative phase gradient, as described in an embodiment of the present invention. Figure 2 For 2025 Mw Location of the epicenter and surrounding geological structure of the 6.9 Philippine Islands earthquake; Figure 3 For 2025 Mw 6.9 Fitting of InSAR data for earthquakes in the Philippine Islands and distribution of residuals. Detailed Implementation

[0011] The technical solutions of the embodiments of the present invention will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of the present invention, and not all embodiments. Based on the embodiments of the present invention, all other embodiments obtained by those skilled in the art without creative effort are within the scope of protection of the present invention.

[0012] To make the above-mentioned objects, features and advantages of the present invention more apparent and understandable, the present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments.

[0013] Example 1 like Figure 1 As shown, this invention provides a method for inverting seismic rupture processes using a combination of InSAR relative phase gradients, comprising: Step S1: Obtain the InSAR relative gradient by performing gradient calculation on the InSAR phase difference; Step S2: Using the relative gradient as the observation data and minimizing the residual of the relative gradient as the objective function, the curvature and arc of the interference fringes are constrained to realize the inversion of the coseismic slip distribution of near-field islands in the subduction zone.

[0014] In one embodiment of the present invention, in step S1, after removing the effects of flat terrain and topographic relief from the island InSAR observations, the observation equation is obtained as follows: Where, on the left side of the equation The number of integer cycles for the interference fringes is unknown. These are InSAR phase observations; the right side of the equation , , These are the phase changes caused by deformation, atmospheric errors, and observation errors, respectively. Choose any point on the island Using this as a reference point, we can obtain any other point. The deformation difference between reference points is: Because the island is small, the atmospheric error between two adjacent points is... It can be approximated as 0. The simplified equation (2) is: As one embodiment of the present invention, in step S2, gradient calculation is performed on both sides of equation (3) to obtain the following equation: In the formula Since it is a constant, If the value is 0, the above formula can be further simplified to: On the left side of the equation Calculate using the following formula: in, For point Nearby Points The relative phase difference, , Points , Distance in the longitude and latitude directions.

[0015] Because the right side of equation (5) That is, any point Deformation gradient and reference point The difference in deformation gradient, and Following a standard normal distribution, after performing gradient calculation Since it still follows a normal distribution, the right side of the equation represents the coseismic deformation gradient and random error. Combining equations (5) and (6), the relative gradient is used as the observation data, and the minimum residual of the relative gradient is used as the objective function to constrain the curvature and arc of the interference fringes.

[0016] After obtaining the coseismic deformation gradient observations at various points on the island, the seismogenic fault plane determined by seismic observations can be discretized into several rectangular or triangular sub-faults. The slip of each sub-fault is determined by four parameters: slip magnitude, slip azimuth, rupture rise time, and rupture velocity. The mapping relationship between the relative coseismic deformation observations at various points on the island and the unit slip of each sub-fault is a Green's function. The sum of the effects of all sub-fault slips on the same surface station represents the overall response of the surface station during the earthquake.

[0017] By combining static GNSS, seismic waves, and other geodetic and geophysical observation data, a multi-source data-constrained seismic rupture inversion method integrating InSAR floating-point observations can be constructed. Static GNSS data can be used with Bernese software to calculate pre- and post-seismic three-dimensional displacements, subtracting linear trends and seasonal variations to obtain coseismic static displacements and observation errors. Seismic wave data is processed by filtering P / SH waves at epicentral distances of 30-90°, applying 0.01-1Hz bandpass filtering and removing instrument response, resampling to 5Hz, and manually picking out P / SH wave arrival times.

[0018] The response of each sub-fault to the observation point is calculated using a uniform elastic half-space Okada dislocation model. The objective function for InSAR observation is to minimize the residual of the relative gradient, as follows: in The phase gradient difference is the predicted phase gradient difference for the sub-fault slip.

[0019] The objective function for GNSS data is , For the observed values, These are the model's predicted values. The observation error is represented by the objective function for seismic wave data. , For low-frequency wavelet coefficient residuals, The residuals are high-frequency cross-correlation. Spatial Laplace smoothing and temporal smoothing are added to suppress non-physical slip oscillations and constrain the continuity of the fracture velocity, resulting in the joint inversion objective function: in , , The weighting coefficients are determined using the L-curve method to find the optimal values.

[0020] Finally, a simulated annealing nonlinear inversion algorithm was used to search and determine the coseismic slip and direction on each sub-fault. First, an initial temperature and cooling coefficient were set, and an initial slip model was randomly generated, including parameters such as sub-fault slip and direction. Then, an iterative search was performed, with each iteration randomly perturbing the parameters of a single sub-fault, and the objective function of the new model was calculated. According to the Metropolis criterion, the new model is accepted: if If the value decreases, accept it directly; if it increases, accept it based on probability. accept, The current temperature is used. The temperature is gradually decreased, and the iteration is repeated until convergence, i.e., the residual stabilizes or the maximum number of iterations is reached.

[0021] After convergence, parameters of each sub-fault are extracted to generate a slip distribution, thereby making full use of near-field island InSAR floating-point observations, strengthening the ability of near-field island InSAR observations to constrain the seismic rupture process, and realizing the inversion of the coseismic slip distribution of earthquakes.

[0022] In this embodiment of the invention, 2025 is used. Mw 6.9 The Philippine Islands earthquake is used as a case study of subduction zone earthquakes. This paper introduces the coseismic slip distribution inversion method based on InSAR relative gradient. In embodiments of the present invention, such as Figure 2 (2025) Mw As shown in the image (Epicenter location and surrounding geological structure of the 6.9 Philippine Islands earthquake), according to the U.S. Geological Survey (USGS), an earthquake occurred in the sea near the Philippine Islands at 13:59:43 (UTC) on September 30, 2025. Mw 6.9 earthquake, Figure 2 It shows the geographical environment around the epicenter of this earthquake.

[0023] InSAR observations are limited by the discontinuity of the observation fringes, making it impossible to accurately constrain near-field coseismic deformation.

[0024] This invention embodiment is based on GAMMA software and uses the bisection method to process 32 scenes of 69 orbits of Sentinel-1 data from September 24, 2025 and October 6, 2025 to obtain InSAR observation stripes generated by seismic deformation.

[0025] In this embodiment of the invention, the earthquake epicenter is located on the seabed, with a lack of continuous land in the surrounding area, and most of the deformation occurs on the Philippine Islands. Due to the obstruction of seawater, the interferometric fringes between the islands are not continuous, and the near-field InSAR interferometric images cannot be directly unwrapped into a single LOS displacement field, which limits the acquisition of the coseismic deformation field. This invention uses SNAPHU software to unwrap the phase of each island, thereby obtaining the relative phase of the InSAR fringes and its gradient observations on each island.

[0026] Figure 3 The figures show the InSAR data fitting and residual distribution in this embodiment of the invention. (a), (b), and (c) represent the phase observation, fitted value, and residual, respectively. The results show that the rupture model has a good fit to InSAR deformation observations. On islands constrained only by InSAR relative gradients, the InSAR interferometric fringes are well recovered, verifying the theoretical correctness of the invention and demonstrating its strong potential application prospects.

[0027] The embodiments described above are merely preferred embodiments of the present invention and are not intended to limit the scope of the present invention. Various modifications and improvements made by those skilled in the art to the technical solutions of the present invention without departing from the spirit of the present invention should fall within the protection scope defined by the claims of the present invention.

Claims

1. A method for inverting seismic rupture processes using joint InSAR relative phase gradient, characterized in that, include: The InSAR relative gradient is obtained by performing gradient calculation on the InSAR phase difference; By using the relative gradient as the observation data and minimizing the residual of the relative gradient as the objective function, the curvature and arc of the interference fringes are constrained to achieve the inversion of the coseismic slip distribution of near-field islands in the subduction zone.

2. The seismic rupture process inversion method based on joint InSAR relative phase gradient as described in claim 1, characterized in that, By removing the flatland effect and topographic undulation influence of the unreferenced InSAR deformation observations in the test area, the mathematical expression of the observations is optimized, and the InSAR relative gradient between the reference point and any point is obtained.

3. The seismic rupture process inversion method based on joint InSAR relative phase gradient as described in claim 2, characterized in that, By combining static GNSS, seismic wave, strong ground motion, and high-frequency GNSS waveform observation data, the Green's function between InSAR observations and sub-fault slip is calculated using dislocation theory. Then, a simulated annealing nonlinear inversion algorithm is used, with the minimum residual of the relative gradient as the objective function, to search and determine the coseismic slip amount and direction on each sub-fault, thereby achieving the inversion of the coseismic slip distribution of earthquakes.

Images