Fault slip model smoothing method based on SAR observation data

Abstract

The invention discloses a fault slip model smoothing method based on SAR observation data, and the method comprises the steps: firstly, obtaining multi-source observation data (including SAR, GNSS and optical images), and carrying out the joint calculation, and generating a high-resolution earth surface deformation field; the method comprises the following steps: combining a remote sensing image, visually interpreting and drawing a fracture trace, then carrying out data downsampling along an inclination angle direction, and finally constructing an irregular fault geometric model by using a TIN method so as to obtain an initial fault plane model. On this basis, a Laplacian operator suitable for triangular dislocation is defined, a slip inversion model is constructed in combination with an elastic dislocation theory, and a smoothing factor is introduced to balance slip details and calculation efficiency. And finally, optimizing a fault inclination angle and a smoothing factor through a grid search method, evaluating fitting errors and slippage roughness under different combinations, and determining optimal parameters. The final model not only comprises a triangular connection relation, but also provides slip amount distribution and fault dip angle information, and an important basis is provided for seismic activity evaluation and disaster response.

Description

Technical Field

[0001] This invention belongs to the field of earthquake monitoring and assessment technology, specifically relating to a method for smoothing fault slip models based on SAR observation data. Background Art

[0002] As a global geological hazard, earthquakes can induce severe casualties and property damage. With the development of geodesy, especially the supplementation of SAR imagery data, it has become a reality to conduct large-scale, long-term, high-precision, non-contact monitoring of surface deformation in seismic zones and areas with high fault activity. A wealth of surface deformation observation data also provides fundamental data for fault model inversion throughout the entire lifecycle of coseismic, inter-seismic, and post-seismic events, which is of great significance for disaster assessment and emergency response.

[0003] Currently, fault model construction is mainly based on "rectangular dislocation" and "irregular dislocation" geometric models. The former divides the fault plane into multiple matrix sub-faults and determines the fault slip distribution through a "two-step" inversion (linear inversion and nonlinear inversion). The latter divides the fault plane into multiple irregular triangles and inverts the slip of each sub-fault plane. Compared with the rectangular dislocation model, this model construction method considers the slip discontinuities and overlaps caused by complex structures such as fault plane steps, bifurcation, and bending, and has significant advantages in complex geometric scenarios, especially strike-slip earthquakes that are often accompanied by long rupture trajectories. However, increasing the number of irregular triangular sub-faults will decrease computational efficiency, while insufficient numbers will prevent the model from accurately depicting the complex rupture trajectories of earthquakes. In addition, the fault dip angle, as one of the important parameters controlling the geometry of the model, plays a crucial role in determining its optimal value for the reliability of slip distribution inversion. Therefore, how to balance fault detail and computational efficiency, and determine the optimal fault geometry, is a problem that needs to be considered in the construction of complex fault models.

[0004] Furthermore, since fault slip distributions typically exhibit a certain degree of smoothness, a smoothing factor is often introduced during the inversion process to smooth the slip gradient between adjacent sub-faults and avoid slip oscillations in the model. However, an excessively large smoothing factor can also cause the model to lose detailed information, especially in complex geometric structures, leading to reduced model reliability. For traditional rectangular dislocation models, some studies have smoothed fault slip distributions by solving for the minimum of the Laplacian second-order smoothing matrix. However, in triangular dislocations, adjacent sub-faults may exhibit non-planarity and bifurcation, requiring the determination of the smoothing factor to consider more complex spatial structural features. Therefore, a method for calculating the smoothing factor in triangular dislocation models that considers the complex spatial relationships of irregular fault models is urgently needed. Summary of the Invention

[0005] In order to solve the technical problems existing in the background art, the present invention aims to provide a smoothing method for fault slip model based on SAR observation data. By constructing a triangular dislocation fault model, the surface rupture trajectory is described more accurately. A grid search algorithm is used to optimize the smoothing factor and the optimal dip angle of the fault, and the optimal fault geometry is determined. A smoothing factor calculation method is proposed to solve the problem of slip distribution oscillation in irregular triangular dislocation models.

