Landslide deformation inversion and motion unit identification method based on dimensionless activity index

Abstract

The application discloses a landslide deformation inversion and motion unit identification method based on a dimensionless activity index, relates to the technical field of landslide terrain inversion, and solves the technical problem of how to realize fine inversion of local activity characteristics in a landslide under the condition of not performing phase unwrapping based on the wrapped phase; the application comprises the following steps: selecting SAR observation data covering a landslide research area, generating a time-series differential interferogram set, and then acquiring the wrapped interferometric phase of a time-series interference pair; preferably selecting an interferogram in the interferogram set, and constructing a global phase model of the landslide area deformation; adopting a sliding window strategy, traversing the landslide coverage area, and calculating a dimensionless activity index time series of the whole image pixel level; constructing a feature vector based on the dimensionless activity index time series, automatically segmenting the landslide area, identifying a plurality of motion units with independent motion characteristics, and revealing the spatiotemporal evolution law of the landslide through time-series analysis; and the application solves the problem of failure of traditional methods in unwrapping difficult areas.

Description

Technical Field

[0001] This invention relates to the field of landslide topographic inversion technology, specifically to a method for landslide deformation inversion and motion unit identification based on a dimensionless activity index. Background Technology

[0002] Landslides are sudden and extremely destructive, making early identification and monitoring crucial. Synthetic Aperture Radar Interferometry (InSAR) is widely used due to its all-weather and high-precision characteristics. However, during the pre-slide stage or accelerated deformation phase of a landslide, a large deformation gradient often exists on the surface, leading to severe decorrelation or phase aliasing in the interferograms. Traditional InSAR techniques (such as SBAS and PS-InSAR) heavily rely on phase unwrapping to obtain deformation values. If unwrapping fails, deformation information will be completely lost or severely biased.

[0003] In complex mountainous environments, the quality of interferometry is difficult to guarantee effectively, and even after phase unwrapping, unwrapping errors may still occur, affecting the accuracy of landslide time-series monitoring. To solve the unwrapping problem, some methods based on direct analysis of the wrapped phase have emerged in existing technologies. For example, by constructing a model phase map representing the overall deformation morphology of the landslide, and matching it with the interferogram in the time series, a global dimensionless activity index can be calculated.

[0004] However, existing activity index methods based on entanglement phases mainly achieve overall activity identification at the whole landslide scale. They are still difficult to characterize the local differential deformation process of different parts inside the landslide without unentanglement. Therefore, it is difficult to automatically identify the motion units with independent evolution laws inside the landslide, and it is also difficult to reveal the initiation sequence, temporal lag relationship and potential driving mechanism between different motion units. Summary of the Invention

[0005] Based on the problems existing in the prior art, this application further addresses the following more detailed technical issues:

[0006] 1. How to expand a single activity index of the entire landslide into a local scale, pixel-level activity intensity inversion result without phase unwrapping;

[0007] 2. How to automatically identify moving units with similar temporal evolution patterns within a landslide based on entangled phase time series, rather than just obtaining conclusions about the overall activity;

[0008] 3. How to further analyze the temporal correlation and time-delay propagation relationship between each unit based on the division of motion units, so as to reveal the spatiotemporal evolution mechanism inside the landslide.

[0009] Therefore, this application proposes a landslide deformation inversion and motion unit identification method based on a dimensionless activity index. The technical problem it solves compared to existing technologies is: how to achieve refined inversion of local activity characteristics within a landslide based on a tangled phase without phase unwrapping, and further automatically identify motion units with independent evolutionary laws, ultimately revealing the spatiotemporal evolution differences and time-delay response relationships between different parts of the landslide. This is mainly achieved through the following scheme:

[0010] A method for landslide deformation inversion and kinematic unit identification based on a dimensionless activity index includes:

[0011] S1: Select SAR observation data covering the landslide study area and extract the time-series wrapped interferometric phase atlas from it;

[0012] S2: Select several interferograms with clear fringes and high coherence quality from the time-series entangled interferometric phase image set, extract the spatial texture of landslide deformation by complex superposition and averaging, and construct a global phase model of landslide deformation.

[0013] S3: Using a sliding window strategy, the landslide-covered area is traversed to establish a local matching relationship between the observed phase of each pixel and the global phase model, and the dimensionless activity index time series at the pixel level of the whole map is calculated.

