A multi-station Beidou displacement time sequence separation and decoupling method, device and medium based on independent component analysis

Abstract

The application relates to the technical field of radio orientation, and discloses a multi-station Beidou displacement time sequence separation and decoupling method, equipment and medium based on independent component analysis, which comprises the following steps: acquiring three-component daily coordinate displacement time sequence of a continuous operating Beidou reference station, eliminating abnormal gross errors and step signals through a standard deviation criterion and bidirectional search, and deducting a physical model correction value to generate a detrended residual time sequence; combining each detrended residual time sequence to construct a multi-station space-time combined observation matrix, and iteratively solving a potential independent source component matrix and a spatial weight coefficient matrix by using a negative entropy maximization criterion; identifying geophysical mechanisms according to the geological distribution form of a fault zone and the characteristics of an external physical environment load model, eliminating independent source components representing common mode errors of a propagation path, and reconstructing signals; and the application breaks the multi-source non-Gaussian signal aliasing bottleneck, and improves the decoupling precision and the solving stability.

Description

Technical Field

[0001] This invention relates to a method, device, and medium for multi-site BeiDou displacement time-series separation and decoupling based on independent component analysis, belonging to the field of radio direction finding technology. Background Technology

[0002] Currently, using a multi-base station network to obtain high-precision coordinate sequences for extracting tectonic deformation features is the mainstream approach. Conventional methods involve processing the observation time series by constructing spatial filters or using principal component analysis, which can deduct some systematic error terms from the known data matrix to achieve preliminary separation of characteristic signals. The displacement time series obtained by multi-station collaborative observation consists of the linear superposition of various physical source signals with different properties. Each physical signal exhibits differences in statistical distribution and spatial correlation. The blind source separation method precisely solves the unmixing matrix by mining the high-order statistical features contained in the observation matrix, under the premise of unknown channel transmission matrix and source signal characteristics, and using the maximization of the independence of unknown components as the criterion, thereby unmixing the multi-source mixed signal.

[0003] To improve decoupling accuracy, simple linear improvement approaches such as adding environmental sensors or adjusting the filter cutoff frequency are commonly used. However, adjusting the cutoff frequency can easily damage low-frequency geological displacement components, and adding sensors can cause interactive interference due to alignment conflicts of heterogeneous sources. This inherent limitation of performance trade-offs prevents traditional open-loop improvement paths from achieving high-purity physical decoupling of overlapping signals while reducing model dependence. The aforementioned improvement paths for overlapping signals are mostly limited to reconstructing the external pier morphology and other mechanical component morphology of monitoring network nodes or expanding the hardware of the control system. They do not fundamentally solve the aliasing interference of high-frequency non-Gaussian noise from the software algorithm dimension of displacement residual spatiotemporal coupling. As a result, the final data calculation is still constrained by the inherent limitations of the open-loop control strategy. For example, Chinese invention patent with publication number CN117807369A... The patent application discloses a GNSS time series denoising method and system based on variational Bayesian independent component analysis. It constructs a variational Bayesian network model for iterative noise learning. The core algorithm implicitly relies on the statistical stationarity of the overall spatiotemporal distribution of the global network. When facing active fault zone conditions, transient slip caused by weak earthquakes or nonlinear steps introduced by receiver hardware updates often exhibit strong local nonstationarity in the time and frequency domains. The variational Bayesian prior smoothing and training parameter iterative update mechanism, in such sudden change scenarios, because the physical action base point cannot offset the impact of nonstationary abrupt changes, is very likely to misjudge the sudden real crustal deformation details in the small band as high-order noise and smooth them out; or, during parameter iterative optimization, the non-ideal fluctuations of the signal cause the coefficients of the unmixing matrix to diverge, thus falling into the local extremum trap.

[0004] Therefore, in order to overcome the bottleneck of signal spatiotemporal decoupling, how to construct the spatiotemporal coupling matrix of the measurement network displacement residual, and achieve high-precision extraction of crustal deformation, periodic load and common mode error by using a high-order statistical blind source separation mechanism without relying on complex external prior models, has become the technical problem to be solved by this invention. Summary of the Invention

[0005] To address the problems in the background art, the technical solution of this invention is as follows: A multi-site BeiDou displacement time series separation and decoupling method based on independent component analysis, comprising the following steps:

[0006] Step S101: Obtain the three-component daily coordinate displacement time series of the continuously running Beidou reference station. Use the standard deviation criterion to remove abnormal gross errors in the three-component daily coordinate displacement time series and use a linear interpolation algorithm to complete the missing data. Based on the step time of the station's spatial coordinate time series, bidirectionally retrieve the three-component daily coordinate displacement time series to correct the step signal. Subtract the solid tide theory correction value, polar tide theory correction value, and ocean tidal load correction value calculated based on the Earth's elastic deformation theory and the regional ocean tidal load Green's function to generate a detrended residual time series.

[0007] Step S102: Combine the detrended residual time series of continuously running BeiDou reference stations to construct a multi-station spatiotemporal joint observation matrix;

[0008] Step S103: Orthogonally whiten the multi-station spatiotemporal joint observation matrix, and combine the negative entropy maximization criterion to carry out fixed-point iterative solution to separate the potential independent source component matrix and the spatial weight coefficient matrix.

[0009] Step S104: Based on the geological distribution morphology of the fault zone, the consistency of spatial coefficients of all stations in the network, and the temporal evolution characteristics of the external physical environment load model, the potential independent source component matrix and the spatial weight coefficient matrix are matched to identify the geophysical mechanism of each independent source component. In the signal reconstruction, the independent source components that characterize the common mode error of the propagation path are removed from the multi-station spatiotemporal joint observation matrix to separate and obtain the reconstructed target geophysical deformation signal.

[0010] Preferably, the method includes the following steps: Step S201, according to the principles of spatial uniform distribution, geological representativeness, low temporal missing rate and spatial independence, select verification stations from each continuously operating Beidou reference station as independent verification nodes, and the independent verification nodes do not participate in the blind source decoupling modeling of the multi-station spatiotemporal joint observation matrix; Step S202, compare the detrended residual time series of the independent verification nodes with the existing velocity field data, and calculate the deformation reconstruction accuracy evaluation index covering the internal residual level and the external deformation rate estimation deviation.

[0011] Preferably, step S101 involves bidirectionally retrieving the three-component daily coordinate displacement time series based on the step time of the station's spatial coordinate time series to correct the step signal, including the following sub-steps: Step S1011, constructing a deformation characteristic basis function containing step terms from known geological seismic examples and equipment replacement step terms; Step S1012, fitting the detrended residual time series using the least squares method to calculate the fitted residual vector for each step time; Step S1013, using a bidirectional sliding window to cyclically scan the time axis in both the forward and reverse directions to retrieve and capture step abrupt change points in the fitted residual vector; Step S1014, performing step offset correction on the detrended residual time series based on the amplitude of the captured step abrupt change points to eliminate the interference of nonlinear local abrupt changes on the iterative solution of non-Gaussian independent sources.

