Wide-range landslide mass displacement early warning method based on GNSS and radar data fusion

Abstract

The invention discloses a large-range landslide mass displacement early warning method based on GNSS and radar data fusion, and the method comprises the steps: obtaining GNSS monitoring sequence data in a landslide region based on the position information of a reference station and a GNSS monitoring station in a GNSS monitoring system, and further obtaining the displacement related information of the landslide region; collecting data through a satellite InSAR, and generating deformation field data of a landslide area by comparing radar images at different time points; performing space-time alignment on the obtained data, and performing two types of data fusion S by adopting an adaptive Kalman filtering method after alignment; decomposing time series data of the fused landslide mass displacement by using wavelet transform, and extracting multi-scale features of a time domain and a frequency domain; the proposed features are combined with historical data of landslide mass displacement, a landslide mass displacement trend prediction model using an LSTM method added with a memory decline factor is input, and future trend prediction of landslide displacement is carried out; setting a multi-level threshold value, judging the change trend of the landslide mass displacement according to the displacement rate, the acceleration and the prediction result, and generating a multi-level early warning signal.

Description

Technical Field

[0001] This invention belongs to the field of landslide prediction and early warning technology, and in particular relates to a large-scale landslide displacement early warning method based on GNSS and radar data fusion. Background Technology

[0002] Landslides are a common geological hazard in mountainous and hilly areas. Their formation process is complex and influenced by various factors, including geological structure, rainfall, earthquakes, and human engineering activities. Once a landslide occurs, it often causes severe loss of life and property and environmental damage. Therefore, real-time monitoring and trend prediction of landslide displacement have significant scientific and social value.

[0003] Current landslide monitoring and prediction methods mainly include GNSS-based detection technology, which can obtain high-precision three-dimensional point displacement information by deploying reference stations and monitoring stations on the landslide body. The measurement accuracy is high and it can be applied to various complex terrain conditions, especially in areas with sparse vegetation cover. The second is satellite InSAR-based detection technology. Satellite InSAR can acquire surface deformation field information over a large area. By comparing SAR images from different time phases, the displacement of the surface in the line of sight direction can be estimated. This method has a wide spatial coverage and can acquire deformation data over a large area, providing planar deformation field information and compensating for the sparseness of GNSS points. However, it is affected by the atmosphere and changes in surface cover, resulting in lower measurement accuracy and a lack of point-level reliability compared to GNSS.

[0004] In summary, current technologies for landslide displacement prediction still suffer from several limitations: difficulty in integrating point and area data; insufficient time-frequency feature analysis, with existing prediction methods largely confined to the time domain, failing to identify potential multi-scale frequency characteristics and resulting in a lack of in-depth understanding of landslide deformation patterns; and limited prediction and early warning accuracy. Therefore, it is necessary to propose a landslide displacement prediction technology that integrates the advantages of GNSS and InSAR, combining time-frequency analysis with deep learning methods, to improve the comprehensiveness of monitoring, the accuracy of prediction, and the timeliness of early warning. Summary of the Invention

[0005] The purpose of this invention is to overcome the shortcomings of existing technologies and propose a large-scale landslide displacement early warning method based on the fusion of GNSS and radar data. This method involves spatiotemporal registration of GNSS point data and satellite InSAR deformation field data, employing an adaptive dual-weighted Kalman filter for data fusion, extracting time-frequency features using wavelet transform, and combining this with an improved LSTM deep learning model incorporating a memory decay factor for trend prediction. This enables accurate prediction and multi-level early warning of landslide displacement.

[0006] Technical Solution: To achieve the above objectives, this invention provides a method for early warning of large-scale landslide displacement based on GNSS and radar data fusion, comprising the following steps:

[0007] S1: Based on the location information of the base station and GNSS monitoring station in the GNSS monitoring system, obtain GNSS monitoring sequence data in the landslide area to further obtain displacement-related information in the landslide area;

[0008] S2: Data is collected via satellite InSAR, and deformation field data of the landslide area is generated by comparing radar images at different time points;

[0009] S3: Perform spatiotemporal registration of the GNSS displacement information obtained in step S1 and the radar deformation field data obtained in step S2, and fuse them using an adaptive dual-weight Kalman filter method;

[0010] S4: The fused data is decomposed into time series data of landslide displacement using wavelet transform to extract multi-scale features in the time and frequency domains.

[0011] S5: Based on the time-frequency features extracted in S4 and combined with historical landslide displacement data, a landslide displacement trend prediction model using the LSTM method with added memory decay factor is constructed to predict the future trend of landslide displacement.

[0012] S6: Set multiple threshold levels, determine the trend of landslide displacement based on displacement rate, acceleration and prediction results, and generate multiple warning signals.

[0013] Furthermore, the specific steps in step S1 for obtaining GNSS monitoring sequence data in the landslide area to further obtain displacement-related information in the landslide area are as follows:

[0014] A1: Convert the coordinates of the GNSS monitored data from latitude, longitude, and elevation to ENU coordinates. Select a reference point (usually the coordinates of a GNSS station) as the origin of the ENU coordinate system. Assuming the latitude and longitude coordinates of the reference point are (), for the latitude and longitude of a point, use the following formula to convert latitude, longitude, and elevation to ECEF coordinates:

[0015] X=(N(φ)+h)cosφcos λ (1)

[0016] Y=(N(φ)+h)cosφsin λ (2)

[0017] Z=((1-e 2 )N(φ)+h)sinφ (3)

[0018] Where φ is latitude, λ is longitude, and h is elevation; is the Earth's radius of curvature; n is the Earth's semi-major axis; e is the Earth's eccentricity; use the same formula to calculate the ECEF coordinates X0, Y0, Z0 of the reference point;

[0019] A2: Calculate the difference between the target point and the reference point, and transform the difference from the ECEF coordinate system to the ENU coordinate system using a rotation matrix. The rotation matrix of the ENU coordinate system is based on the latitude and longitude of the reference point. The coordinates of the target point in the ECEF coordinate system are (X, Y, Z), and the ECEF coordinates of the reference point are X0, Y0, Z0. Then the local ENU coordinates are calculated as follows:

[0020]

[0021] The rotation matrix R(φ0, λ0) is used to convert ECEF coordinates to ENU coordinates, as shown in the following formula:

[0022]

[0023] The target point's east, north, and altitude in the ENU coordinate system can be obtained;

[0024] A3: After obtaining the ENU coordinates, we can perform displacement calculations. Assuming multiple observations across consecutive time steps, taking k1 and k2 as examples, the ENU coordinates at each time step are: E(k1), N(k1), U(k1) correspond to the ENU coordinates at time k1, and E(k2), N(k2), U(k2) correspond to the ENU coordinates at time k2. Δd GNSS,k The displacement difference relative to the previous moment is calculated using the following formula:

[0025]