[0014] S4: Construct feature vectors based on dimensionless activity index time series, use machine learning clustering algorithms to automatically segment the landslide area, identify several motion units with independent motion characteristics, and extract the average activity index sequence of each motion unit. Through time series analysis, reveal the spatiotemporal evolution law of different parts of the landslide.

[0015] Further, S1 includes:

[0016] S11: Select SAR observation impact data covering the target landslide area, select an image with complete coverage and no obvious atmospheric noise as the main image, and use the remaining images as secondary images. After registering the main image and secondary images, select interferometric pairs with good spatiotemporal coherence.

[0017] S12: Combine the selected interferometric pairs, perform conjugate multiplication of the master and slave images to generate the original interferogram, remove the topographic phase and atmospheric phase, and retain the entangled phase containing deformation information;

[0018] S13: Filter the wrapped phase and then geocode it to convert it into a wrapped phase result in a geographic coordinate system.

[0019] Further, S2 includes:

[0020] S21: Based on the average coherence coefficient diagram and manual visual interpretation, select n interferograms as the sample set for constructing the global phase model of deformation;

[0021] S22: The sample set is superimposed and averaged in the complex domain to obtain the complex superposition result, which cancels out or suppresses random atmospheric noise and enhances the phase consistency of landslide deformation characteristics.

[0022] S23: Extract the phase angle of the complex superposition result to construct a global phase model of landslide deformation.

[0023] Further, S3 includes:

[0024] S31: For any pixel within the target landslide area, compare the observed phase of the pixel with the global phase model, and find an optimal scaling factor to maximize the matching degree between the two, which serves as the dimensionless activity index of the pixel; obtain the time series of the activity index of the pixel by traversing all interferograms.

[0025] S32: Standardize the activity index time series of each pixel to obtain the standardized feature vector of the pixel.

[0026] Further, S4 includes:

[0027] S41: Perform cluster analysis on the standardized feature vectors of all pixels to determine the optimal number of clusters K, group pixels with similar temporal evolution patterns into one class, and generate a cluster label distribution map; then perform morphological opening and closing operations on the cluster label distribution map to remove small fragmented patches, and finally automatically divide the landslide area into K landslide motion units with independent motion characteristics.

[0028] S42: Calculate the spatial arithmetic mean of the dimensionless activity index sequence of all pixels in the landslide movement unit to obtain the average activity index time series of the unit, which is used as the evolution time series of the landslide movement unit;

[0029] S43: Combining an external digital elevation model, use GIS software to analyze the landslide area, obtain the slope and aspect of each landslide movement unit as spatial location attributes; calculate the cross-correlation function between the time series of the average activity index of the landslide movement units, analyze the temporal relationship between different movement units, and determine whether there is synchronicity or time lag in the activity patterns and evolution laws; then extract the time delay corresponding to the maximum value of the cross-correlation coefficient of the cross-correlation function, and determine the spatiotemporal evolution law and driving mechanism of the landslide based on the time delay.

[0030] The beneficial effects of this invention include:

[0031] 1. This invention is based on landslide activity information extraction using entangled phases, without involving unentanglement operations. By establishing complex correlation matching between the observed phase and the model phase, it can still accurately extract landslide activity information under harsh conditions of low coherence and high gradient, solving the problem of traditional methods failing in areas where unentanglement is difficult.

[0032] 2. This invention introduces a motion unit zoning strategy for landslide morphology. By inverting the independent activity index sequences of different landslide regions, it can keenly capture the time difference (temporal lag relationship) between the initiation and acceleration of different parts of the landslide, thereby effectively determining the driving mechanism of the landslide and providing a deep physical basis for the analysis of the causes of landslide disasters and the formulation of treatment plans.

[0033] 3. This invention, through a deformation model of the entangled phase, introduces a dimensionless activity index that can intuitively transform complex interference fringe density into a numerical value reflecting the relative activity level of landslides. This index not only avoids the uncertainties of absolute deformation-level inversion but is also extremely sensitive to the acceleration precursors of landslides (fringe densification). By monitoring the temporal abrupt changes in this index, landslide instability trends can be detected earlier than traditional displacement time series, demonstrating significant early warning value. Attached Figure Description

[0034] Figure 1 This is a flowchart of the landslide deformation inversion and motion unit identification method based on the dimensionless activity index involved in the embodiments of this application.