[0012] Preferably, step S101, which deducts the solid tide theory correction, polar tide theory correction, and ocean tidal load correction, includes the following sub-steps: Step S1015, based on the Earth's elastic deformation theory and the regional ocean tidal load Green's function, calculates the solid tide theory correction, polar tide theory correction, and ocean tidal load correction for the spatial location of each continuously operating BeiDou reference station; Step S1016, deducts the calculated solid tide theory correction, polar tide theory correction, and ocean tidal load correction from the detrended residual time series to eliminate the interference of external periodic deformation fluctuations.

[0013] Preferably, step S103, which iteratively decomposes the multi-station spatiotemporal joint observation matrix to independently decouple non-Gaussian independent sources, includes the following sub-steps: Step S1031, using the negative entropy maximization criterion as the non-Gaussianity measurement objective function, orthogonally whitens the multi-station spatiotemporal joint observation matrix to achieve spatial decoupling; Step S1032, using a fixed-point iterative algorithm to update the separation matrix coefficients until the mutual information of the potential independent source components in the potential independent source component matrices of two adjacent iterations satisfies the convergence condition, and outputs the statistically independent potential independent source component matrices and the corresponding spatial weight coefficient matrices.

[0014] Preferably, in step S104, based on the geological distribution morphology of the fault zone, the consistency of spatial coefficients across the entire network of stations, and the temporal evolution characteristics of the external physical environmental load model, the potential independent source component matrix and the spatial weight coefficient matrix are matched to identify the geophysical mechanism of each independent source component. This includes the following sub-steps: Step S1041, in the one-dimensional time dimension, calculate the temporal correlation coefficient between each independent source component and the existing atmospheric load time series model, hydrological load time series model, and asymmetric seasonal load model to identify the environmental load components; Step S1042 In the two-dimensional spatial dimension, the consistency of the spatial coefficients of all stations in the network corresponding to each independent source component in the spatial weight coefficient matrix is ​​analyzed. If the spatial weights show a uniform distribution across the network, they are determined to be regional common mode error components. If the spatial weights are limited to the edge of the fault zone, they are determined to be crustal tectonic deformation components. In step S1043, in the three-dimensional geological mechanism dimension, the distribution pattern of stations with high spatial weights is matched with the geometrical spatial overlap of the known geological distribution pattern of the fault zone. The spatial trajectory fit is calculated, and the tectonic deformation components reflecting fault zone closure and post-earthquake relaxation are identified.

[0015] Preferably, in step S104, the process of removing independent source components representing propagation path common-mode errors from the multi-station spatiotemporal joint observation matrix during signal reconstruction, and separating to obtain the reconstructed target geophysical deformation signal, includes the following sub-steps: Step S1044, performing matrix multiplication on the spatial weight coefficients of the corresponding columns and the independent source components determined to be regional common-mode error components, and inverting and reconstructing the propagation path common-mode error matrix caused by the superposition of satellite orbit errors, satellite clock errors, and regional atmospheric delays; Step S1045, subtracting the propagation path common-mode error matrix and the environmental load component matrix from the multi-station spatiotemporal joint observation matrix, and inversely transforming to reconstruct the remaining tectonic deformation components, to obtain the reconstructed target geophysical deformation signal after removing path common-mode errors and environmental interference.

[0016] Preferably, after obtaining the reconstructed target geophysical deformation signal in step S104, the following follow-up main steps are included: Step S801, the sliding window time difference calculation of the obtained reconstructed target geophysical deformation signal is performed, the horizontal and vertical movement rates of the crust are extracted under linear least squares fitting, and the first time derivative of the movement rate is obtained to obtain the tectonic deformation acceleration index; Step S802, a tectonic deformation anomaly early warning mechanism is established. When the tectonic deformation acceleration index continuously exceeds three times the existing system reference noise variance envelope, a deformation anomaly early warning signal reflecting the accelerated accumulation of strain in the fault zone and the crustal steady-state imbalance is output. In step S101, continuously operating Beidou reference stations with an annual average time series data missing rate of less than 5% are selected as input sources, and the spatial distribution density of each continuously operating Beidou reference station is not less than one station per 10,000 km², and the total time series length of the three-component daily coordinate displacement time series is not less than 1826 days, so as to ensure the convergence stability of the multi-station spatiotemporal joint observation matrix during blind source separation iterative solution.

[0017] A multi-site BeiDou displacement time-series separation and decoupling device based on independent component analysis includes: a memory for storing computer program instructions; and a processor for executing the computer program instructions to implement the steps of a multi-site BeiDou displacement time-series separation and decoupling method based on independent component analysis.

[0018] A multi-site BeiDou displacement timing separation and decoupling medium based on independent component analysis is disclosed. The multi-site BeiDou displacement timing separation and decoupling medium based on independent component analysis stores a computer program. When the computer program is executed by a processor, it implements the steps of a multi-site BeiDou displacement timing separation and decoupling method based on independent component analysis.

[0019] Compared with the prior art, the beneficial effects of the present invention are:

[0020] 1. In the multi-site BeiDou displacement time series separation, by uniformly constructing the three-dimensional coordinate time series of multiple stations into a detrended normalized observation matrix, and using a fixed-point iterative algorithm based on maximizing negative entropy to optimize the unmixing matrix, the independent non-Gaussian components can be separated based on the statistical independence characteristics of each physical source signal. This multi-site spatiotemporal collaborative processing path utilizes the spatial correlation within the reference network to break the constraint that traditional single-site filtering cannot distinguish overlapping signals in the frequency domain, enabling the decoupling of tectonic deformation signals with similar periodic characteristics and seasonal environmental load components in terms of spatial distribution and evolution, effectively avoiding the signal aliasing problem when separating non-Gaussian sources using conventional second-order statistical methods.

[0021] 2. A three-dimensional identification system is constructed based on temporal evolution characteristics, spatial weight distribution, and geophysical mechanism matching. By comparing independent components with the activity distribution pattern of fault zones, the consistency of spatial coefficients of all stations in the network, and the prediction patterns of various external environmental load models, it can identify the physical mechanism of each separated component. This mechanism transforms the traditional black-box data filtering into a signal judgment process with a definite physical causal relationship. In this way, it can accurately remove the common mode error of propagation path caused by satellite orbit, satellite clock error, and regional atmospheric delay from the self-observation matrix during the reconstruction stage, ensuring the physical interpretability of the final reconstructed target signal.

[0022] 3. By employing the three-standard-deviation criterion to eliminate isolated outliers and using linear interpolation to complete the sequence, combined with a bidirectional ergonomic scanning mechanism of least-squares fitting residuals of deformed characteristic functions, the step signal is accurately corrected. At the same time, the theoretical correction values ​​of solid tides, polar tides, and ocean tide models are deducted, which can provide a clean detrended residual sequence with the system deformation term removed for subsequent blind source separation. This multi-task sequence full-link preprocessing mechanism is deeply coupled with the construction of a multi-station spatiotemporal joint observation matrix, which effectively avoids the mathematical interference of step abrupt changes and nonlinear local high-frequency noise on the iterative solution process of non-Gaussian independent sources, prevents the calculation convergence from getting trapped in local extrema, and ensures the solution stability of the system when processing large-scale reference networks in batches. Attached Figure Description

[0023] Figure 1 This is a flowchart of the multi-site temporal blind source decomposition and target deformation separation decoupling method of the present invention;