[0026] Furthermore, the method for generating the deformation field of the landslide area using satellite InSAR in step S2 is as follows:

[0027] B1: Radar imagery acquired via satellite platforms typically contains phase information of radar waves, reflecting surface deformation at different times. In each SAR image, the phase information φ is related to the radar wave propagation path R of the target, and can be expressed as:

[0028]

[0029] Where φ is the phase, R is the path length of the radar wave propagation, and λ is the radar wavelength;

[0030] B2: Using differential interferometry, the interferogram between two SAR images is calculated. Differential interferometry can reveal minute surface deformations, such as the phase change φ at a given time t. ps The relationship between (k) and displacement d(k):

[0031]

[0032] The displacement difference is calculated using the phase difference, according to the following formula:

[0033] Δφ=φ2-φ1 (9)

[0034] Where Δφ is the differential phase of the two images, and φ2 and φ1 are the phases of the two images at times k2 and k1, respectively. According to formula (8), the relationship between the phase difference and the displacement difference can be obtained:

[0035]

[0036] The displacement difference Δd can be calculated from this formula (10):

[0037]

[0038] B3: After calculating the phase difference at multiple time points and obtaining the corresponding displacement difference, a projection transformation is performed to generate the deformation field of the surface displacement. In order to obtain deformation data consistent with the actual landslide direction, the displacement data needs to be projected from the radar line-of-sight direction to the downslope direction. The calculation formula is as follows:

[0039] Δd downslop =Δd·cosθ (12)

[0040] Where Δd is the displacement difference calculated from the phase difference, Δd downslop It is the transformed displacement, and θ is the angle between the LOS direction and the direction perpendicular to the ground.

[0041] Further, in step S3, the GNSS data and the obtained deformation field data are spatiotemporally registered, and an adaptive dual-weight Kalman filter method is used for data fusion. The specific method is as follows:

[0042] C1: Because GNSS data and InSAR data have different time intervals, GNSS data is typically acquired at a higher frequency (e.g., hourly or hourly), while InSAR data is acquired at a lower frequency (e.g., every few days). Therefore, to ensure that the data from both are consistent in time, interpolation processing is required for the InSAR data. Here, the GNSS data is denoted as d. GNSS (k) = [E(k), N(k), U(k)], while InSAR data is represented as Δd downslop (k), where k is the time point. For interpolation, we interpolate the InSAR data based on the time point k of the GNSS data. Assume the InSAR data at time points k1 and k2 is denoted as Δd. downslop (k1), Δddownslop (k2), to calculate the InSAR data corresponding to time point k, we can use the following interpolation formula:

[0043]

[0044] C2: To ensure spatial consistency between GNSS and radar data, they need to be mapped to the same coordinate system, spatially aligning GNSS and InSAR data. Since InSAR data is measured based on the satellite's line-of-sight (LOS), and we want to compare it with GNSS data in the same coordinate system, the displacement in the LOS direction needs to be converted to displacement in the ENU coordinate system. Let the angular relationship between the satellite and the target point be p. Then, a rotation matrix can be used to convert the LOS displacement to displacement in the ENU coordinate system. Assume the LOS displacement in the InSAR data is Δd. LOS It can be converted into displacement in the ENU coordinate system using the following formula:

[0045]

[0046] Where R(β) is the rotation matrix, which converts the LOS displacement into three components in the ENU coordinate system, namely the displacements in the east, north, and sky directions.

[0047]

[0048] GNSS data is already represented in the ENU coordinate system, and the transformed InSAR data is also represented in the ENU coordinate system, so their spatial positions are already aligned;

[0049] C3: After spatiotemporal alignment, to obtain the displacement data of GNSS data and satellite InSAR data at the same temporal and spatial location at each time step, an adaptive dual-weight Kalman filter method is used to fuse the GNSS data and satellite InSAR data. In the time dimension, a time weighting coefficient w is introduced. t The proportion of GNSS and InSAR data in the Kalman update is dynamically adjusted according to different observation time intervals; in the spatial dimension, a spatial weighting coefficient w is introduced. s The weights are dynamically allocated based on the complexity of the terrain. In areas with broken or heavily obstructed mountains, the weight of GNSS points is increased, while in areas with relatively flat terrain and good visibility, the weight of InSAR networks is increased.

[0050] Using GNSS and satellite InSAR data as measurement inputs, and calculating the Kalman gain using their measurement noise covariance, the state estimate is then updated based on their respective gain weights, as follows:

[0051] Update GNSS data:

[0052]

[0053] This represents the updated GNSS state estimate, and this represents the displacement estimate at time k. This represents the GNSS state estimate at the previous time k-1, i.e., the displacement estimate at the previous time. Since the displacement data we obtained before this time is the displacement difference with the previous time, to obtain the displacement estimate at a certain time, we need to give the initial reference position displacement as d0, as shown below:

[0054] in For the spatiotemporal integration weights of GNSS, The time weighting is applied; the smaller the observation interval, the greater the weight. For spatial weighting, the more complex the terrain, the greater the GNSS weight; K GNSS,k Δd represents the GNSS Kalman gain, indicating the degree of response to observation errors in the Kalman filter. A larger gain means a greater weight is given to the observation data during updates. GNSS,k H represents the displacement difference from the previous moment; GNSS It is the GNSS observation matrix, used to map state estimates to the observation space, and is generally in the form of an identity matrix;

[0055] Update InSAR data:

[0056]

[0057] This represents the updated InSAR state estimate, and this represents the displacement estimate at time k. This represents the InSAR state estimate at the previous time k-1, i.e., the displacement estimate at the previous time. The displacement calculation method for InSAR data at a certain time is the same as that for GNSS data described above; where... For the spatiotemporal integrated weights of InSAR, The time weighting is applied; the smaller the observation interval, the greater the weight. For spatial weights, the flatter the terrain, the greater the InSAR weight; here... θ represents the slope, ε represents the coherence of the InSAR data, typically ε∈(0,1); K InSAR,k Δd represents the InSAR Kalman gain, indicating the degree of response to observation errors in the Kalman filter. A larger gain means a greater weight is given to the observation data during updates. InSAR，k H represents the displacement difference from the previous moment; InSARThis represents the observation matrix, used to map the state estimate to the observation space; it is generally in the form of an identity matrix.

[0058] Finally, by combining the estimation results from both methods, the fused state estimate is obtained:

[0059]

[0060] Represents the fused state at time k, and α represents the displacement estimate after fusing GNSS and InSAR data; g and α i These are weighting factors for GNSS and InSAR data, respectively, controlling the contribution of each data source to the final estimate; This represents the updated GNSS state estimate. This represents the updated InSAR state estimate; Where R GNSS It is the GNSS observation noise covariance matrix, P k|k-1 The prediction error covariance matrix of the previous time step; Where R InSAR It is the InSAR observation noise covariance matrix; K InSAR，k Similarly;