[0006] To solve the technical problem, the technical solution of the present invention is as follows:

[0007] A method for smoothing fault slip models based on SAR observation data, the method comprising:

[0008] S1: Acquire surface deformation observation data and remote sensing-aided information, and jointly calculate to generate a high-resolution surface deformation field; combine the deformation field and remote sensing image visual interpretation to draw rupture traces; perform data downsampling along the dip direction; and use the TIN method to construct an irregular fault geometry based on the sampling points to obtain an initial fault plane model.

[0009] S2: Based on the initial fault geometry model, node positions and triangular patch connections, determine the spatial relationships between sub-faults, construct the Laplacian operator and the smoothing matrix, and obtain the Laplacian operator matrix.

[0010] S3: Input surface deformation observation data, initial fault geometric model and Laplacian operator matrix, construct Green's function matrix, introduce smoothing factor, perform slip inversion, and output initial slip distribution model, including the slip amount of each sub-fault segment;

[0011] S4: Perform a grid search on the initial slip distribution model, the fault dip angle to be optimized, and the range and step size of the smoothing factor. Repeat the modeling process from S1 to S3 for different combinations of dip angles and smoothing factors. Then, evaluate the model and compare the fitting error and slip roughness under different combinations. When increasing the roughness no longer significantly reduces the error, determine the optimal parameters. Finally, extract the optimal solution. Using slip roughness and fitting error as the standard, output the optimal dip angle and smoothing factor, which yields the optimal geometric fault model, containing triangular connection relationships, slip distribution, and complete dip angle information.

[0012] Furthermore, step S1 includes:

[0013] S101: Obtain Data

[0014] Collect surface deformation observation data and remote sensing auxiliary information. The surface deformation observation data includes SAR, optical imagery, and GNSS; the remote sensing auxiliary information includes remote sensing imagery of fault rupture lines and existing seismic zone data.

[0015] S102: Joint Solution

[0016] By using data fusion and computation techniques, multi-source observation data are comprehensively processed to generate a high-resolution surface deformation field;

[0017] S103: Drawing the fracture trace

[0018] By combining the generated surface deformation field with remote sensing imagery, rupture traces are visually interpreted and mapped to show the rupture characteristics and locations on the surface.

[0019] S104: Data downsampling

[0020] Data downsampling is performed along the defined rupture trajectory, along the dip direction of the fault, and the number of sampling points is reduced in areas of greater depth to improve computational efficiency;

[0021] S105: Constructing fault geometry

[0022] Based on the downsampled sampling points, the irregular fault geometry is constructed using the TIN method, and the initial fault plane model is finally output, including node positions, triangular connection relationships and dip direction.

[0023] Furthermore, in step S2, a Laplacian computation matrix suitable for the triangular dislocation model is established:

[0024]

[0025] In the formula, ξ is the Laplacian operator; p0 is the slip on the central sub-fault patch0; p i The slip of the sub-fault adjacent to patch0; l i α is the distance between adjacent sub-faults; i It is the dihedral angle between adjacent sub-faults.

[0026] Furthermore, step S3 includes:

[0027] S301: Based on the theory of elastic dislocations, establish the Green's function matrix G = [G...] that connects earthquake source parameters and surface displacements. s G d ]:

[0028] G s =f(x,y,x0,y0,depth,strike,dip,length,width,s s ,0)

[0029] G d=f(x,y,x0,y0,depth,strike,dip,length,width,0,s d )

[0030] In the formula, G s , G d These are the surface displacement vectors under a unit slip disturbance in the strike and dip directions, respectively. The remaining parameters are the surface deformation observation location (x, y), fault location (x0, y0), depth, strike, dip, length, and width. s and s d The value is usually a constant of 1.0. In addition, the InSAR deformation observation matrix D is a one-dimensional displacement field in the LOS direction. Before calculating the Green's function matrix, the LOS deformation needs to be projected onto the three-dimensional displacement field of the Earth's surface.