[0024] Figure 2 This is a state machine diagram of the dynamic step size control during the independent component analysis calculation process of this invention.

[0025] The objectives, features, and advantages of this invention will be further explained in conjunction with the embodiments and with reference to the accompanying drawings. Detailed Implementation

[0026] The technical solutions of the embodiments of this application will be clearly described below with reference to the accompanying drawings. Obviously, the described embodiments are only some embodiments of this application, not all embodiments. All other embodiments obtained by those skilled in the art based on the embodiments of this application are within the scope of protection of this application.

[0027] A multi-site BeiDou displacement time series separation and decoupling method based on independent component analysis, for use in the field of radio direction finding technology, includes the following steps:

[0028] Step S101: Obtain the three-component daily coordinate displacement time series of the continuously running Beidou reference station. Use the standard deviation criterion to remove abnormal gross errors in the three-component daily coordinate displacement time series and use a linear interpolation algorithm to complete the missing data. Based on the step time of the station's spatial coordinate time series, bidirectionally retrieve the three-component daily coordinate displacement time series to correct the step signal. Subtract the solid tide theory correction value, polar tide theory correction value, and ocean tidal load correction value calculated based on the Earth's elastic deformation theory and the regional ocean tidal load Green's function to generate a detrended residual time series.

[0029] Step S102: Combine the detrended residual time series of continuously running BeiDou reference stations to construct a multi-station spatiotemporal joint observation matrix;

[0030] Step S103: Orthogonally whiten the multi-station spatiotemporal joint observation matrix, and combine the negative entropy maximization criterion to carry out fixed-point iterative solution to separate the potential independent source component matrix and the spatial weight coefficient matrix.

[0031] Step S104: Based on the geological distribution morphology of the fault zone, the consistency of spatial coefficients of all stations in the network, and the temporal evolution characteristics of the external physical environment load model, the potential independent source component matrix and the spatial weight coefficient matrix are matched to identify the geophysical mechanism of each independent source component. In the signal reconstruction, the independent source components that characterize the common mode error of the propagation path are removed from the multi-station spatiotemporal joint observation matrix to separate and obtain the reconstructed target geophysical deformation signal.

[0032] Preferably, the method includes the following steps: Step S201, according to the principles of spatial uniform distribution, geological representativeness, low temporal missing rate and spatial independence, select verification stations from each continuously operating Beidou reference station as independent verification nodes, and the independent verification nodes do not participate in the blind source decoupling modeling of the multi-station spatiotemporal joint observation matrix; Step S202, compare the detrended residual time series of the independent verification nodes with the existing velocity field data, and calculate the deformation reconstruction accuracy evaluation index covering the internal residual level and the external deformation rate estimation deviation.

[0033] Preferably, step S101 involves bidirectionally retrieving the three-component daily coordinate displacement time series based on the step time of the station's spatial coordinate time series to correct the step signal, including the following sub-steps: Step S1011, constructing a deformation characteristic basis function containing step terms from known geological seismic examples and equipment replacement step terms; Step S1012, fitting the detrended residual time series using the least squares method to calculate the fitted residual vector for each step time; Step S1013, using a bidirectional sliding window to cyclically scan the time axis in both the forward and reverse directions to retrieve and capture step abrupt change points in the fitted residual vector; Step S1014: performing step offset correction on the detrended residual time series based on the amplitude of the captured step abrupt change points to eliminate the interference of nonlinear local abrupt changes on the iterative solution of non-Gaussian independent sources.

[0034] Preferably, step S101, which deducts the solid tide theory correction, polar tide theory correction, and ocean tidal load correction, includes the following sub-steps: Step S1015, based on the Earth's elastic deformation theory and the regional ocean tidal load Green's function, calculates the solid tide theory correction, polar tide theory correction, and ocean tidal load correction for the spatial location of each continuously operating BeiDou reference station; Step S1016, deducts the calculated solid tide theory correction, polar tide theory correction, and ocean tidal load correction from the detrended residual time series to eliminate the interference of external periodic deformation fluctuations.

[0035] Preferably, step S103, which iteratively decomposes the multi-station spatiotemporal joint observation matrix to independently decouple non-Gaussian independent sources, includes the following sub-steps: Step S1031, using the negative entropy maximization criterion as the non-Gaussianity measurement objective function, orthogonally whitens the multi-station spatiotemporal joint observation matrix to achieve spatial decoupling; Step S1032, using a fixed-point iterative algorithm to update the separation matrix coefficients until the mutual information of the potential independent source components in the potential independent source component matrices of two adjacent iterations satisfies the convergence condition, and outputs the statistically independent potential independent source component matrices and the corresponding spatial weight coefficient matrices.

[0036] Preferably, in step S104, based on the geological distribution morphology of the fault zone, the consistency of spatial coefficients across the entire network of stations, and the temporal evolution characteristics of the external physical environmental load model, the potential independent source component matrix and the spatial weight coefficient matrix are matched to identify the geophysical mechanism of each independent source component. This includes the following sub-steps: Step S1041, in the one-dimensional time dimension, calculate the temporal correlation coefficient between each independent source component and the existing atmospheric load time series model, hydrological load time series model, and asymmetric seasonal load model to identify the environmental load components; Step S1042 In the two-dimensional spatial dimension, the consistency of the spatial coefficients of all stations in the network corresponding to each independent source component in the spatial weight coefficient matrix is ​​analyzed. If the spatial weights show a uniform distribution across the network, they are determined to be regional common mode error components. If the spatial weights are limited to the edge of the fault zone, they are determined to be crustal tectonic deformation components. In step S1043, in the three-dimensional geological mechanism dimension, the distribution pattern of stations with high spatial weights is matched with the geometrical spatial overlap of the known geological distribution pattern of the fault zone. The spatial trajectory fit is calculated, and the tectonic deformation components reflecting fault zone closure and post-earthquake relaxation are identified.

[0037] Preferably, in step S104, the process of removing independent source components representing propagation path common-mode errors from the multi-station spatiotemporal joint observation matrix during signal reconstruction, and separating to obtain the reconstructed target geophysical deformation signal, includes the following sub-steps: Step S1044, performing matrix multiplication on the spatial weight coefficients of the corresponding columns and the independent source components determined to be regional common-mode error components, and inverting and reconstructing the propagation path common-mode error matrix caused by the superposition of satellite orbit errors, satellite clock errors, and regional atmospheric delays; Step S1045, subtracting the propagation path common-mode error matrix and the environmental load component matrix from the multi-station spatiotemporal joint observation matrix, and inversely transforming to reconstruct the remaining tectonic deformation components, to obtain the reconstructed target geophysical deformation signal after removing path common-mode errors and environmental interference.