[0061] Furthermore, the data fused in step S4 is decomposed into time-series landslide displacement data using wavelet transform to extract multi-scale features in the time and frequency domains. In the time-domain features, not only the mean and standard deviation are extracted, but the first derivative (velocity) and second derivative (acceleration) of the displacement sequence are also calculated to characterize the dynamic evolution of landslide displacement. The displacement rate reflects the speed of displacement change, and the displacement acceleration reflects the trend of displacement rate change. When the acceleration is consistently positive and gradually increases, it indicates that the landslide has entered an accelerated deformation stage, a precursor to a significantly increased probability of instability. Through joint analysis of the velocity and acceleration features with the original displacement features, the entire evolution of the landslide from stability to acceleration and finally to the point of imminent landslide can be more comprehensively revealed. The steps are as follows:

[0062] D1: Based on the landslide displacement time series, calculate the velocity and acceleration time series of the landslide. The formulas are as follows: displacement acceleration

[0063] D2: Morlet wavelets are signals with obvious frequency characteristics and are often used for frequency analysis of landslide data. Therefore, Morlet wavelets are used as the mother wavelet to decompose the signal. Wavelet transform is performed on the time series data x(t) of landslide displacement. Here, continuous wavelet transform (CWT) is selected to decompose the signal.

[0064] The Morlet wavelet is defined as:

[0065]

[0066] Where ψ(t) is the time-domain representation of the Morlet wavelet, describing the wavelet's shape in the time domain, and is used to generate the multi-scale transformation of the signal; w0 is the center frequency, the value of which depends on the frequency range of the signal we want to analyze.

[0067] The formula for Continuous Wavelet Transform (CWT):

[0068]

[0069] Where x(t) represents the time series data of the landslide displacement, i.e., the aforementioned... ψ * W(a,b) is the complex conjugate of the mother wavelet. 'a' is the scale factor, which determines the time-frequency resolution of the wavelet function. A larger 'a' corresponds to a wider wavelet (for analyzing low-frequency components), and a smaller 'a' corresponds to a narrower wavelet (for analyzing high-frequency components). 'b' is the translation factor, which controls the position of the wavelet function on the time axis, indicating that the time of the wavelet is unknown. This is used to analyze the local features of the signal at different time points. W(a,b) is the result of the wavelet transform, representing the characteristics of the signal at scale 'a' and time position 'b'.

[0070] D3: Through multiple wavelet transforms, the landslide displacement data can be ultimately decomposed into low-frequency and high-frequency components at multiple scales, respectively represented by C. a C d The low-frequency and high-frequency coefficients are used to represent the low-frequency components. The low-frequency components can be used to extract the time-domain features of the data, including the mean and standard deviation. The high-frequency components can help analyze the local changes in landslide displacement.

[0071] By performing multi-scale decomposition of landslide displacement data using CWT, frequency components of different frequency bands are obtained. By analyzing the frequency range of each decomposition layer and extracting the frequency characteristics of landslide displacement data at different scales, its frequency domain characteristics, instantaneous frequency, dominant frequency, and power spectral density can be obtained.

[0072] Furthermore, the landslide displacement trend prediction model in step S5 introduces a memory decay factor γ in the forget gate of the LSTM network. t It is used to dynamically adjust memory retention time and predict trends based on historical data and time-frequency characteristics;

[0073] E1: After cleaning and preprocessing the data fused in the previous step, and using the obtained time-frequency features, construct the input sequence, determine the input time window of the model, and divide the training data and target prediction data according to the fixed-length time window to construct the training set and test set;

[0074] E2: An improvement on the forget gate based on the LSTM network, introducing a memory decay factor γ. t This factor adaptively adjusts with time and the degree of data fluctuation. When the landslide displacement fluctuates drastically at a certain stage, γ t A larger value of γ accelerates the forgetting gate, thereby suppressing the interference of sudden abnormalities on long-term memory; when the landslide displacement is in a relatively stable state, γ... t With a smaller value, the forgetting gate extends memory retention time to capture long-term dependencies; here, γ... t It can be adaptively adjusted based on the sliding window standard deviation, which is as follows:

[0075]

[0076] Where W is the window size, x is the time series, μ is the mean of the data within the sliding window, and m is a small constant used to control the magnitude of the decay factor change;

[0077] E3: In the improved network structure, the formula for calculating the forget gate is updated as follows:

[0078] f t =σ(W f .[h t-1 x t ]+b f )·γ t (twenty three)

[0079] Among them, f t The output of the forget gate is σ, where σ is the sigmoid function and h is the output of the forget gate. t-1 x is the hidden state from the previous time step. t For the current input, W f and b f These represent the weights and bias parameters of the forget gate, γ. t An adaptive memory decay factor;

[0080] E4: After constructing the training and testing sets, the improved model is trained. After training, the improved LSTM network can adaptively adjust the "memory half-life" according to the different period characteristics of the landslide displacement data. It can respond quickly during violent fluctuations and maintain long-term dependence during stable phases, thus outputting displacement trend prediction results that are more consistent with reality. Through further analysis of the prediction results, more accurate displacement change trends can be provided for the landslide early warning system, thereby providing early warning of the occurrence of landslide disasters.

[0081] Furthermore, in step S6, different thresholds are set to generate multi-level early warning signals. These multi-level thresholds include not only displacement thresholds but also displacement rate thresholds and displacement acceleration thresholds: when the displacement is within a safe range but the rate or acceleration increases significantly, but has not yet reached a dangerous level, it is determined to be an early warning state; when the displacement approaches a dangerous level and the rate continues to increase, it is determined to be an enhanced early warning state; when both the displacement rate and acceleration exceed the thresholds, it indicates that the landslide has entered the pre-slip stage, is determined to be a warning state, and the highest level of early warning is issued. By introducing rate and acceleration indicators, this system can identify instability trends in advance before the cumulative displacement reaches the dangerous threshold, thereby achieving more sensitive and timely early warnings.

[0082] Beneficial effects: Compared with the prior art, the present invention has the following advantages:

[0083] 1. This invention overcomes the shortcomings of both by fusing high-precision GNSS location data and large-scale deformation field data from satellite InSAR. GNSS data can provide accurate location displacement information, while InSAR can cover a wide area and provide comprehensive deformation field information. By fusing the data through spatiotemporal registration and adaptive dual-weight Kalman filtering, the accuracy of landslide displacement prediction is improved, providing more comprehensive and reliable monitoring data for landslide monitoring.

[0084] 2. This invention employs wavelet transform to perform multi-scale time-frequency analysis on landslide displacement time-series data, effectively extracting multi-level features in both the time and frequency domains. Through wavelet transform, it can reveal periodic changes and local nonlinear deformation characteristics in landslide displacement, enhancing the predictive ability of landslide displacement trends. Compared to traditional methods that rely solely on time-domain data, this invention can more comprehensively capture the patterns of landslide deformation, thereby improving prediction accuracy.