[0035] Figure 2 This is a schematic diagram of the global deformation model of a landslide involved in an embodiment of this application.

[0036] Figure 3 This is a time series evolution diagram of the dimensionless activity index of different motion units of the landslide involved in the embodiments of this application. Detailed Implementation

[0037] To make the objectives, technical solutions, and advantages of the embodiments of this application clearer, the technical solutions of the embodiments of this application will be clearly and completely described below with reference to the accompanying drawings. Obviously, the described embodiments are only a part of the embodiments of this application, and not all of the embodiments. Therefore, the following detailed description of the embodiments of this application provided in the accompanying drawings is not intended to limit the scope of the claimed application, but merely represents selected embodiments of this application. All other embodiments obtained by those skilled in the art based on the embodiments of this application without creative effort are within the scope of protection of this application.

[0038] Example 1

[0039] The following is in conjunction with the appendix Figure 1 Specific embodiments of the present invention will be described in detail;

[0040] A method for landslide deformation inversion and kinematic unit identification based on a dimensionless activity index includes the following steps:

[0041] S1: Select Sentinel-1 SAR observation data covering the landslide study area, generate a time-series differential interferometric atlas based on the short baseline principle of 48-day time baseline and 150m spatial baseline threshold, and then use DInSAR technology to obtain the entangled interferometric phase of the time-series interferometric pair.

[0042] S2: Select several interferograms with clear fringes and high coherence quality from the interferogram set, extract the spatial texture of landslide deformation by complex superposition and averaging, and construct a global phase model of landslide deformation.

[0043] S3: Using a sliding window strategy, the landslide-covered area is traversed to establish a local matching relationship between the observed phase of each pixel and the global phase model, and the dimensionless activity index time series at the pixel level of the whole map is calculated.

[0044] S4: Construct feature vectors based on dimensionless activity index time series, use machine learning clustering algorithms to automatically segment the landslide area, identify several kinematic units with independent motion characteristics, and extract the average activity index sequence of each kinematic unit. Through time series analysis, reveal the spatiotemporal evolution law of different parts of the landslide.

[0045] In one embodiment, the specific steps of S1 include:

[0046] Step S11: Select Sentinel-1 SAR single-look complex imagery covering the target landslide area, spanning from May 2021 to May 2024, and introduce corresponding precise orbital ephemeris data to correct orbital errors. Select one image with complete coverage and no significant atmospheric noise as the master image, and use the remaining images as slave images. Using precise orbital data and DEM assistance, register all slave images to the geometric space of the master image with high precision, achieving sub-pixel accuracy. Subsequently, according to the SBAS principle, set a 48-day temporal baseline threshold and a 150-meter spatial vertical baseline threshold to select interferometric pair combinations with good spatiotemporal coherence.

[0047] S12: For the selected interferometric pairs, the master and slave images are conjugately multiplied to generate the original interferogram. An external high-precision digital elevation model (DEM) is introduced to simulate the terrain phase and subtract it from the original interferogram. Then, external GACOS data is used to remove the atmospheric phase to retain the entangled phase containing deformation information.

[0048] S13: The entangled phase obtained in S12 is filtered using an adaptive Goldstein filtering algorithm. Utilizing a range-Doppler model, the filtered entangled phase in the radar coordinate system is geocoded and converted into an entangled phase result in the geographic coordinate system. The filtering process includes setting appropriate filtering windows and filtering intensity parameters to reduce phase noise while maintaining the clarity of the landslide deformation stripe edges.

[0049] In one embodiment, the specific steps of S2 include:

[0050] S21: Based on the average coherence coefficient map and manual visual interpretation, 10 interferograms with clear, continuous deformation fringes and high coherence are selected from the entanglement phase results in the geographic coordinate system generated in S13 as a sample set for constructing the global deformation phase model.

[0051] When calculating the entangled interference phase map, the coherence coefficient map of each interference pair is also obtained simultaneously. The average coherence coefficient map can be obtained by averaging all the coherence coefficient maps.

[0052] S22: To avoid phase jump errors introduced by directly arithmetically averaging the entangled phases, the 10 interferograms selected in S21 are superimposed and averaged in the complex domain. Through complex superposition, random atmospheric noise is mutually canceled or suppressed, while the dominant landslide deformation characteristics are effectively enhanced due to phase consistency. The complex superposition formula is constructed as follows:

[0053]

[0054] In the formula, The result is a complex superposition, where j is the imaginary unit, x is the column coordinate of the range direction in the interferogram, and y is the row coordinate of the orientation direction in the interferogram. The observation phase of the nth interferogram, where obs represents the observation value and n represents the nth scene; This represents the observed interference phase at pixel coordinates (x, y) of the nth interferogram.

[0055] S23: Extract the phase angle from the complex superposition result obtained in S22 to construct a global phase model of deformation in the landslide area. This model represents the most representative spatial deformation morphology of the landslide during the observation period, specifically as follows: Figure 2 As shown, the color stripes used are Jet color stripes, and the mapping unit from blue to red is [-3.14, 3.14], which is the phase in radians. The change in color is the change in phase in radians.

[0056] In one embodiment, step S3 specifically includes:

[0057] S31: Define a W×W sliding window, where W=9. For any pixel (x,y) within the target landslide area, construct a local window centered on it. Assume that within this local window, the observed phase is a linear scaling of the global phase model. Compare the observed phase of each pixel with the global phase model, and find an optimal scaling factor to maximize the matching degree between the two. This scaling factor serves as a dimensionless activity index for each pixel, reflecting the changes of that pixel at different time points. Finally, an activity index time series is generated, providing data support for subsequent motion unit analysis and the revelation of spatiotemporal evolution laws.

[0058] For the t-th interferogram, the dimensionless activity index a(x,y,t) at that pixel is calculated, and the objective function is constructed as follows:

[0059]

[0060] In the formula, This represents the estimated value to be solved. Representative parameter optimization operator, Let p be the set of pixels within a sliding window centered at (x, y), where p represents the index of any pixel within the set. This represents the entangled observation phase value at pixel p in the t-th interferogram. This indicates the phase entanglement operator. This represents the model phase value at pixel p in the global phase model.

[0061] By traversing all temporal interferograms (t=1...M), the activity exponent time series vector of this pixel is obtained:

[0062]

[0063] In the formula, Let represent the dimensionless activity index at coordinates (x, y) of the t-th interferogram, and M represent the total number of time-series interferograms.

[0064] S32: Activity exponential time series vector of each pixel obtained from S31 To perform Z-score normalization, first calculate the arithmetic mean of the time series of that pixel. and standard deviation The original activity index sequence was then converted into a standardized feature vector using the standard deviation formula. This ensures that the standardized feature vectors of all processed pixels have zero mean and unit variance. The specific calculation expression is as follows:

[0065]

[0066] In the formula, The original dimensionless activity index, These are the standardized eigenvalues. It is a tiny constant.

[0067] In another embodiment, the specific steps of S4 are as follows:

[0068] S41: Utilize the K-Means clustering algorithm from machine learning to standardize the feature vectors of all pixels. Cluster analysis is performed by first determining the optimal number of clusters K based on the elbow rule or silhouette coefficient, grouping pixels with similar temporal evolution patterns into one class, and generating a cluster label distribution map. Then, morphological opening and closing operations are performed on the cluster label distribution map to remove fragmented patches, ultimately automatically segmenting the landslide area into K landslide motion units with independent motion characteristics.

[0069] S42: Based on the landslide motion units divided in S41, the evolution time sequence of each unit is extracted.

[0070] For each motion unit Calculate the original dimensionless activity index of all pixels within the cell range in S31. The spatial arithmetic mean is used to obtain the average activity index time series of this unit. The calculation formula is:

[0071]

[0072] In the formula, This represents the total number of pixels within the k-th motion unit.

[0073] S43: Using an external digital elevation model (DEM), GIS software was employed to analyze the landslide area. First, the DEM provided elevation information for the area. Based on this elevation data, the slope and aspect of each movement unit were calculated as spatial location attributes for each landslide movement unit, describing its geospatial characteristics. The average activity index time series for different landslide movement units was then calculated. The cross-correlation function between different landslide movement units is used to analyze the temporal relationship between them, i.e., whether their activity patterns and evolution laws exhibit synchronicity or temporal lag. The time delay corresponding to the maximum value of the cross-correlation coefficient is extracted. Based on the time delay, the spatiotemporal evolution law and driving mechanism of the landslide are determined. Specific results are as follows: Figure 3 As shown.