[0031] S302: In the classical elastic dislocation theory model GX=D, a Laplacian operator and smoothing factor applicable to triangular dislocations are introduced, and the slip of each sub-fault is estimated based on the least squares criterion.

[0032] S303: The constructed slip model is as follows:

[0033] [D 0] T =[G ε 2 ξ] T p

[0034] In the formula, D is the surface deformation matrix, G is the Green's function, ε is the smoothing factor, ξ is the Laplacian operator, and p is the fault slip matrix.

[0035] Furthermore, in step S4, considering the one-to-one correspondence between slip roughness and smoothing factor, it can be used as the selection index for the optimal value. The specific solution method is as follows:

[0036]

[0037] In the formula, n is the number of sub-fault segments in the fault model, r is the total roughness of the fault model, i.e., the average slip gradient, usually expressed in cm / km or m / km. i It is a weighted slip distribution matrix with smoothing applied.

[0038] Compared with the prior art, the advantages of the present invention are:

[0039] Based on high-resolution geodetic data, this paper realizes the inversion and slip modeling of faults with complex geometries (including bifurcation, step regions, and bending features), providing a more accurate description of surface rupture trajectories. A novel, computationally simple, and widely applicable method for calculating the smoothing factor is proposed. This method effectively handles the dislocation gradient problem between irregular triangular dislocation sub-faults, avoiding unnecessary oscillations in the slip distribution in the model. A grid search algorithm is implemented to optimize the smoothing factor and the optimal dip angle of the fault, determining the optimal fault geometry and achieving accurate inversion of the fault slip distribution. Attached Figure Description

[0040] Figure 1 A technical roadmap for a fault slip model smoothing method based on SAR observation data according to the present invention;

[0041] Figure 2 Map showing the tectonic background and surface fracture identification results of the Maduo earthquake;

[0042] Figure 3 Optimal tilt angle and smoothing factor search results chart;

[0043] Figure 4 , three 3D optimal slip model diagram;

[0044] Figure 5 Two-dimensional planar projection of the sliding model. Detailed Implementation

[0045] The specific implementation of the present invention is described below with reference to embodiments:

[0046] It should be noted that the structures, proportions, sizes, etc. shown in this specification are only used to complement the content disclosed in the specification for those skilled in the art to understand and read, and are not intended to limit the conditions under which the present invention can be implemented. Any modifications to the structure, changes in the proportions, or adjustments to the size, without affecting the effects and objectives that the present invention can produce, should still fall within the scope of the technical content disclosed in the present invention.

[0047] Furthermore, the terms such as "upper," "lower," "left," "right," "middle," and "one" used in this specification are merely for clarity of description and are not intended to limit the scope of the invention. Any changes or adjustments to their relative relationships, without substantially altering the technical content, should also be considered within the scope of the invention.

[0048] Key terms and abbreviations:

[0049] SAR: Synthetic Aperture Radar

[0050] GNSS: Global Navigation Satellite System

[0051] TIN: Triangulated Irregular Network

[0052] LOS: Line of sight

[0053] Example 1:

[0054] like Figure 1 As shown in the figure, this embodiment provides a method for constructing an irregular geometric fault model based on surface deformation observation data, which specifically includes:

[0055] The first step is to construct a fault geometry that fully reflects the complexity of surface rupture. ① First, geodetic data such as SAR, optical, and GNSS are processed to obtain a large-scale, high spatial resolution surface deformation field. Surface rupture trajectories are then visually interpreted based on the surface deformation field and remote sensing data. ② Next, to improve computational efficiency, the surface rupture trajectory is downsampled along the dip angle (the optimal dip angle of the fault will be determined in step 4), so that the number of sampling points decreases with increasing depth. ③ Finally, the TIN method is used to connect the sampling points to construct an irregular fault plane.