[0085] 3. This invention employs LSTM (Long Short-Term Memory) network in predicting landslide displacement trends. This method can effectively process time series data and capture long-term dependencies in the data. In addition, by introducing a memory decay factor into the forget gate of the LSTM network, the model can automatically adjust the memory retention time according to the degree of data fluctuation, thereby more accurately identifying the changing trend of landslide displacement. Compared with traditional regression models, LSTM network can more accurately predict complex nonlinear time series data and avoid overfitting problems.

[0086] 4. By introducing dynamic characteristics such as displacement rate and acceleration, this invention can identify the stage of landslide evolution from slow creep to accelerated deformation, overcoming the shortcomings of traditional methods that rely solely on displacement thresholds. When the displacement rate and acceleration of the landslide show a continuous increasing trend, the system can issue an early warning even if the cumulative displacement has not yet exceeded the danger threshold, thereby improving the sensitivity and reliability of the warning. The multi-level early warning mechanism of this invention can provide different levels of warning when the landslide is at different stages, thereby providing decision-makers with more accurate and timely disaster prevention information, helping to prevent disasters and reduce their losses. Attached Figure Description

[0087] Figure 1 This is a flowchart of the method of the present invention;

[0088] Figure 2 This is a network structure diagram of the method of the present invention; Detailed Implementation

[0089] The present invention will be further illustrated below with reference to the accompanying drawings and specific embodiments. It should be understood that these embodiments are for illustrative purposes only and are not intended to limit the scope of the invention. After reading this invention, any modifications of the invention in various equivalent forms by those skilled in the art will fall within the scope defined by the appended claims.

[0090] S1: Based on the location information of the base station and GNSS monitoring station in the GNSS monitoring system, obtain GNSS monitoring sequence data in the landslide area to further obtain displacement-related information in the landslide area;

[0091] S2: Data is collected via satellite InSAR, and deformation field data of the landslide area is generated by comparing radar images at different time points;

[0092] S3: Perform spatiotemporal registration of the GNSS displacement information obtained in step S1 and the radar deformation field data obtained in step S2, and fuse them using an adaptive dual-weight Kalman filter method;

[0093] S4: The fused data is decomposed into time series data of landslide displacement using wavelet transform to extract multi-scale features in the time and frequency domains.

[0094] S5: Based on the time-frequency features extracted in S4 and combined with historical landslide displacement data, a landslide displacement trend prediction model using the LSTM method with added memory decay factor is constructed to predict the future trend of landslide displacement.

[0095] S6: Set multiple threshold levels, determine the trend of landslide displacement based on displacement rate, acceleration and prediction results, and generate multiple warning signals.

[0096] 2. The method for early warning of large-scale landslide displacement based on GNSS and radar data fusion according to claim 1, characterized in that the specific steps in step S1 of obtaining GNSS monitoring sequence data in the landslide area to further obtain relevant information on landslide area displacement are as follows:

[0097] A1: Convert the coordinates of the data monitored by GNSS from latitude, longitude, and elevation to ENU coordinates. Select a reference point (usually the coordinates of a GNSS station) as the origin of the ENU coordinate system. Assuming the latitude and longitude coordinates of the reference point are (φ0, λ0, h0), for the latitude and longitude (φ, λ, h) of a point, use the following formula to convert latitude, longitude, and elevation to ECEF coordinates:

[0098] X=(N(φ)+h)cosφcos λ (1)

[0099] Y=(N(φ)+h)cosφsin λ (2)

[0100] Z=((1-e 2 )N(φ)+h)sinφ (3)

[0101] Where φ is latitude, λ is longitude, and h is elevation; Here, n is the Earth's radius of curvature, e is the Earth's semi-major axis, and e is the Earth's eccentricity. The same formula is used to calculate the ECEF coordinates (X0, Y0, Z0) of the reference point.

[0102] A2: Calculate the difference between the target point and the reference point, and transform the difference from the ECEF coordinate system to the ENU coordinate system using a rotation matrix. The rotation matrix of the ENU coordinate system is based on the latitude and longitude of the reference point. The coordinates of the target point in the ECEF coordinate system are (X, Y, Z), and the ECEF coordinates of the reference point are X0, Y0, Z0. Then the local ENU coordinates are calculated as follows:

[0103]

[0104] The rotation matrix R(φ0, λ0) is used to convert ECEF coordinates to ENU coordinates, as shown in the following formula:

[0105]

[0106] The target point's east, north, and altitude in the ENU coordinate system can be obtained;

[0107] A3: After obtaining the ENU coordinates, we can perform displacement calculations. Assuming multiple observations across consecutive time steps, taking k1 and k2 as examples, the ENU coordinates at each time step are: E(k1), N(k1), U(k1) correspond to the ENU coordinates at time k1, and E(k2), N(k2), U(k2) correspond to the ENU coordinates at time k2. Δd GNSS,k The displacement difference relative to the previous moment is calculated using the following formula:

[0108]

[0109] 3. The method for early warning of large-scale landslide displacement based on GNSS and radar data fusion according to claim 1, characterized in that the method for generating the deformation field of the landslide area through satellite InSAR in step S2 is as follows:

[0110] B1: Radar imagery acquired via satellite platforms typically contains phase information of radar waves, reflecting surface deformation at different times. In each SAR image, the phase information φ is related to the radar wave propagation path R of the target, and can be expressed as:

[0111]

[0112] Where φ is the phase, R is the path length of the radar wave, and λ is the radar wavelength.

[0113] B2: Using differential interferometry, the interferogram between two SAR images is calculated. Differential interferometry can reveal minute surface deformations, such as the phase change φ at a given time t. ps The relationship between (k) and displacement d(k):

[0114]

[0115] The displacement difference is calculated using the phase difference, according to the following formula:

[0116] Δφ=φ2-φ1 (9)

[0117] Where Δφ is the differential phase of the two images, and φ2 and φ1 are the phases of the two images at times k2 and k1, respectively. According to formula (8), the relationship between the phase difference and the displacement difference can be obtained:

[0118]

[0119] The displacement difference Δd can be calculated from this formula (10):

[0120]

[0121] B3: After calculating the phase difference at multiple time points and obtaining the corresponding displacement difference, a projection transformation is performed to generate the deformation field of the surface displacement. In order to obtain deformation data consistent with the actual landslide direction, the displacement data needs to be projected from the radar line-of-sight direction to the downslope direction. The calculation formula is as follows:

[0122] Δd downslop =Δd·cos θ (12)

[0123] Where Δd is the displacement difference calculated from the phase difference, Δd downslop It is the transformed displacement, and θ is the angle between the LOS direction and the direction perpendicular to the ground.