[0038] Preferably, after obtaining the reconstructed target geophysical deformation signal in step S104, the following follow-up main steps are included: Step S801, the sliding window time difference calculation of the obtained reconstructed target geophysical deformation signal is performed, the horizontal and vertical movement rates of the crust are extracted under linear least squares fitting, and the first time derivative of the movement rate is obtained to obtain the tectonic deformation acceleration index; Step S802, a tectonic deformation anomaly early warning mechanism is established. When the tectonic deformation acceleration index continuously exceeds three times the existing system reference noise variance envelope, a deformation anomaly early warning signal reflecting the accelerated accumulation of strain in the fault zone and the crustal steady-state imbalance is output. In step S101, continuously operating Beidou reference stations with an annual average time series data missing rate of less than 5% are selected as input sources, and the spatial distribution density of each continuously operating Beidou reference station is not less than one station per 10,000 km², and the total time series length of the three-component daily coordinate displacement time series is not less than 1826 days, so as to ensure the convergence stability of the multi-station spatiotemporal joint observation matrix during blind source separation iterative solution.

[0039] A multi-site BeiDou displacement time-series separation and decoupling device based on independent component analysis includes: a memory for storing computer program instructions; and a processor for executing the computer program instructions to implement the steps of a multi-site BeiDou displacement time-series separation and decoupling method based on independent component analysis.

[0040] A multi-site BeiDou displacement timing separation and decoupling medium based on independent component analysis is disclosed. The multi-site BeiDou displacement timing separation and decoupling medium based on independent component analysis stores a computer program. When the computer program is executed by a processor, it implements the steps of a multi-site BeiDou displacement timing separation and decoupling method based on independent component analysis.

[0041] Example 1: When a radio directional monitoring system faces the scenario of continuous BeiDou displacement time series analysis of multiple observation stations in an active fault zone, the daily coordinate displacement time series of the three components (eastward, northward, and zenith) acquired are mixed with common-mode errors caused by continuous crustal deformation, satellite orbit errors, satellite clock errors, and regional atmospheric delays, as well as seasonal periodic fluctuations caused by atmospheric pressure, surface hydrology, and non-tidal ocean loads. These non-Gaussian source signals with different properties overlap in the time and frequency domains, causing conventional model-driven fitting techniques to be unable to adapt to blind separation scenarios with unknown coupling interference due to their reliance on the accuracy of prior geophysical models. The time-frequency filtering process of a single station is prone to synchronously erroneously removing target structural deformation signals and losing effective information while removing high-frequency random noise, because it ignores the spatial coordination characteristics between the station nodes of the network.

[0042] The radio direction-of-care monitoring system selects continuously operating BeiDou reference stations with an annual average time-series data missing rate of less than 5% and a spatial distribution density of no less than one station per 10,000 square kilometers. It collects daily coordinate displacement time-series sequences of the eastward, northward, and zenith directions for a total length of no less than 1826 days. For isolated outliers in the daily coordinate displacement time-series sequences, their mean and standard deviation are calculated. When the absolute deviation of an observed value from the mean at a certain moment is greater than three times the standard deviation, the observed value is determined to be an outlier and is discarded. A linear interpolation algorithm is used to complete the missing data to ensure temporal continuity. A system containing the daily coordinate displacement time-series sequences of the three components is constructed. The deformed characteristic functions of the constant term, velocity term, annual cycle term, semi-annual cycle term, and seasonal cycle term are fitted using the least squares method to obtain the fitting residual vector after removing gross errors. The fitting residual vector is scanned using a double traversal method of forward and reverse time order. The step signals in the daily coordinate displacement time series of the three components are corrected based on the amplitude of the captured step abrupt change points to eliminate nonlinear local abrupt change interference. In the actual construction of the deformed characteristic function as the fitting basis, it is specifically composed of a linear combination of constants, multi-band harmonic components, and a unit step function superimposed. Specifically, the time independent variable of this deformed characteristic function is in days, and the linear part includes... The system includes a constant term and a velocity term multiplied by the time variable to fit long-term crustal tectonic movement trends; the periodic fluctuation component includes annual and semi-annual periodic terms, respectively composed of sine and cosine combinations with a base period of 365.25 days and sine and cosine combinations with a base period of 182.62 days, to capture steady-state periodic fluctuations caused by seasonal loads; the abrupt step component, based on the seismic occurrence time of previously recorded geological earthquakes and the equipment replacement time determined by changes in station measurement and receiver hardware, introduces a deformation characteristic basis function composed of multiple unit step functions with fixed time starting points. The step function takes a value of 0 before its corresponding step moment and a value of 1 at and after its corresponding step moment. Each step term is preceded by an amplitude coefficient to be solved. The displacement jump amplitude caused by each known interference source can be directly and accurately calculated through global least squares fitting. The solid tide theory correction value calculated by the elastic Earth tide model, the polar tide theory correction value calculated by the polar motion tide model, and the ocean tidal load correction value calculated by the global tide model are deducted from the three-component daily coordinate displacement time series to eliminate known systematic errors. The linear trend term and the annual and semi-annual periodic terms are removed by least squares fitting to generate a detrended residual time series.

[0043] The radio directional monitoring system combines the detrended residual time series of each station to construct an observation matrix consisting of the number of reference stations multiplied by the time series length dimension. Each row of the observation matrix corresponds to the residual sequence of one reference station, and each column corresponds to the observation values ​​of all reference stations at the same time. The observation matrix is ​​standardized with zero mean and unit variance to eliminate spatial differences in time series magnitude between different stations. The standardized observation matrix is ​​then orthogonally whitened to obtain a whitening matrix, which eliminates the linear correlation between signal components and reduces the computational complexity of iteration. Negative entropy is used as the criterion for non-Gaussianity and statistical independence, and maximizing the negative entropy of each component is used as the objective function for optimization. The separation matrix coefficients are dynamically updated through a fixed-point iterative algorithm until the absolute value of the change in the separation matrix coefficients between two adjacent iterations is less than 0.00001. The system then separates and outputs the potential independent source component matrix and the corresponding spatial weight coefficient matrix representing the spatial response weights of each independent source signal at different stations.