[0074] Cross-correlation functions help determine the interactions between two moving units by calculating the similarity of their activity index time series. For example, if the activity index time series of two moving units are highly correlated, it may mean that they influence each other or have similar deformation characteristics during a landslide. Spatial location attributes (such as elevation, slope, and aspect) describe the physical location of each moving unit in geographic space. Spatial attributes have a direct impact on landslide evolution because different topographic features affect the landslide occurrence mechanisms and evolution patterns in different regions. Therefore, by combining the results of spatial location attributes and cross-correlation functions, it is possible to determine whether the temporal evolution of different moving units is related to their spatial location. For example, moving units located at similar slopes and aspects may exhibit similar activity patterns, which can be confirmed by cross-correlation functions.

[0075] The cross-correlation coefficient is extracted from the cross-correlation function, representing the similarity between two time series at different time delays. The cross-correlation function calculates the relationship between the similarity between two time series and the time delay. The time delay between the two time series, i.e., their best temporal match, is obtained through the time point corresponding to the maximum cross-correlation coefficient. The larger the cross-correlation coefficient, the higher the similarity between the two time series at that time delay. The time point corresponding to the maximum cross-correlation coefficient represents the optimal temporal match point between the two motion units, i.e., the best moment when their activity patterns match each other temporally. The time delay is the time difference between the two motion units at this maximum cross-correlation coefficient moment. It reflects whether one motion unit lags behind or precedes the other in time. In this way, the relationship between the cross-correlation coefficient and the time delay indicates whether the activity patterns of different motion units are synchronized, or whether there is a temporal lag relationship between them. Based on the analysis of the temporal evolution trend of the deformation time delay of different motion units, the correlation of different spatial locations of the entire landslide is analyzed, thereby deducing the type of driven landslide.

[0076] Determining the spatiotemporal evolution and driving mechanism of landslides based on time delay is achieved by analyzing the time differences between the activity index time series of different moving units. The time delay reflects the activity sequence of each moving unit. Landslides typically exhibit two modes: a push-driven mode (the leading edge slides first, followed by the trailing edge) and a traction-driven mode (the trailing edge slides first, driving the leading edge to slide). For example, as... Figure 3As shown, the diagram illustrates five motion units divided into zones from the rear edge to the front edge of a landslide. The red line represents the rear edge region (Zone 1), while the blue line represents the front edge region (Zone 5). The diagram reveals that the activity index at the rear edge is earlier and larger than that at the front edge, indicating that the rear edge slides first, driving the front edge to slide—a traction-driven pattern. Based on the time delay of each motion unit, it is possible to determine which areas slide first and which slide later, thus revealing the landslide type and its sequence of activity. This approach allows for a clear understanding of the spatiotemporal evolution of landslides and the inference of their driving mechanisms, such as how external factors like gravity, precipitation, or earthquakes affect different parts of the landslide.

[0077] Because landslides have various driving modes, including shove landslides, traction landslides, and translational landslides, the deformation characteristics of different zones within the landslide will differ depending on the driving mode. For example, in shove landslides, the rear edge slides first, relying on gravity to push the central and front areas of the landslide, ultimately triggering a complete landslide. The deformation pattern shows that the rear edge slides before the front edge. Traction landslides are the opposite of shove landslides, where the front edge slides first, pulling the central and rear sections of the landslide body backward, gradually causing instability and ultimately forming a complete landslide. Translational landslides involve approximately uniform deformation time across the entire landslide, with the entire landslide sliding at roughly the same pace. Therefore, based on the methods used in this study, combined with dimensionless activity index monitoring of landslide zones, the entire landslide was divided into motion units (motion units 1-5) from the rear edge to the front edge. Dimensionless activity indices were then calculated for each unit, and their deformation trends over time were statistically analyzed. For example... Figure 3 It can be seen that the landslide was relatively stable before January 2023. In February 2023, all the moving units of the landslide showed an exponential increase at the same time, indicating that a small deformation had occurred in this part of the landslide. Until October 2023, the exponential changes of the moving units of the landslide showed a time difference. It can be seen that the rear edge started to increase the exponential earliest, followed by other moving units, and finally the moving units at the front edge of the landslide. Based on the sliding time delay of the landslide and the landslide driving mode theory, it can be inferred that the landslide is a push-type landslide, and the landslide began to show significant deformation in October 2023, which continued until May 2024.