[0124] 4. A method for early warning of large-scale landslide displacement based on GNSS and radar data fusion according to claim 1, characterized in that, in step S3, the GNSS data and the obtained deformation field data are spatiotemporally registered, and an adaptive dual-weight Kalman filter method is used for data fusion, the specific method being:

[0125] C1: Because GNSS data and InSAR data have different time intervals, GNSS data is typically acquired at a higher frequency (e.g., hourly or hourly), while InSAR data is acquired at a lower frequency (e.g., every few days). Therefore, to ensure that the data from both are consistent in time, interpolation processing is required for the InSAR data. Here, the GNSS data is denoted as d. GNSS (k) = [E(k), N(k), U(k)], while InSAR data is represented as Δd downslop (k), where k is the time point. For interpolation, we interpolate the InSAR data based on the time point k of the GNSS data. Assume the InSAR data at time points k1 and k2 is denoted as Δd. downslop (k1), Δd downslop (k2), to calculate the InSAR data corresponding to time point k, we can use the following interpolation formula:

[0126]

[0127] C2: To ensure spatial consistency between GNSS and radar data, they need to be mapped to the same coordinate system, spatially aligning GNSS and InSAR data. Since InSAR data is measured based on the satellite's line-of-sight (LOS), and we want to compare it with GNSS data in the same coordinate system, the displacement in the LOS direction needs to be converted to displacement in the ENU coordinate system. Let the angular relationship between the satellite and the target point be β. Then, a rotation matrix can be used to convert the LOS displacement to displacement in the ENU coordinate system. Assume the LOS displacement Δd in the InSAR data... LOS It can be converted into displacement in the ENU coordinate system using the following formula:

[0128]

[0129] Where R(β) is the rotation matrix, which converts the LOS displacement into three components in the ENU coordinate system, namely the displacements in the east, north, and sky directions.

[0130]

[0131] GNSS data is already represented in the ENU coordinate system, and the transformed InSAR data is also represented in the ENU coordinate system, so their spatial positions are already aligned;

[0132] C3: After spatiotemporal alignment, to obtain the displacement data of GNSS data and satellite InSAR data at the same temporal and spatial location at each time step, an adaptive dual-weight Kalman filter method is used to fuse the GNSS data and satellite InSAR data. In the time dimension, a time weighting coefficient w is introduced. t The proportion of GNSS and InSAR data in the Kalman update is dynamically adjusted according to different observation time intervals; in the spatial dimension, a spatial weighting coefficient w is introduced. s The weights are dynamically allocated based on the complexity of the terrain. In areas with broken or heavily obstructed mountains, the weight of GNSS points is increased, while in areas with relatively flat terrain and good visibility, the weight of InSAR networks is increased.

[0133] Using GNSS and satellite InSAR data as measurement inputs, and calculating the Kalman gain using their measurement noise covariance, the state estimate is then updated based on their respective gain weights, as follows:

[0134] Update GNSS data:

[0135]

[0136] This represents the updated GNSS state estimate, and this represents the displacement estimate at time k. This represents the GNSS state estimate at the previous time k-1, i.e., the displacement estimate at the previous time. Since the displacement data we obtained before this time is the displacement difference with the previous time, to obtain the displacement estimate at a certain time, we need to give the initial reference position displacement as d0, as shown below:

[0137] in For the spatiotemporal integration weights of GNSS, The time weighting is applied; the smaller the observation interval, the greater the weight. For spatial weighting, the more complex the terrain, the greater the GNSS weight; K GNSS,k Δd represents the GNSS Kalman gain, indicating the degree of response to observation errors in the Kalman filter. A larger gain means a greater weight is given to the observation data during updates. GNSS,k H represents the displacement difference from the previous moment; GNSS It is the GNSS observation matrix, used to map state estimates to the observation space, and is generally in the form of an identity matrix;

[0138] Update InSAR data:

[0139]

[0140] This represents the updated InSAR state estimate, and this represents the displacement estimate at time k. This represents the InSAR state estimate at the previous time k-1, i.e., the displacement estimate at the previous time. The displacement calculation method for InSAR data at a certain time is the same as that for GNSS data described above; where... For the spatiotemporal integrated weights of InSAR, The time weighting is applied; the smaller the observation interval, the greater the weight. For spatial weights, the flatter the terrain, the greater the InSAR weight; here... θ represents the slope, ε represents the coherence of the InSAR data, typically ε∈(0,1); K InSAR,k Δd represents the InSAR Kalman gain, indicating the degree of response to observation errors in the Kalman filter. A larger gain means a greater weight is given to the observation data during updates. InSAR,k H represents the displacement difference from the previous moment; InSAR This represents the observation matrix, used to map the state estimate to the observation space; it is generally in the form of an identity matrix.

[0141] Finally, by combining the estimation results from both methods, the fused state estimate is obtained:

[0142]

[0143] Represents the fused state at time k, and α represents the displacement estimate after fusing GNSS and InSAR data; g and α i These are weighting factors for GNSS and InSAR data, respectively, controlling the contribution of each data source to the final estimate; This represents the updated GNSS state estimate. This represents the updated InSAR state estimate; Where R GNSs It is the GNSS observation noise covariance matrix, P k|k-1 The prediction error covariance matrix of the previous time step; Where R InSAR It is the InSAR observation noise covariance matrix; K InSAR,k Similarly;

[0144] 5. A large-scale landslide displacement early warning method based on GNSS and radar data fusion according to claim 1, characterized in that, in step S4, the fused data is decomposed into time-series landslide displacement data using wavelet transform to extract multi-scale features in the time and frequency domains; in the time-domain features, not only the mean and standard deviation are extracted, but also the first derivative (velocity) and second derivative (acceleration) of the displacement sequence are calculated to characterize the dynamic evolution process of landslide displacement. The displacement rate reflects the speed of displacement change, and the displacement acceleration reflects the trend of displacement rate change. When the acceleration is consistently positive and gradually increases, it indicates that the landslide has entered an accelerated deformation stage, which is a precursor to a significantly increased probability of instability. Through joint analysis of the velocity and acceleration features with the original displacement features, the entire evolutionary pattern of the landslide from stability to acceleration to imminent landslide can be more comprehensively revealed. The steps are as follows:

[0145] D1: Based on the landslide displacement time series, calculate the velocity and acceleration time series of the landslide. The formulas are as follows: displacement acceleration

[0146] D2: Morlet wavelets are signals with obvious frequency characteristics and are often used for frequency analysis of landslide data. Therefore, Morlet wavelets are used as the mother wavelet to decompose the signal. Wavelet transform is performed on the time series data x(t) of landslide displacement. Here, continuous wavelet transform (CWT) is selected to decompose the signal.

[0147] The Morlet wavelet is defined as:

[0148]