[0044] The radio directional monitoring system constructs a three-dimensional identification system encompassing temporal evolution characteristics, spatial weight distribution, and geophysical mechanism matching. In the one-dimensional time dimension, it calculates the temporal correlation coefficient between the time series of each independent source component and existing atmospheric, hydrological, and non-tidal ocean load time series models to identify environmental load components. In the two-dimensional spatial dimension, it analyzes the consistency of spatial coefficients across the entire network corresponding to each independent source component in the spatial weight coefficient matrix. When the spatial weights exhibit a uniform distribution across the entire network, the independent source components are identified as propagation path common-mode error components caused by the superposition of satellite orbital errors, satellite clock errors, and regional atmospheric delays. In practical processing, to accurately map the spatial weight coefficient matrix, which belongs to multidimensional continuous values... To label discrete physical mechanisms, this invention introduces spatial variation coefficient and spatial clustering degree as dimensionality reduction criteria. Specifically, for any independent source component, its corresponding spatial weight coefficient vector is extracted. The dimension of this vector is completely equivalent to the total number of reference stations participating in blind source decoupling in the measurement network. The average and standard deviation of the absolute values ​​of the weights of all stations in this weight coefficient vector are calculated, and the standard deviation is divided by the average to obtain the dimensionless spatial variation coefficient. When the spatial variation coefficient is less than or equal to 0.35, it indicates that the weight fluctuation of the independent component in geographic space is weak, satisfying the physical characteristics of uniform distribution throughout the network. Based on this, the system determines it as a regional common mode error component. Conversely, when the spatial variation coefficient is greater than 0.At 35 o'clock, the system further calculates the distance from each station to the nearest fault edge and extracts the sum of the absolute weights of all fault edge stations with a distance of less than 25 kilometers. This sum is then divided by the sum of the absolute weights of all stations in the entire network to derive the spatial clustering index. When the spatial clustering index is greater than or equal to 65%, it indicates that the energy of this component is highly confined to the fault edge, satisfying the physical characteristics of crustal tectonic deformation. Based on this, the system accurately maps and determines it as a crustal tectonic deformation component. When the distribution pattern of high-weight stations overlaps geometrically with the geological distribution pattern of active fault zones in the study area, the spatial trajectory fit is calculated in the three-dimensional geological mechanism dimension to delineate the tectonic deformation target components reflecting fault zone closure and post-earthquake relaxation. When matching in the three-dimensional geological mechanism dimension, the quantitative calculation of the spatial trajectory fit adopts a point-line geometric matching algorithm based on buffer layered projection. The specific execution steps are: digitizing the known active fault zone distribution pattern into... A spatial polyline segment, formed by connecting discrete geographic coordinates, is used to construct a tectonic influence buffer zone extending horizontally 25 kilometers to both sides of the polyline segment. Next, the spatial weight coefficient vectors of the corresponding independent source components are extracted, and the average and standard deviation of the absolute weights of all stations in the network are calculated. Stations with spatial weights greater than the average plus twice the standard deviation are marked as high-value response stations, forming a discrete point set. Then, the shortest geometric vertical distance from each high-value response station to the spatial polyline segment of the fault zone is calculated, and the number of high-value response stations whose shortest geometric vertical distance falls within the 25-kilometer buffer zone is counted. Finally, the number of stations falling within the buffer zone is divided by the total number of high-value response stations to obtain a percentage value as the spatial trajectory fit. When this fit percentage value is greater than or equal to 70%, the independent component is determined to be highly consistent with the geological fault zone in three-dimensional spatial distribution, thus identifying the tectonic deformation target component that accurately reflects fault zone closure and post-earthquake relaxation.

[0045] The radio direction finding system calculates the product of the spatial weight coefficients of the corresponding columns and the independent source components identified as regional common-mode error components to reconstruct the propagation path common-mode error matrix. This propagation path common-mode error matrix and the environmental load component matrix are subtracted from the self-observation matrix. An inverse transformation is then performed to reconstruct the remaining tectonic deformation components, and the constant and velocity terms removed in the preprocessing stage are superimposed. The reconstructed target geophysical deformation signal is then obtained separately. In the physical logic of this invention, the propagation path common-mode error matrix is ​​a spatial projection representation of the unmodeled residual errors in the range and phase domains onto the position domain. Specifically, satellite orbital residuals, satellite clock errors, and wide-area atmospheric delays exhibit high spatial coherence within short-range or regional measurement networks. When these range-domain errors are projected onto the station location domain through the time-varying observation geometric matrix of each station, they induce coordinated fluctuations with the same frequency, phase, and extensive spatial correlation in the three-dimensional coordinate displacement residual sequence. Since all reference stations in the network have highly similar sensitivity matrices to satellite errors in the same line-of-sight direction, independent component analysis utilizes this high-order statistical independence to capture the overall consistency of this geometric projection relationship through the spatial weighting coefficient matrix without relying on the original carrier phase in the range domain. Therefore, the common-mode error matrix reconstructed by inverting the matrix product of the spatial weighting coefficients and the corresponding independent source components mathematically and physically represents the projection of the propagation path physical error onto the three-dimensional coordinate space. The total comprehensive response can be effectively offset by directly performing a subtraction operation in the location domain, thus eliminating the nonlinear geometric deformation interference caused by the line-of-sight propagation path error on the location coordinate time series. Before constructing the observation matrix, based on the principles of spatial uniformity, geological representativeness, low temporal missing rate, and spatial independence, 15% to 25% of all continuously operating BeiDou reference stations, with a number of three or more, are selected as independent verification nodes. The independent verification nodes do not participate in the independent component analysis modeling process. The radio direction-of-sight monitoring system directly compares the detrended residual time series of the independent verification nodes with the publicly released velocity field data to calculate the estimation of internal residual level and external deformation rate. The dual-constraint reconstructed accuracy evaluation index of the bias is achieved by the radio directional monitoring system, which reduces the root mean square residual of the vertical displacement time series of BeiDou from 8.2 mm to 4.5 mm, a 45% reduction. At the same time, it reduces the system bias of the horizontal deformation rate estimation from 2.1 mm per year to 0.8 mm per year, a 62% reduction. The spatial correlation constraint and the high-order statistical blind source separation mechanism within the observation benchmark network work together to separate the frequency domain overlapping signals in the absence of an external prior geophysical model. The separation error between the deformation signal and aliasing interference is controlled through a data-driven processing flow, and the reconstructed target geophysical deformation signal with spatial correspondence and temporal evolution law is output.

[0046] Example 2: When the current radio direction finding system processes the coordinate deformation observation time series generated by the continuously operating BeiDou reference network in the active fault zone region, a comparative test system is established, consisting of multiple spatially overlapping interference gradients and parameter constraint boundaries. The original input data required by this test system comes from a regional deformation monitoring network consisting of 25 continuously operating BeiDou reference stations. The spatial distribution density of the monitoring network is one station per 10,000 square kilometers. The total data acquisition length is selected as 1826 days, and the data sampling frequency is fixed at daily acquisition. The test system actively superimposes Gaussian white noise with a signal-to-noise ratio of 20dB and power frequency interference harmonics with a frequency of 50Hz into the data source to simulate the non-ideal multi-source interference field in the real operating environment, which is composed of receiver multipath effect, radiation from peripheral radio equipment, and interweaving of higher-order ionospheric delay.

[0047] In establishing the aforementioned original input data specifications, the control logic for the upper limit of the missing data rate of the input station and the iterative convergence threshold of the unmixing matrix in the system parameters depends on the nonlinear constraint relationship between the geometric integrity of the observation network and the local load noise power. The setting of the upper limit of the missing data rate of the input station balances the completeness of the reference station for spatial common mode error extraction and the pseudo-periodic high-frequency noise power introduced by continuous data stream interpolation. The setting of the iterative convergence threshold balances the computational load of the independent source signal separation algorithm and the convergence accuracy of weak crustal deformation feature extraction. When the environmental load noise power in the area where the monitoring network is located increases or the reference stability of the station decreases, in order to ensure the identification accuracy of deformation signals within the wavelet band and avoid the risk of noise overfitting, the iterative convergence threshold moves towards the high-precision range. The system determines the upper limit of the missing data rate of the input station to 5% and the iterative convergence threshold of the unmixing matrix in independent component analysis calculation to 0.00001 through the judgment rules. This parameter configuration constitutes a closed calibration instance for the current wide-area high-density monitoring scenario of geological fault zones.