[0078] The embodiments described above merely illustrate specific implementation methods of this application, and while the descriptions are detailed and specific, they should not be construed as limiting the scope of protection of this application. It should be noted that those skilled in the art can make various modifications and improvements without departing from the concept of the technical solution of this application, and these modifications and improvements all fall within the scope of protection of this application.

Claims

1. A method for landslide deformation inversion and motion unit identification based on a dimensionless activity index, characterized in that, Includes the following steps: S1: Select SAR observation data covering the landslide study area and extract the time-series wrapped interferometric phase atlas from it; S2: Select several interferograms with clear fringes and high coherence quality from the time-series entangled interferometric phase image set, extract the spatial texture of landslide deformation by complex superposition and averaging, and construct a global phase model of landslide deformation. S3: Using a sliding window strategy, the landslide-covered area is traversed to establish a local matching relationship between the observed phase of each pixel and the global phase model, and the dimensionless activity index time series at the pixel level of the whole map is calculated. S4: Construct feature vectors based on dimensionless activity index time series, use machine learning clustering algorithms to automatically segment the landslide area, identify several motion units with independent motion characteristics, and extract the average activity index sequence of each motion unit. Through time series analysis, reveal the spatiotemporal evolution law of different parts of the landslide.

2. The landslide deformation inversion and motion unit identification method based on dimensionless activity index according to claim 1, characterized in that, S1 includes: S11: Select SAR observation impact data covering the target landslide area, select an image with complete coverage and no obvious atmospheric noise as the main image, and use the remaining images as secondary images. After registering the main image and secondary images, select interferometric pairs with good spatiotemporal coherence. S12: Combine the selected interferometric pairs, perform conjugate multiplication of the master and slave images to generate the original interferogram, remove the topographic phase and atmospheric phase, and retain the entangled phase containing deformation information; S13: The entangled phase is filtered and then geocoded to convert it into a entangled phase result in a geographic coordinate system to obtain a time-series entangled interferometric phase atlas.

3. The landslide deformation inversion and motion unit identification method based on dimensionless activity index according to claim 1, characterized in that, S2 includes: S21: Based on the average coherence coefficient map and manual visual interpretation, n time-series winding interferometric phase maps are selected from the time-series winding interferometric phase map set as a sample set for constructing the global phase model of deformation; S22: The sample set is superimposed and averaged in the complex domain to obtain the complex superposition result, which cancels out or suppresses random atmospheric noise and enhances the phase consistency of landslide deformation characteristics. S23: Extract the phase angle of the complex superposition result to construct a global phase model of landslide deformation.

4. The landslide deformation inversion and motion unit identification method based on dimensionless activity index according to claim 1, characterized in that, S3 includes: S31: For any pixel within the target landslide area, compare the observed phase of the pixel with the global phase model, and find an optimal scaling factor to maximize the matching degree between the two, which serves as the dimensionless activity index of the pixel; obtain the time series of the activity index of the pixel by traversing all interferograms. S32: Standardize the activity index time series of each pixel to obtain the standardized feature vector of the pixel.

5. The landslide deformation inversion and motion unit identification method based on dimensionless activity index according to claim 1, characterized in that, S4 includes: S41: Perform cluster analysis on the standardized feature vectors of all pixels to determine the optimal number of clusters K, group pixels with similar temporal evolution patterns into one class, and generate a cluster label distribution map; then perform morphological opening and closing operations on the cluster label distribution map to remove small fragmented patches, and finally automatically divide the landslide area into K landslide motion units with independent motion characteristics. S42: Calculate the spatial arithmetic mean of the dimensionless activity index sequence of all pixels in the landslide movement unit to obtain the average activity index time series of the unit, which is used as the evolution time series of the landslide movement unit; S43: Combining an external digital elevation model, use GIS software to analyze the landslide area, obtain the slope and aspect of each landslide movement unit as spatial location attributes; calculate the cross-correlation function between the time series of the average activity index of the landslide movement units, analyze the temporal relationship between different movement units, and determine whether there is synchronicity or time lag in the activity patterns and evolution laws; then extract the time delay corresponding to the maximum value of the cross-correlation coefficient of the cross-correlation function, and determine the spatiotemporal evolution law and driving mechanism of the landslide based on the time delay.

Images