[0149] Where ψ(t) is the time-domain representation of the Morlet wavelet, describing the wavelet's shape in the time domain, and is used to generate the multi-scale transformation of the signal; w0 is the center frequency, the value of which depends on the frequency range of the signal we want to analyze.

[0150] The formula for Continuous Wavelet Transform (CWT):

[0151]

[0152] Where x(t) represents the time series data of the landslide displacement, as mentioned above. ψ * W(a,b) is the complex conjugate of the mother wavelet. 'a' is the scale factor, which determines the time-frequency resolution of the wavelet function. A larger 'a' corresponds to a wider wavelet (for analyzing low-frequency components), and a smaller 'a' corresponds to a narrower wavelet (for analyzing high-frequency components). 'b' is the translation factor, which controls the position of the wavelet function on the time axis, indicating that the time of the wavelet is unknown. This is used to analyze the local features of the signal at different time points. W(a,b) is the result of the wavelet transform, representing the characteristics of the signal at scale 'a' and time position 'b'.

[0153] D3: Through multiple wavelet transforms, the landslide displacement data can be ultimately decomposed into low-frequency and high-frequency components at multiple scales, respectively represented by C. a C d The low-frequency coefficient and high-frequency coefficient are used to represent the low-frequency component. The low-frequency component can be used to extract the characteristics of the data in the time domain, including the mean and standard deviation. The high-frequency component can help analyze the local changes in landslide displacement.

[0154] By performing multi-scale decomposition of landslide displacement data using CWT, frequency components of different frequency bands are obtained. By analyzing the frequency range of each decomposition layer and extracting the frequency characteristics of landslide displacement data at different scales, its frequency domain characteristics, instantaneous frequency, dominant frequency, and power spectral density can be obtained.

[0155] 6. A method for early warning of large-scale landslide displacement based on GNSS and radar data fusion according to claim 1, characterized in that the landslide displacement trend prediction model in step S5 introduces a memory decay factor γ in the forget gate of the LSTM network. t It is used to dynamically adjust memory retention time and predict trends based on historical data and time-frequency characteristics;

[0156] E1: After cleaning and preprocessing the data fused in the previous step, and using the obtained time-frequency features, construct the input sequence, determine the input time window of the model, and divide the training data and target prediction data according to the fixed-length time window to construct the training set and test set;

[0157] E2: An improvement on the forget gate based on the LSTM network, introducing a memory decay factor γ. t This factor adaptively adjusts with time and the degree of data fluctuation. When the landslide displacement fluctuates drastically at a certain stage, γ t A larger value of γ accelerates the forgetting gate, thereby suppressing the interference of sudden abnormalities on long-term memory; when the landslide displacement is in a relatively stable state, γ... t With a smaller value, the forgetting gate extends memory retention time to capture long-term dependencies; here, γ... t It can be adaptively adjusted based on the sliding window standard deviation, which is as follows:

[0158]

[0159] Where W is the window size, x is the time series, μ is the mean of the data within the sliding window, and m is a small constant used to control the magnitude of the decay factor change;

[0160] E3: In the improved network structure, the formula for calculating the forget gate is updated as follows:

[0161] f t =σ(W f ·[h t-l x t ]+b f )·γ t (twenty three)

[0162] Among them, f t The output of the forget gate is σ, where σ is the sigmoid function and h is the output of the forget gate. t-1 x is the hidden state from the previous time step. t For the current input, W f and b f These represent the weights and bias parameters of the forget gate, γ. t An adaptive memory decay factor;

[0163] E4: After constructing the training and testing sets, the improved model is trained. After training, the improved LSTM network can adaptively adjust the "memory half-life" according to the different period characteristics of the landslide displacement data. It can respond quickly during violent fluctuations and maintain long-term dependence during stable phases, thus outputting displacement trend prediction results that are more consistent with reality. Through further analysis of the prediction results, more accurate displacement change trends can be provided for the landslide early warning system, thereby providing early warning of the occurrence of landslide disasters.

[0164] 7. A large-scale landslide displacement early warning method based on GNSS and radar data fusion according to claim 1, characterized in that, in step S6, different thresholds are set to generate multi-level early warning signals. The multi-level thresholds include not only displacement thresholds but also displacement rate thresholds and displacement acceleration thresholds: when the displacement is within a safe range but the rate or acceleration increases significantly, but has not yet reached a dangerous level, it is determined to be in a normal state; when the displacement is close to a dangerous level and the rate continues to increase, it is determined to be in a state of concern; when the displacement rate and acceleration both exceed the thresholds, it indicates that the landslide has entered the pre-slide stage and is determined to be in an alarm state. By introducing rate and acceleration indicators, this system can identify instability trends in advance before the cumulative displacement reaches the dangerous threshold, thereby achieving more sensitive and timely early warning.

Claims

1. A method for early warning of large-scale landslide displacement based on GNSS and radar data fusion, characterized in that, Includes the following steps: S1: Based on the location information of the base station and GNSS monitoring station in the GNSS monitoring system, obtain GNSS monitoring sequence data in the landslide area to further obtain displacement-related information in the landslide area; S2: Data is collected via satellite InSAR, and deformation field data of the landslide area is generated by comparing radar images at different time points; S3: Perform spatiotemporal registration of the GNSS displacement information obtained in step S1 and the radar deformation field data obtained in step S2, and fuse them using an adaptive dual-weight Kalman filter method; S4: The fused data is decomposed into time series data of landslide displacement using wavelet transform to extract multi-scale features in the time and frequency domains. S5: Based on the time-frequency features extracted in S4 and combined with historical landslide displacement data, a landslide displacement trend prediction model using the LSTM method with added memory decay factor is constructed to predict the future trend of landslide displacement. S6: Set multiple threshold levels, determine the trend of landslide displacement based on displacement rate, acceleration and prediction results, and generate multiple warning signals.

2. The method for early warning of large-scale landslide displacement based on GNSS and radar data fusion as described in claim 1, characterized in that, The specific steps in step S1 to obtain GNSS monitoring sequence data in the landslide area to further obtain displacement-related information of the landslide area are as follows: AI: Convert the coordinates of GNSS monitored data from latitude, longitude, and elevation to ENU coordinates. Select a reference point (usually the coordinates of a GNSS station) as the origin of the ENU coordinate system. Assuming the latitude and longitude coordinates of the reference point are (), for the latitude and longitude of a point, use the following formula to convert latitude, longitude, and elevation to ECEF coordinates: X=(N(φ)+h)cosφcosλ (1) Y=(N(φ)+h)cosφsinλ (2) Z=((1-e 2 )N(φ)+h)sinφ (3) Where φ is latitude, λ is longitude, and h is elevation; is the Earth's radius of curvature; n is the Earth's semi-major axis; e is the Earth's eccentricity; use the same formula to calculate the ECEF coordinates of the reference point X0, Y0, Z0; A2: Calculate the difference between the target point and the reference point, and transform the difference from the ECEF coordinate system to the ENU coordinate system using a rotation matrix. The rotation matrix of the ENU coordinate system is based on the latitude and longitude of the reference point. The coordinates of the target point in the ECEF coordinate system are (X, Y, Z), and the ECEF coordinates of the reference point are [missing information]. Then, the local ENU coordinates are calculated as follows: The rotation matrix R(φ0, λ0) is used to convert ECEF coordinates to ENU coordinates, as shown in the following formula: The target point's east, north, and altitude in the ENU coordinate system can be obtained; A3: After obtaining the ENU coordinates, we can perform displacement calculations. Assuming multiple observations across consecutive time steps, taking k1 and k2 as examples, the ENU coordinates at each time step are: E(k1), N(k1), U(k1) correspond to the ENU coordinates at time k1, and E(k2), N(k2), U(k2) correspond to the ENU coordinates at time k2. Δd GNSS，k The displacement difference relative to the previous moment is calculated using the following formula:

3. The method for early warning of large-scale landslide displacement based on GNSS and radar data fusion according to claim 1, characterized in that, The method for generating the landslide area deformation field using satellite InSAR in step S2 is as follows: B1: Radar imagery acquired via satellite platforms typically contains phase information of radar waves, reflecting surface deformation at different times. In each SAR image, the phase information φ is related to the radar wave propagation path R of the target, and can be expressed as: Where φ is the phase, R is the path length of the radar wave propagation, and λ is the radar wavelength; B2: Using differential interferometry, the interferogram between two SAR images is calculated. Differential interferometry can reveal minute surface deformations, such as the phase change φ at a given time t. ps The relationship between (k) and displacement d(k): The displacement difference is calculated using the phase difference, according to the following formula: Δφ=φ2-φ1 (9) Where Δφ is the differential phase of the two images, and φ2 and φ1 are the phases of the two images at times k2 and k1, respectively. According to formula (8), the relationship between the phase difference and the displacement difference can be obtained: The displacement difference Δd can be calculated from this formula (10): B3: After calculating the phase difference at multiple time points and obtaining the corresponding displacement difference, a projection transformation is performed to generate the deformation field of the surface displacement. In order to obtain deformation data consistent with the actual landslide direction, the displacement data needs to be projected from the radar line-of-sight direction to the downslope direction. The calculation formula is as follows: Δd downslop =Δd·cosθ (12) Where Δd is the displacement difference calculated from the phase difference, Δd downslop It is the transformed displacement, and θ is the angle between the LOS direction and the direction perpendicular to the ground.

4. The method for early warning of large-scale landslide displacement based on GNSS and radar data fusion according to claim 1, characterized in that, Step S3 involves spatiotemporal registration of the GNSS data and the obtained deformation field data, followed by data fusion using an adaptive dual-weight Kalman filter method. The specific method is as follows: C1: Because GNSS data and InSAR data have different time intervals, GNSS data is typically acquired at a higher frequency (e.g., hourly or hourly), while InSAR data is acquired at a lower frequency (e.g., every few days). Therefore, to ensure that the data from both are consistent in time, interpolation processing is required for the InSAR data. Here, the GNSS data is denoted as d. GNSS (k) = [E(k), N(k), U(k)], while InSAR data is represented as Δd downslop (k), where k is the time point. For interpolation, we interpolate the InSAR data based on the time point k of the GNSS data. Assume the InSAR data at time points k1 and k2 is denoted as Δd. downslop (k1), Δd downslop (k2), to calculate the InSAR data corresponding to time point k, we can use the following interpolation formula: C2: To ensure spatial consistency between GNSS and radar data, they need to be mapped to the same coordinate system, spatially aligning GNSS and InSAR data. Since InSAR data is measured based on the satellite's line-of-sight (LOS), and we want to compare it with GNSS data in the same coordinate system, the displacement in the LOS direction needs to be converted to displacement in the ENU coordinate system. Let the angular relationship between the satellite and the target point be β. Then, a rotation matrix can be used to convert the LOS displacement to displacement in the ENU coordinate system. Assume the LOS displacement Δd in the InSAR data... LOS It can be converted into displacement in the ENU coordinate system using the following formula: Where R(β) is the rotation matrix, which converts the LOS displacement into three components in the ENU coordinate system, namely the displacements in the east, north, and sky directions. GNSS data is already represented in the ENU coordinate system, and the transformed InSAR data is also represented in the ENU coordinate system, so their spatial positions are already aligned; C3: After spatiotemporal alignment, to obtain the displacement data of GNSS data and satellite InSAR data at the same temporal and spatial location at each time step, an adaptive dual-weight Kalman filter method is used to fuse the GNSS data and satellite InSAR data. In the time dimension, a time weighting coefficient w is introduced. t The proportion of GNSS and InSAR data in the Kalman update is dynamically adjusted according to different observation time intervals; in the spatial dimension, a spatial weighting coefficient w is introduced. s The weights are dynamically allocated based on the complexity of the terrain. In areas with broken or heavily obstructed mountains, the weight of GNSS points is increased, while in areas with relatively flat terrain and good visibility, the weight of InSAR networks is increased. Using GNSS and satellite InSAR data as measurement inputs, and calculating the Kalman gain using their measurement noise covariance, the state estimate is then updated based on their respective gain weights, as follows: Update GNSS data: This represents the updated GNSS state estimate, and this represents the displacement estimate at time k. This represents the GNSS state estimate at the previous time k-1, i.e., the displacement estimate at the previous time. Since the displacement data we obtained before this time is the displacement difference with the previous time, to obtain the displacement estimate at a certain time, we need to give the initial reference position displacement as d0, as shown below: in For the spatiotemporal integration weights of GNSS, The time weighting is applied; the smaller the observation interval, the greater the weight. For spatial weighting, the more complex the terrain, the greater the GNSS weight; K GNSS,k Δd represents the GNSS Kalman gain, indicating the degree of response to observation errors in the Kalman filter. A larger gain means a greater weight is given to the observation data during updates. GNSS，k H represents the displacement difference from the previous moment; GNSS It is the GNSS observation matrix, used to map state estimates to the observation space, and is generally in the form of an identity matrix; Update InSAR data: This represents the updated InSAR state estimate, and this represents the displacement estimate at time k. This represents the InSAR state estimate at the previous time k-1, i.e., the displacement estimate at the previous time. The displacement calculation method for InSAR data at a certain time is the same as that for GNSS data described above; where... For the spatiotemporal integrated weights of InSAR, The time weighting is applied; the smaller the observation interval, the greater the weight. For spatial weights, the flatter the terrain, the greater the InSAR weight; here... θ represents the slope, ε represents the coherence of the InSAR data, typically ε∈(0,1); K InSAR,k Δd represents the InSAR Kalman gain, indicating the degree of response to observation errors in the Kalman filter. A larger gain means a greater weight is given to the observation data during updates. InSAR,k H represents the displacement difference from the previous moment; InSAR This represents the observation matrix, used to map the state estimate to the observation space; it is generally in the form of an identity matrix. Finally, by combining the estimation results from both methods, the fused state estimate is obtained: Represents the fused state at time k, and α represents the displacement estimate after fusing GNSS and InSAR data; g and α i These are weighting factors for GNSS and InSAR data, respectively, controlling the contribution of each data source to the final estimate; This represents the updated GNSS state estimate. This represents the updated InSAR state estimate; Where R GNSS It is the GNSS observation noise covariance matrix, P k|k-1 The prediction error covariance matrix of the previous time step; Where R InSAR It is the InSAR observation noise covariance matrix; K InSAR,k Similarly.