[0048] To verify the practical effect of the multi-feature collaborative separation mechanism and boundary constraints, the test system divided the collected detrended residual time series into four independent test groups. The test group of this invention includes complete bidirectional step mutation correction, three-dimensional geophysical mechanism identification, and fixed-point iterative algorithm based on maximizing negative entropy. The first control group uses the principal component analysis stack method instead of the independent component analysis algorithm as a reference for conventional second-order statistical separation methods. The second control group serves as a partially missing control group, and the step mutation correction process of forward and reverse time double traversal scanning is removed from the architecture of the test group of this invention. The third control group serves as an out-of-range control group, and the total length of input data collection is shortened to 365 days, while the missing rate of input station data is relaxed to 12%. Each test group uniformly faces the same multi-site BeiDou displacement time series containing Gaussian white noise and power frequency interference harmonics, thereby isolating and verifying the contribution of each process in the signal purification process.

[0049] In the parameter sensitivity optimization experiment spanning data acquisition duration, the test system observed the performance inflection point change law of nonlinear effects and deformation extraction accuracy. When the total acquisition length of the input data is less than 730 days, the system's separation resolution for seasonal environmental loads and path common mode errors decreases because the time span cannot cover more than two complete annual and semi-annual load fluctuations. This leads to an increasing trend in the deformation rate estimation deviation of the zenith displacement time series, indicating that the time domain length has not reached the minimum cutoff period limit. When the total acquisition length is extended from 1095 days to 1826 days, the reconstruction accuracy of the target geophysical deformation signal shows an upward trend, and the root mean square error in the vertical direction decreases accordingly. However, when the total acquisition length continues to increase and exceeds 1826 days or even reaches 3650 days, the data shows a saturation effect, and the rate of decrease in the root mean square error slows down and tends to flatten. This indicates that the signal gain brought by increasing the data length is close to the noise suppression limit of the algorithm itself, and further extending the sampling period cannot bring a performance improvement commensurate with the measurement benefits. This confirms that the 1826-day data length limit constitutes a direct experimental basis for determining the boundary of the numerical range.

[0050] The final calculation output data of each test group confirmed that the unprocessed first control group, due to its strong reliance on the Gaussian distribution assumption and inability to distinguish overlapping frequency domain signals, had residual load interference in its zenith displacement time series reconstruction results. This resulted in the root mean square error of the calculated vertical displacement time series residual remaining at 8.2 mm, and the systematic bias of the horizontal deformation rate estimation reaching 2.1 mm per year. The experimental group of this invention utilized the synergistic effect of spatial correlation constraints within the reference network and the separation of higher-order statistical blind sources to reduce the root mean square error of the vertical displacement time series residual to 4.5 mm, improving the reconstruction accuracy by 45%, while simultaneously reducing the systematic bias of the horizontal deformation rate estimation to 0. With a deformation rate of 8 mm per year, the system bias is reduced by 62%. In contrast, the second control group, which eliminated step correction, failed to eliminate the interference of nonlinear local mutations, and its deformation rate estimation bias rose back to 1.4 mm per year. The third control group, which exceeded the parameter working window, had its vertical residual root mean square error increase to 7.8 mm due to the limited sampling length. The above quantitative comparison results empirically demonstrate the nonlinear error hedging effectiveness generated by the combination of multi-site spatiotemporal joint observation matrix construction and negative entropy maximization blind source decoupling process. Under the condition of no external prior geophysical model, it achieves the isolation and independent extraction of crustal tectonic deformation signals, seasonal environmental loads and path common mode errors.

[0051] Example 3: This example combines Figures 1 to 2 This document describes a method, equipment, and medium for multi-site BeiDou displacement time-separation and decoupling based on independent component analysis, such as... Figure 1As shown, the multi-site BeiDou displacement time series is input to the standardization process through gross error automatic detection and completion, step signal identification and correction, and deterministic system error correction. After standardization, the standardized observation matrix is ​​output. The standardized observation matrix is ​​input to the multi-site time series blind source decomposition. The multi-site time series blind source decomposition outputs the spatial response weight coefficient and independent component time series sequence. The spatial response weight coefficient is output to the time series evolution characteristics, spatial weight distribution, and geophysical mechanism matching. The independent component time series sequence is output to the time series evolution characteristics, spatial weight distribution, and geophysical mechanism matching. The time series evolution characteristics, spatial weight distribution, and geophysical mechanism matching are input to the target signal reconstruction and clean time series extraction. The target signal reconstruction and clean time series extraction output the target deformation time series.

[0052] like Figure 2 As shown, the initial weight vector construction state is input to the variance calculation state by reading the multidimensional column vector. When the variance calculation state is greater than the variance tolerance threshold for three consecutive generations, it is input to the step factor downward adjustment state. The step factor downward adjustment state is input to the orthogonal projection superposition state by reducing the step factor by 50%. When the variance calculation state is less than the variance tolerance threshold for three consecutive generations, it is input to the step factor upward adjustment state. The step factor upward adjustment state is input to the orthogonal projection superposition state by increasing the step factor by 20%. When the absolute value of the variation is greater than 0.00001, the orthogonal projection superposition state is returned to the variance calculation state. When the deviation from the average value within the sliding time window exceeds three times the standard deviation, the orthogonal projection superposition state is input to the calibration loop recalibration state.

[0053] Example 4: When the radio directional monitoring system faces the distributed collaborative analysis of the coordinate displacement residual time series continuously generated by a high-density monitoring network composed of 25 continuously operating Beidou reference stations, the waveform characteristics of non-Gaussian implicit negative event records caused by local geological environment disturbances at the stations, seasonal groundwater extraction, and receiver antenna hardware updates are similar to the actual instantaneous tectonic deformation inside the Earth's crust. This causes the back-end processing unit to easily fall into local extremum traps when running a fixed-point iterative algorithm based on maximizing negative entropy to optimize the unmixing matrix. As a result, the separated potential independent source components will experience spatial phase confusion and misjudgment of physical properties. Due to the lack of construction details for surface-level program action control and parameter adaptive convergence criteria, the computation flow is prone to implementation obstacles such as iterative dead loops or weight coefficient divergence when facing the continuous impact of high-order non-Gaussian aliasing noise.

[0054] To address the challenges of computational flow getting stuck in an infinite loop and coefficient divergence, the system utilizes a distributed data bus architecture and a processing unit containing a double-precision floating-point storage array and a dedicated logic operation core. It introduces a deterministic procedure for implementing surface-level program action control and dynamic calibration of convergence thresholds for implicit local mutation events.