[0056] The second step is to define a Laplacian operator computational model applicable to triangular dislocations. Based on the geometric relationships between sub-faults and their slip amounts, a Laplacian computational matrix suitable for the triangular dislocation model is established.

[0057]

[0058] In the formula, ξ is the Laplacian operator; p0 is the slip on the central sub-fault patch0; p i The slip of the sub-fault adjacent to patch0; l i α is the distance between adjacent sub-faults; i This refers to the dihedral angle between adjacent sub-faults (see Figure b below for details).

[0059] Step 3 involves introducing the Laplacian operator and smoothing factor to construct a fault slip model. The process is as follows: ① First, based on elastic dislocation theory, establish the Green's function matrix G = [G...]. s G d ]:

[0060] G s =f(x,y,x0,y0,depth,strike,dip,length,width,ss ,0)

[0061] G d =f(x,y,x0,y0,depth,strike,dip,length,width,0,s d )

[0062] In the formula, G s , G d These are the surface displacement vectors under a unit slip disturbance in the strike and dip directions, respectively. The remaining parameters are the surface deformation observation location (x, y), fault location (x0, y0), depth, strike, dip, length, and width. s and s d Typically, it is a constant of 1.0. Furthermore, the InSAR deformation observation matrix D represents a one-dimensional displacement field in the LOS direction; before calculating the Green's function matrix, the LOS deformation needs to be projected onto the three-dimensional displacement field at the Earth's surface. ② Then, in the classical elastic dislocation theory model GX=D, a Laplacian operator and smoothing factor suitable for triangular dislocations are introduced, and the slip of each sub-fault is estimated based on least-squares criterion. ③ Finally, the constructed slip model is:

[0063] [D 0] T =[G ε 2 ξ] T p

[0064] In the formula, D is the surface deformation matrix, G is the Green's function, ε is the smoothing factor, ξ is the Laplacian operator, and p is the fault slip matrix.

[0065] Step 4 is to determine the optimal dip angle and smoothing factor of the fault. Although Step 1 constructed the fault geometry model, it only determined the connection mode and strike of the triangular dislocations on the fault plane, and could not determine the optimal dip angle. To determine the optimal fault parameters, a grid search method was introduced. Appropriate ranges and step sizes were set for the fault dip angle and smoothing factor, and steps 1 to 3 were repeated to test the fault distribution and data fitting error under different combinations of dip angles and smoothing factors. When increasing the slip roughness no longer significantly reduces the fitting error (e.g....), Figure 1 c) Determine the optimal values ​​for the tilt angle and smoothing factor. Since there is a one-to-one correspondence between slip roughness and the smoothing factor, they can be used as indicators for selecting the optimal values. The specific solution method is as follows:

[0066]

[0067] In the formula, n is the number of sub-fault segments in the fault model, r is the total roughness of the fault model, i.e., the average slip gradient, usually expressed in cm / km or m / km. iIt is a weighted slip distribution matrix with smoothing applied.

[0068] Example 2:

[0069] Take the 2021 Madoi earthquake as an example:

[0070] Figure 2 This is the tectonic setting of the 2021 Mw 7.4 Mado strike-slip earthquake event, and its surface rupture exhibits complex geometric structures such as bifurcation and bending. Figure 3 This paper presents the optimal dip angle and smoothing factor of a fault, obtained using a grid search method. To improve computational efficiency, sampling points are downsampled along the fault dip. Figure 4 The constructed fault model shows that the size of the irregular triangular dislocation fault gradually increases with depth, thus maximizing computational efficiency and ensuring the reliability of deformation inversion. Furthermore, by testing multiple sets of slip models estimated with fault dip angles and smoothing factors, the optimal value was determined when the data fitting error no longer decreased with increasing slip roughness. Figure 3 As can be seen, the optimal sliding factor coincides with the trade-off position of the data fitting error, and the optimal fault dip angle corresponds to the point where the fitting error is minimized. Figure 4 A fault slip model for seismogenic faults was established under the optimal dip angle and smoothing factor. The model shows that the seismic slip distribution is concentrated in the shallow part of the fault, and the slip between adjacent sub-faults exhibits good continuity, indicating that the method for determining the smoothing factor is highly feasible. Figure 5 The figure shows the relationship between the fault slip model projected onto a two-dimensional plane and the spatial distribution of aftershocks. As can be seen from the figure, aftershocks are mainly distributed at deeper points or along the sides of faults with less slip, exhibiting a strong spatial correlation with the slip distribution. This result further validates that the proposed method for determining the smoothing factor and fault geometric parameters has a good explanatory effect on fault slip in complex geometries and the mechanism of earthquake generation.

[0071] It will be understood by those skilled in the art that embodiments of the present invention may be provided as methods, systems, or computer program products. Thus, the present invention may take the form of an entirely hardware embodiment, an entirely software embodiment, or an embodiment combining software and hardware. Furthermore, the present invention may take the form of a computer program product implemented on one or more computer-usable storage media (including but not limited to magnetic disk storage, CD-ROM, optical storage, etc.) containing computer-usable program code.

[0072] The present invention is described with reference to flowcharts and / or block diagrams of methods, devices (systems), and computer program products according to embodiments of the present invention. It should be understood that each process and / or block in the flowcharts and / or block diagrams, as well as combinations of processes and / or blocks in the flowcharts and / or block diagrams, can be implemented by computer program instructions. These computer program instructions can be provided to a processor of a general-purpose computer, a special-purpose computer, an embedded processor, or other programmable data processing device to produce a machine, so that the instructions executed by the processor of the computer or other programmable data processing device generate instructions for implementing the processes in the flowcharts and / or block diagrams. Figure 1 a process or multiple processes and / or boxes Figure 1 A device that provides the functions specified in one or more boxes.

[0073] These computer program instructions may also be stored in a computer readable memory that can direct a computer or other programmable data processing device to work in a specific manner, so that the instructions stored in the computer readable memory produce an article of manufacture comprising an instruction device, which implements the process Figure 1 a process or multiple processes and / or boxes Figure 1 The function specified in one or more boxes.

[0074] These computer program instructions may also be loaded onto a computer or other programmable data processing equipment to cause a series of operational steps to be performed on the computer or other programmable equipment to produce a computer-implemented process, thereby providing instructions that execute on the computer or other programmable equipment for implementing the process. Figure 1 a process or multiple processes and / or boxes Figure 1 The steps of the function specified in one or more boxes.

[0075] The preferred embodiments of the present invention are described in detail above, but the present invention is not limited to the above embodiments. Various changes can be made within the knowledge of ordinary technicians in this field without departing from the scope of the present invention.

[0076] Many other changes and modifications can be made without departing from the concept and scope of this invention. It should be understood that this invention is not limited to the specific embodiments, and the scope of this invention is defined by the appended claims.

Claims

1. A method for smoothing fault slip models based on SAR observation data, characterized in that, The method includes: S1: Acquire surface deformation observation data and remote sensing-aided information, and jointly calculate to generate a high-resolution surface deformation field; combine the deformation field and remote sensing image visual interpretation to draw rupture traces; perform data downsampling along the dip direction; and use the TIN method to construct an irregular fault geometry based on the sampling points to obtain an initial fault plane model. S2: Based on the initial fault geometry model, node positions and triangular patch connections, determine the spatial relationships between sub-faults, construct the Laplacian operator and the smoothing matrix, and obtain the Laplacian operator matrix. S3: Input surface deformation observation data, initial fault geometric model and Laplacian operator matrix, construct Green's function matrix, introduce smoothing factor, perform slip inversion, and output initial slip distribution model, including the slip amount of each sub-fault segment; S4: Perform a grid search on the initial slip distribution model, the fault dip angle to be optimized, and the range and step size of the smoothing factor. Repeat the modeling process from S1 to S3 for different combinations of dip angles and smoothing factors. Then, evaluate the model and compare the fitting error and slip roughness under different combinations. When increasing the roughness no longer significantly reduces the error, determine the optimal parameters. Finally, extract the optimal solution. Using slip roughness and fitting error as the standard, output the optimal dip angle and smoothing factor, which yields the optimal geometric fault model, containing triangular connection relationships, slip distribution, and complete dip angle information.