5. The method for early warning of large-scale landslide displacement based on GNSS and radar data fusion according to claim 1, characterized in that, The data fused in step S4 is used to decompose the time-series data of landslide displacement using wavelet transform, extracting multi-scale features in the time and frequency domains. In the time domain features, not only are the mean and standard deviation extracted, but the first derivative (velocity) and second derivative (acceleration) of the displacement sequence are further calculated to characterize the dynamic evolution of landslide displacement. The displacement rate reflects the speed of displacement change, and the displacement acceleration reflects the trend of the displacement rate. When the acceleration is consistently positive and gradually increases, it indicates that the landslide has entered an accelerated deformation stage, a precursor to a significantly increased probability of instability. Through joint analysis of the velocity and acceleration features with the original displacement features, the entire process of landslide evolution from stability to acceleration and finally to the point of imminent landslide can be more comprehensively revealed. The steps are as follows: D1: Based on the landslide displacement time series, calculate the velocity and acceleration time series of the landslide. The formulas are as follows: displacement acceleration D2: Morlet wavelets are signals with obvious frequency characteristics and are often used for frequency analysis of landslide data. Therefore, Morlet wavelets are used as the mother wavelet to decompose the basic function of the signal. Wavelet transform is performed on the time series data x(t) of landslide displacement. Here, continuous wavelet transform (CWT) is selected to decompose the signal. The Morlet wavelet is defined as: Where ψ(t) is the time-domain representation of the Morlet wavelet, describing the wavelet's shape in the time domain, and is used to generate the multi-scale transformation of the signal; w0 is the center frequency, the value of which depends on the frequency range of the signal we want to analyze. The formula for Continuous Wavelet Transform (CWT): Where x(t) represents the time series data of the landslide displacement, as mentioned above. ψ * W(a,b) is the complex conjugate of the mother wavelet. 'a' is the scale factor, which determines the time-frequency resolution of the wavelet function. A larger 'a' corresponds to a wider wavelet (for analyzing low-frequency components), and a smaller 'a' corresponds to a narrower wavelet (for analyzing high-frequency components). 'b' is the translation factor, which controls the position of the wavelet function on the time axis, indicating that the time of the wavelet is unknown. This is used to analyze the local features of the signal at different time points. W(a,b) is the result of the wavelet transform, representing the characteristics of the signal at scale 'a' and time position 'b'. D3: Through multiple wavelet transforms, the landslide displacement data can be ultimately decomposed into low-frequency and high-frequency components at multiple scales, respectively represented by C. a C d The low-frequency coefficient and high-frequency coefficient are used to represent the low-frequency component. The low-frequency component can be used to extract the characteristics of the data in the time domain, such as the mean and standard deviation. The high-frequency component can help analyze the local changes in landslide displacement. By performing multi-scale decomposition of landslide displacement data using CWT, frequency components of different frequency bands are obtained. By analyzing the frequency range of each decomposition layer and extracting the frequency characteristics of landslide displacement data at different scales, its frequency domain characteristics, instantaneous frequency, dominant frequency, and power spectral density can be obtained.

6. The method for early warning of large-scale landslide displacement based on GNSS and radar data fusion according to claim 1, characterized in that, The landslide displacement trend prediction model in step S5 introduces a memory decay factor γ into the forget gate of the LSTM network. t It is used to dynamically adjust memory retention time and predict trends based on historical data and time-frequency characteristics; E1: After cleaning and preprocessing the data fused in the previous step, and using the obtained time-frequency features, construct the input sequence, determine the input time window of the model, and divide the training data and target prediction data according to the fixed-length time window to construct the training set and test set; E2: An improvement on the forget gate based on the LSTM network, introducing a memory decay factor γ. t This factor adaptively adjusts with time and the degree of data fluctuation. When the landslide displacement fluctuates drastically at a certain stage, γ t A larger value of γ accelerates the forgetting gate, thereby suppressing the interference of sudden abnormalities on long-term memory; when the landslide displacement is in a relatively stable state, γ... t With a smaller value, the forgetting gate extends memory retention time to capture long-term dependencies; here, γ... t It can be adaptively adjusted based on the sliding window standard deviation, which is as follows: Where W is the window size, x is the time series, μ is the mean of the data within the sliding window, and m is a small constant used to control the magnitude of the change in the decay factor; E3: In the improved network structure, the formula for calculating the forget gate is updated as follows: f t =σ(W f ·[h t-1 ,x t ]+b f )·c t (23) Among them, f t The output of the forget gate is σ, where σ is the sigmoid function and h is the output of the forget gate. t-1 x is the hidden state from the previous time step. t For the current input, W f and b f These represent the weights and bias parameters of the forget gate, γ. t An adaptive memory decay factor; E4: After constructing the training and testing sets, the improved model is trained. After training, the improved LSTM network can adaptively adjust the "memory half-life" according to the different period characteristics of the landslide displacement data. It can respond quickly during violent fluctuations and maintain long-term dependence during stable phases, thus outputting displacement trend prediction results that are more consistent with reality. Through further analysis of the prediction results, more accurate displacement change trends can be provided for the landslide early warning system, thereby providing early warning of the occurrence of landslide disasters.

7. The method for early warning of large-scale landslide displacement based on GNSS and radar data fusion according to claim 1, characterized in that, In step S6, different thresholds are set to generate multi-level early warning signals. These multi-level thresholds include not only displacement thresholds but also displacement rate thresholds and displacement acceleration thresholds: when the displacement is within a safe range but the rate or acceleration increases significantly, but has not yet reached a dangerous level, it is determined to be in a normal state; when the displacement is close to a dangerous level and the rate continues to increase, it is determined to be in a state of concern; when both the displacement rate and acceleration exceed the thresholds, it indicates that the landslide has entered the pre-slip stage and is determined to be in an alarm state. By introducing rate and acceleration indicators, this system can identify instability trends in advance before the cumulative displacement reaches the dangerous threshold, thereby achieving more sensitive and timely early warning.

Images