[0055] In the initial state definition phase, the processing unit reads the multi-dimensional column vectors segmented from the previous de-trending and normalized observation matrix from the preset address space of the cache. The core object of the multi-dimensional column vectors is the three-dimensional displacement residual data stream containing implicit local mutations at each station. The functional specifications of the enabled environment are set so that the floating-point operation precision of the dedicated logic operation core is not less than 0.000001, and the sampling frequency is maintained at once per day. To eliminate computational uncertainty, the processing unit executes the following peripheral control steps: the processor performs a linear combination product of the separation matrix coefficients of the current generation and the whitening matrix to construct the initial weight vector, uses a nonlinear mapping function to activate the column vectors under the action of the initial weight vector, and calculates the current generation and the previous generation. To address the variance of the weight vector across iterations, the system introduces a step-size calibration loop that includes a bidirectional switching mechanism for uplink and downlink adjustment. During each weight vector update, the variance generated in the previous step is compared in real-time with a preset tolerance threshold. If the variance exceeds the tolerance threshold for three consecutive iterations, it indicates the computational flow is in a divergent oscillation region. The processor triggers a downlink adjustment instruction, reducing the current step-size factor by 50% to suppress nonlinear jumps in the computational loop. Conversely, if the variance is less than the tolerance threshold for three consecutive iterations, it indicates the computational flow has entered a smooth approximation trajectory. The processor triggers an uplink adjustment instruction, increasing the current step-size factor by 20% to accelerate the approximation. The near-maximum negative entropy independent state, and the aforementioned selection of the number of consecutive iterations for the variable variance, as well as the proportional settings of decreasing the step size factor by 50% and increasing it by 20%, are determined based on the nonlinear convergence boundary test of the negative entropy objective function under strong non-Gaussian noise impact in independent component analysis. Specifically, when the number of iterations within the processor is less than 3, the system is easily mistriggered by small fluctuations in random noise, leading to frequent oscillations in the step size and an inability to accurately identify the true divergence trend. If the number of consecutive iterations is greater than 3, for example, reaching 5 or more, the system stays in the divergence region for too long, which can easily cause the weight coefficients of the unmixing matrix to diverge. Therefore, 3 iterations constitute a critical boundary that balances noise robustness and divergence prevention timeliness. Similarly, if the step size reduction ratio is less than 50%, its deceleration effect is insufficient to counteract the risk of gradient explosion caused by non-Gaussian pulses when facing strong local mutations such as hardware antenna replacement. If the reduction ratio exceeds 50%, it will cause the step size to become excessively blunted, causing the algorithm to stall prematurely and fall into a local extremum trap. On a smooth track, if the step size increase ratio is higher than 20%, it is easy for the computation flow to directly cross the peak of the maximum negative entropy of the global optimal solution and cause overshoot. If it is lower than 20%, it will cause the approximation speed to be too slow, resulting in redundant consumption of distributed computing resources. The processor will orthogonally project and superimpose the updated step size factor with the weight vector of the current generation and write it back to the weight register to complete the single iteration loop.

[0056] To verify the reproducibility of the above calibration procedure, under the condition that the multidimensional column vector contains three non-Gaussian implicit local mutation pulses caused by receiver antenna updates, the variance of the variation calculated by the processor in the 12th generation is: This value is greater than the preset tolerance threshold. After the dedicated logic operation core continuously scanned and maintained the out-of-state state in the 13th and 14th generations, it initiated a downlink adjustment action in the 15th generation, changing the current step size factor from... Adjust directly to This causes the variance of the 16th generation to converge to This avoids the deadlock trap. At the same time, to ensure the long-term operational stability of the dynamic control component, the system is configured with a sliding time window timeliness guarantee mechanism. The sliding time window length is set at 90 days to eliminate outdated historical deformation baselines. When the estimated value of the deformation residual in the latest cycle deviates from the average value within the sliding time window by more than 3 times the standard deviation, the recalibration of the calibration loop is automatically triggered.

[0057] By employing the aforementioned working mode encompassing dynamic bidirectional conversion of the step size factor and the timeliness guarantee mechanism for the historical deformation baseline via a sliding time window, the system eliminates coefficient divergence during the calculation process. When abrupt disturbances caused by non-Gaussian implicit negative events are injected again into the top-direction displacement residual data stream, the dedicated logic operation core stably controls the absolute value of the change in the separation matrix coefficients between adjacent iterations within a specified range. Within this range, the potential independent source component matrix is ​​separated, thereby stabilizing the reconstruction index of the vertical displacement time series within 4.5 mm, which confirms the ability of this method to solve specific environmental boundary constraints and the stability of engineering reproduction.

[0058] Example 5: In the deformation monitoring scenario of the BeiDou continuously operating reference station network in the active fault zone, when the transient displacement caused by weak seismic signals in the regional crust overlaps with the sudden coordinate jump caused by the update of the station receiver antenna hardware, the existing demixing process lacks a differentiated reconstruction mechanism for the implicit mutation event record, which leads to coefficient divergence and getting trapped in a local optimization trap during the calculation process.

[0059] The radio direction finding system deploys a data processing unit comprising a double-precision floating-point storage array and a dedicated logic operation core. It performs deterministic calibration and dynamic anomaly reconstruction on the raw displacement time-series data. The processor loads the daily displacement residual matrix of the three components of the reference station network into a dynamic memory region within a preset address space. The floating-point precision of the dedicated logic operation core is set to 0.000001, and the sampling period is synchronized once daily. The processor constructs a step-size control loop with uplink and downlink adjustment conversion mechanisms by comparing the variance of the separation matrix coefficients in real time. Specifically, when the variance exceeds the preset tolerance threshold of 0.001 for three consecutive iterations, the processor automatically reduces the current step-size factor by 50%, limiting the step jump amplitude of the computational weight vector. When the variance remains below the 0.001 tolerance threshold for three consecutive iterations... At this time, the processor increases the step size factor by 20% to accelerate convergence to the optimal solution space. For non-Gaussian step signals in the observation matrix caused by receiver hardware mutations, the processor performs bidirectional time-series scanning, marks segments with variance exceeding the threshold as implicit local mutation events, uses adjacent time window data for interpolation coverage, and uses this corrected displacement residual as the underlying reference data for constructing a multi-station spatiotemporal joint observation matrix, thereby achieving physical state stability during algorithm iteration. This not only eliminates the risk of the calculation process getting stuck in an infinite loop, but also locks the deformation baseline update cycle at 90 days through the sliding time window mechanism, maintaining the high-precision decoupling stability of the BeiDou reference station network under complex crustal deformation conditions. According to actual measurements, this deterministic procedure stably controls the root mean square error of the vertical residual within 4.5 mm, ensuring the reliability of observation data under implicit mutation interference scenarios.

[0060] It will be apparent to those skilled in the art that the present invention is not limited to the details of the exemplary embodiments described above, and that the present invention can be implemented in other specific forms without departing from the spirit or essential characteristics of the present invention.

[0061] Finally, it should be noted that the above embodiments are only used to illustrate the technical solutions of the present invention and are not intended to limit it. Although the present invention has been described in detail with reference to preferred embodiments, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical solutions of the present invention without departing from the spirit and scope of the technical solutions of the present invention.

Claims