2. The method for smoothing a fault slip model based on SAR observation data according to claim 1, characterized in that, Step S1 includes: S101: Obtain Data Collect surface deformation observation data and remote sensing auxiliary information. The surface deformation observation data includes SAR, optical imagery, and GNSS; the remote sensing auxiliary information includes remote sensing imagery of fault rupture lines and existing seismic zone data. S102: Joint Solution By using data fusion and computation techniques, multi-source observation data are comprehensively processed to generate a high-resolution surface deformation field; S103: Drawing the fracture trace By combining the generated surface deformation field with remote sensing imagery, rupture traces are visually interpreted and mapped to show the rupture characteristics and locations on the surface. S104: Data downsampling Data downsampling is performed along the defined rupture trajectory, along the dip direction of the fault, and the number of sampling points is reduced in areas of greater depth to improve computational efficiency; S105: Constructing fault geometry Based on the downsampled sampling points, the irregular fault geometry is constructed using the TIN method, and the initial fault plane model is finally output, including node positions, triangular connection relationships and dip direction.

3. The method for smoothing a fault slip model based on SAR observation data according to claim 1, characterized in that, In step S2, a Laplacian computation matrix suitable for the triangular dislocation model is established: In the formula, ξ is the Laplacian operator; p0 is the slip on the central sub-fault patch0; p i The slip of the sub-fault adjacent to patch0; l i α is the distance between adjacent sub-faults; i It is the dihedral angle between adjacent sub-faults.

4. The method for smoothing a fault slip model based on SAR observation data according to claim 1, characterized in that, Step S3 includes: S301: Based on the theory of elastic dislocations, establish the Green's function matrix G = [G...] that connects earthquake source parameters and surface displacements. s G d ]: G s =f(x,y,x0,y0,depth,strike,dip,length,width,s s ,0) G d =f(x,y,x0,y0,depth,strike,dip,length,width,0,s d ) In the formula, G s , G d These are the surface displacement vectors under a unit slip disturbance in the strike and dip directions, respectively. The remaining parameters are the surface deformation observation location (x, y), fault location (x0, y0), depth, strike, dip, length, and width. s and s d The value is usually a constant of 1.

0. In addition, the InSAR deformation observation matrix D is a one-dimensional displacement field in the LOS direction. Before calculating the Green's function matrix, the LOS deformation needs to be projected onto the three-dimensional displacement field of the Earth's surface. S302: In the classical elastic dislocation theory model GX=D, a Laplacian operator and smoothing factor applicable to triangular dislocations are introduced, and the slip of each sub-fault is estimated based on the least squares criterion. S303: The constructed slip model is as follows: [D 0] T =[G e 2 [ξ] T p In the formula, D is the surface deformation matrix, G is the Green's function, ε is the smoothing factor, ξ is the Laplacian operator, and p is the fault slip matrix.

5. The method for smoothing a fault slip model based on SAR observation data according to claim 1, characterized in that, In step S4, considering the one-to-one correspondence between slip roughness and smoothness factor, it can be used as the selection index for the optimal value. The specific solution method is as follows: In the formula, n is the number of sub-fault segments in the fault model, r is the total roughness of the fault model, i.e., the average slip gradient, usually expressed in cm / km or m / km. i It is a weighted slip distribution matrix with smoothing applied.

Images