1. A method for separating and decoupling the time series of BeiDou displacement data at multiple stations based on independent component analysis, characterized in that, Includes the following steps: Step S101: Obtain the three-component daily coordinate displacement time series of continuously operating BeiDou reference stations. Select continuously operating BeiDou reference stations with an annual average time series data missing rate of less than 5% as input sources. The spatial distribution density of each continuously operating BeiDou reference station is not less than one station per 10,000 km², and the total time series length of the three-component daily coordinate displacement time series is not less than 1826 days to ensure the convergence stability of the multi-station spatiotemporal joint observation matrix during blind source separation iterative solution. Use the standard deviation criterion to remove abnormal gross errors in the three-component daily coordinate displacement time series and use a linear interpolation algorithm to complete the missing data. Scan the three-component daily coordinate displacement time series in both forward and reverse order according to the step time of the station spatial coordinate time series to correct the step signal. Subtract the solid tide theory correction value and polar tide theory correction value calculated based on the Earth elastic deformation theory, as well as the ocean tidal load correction value calculated based on the regional ocean tidal load Green's function, to generate a detrended residual time series. Step S102: Combine the detrended residual time series of continuously running BeiDou reference stations to construct a multi-station spatiotemporal joint observation matrix; Step S103: Orthogonally whiten the multi-station spatiotemporal joint observation matrix, and combine the negative entropy maximization criterion to carry out fixed-point iterative solution to separate the potential independent source component matrix and the spatial weight coefficient matrix. Step S104: Based on the geological distribution morphology of the fault zone, the consistency of spatial coefficients of all stations in the network, and the temporal evolution characteristics of the external physical environment load model, the potential independent source component matrix and the spatial weight coefficient matrix are matched to identify the geophysical mechanism of each independent source component. In the signal reconstruction, the independent source components that characterize the common mode error of the propagation path are removed from the multi-station spatiotemporal joint observation matrix to separate and obtain the geophysical deformation signal of the reconstruction target. The process includes the following sub-steps: Step S1041, in the one-dimensional time dimension, calculate the temporal correlation coefficient between each independent source component and the existing atmospheric load time series model, hydrological load time series model, and asymmetric seasonal load model to identify environmental load components; Step S1042, in the two-dimensional spatial dimension, analyze the consistency of the spatial coefficients of all stations in the spatial weight coefficient matrix corresponding to each independent source component. If the spatial weights show a uniform distribution across the entire network, they are determined to be regional common mode error components. If the spatial weights are limited to the edge of the fault zone, they are determined to be crustal tectonic deformation components; Step S1043, in the three-dimensional geological mechanism dimension, match the distribution pattern of stations with high spatial weights with the geometrical spatial overlap of the known geological distribution pattern of the fault zone, calculate the spatial trajectory fit, and delineate the tectonic deformation components reflecting fault zone closure and post-earthquake relaxation.

2. The method for multi-site BeiDou displacement time series separation and decoupling based on independent component analysis according to claim 1, characterized in that, Includes the following steps: Step S201: Based on the principles of uniform spatial distribution, geological representativeness, low temporal missing rate, and spatial independence, verification stations are selected from each continuously operating BeiDou reference station as independent verification nodes. The independent verification nodes do not participate in the blind source decoupling modeling of the multi-station spatiotemporal joint observation matrix. Step S202: Compare the detrended residual time series of the independent verification nodes with the existing velocity field data, and calculate the deformation reconstruction accuracy evaluation index that covers the internal residual level and the external deformation rate estimation bias.

3. The method for multi-site BeiDou displacement time series separation and decoupling based on independent component analysis according to claim 1, characterized in that, Step S101 involves a double-cross scan of the three-component daily coordinate displacement time series based on the step time of the station's spatial coordinate time series to correct the step signal. This includes the following sub-steps: Step S1011: Constructing a deformation characteristic basis function that includes step terms from known geological earthquake examples and equipment replacement step terms; Step S1012: Fitting the detrended residual time series using the least squares method to calculate the fitted residual vector for each step time; Step S1013: Using a bidirectional sliding window to cyclically scan the time axis in both forward and reverse directions to retrieve and capture step abrupt change points in the fitted residual vector; Step S1014: Performing step offset correction on the detrended residual time series based on the amplitude of the captured step abrupt change points to eliminate the interference of nonlinear local abrupt changes on the iterative solution of non-Gaussian independent sources.

4. The method for multi-site BeiDou displacement time series separation and decoupling based on independent component analysis according to claim 1, characterized in that, Step S101 involves deducting the theoretical corrections for solid tides, polar tides, and ocean tidal loads, including the following sub-steps: Step S1015: Calculate the theoretical corrections for solid tides and polar tides for the corresponding spatial locations of each continuously operating BeiDou reference station based on the Earth's elastic deformation theory, as well as the ocean tidal load correction calculated based on the regional ocean tidal load Green's function; Step S1016: Deduct the calculated theoretical corrections for solid tides, polar tides, and ocean tidal loads from the detrended residual time series to eliminate interference from external periodic deformation fluctuations.

5. The method for multi-site BeiDou displacement time series separation and decoupling based on independent component analysis according to claim 1, characterized in that, Step S103 involves iteratively decomposing the multi-station spatiotemporal joint observation matrix to independently decouple non-Gaussian independent sources. This includes the following sub-steps: Step S1031: Using the negative entropy maximization criterion as the objective function for non-Gaussianity measurement, the multi-station spatiotemporal joint observation matrix is ​​orthogonally whitened to achieve spatial decoupling; Step S1032: The coefficients of the separation matrix are updated using a fixed-point iterative algorithm until the mutual information of the potential independent source components in the potential independent source component matrices of two adjacent iterations satisfies the convergence condition, and the statistically independent potential independent source component matrices and their corresponding spatial weight coefficient matrices are output.

6. The method for multi-site BeiDou displacement time series separation and decoupling based on independent component analysis according to claim 1, characterized in that, In step S104, the independent source components representing the common-mode error of the propagation path are removed from the multi-station spatiotemporal joint observation matrix during signal reconstruction, and the reconstructed target geophysical deformation signal is obtained. This includes the following sub-steps: Step S1044, performing matrix multiplication on the spatial weight coefficients of the corresponding columns and the independent source components determined to be regional common-mode error components, and reconstructing the propagation path common-mode error matrix caused by the superposition of satellite orbit error, satellite clock error and regional atmospheric delay; Step S1045, subtracting the propagation path common-mode error matrix and the environmental load component matrix from the multi-station spatiotemporal joint observation matrix, and reconstructing the remaining tectonic deformation components by inverse transformation, to obtain the reconstructed target geophysical deformation signal after removing the path common-mode error and environmental interference.

7. The method for multi-site BeiDou displacement time series separation and decoupling based on independent component analysis according to claim 1, characterized in that, After obtaining the reconstructed target geophysical deformation signal in step S104, the following sub-steps are also included: Step S801, the sliding window time difference calculation is performed on the obtained reconstructed target geophysical deformation signal, the horizontal and vertical movement rates of the crust are extracted under linear least squares fitting, and the first time derivative of the movement rate is obtained to obtain the tectonic deformation acceleration index; Step S802, a tectonic deformation anomaly early warning mechanism is established, and when the tectonic deformation acceleration index continues to exceed three times the system reference noise variance envelope, a deformation anomaly early warning signal reflecting the accelerated accumulation of strain in the fault zone and the crustal steady-state imbalance is output.

8. A multi-site BeiDou displacement time series separation and decoupling device based on independent component analysis, characterized in that, include: Memory is used to store computer program instructions; A processor for executing computer program instructions to implement the steps of the method as described in claim 1.

9. A multi-site BeiDou displacement time series separation and decoupling medium based on independent component analysis, characterized in that, The multi-site BeiDou displacement timing separation and decoupling medium based on independent component analysis stores a computer program, which, when executed by a processor, implements the steps of the method described in claim 1.

Images