Abnormal positioning data cleaning method based on multi-base station positioning

Abstract

The application belongs to the technical field of data cleaning, and particularly relates to an abnormal positioning data cleaning method based on multi-base station positioning, comprising the following steps: acquiring original positioning data collected by multiple base stations within a preset time window and measurement noise standard deviation corresponding to each positioning attribute of a terminal to be positioned; calculating a time gradient sequence of each positioning attribute for continuous data points within the time window, and constructing a cross-attribute time sequence gradient correlation matrix; solving an inter-attribute conflict index based on the cross-attribute time sequence gradient correlation matrix, and extracting a dominant abnormal direction vector; and mapping the inter-attribute conflict index to a scalar abnormal amplitude through a probability model trained based on historical normal data. The application constructs a cross-attribute time sequence gradient matrix to identify positioning abnormalities, constructs an abnormal covariance, and corrects filtering parameters through a Jacobian matrix to adaptively suppress abnormal interference, thereby significantly improving positioning accuracy and system stability.

Description

Technical Field

[0001] This invention belongs to the field of data cleaning technology, specifically relating to a method for cleaning abnormal location data based on multi-base station positioning. Background Technology

[0002] Among numerous positioning technologies, multi-base station collaborative positioning methods are widely used in cellular network positioning, indoor positioning, and vehicle-to-everything (V2X) applications due to their advantages such as wide coverage and no need for additional terminal hardware modifications. This method calculates the terminal's location by comprehensively utilizing measurement information from multiple base stations, such as signal arrival time, time difference of arrival, and angle of arrival. In complex environments, the propagation path of wireless signals is easily interfered with by factors such as obstacles, reflection, and diffraction, resulting in non-line-of-sight propagation effects. This effect can cause deviations in the measurement values ​​of some positioning attributes, resulting in abnormal data points or outliers. If these are used directly for location calculation without processing, it will severely reduce the accuracy and stability of the positioning system, and may even cause drastic changes in the positioning trajectory.

[0003] Chinese invention patent CN114363949B discloses an abnormal data detection method for a UWB positioning system. The method involves: S1, deploying N UWB positioning base stations in the environment to be located. The positions of the base stations and reference points are calibrated using a total station. Distance is measured by interacting with the base stations using UWB positioning tags in the environment. The sampling frequency is set to 5Hz, and the distance from the tags to each base station is recorded. A deep learning RNN network model is used to train the detection model based on historical data and data from the next moment. Simultaneously, data prediction is performed based on real-time data and compared with the actual values, thereby achieving dynamic and real-time abnormal data detection.

[0004] Kalman filters and their variants are commonly used tools for fusing time-series data and smoothing positioning trajectories. Existing techniques often employ chi-square tests of residuals, random sampling consensus algorithms, or threshold decision methods to identify and eliminate potentially anomalous measurement data points. Thresholding or statistical tests are sensitive to environmental changes, difficult to set, and prone to missed detections or false positives. Completely eliminating identified data points is a crude approach, failing to consider that in a multi-attribute measurement, only some attributes may be anomalous while others remain reliable, resulting in information loss. Most existing methods treat each measurement in isolation or evaluate them as a whole, failing to deeply explore and utilize the inherent correlation and consistency that different positioning attributes should possess in temporal changes. When anomalies occur, they cannot accurately diagnose the source, direction, and specific form of the anomaly in the multi-dimensional attribute space, causing the filter to be unable to accurately suppress the impact of the anomaly. Summary of the Invention

[0005] This invention provides an abnormal location data cleaning method based on multi-base station positioning, in order to solve the technical problem that the existing technology fails to deeply explore and utilize the inherent correlation and consistency that different location attributes should have in the temporal changes, and cannot accurately diagnose the source, direction and specific form of the anomaly in the multi-dimensional attribute space.

[0006] This invention provides a method for cleaning abnormal location data based on multi-base station positioning, comprising the following steps: S1, acquire the original positioning data of the terminal to be positioned collected by multiple base stations within a preset time window and the standard deviation of measurement noise corresponding to each positioning attribute; S2, for continuous data points within the time window, calculate the temporal gradient sequence of each localization attribute and construct a cross-attribute temporal gradient correlation matrix; S3, based on the cross-attribute temporal gradient correlation matrix, solve for the inter-attribute conflict index and extract the dominant anomaly direction vector; S4 maps the inter-attribute conflict index to a scalar anomaly magnitude using a probability model trained on historical normal data. S5. Based on the scalar anomaly amplitude, the dominant anomaly direction vector, and the standard deviation of measurement noise corresponding to each positioning attribute, an anomaly covariance matrix with physical dimensions is constructed in the positioning attribute space. S6, obtain the positioning solution function used to solve the original positioning data to calculate the measurement position, and calculate the Jacobian matrix of the positioning solution function at the current working point; S7. Using the Jacobian matrix, the anomaly covariance matrix with physical dimensions in the attribute space is projected onto the location space to generate the measurement noise adjustment matrix in the location space. S8. Add the measurement noise adjustment matrix to the measurement noise covariance matrix preset by the Kalman filter to obtain the updated measurement noise covariance matrix. S9 employs a Kalman filter with an updated measurement noise covariance matrix to update the state of the measurement position, thereby reducing the impact of abnormal data on the positioning results.

[0007] Its effect is as follows: by constructing a cross-attribute temporal gradient correlation matrix and mapping it to an anomaly covariance matrix with physical dimensions, it solves the problems of information loss caused by simply removing outliers in complex non-line-of-sight environments and sensitivity to environmental changes in traditional methods. It achieves adaptive adjustment of the filtering gain according to the degree of data conflict without discarding the original observations, which significantly improves the positioning continuity and dynamic accuracy in strong interference environments.

[0008] Furthermore, for continuous data points within the time window, the temporal gradient sequence of each localization attribute is calculated, and a cross-attribute temporal gradient correlation matrix is ​​constructed, including: For each positioning attribute at the current moment, the backward difference method is used to calculate the time gradient, and the time gradient sequence of each attribute is obtained. For the original temporal gradient sequence of each attribute, the gradient value at each time step is divided by the standard deviation of the measurement noise corresponding to that attribute to obtain the normalized gradient. Set a sliding time window of length W, and extract the normalized gradient sequence fragments of the current time step and the previous W-1 consecutive time steps. Within the sliding time window, calculate the Pearson correlation coefficient between normalized gradient sequence segments corresponding to any two localization attributes; The correlation coefficients of all attribute pairs are combined into an N×N symmetric matrix, where N is the number of categories of location attributes. This matrix is ​​the cross-attribute temporal gradient correlation matrix.

[0009] Its effect is that by calculating the Pearson correlation coefficient within the sliding time window, it achieves a fine quantification of the consistency of temporal changes among different positioning attributes, and can effectively capture potential conflicts between attributes caused by signal blockage or multipath effects during multi-base station cooperative positioning, providing a reliable spatiotemporal feature basis for accurately diagnosing the source and direction of anomalies.

[0010] Furthermore, based on the cross-attribute temporal gradient correlation matrix, the inter-attribute conflict index is obtained, and the dominant anomaly direction vector is extracted, including: Eigenvalue decomposition is performed on the cross-attribute temporal gradient correlation matrix to obtain N eigenvalues ​​and their corresponding eigenvectors or feature spaces. The feature value with the smallest numerical value is selected as the inter-attribute conflict index; Extract any normalized feature vector in the feature space corresponding to the minimum feature value, and use it as the dominant anomaly direction vector.

[0011] Its effect is that by using eigenvalue decomposition technology to extract the smallest eigenvalue and its corresponding eigenvector in the correlation matrix, it can accurately extract the dimension with the most intense conflict between attributes from the mathematical level, thereby accurately locking the direction of the dominant anomaly in the multidimensional signal space. Compared with the conventional statistical threshold method, it can more sensitively identify subtle inconsistencies and fluctuations between multiple signals.

[0012] Furthermore, the inter-attribute conflict index is mapped to a scalar anomaly magnitude using a probabilistic model trained on historical normal data, including: Using historical normal location datasets, a series of inter-attribute conflict indices are calculated, and the probability distribution of the inter-attribute conflict indices is fitted using the maximum likelihood estimation method to obtain a gamma distribution model. For the currently calculated inter-attribute conflict index, the corresponding anomaly confidence level is calculated using the cumulative distribution function of the gamma distribution model; The anomaly confidence level is input into a preset nonlinear mapping function to calculate the scalar anomaly magnitude. The nonlinear mapping function outputs zero when the input is zero.

[0013] Its effect is that, by using a gamma distribution model based on historical data and a nonlinear mapping function, the abstract degree of attribute conflict is transformed into a statistically significant scalar anomaly amplitude, which realizes the scientific classification and quantification of the severity of anomalies and ensures that the system has extremely high robustness and stability in anomaly judgment under different electromagnetic environments and signal-to-noise ratio conditions.

[0014] Furthermore, based on the scalar anomaly amplitude, the dominant anomaly direction vector, and the standard deviation of measurement noise corresponding to each positioning attribute, an anomaly covariance matrix with physical dimensions in the positioning attribute space is constructed, including: Obtain the dominant anomaly direction vector and its transpose; The outer product of the dominant anomaly direction vector and its transpose is performed to obtain an N×N direction matrix. Scalar multiplication of the scalar anomaly magnitude and the direction matrix yields a dimensionless covariance matrix; A diagonal standard deviation matrix is ​​constructed using the measurement noise standard deviation corresponding to each positioning attribute. The diagonal standard deviation matrix, the dimensionless covariance matrix, and the diagonal standard deviation matrix are then multiplied sequentially to obtain an anomaly covariance matrix with physical dimensions.

[0015] Its effect is that by constructing an anomaly covariance matrix with physical dimensions through the composite operation of the direction matrix and the diagonal standard deviation matrix, the abstract statistical anomaly can be directly matched with the actual physical range of the sensor, solving the problem of the inconsistency between the mathematical model indicators and the physical measurement dimensions, and providing a direct and usable physical basis for the accurate correction of the filter closed-loop parameters.

[0016] Furthermore, using the Jacobian matrix, the anomaly covariance matrix with physical dimensions in the attribute space is projected onto the location space to generate a measurement noise adjustment matrix in the location space, including: The positioning solution function is determined to be a function that maps N positioning attributes to m-dimensional position coordinates, where m is the dimension of the position space; The first-order partial derivatives of the positioning solution function with respect to each input positioning attribute are obtained, and the partial derivatives are arranged into an m×N matrix to obtain the Jacobian matrix. The measurement noise adjustment matrix is ​​obtained by performing matrix multiplication operations sequentially on the Jacobian matrix, the anomaly covariance matrix with physical dimensions, and the transpose of the Jacobian matrix.

[0017] Furthermore, the measurement noise adjustment matrix is ​​added to the measurement noise covariance matrix preset by the Kalman filter to obtain the updated measurement noise covariance matrix, which includes: Obtain the preset measurement noise covariance matrix at the current time from the Kalman filter; The measurement noise adjustment matrix is ​​added to the preset measurement noise covariance matrix to obtain the updated measurement noise covariance matrix. In the measurement update step of the Kalman filter, the Kalman gain is calculated and the state estimate is updated using the updated measurement noise covariance matrix.

[0018] Furthermore, the raw positioning data includes the time difference of arrival, angle of arrival, and reference signal received power.

[0019] Furthermore, the Jacobi iteration method or QR algorithm is used to perform eigenvalue decomposition on the cross-attribute temporal gradient correlation matrix.

[0020] Furthermore, the Jacobian matrix is ​​a matrix composed of all the first-order partial derivatives of a multivariate vector-valued function, used to describe the local linear change relationship of the output vector caused by a small change in the input variable at a specific point in the localization solution function.

[0021] The beneficial effects are as follows: This invention constructs a cross-attribute temporal gradient correlation matrix representing the consistency of temporal changes in multiple positioning attributes, and performs eigenvalue decomposition on the matrix. This enables the identification and representation of anomalies from the perspective of inherent conflicts between attributes, determining not only the severity of anomalies but also identifying the main sources of anomalies. An anomaly covariance matrix with clear physical meaning in the attribute space is constructed, and the covariance matrix is ​​projected onto the location space using the Jacobian matrix. This is used to adjust the measurement noise covariance matrix in the Kalman filter, allowing the filtering process to adjust the confidence level of the current measurement according to the degree of data anomaly, thereby suppressing the impact of abnormal measurement data on the positioning results and improving the overall stability and positioning accuracy of the positioning system. Attached Figure Description

[0022] Figure 1 This is a flowchart of an abnormal location data cleaning method based on multi-base station positioning. Detailed Implementation

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

[0024] An embodiment of the abnormal location data cleaning method based on multi-base station positioning provided by the present invention: like Figure 1 As shown, the abnormal location data cleaning method based on multi-base station positioning includes the following steps: S1, obtain the original positioning data and attribute noise standard deviation of multiple base stations.

[0025] The system acquires the raw positioning data of the terminal to be located from multiple base stations within a preset time window, as well as the standard deviation of the measurement noise corresponding to each positioning attribute.

[0026] Set a time sliding window of fixed length W; at each sampling time t, obtain the original positioning data vector of the terminal to be located from a multi-base station network such as 5G or LTE. This vector contains at least two different types of positioning attributes, such as Time Difference of Arrival (TDOA), Angle of Arrival (AOA), and Reference Received Power (RSRP). Simultaneously, based on the base station equipment manual or system calibration results, a measurement noise standard deviation vector with the same dimension as the positioning attribute vector is obtained. Each element Let be the standard deviation of the measurement noise for the i-th positioning attribute.

[0027] S2, calculate the temporal gradient of the attribute and construct the cross-attribute gradient correlation matrix.

[0028] For continuous data points within a time window, calculate the temporal gradient sequence for each location attribute and construct a cross-attribute temporal gradient correlation matrix.

[0029] For each localization attribute sequence within the time window, the temporal gradient is calculated using the first-order backward difference method. That is, the gradient value of the i-th attribute at time t is equal to the attribute value at time t minus the attribute value at time t-1, and then divided by the sampling time interval. To eliminate differences in dimensions and scales among different positioning attributes, the temporal gradient sequence of each attribute is normalized. The gradient value of the i-th attribute at time t is divided by the corresponding standard deviation of the measurement noise to obtain the normalized gradient. The normalized gradient values ​​of all N attributes at all times within the time window are used to construct an N-row, W-1-column normalized gradient matrix G. Based on this matrix G, the Pearson product-moment correlation coefficient algorithm is used to calculate the correlation between rows of the matrix, generating an N-row, N-column cross-attribute temporal gradient correlation matrix. The element in the i-th row and j-th column of the matrix represents the correlation coefficient between the normalized temporal gradient sequences of the i-th attribute and the j-th attribute.

[0030] In an optional embodiment, for continuous data points within a time window, the temporal gradient sequence of each localization attribute is calculated, and a cross-attribute temporal gradient correlation matrix is ​​constructed, including: For each positioning attribute at the current moment, the backward difference method is used to calculate the time gradient, and the time gradient sequence of each attribute is obtained. For the original temporal gradient sequence of each attribute, the gradient value at each time step is divided by the standard deviation of the measurement noise corresponding to that attribute to obtain the normalized gradient. Set a sliding time window of length W, and extract the normalized gradient sequence fragments of the current time step and the previous W-1 consecutive time steps. Within the sliding time window, calculate the Pearson correlation coefficient between normalized gradient sequence segments corresponding to any two localization attributes; The correlation coefficients of all attribute pairs are combined into an N×N symmetric matrix, where N is the number of categories of location attributes. This matrix is ​​the cross-attribute temporal gradient correlation matrix.

[0031] There are N positioning attributes, for example, N=3, corresponding to Time of Arrival (TOA), Time Difference of Arrival (TDOA), and Angle of Arrival (AOA). Assume the data sampling time interval is... The preferred value is 0.1 seconds. At each discrete time point... For the i-th positioning attribute time gradient At any moment The value is calculated using the backward difference formula: The normalized gradient is Set a sliding time window of length W, with an optimal range of 10 to 30, for example, choosing W=20. At the current moment... Extract each attribute i within the time window The normalized gradient sequence segment within is denoted as Calculate the gradient sequence fragments for any two attributes i and j. and Pearson correlation coefficient between The calculation formula is: Where cov represents covariance. , These are normalized gradient sequence fragments. , The sample standard deviation. Construct an N×N cross-attribute time-series gradient correlation matrix from all calculated correlation coefficients. The element in the i-th row and j-th column of the matrix That is, the corresponding correlation coefficient. Because and =1, this matrix is ​​a symmetric matrix with diagonal elements of 1.

[0032] S3, solve for the attribute conflict index and extract the anomaly dominance direction vector.

[0033] Based on the cross-attribute temporal gradient correlation matrix, the inter-attribute conflict index is obtained, and the dominant anomaly direction vector is extracted.

[0034] Use eigenvalue decomposition algorithms, such as the Jacobi iteration method or the QR algorithm, to analyze the cross-attribute temporal gradient correlation matrix. Decompose to obtain N eigenvalues , , , and the corresponding N feature vectors , , , Sort all eigenvalues ​​and select the eigenvalue with the smallest value. As an index of conflict between attributes; at the same time, extracting and corresponding feature vector , serving as the dominant anomaly direction vector.

[0035] In an optional embodiment, based on the cross-attribute temporal gradient correlation matrix, the inter-attribute conflict index is obtained, and the dominant anomaly direction vector is extracted, including: Eigenvalue decomposition is performed on the cross-attribute temporal gradient correlation matrix to obtain N eigenvalues ​​and their corresponding eigenvectors or feature spaces. The feature value with the smallest numerical value is selected as the inter-attribute conflict index; Extract any normalized feature vector in the feature space corresponding to the minimum feature value, and use it as the dominant anomaly direction vector.

[0036] Obtain the N×N cross-attribute temporal gradient correlation matrix Then, the matrix is ​​subjected to standard eigenvalue decomposition, and the equation is solved. =0, where I is the identity matrix, and det is the abbreviation for determinant in mathematics, yielding N real eigenvalues. , , , And solve the equation =0 yields the values ​​of each eigenvalue. corresponding feature vector .because It is a positive semi-definite matrix, and all eigenvalues ​​are non-negative. The eigenvalues ​​are then sorted in descending order of magnitude. Select the smallest eigenvalue. The inter-attribute conflict index serves as a representation of the degree of anomaly in the current data point. In an ideal, conflict-free scenario, all attributes exhibit completely consistent temporal changes. The rank of the matrix is ​​1, at which point the smallest eigenvalue is... Theoretically, it is 0; however, when there is inconsistency between attributes, The value will be greater than 0; the magnitude of the eigenvalue indicates the severity of the conflict, for example, a value between 0.01 and 0.5. Simultaneously, the minimum eigenvalue is extracted. corresponding feature vector The vector is then normalized using the L2 norm to obtain the dominant anomaly direction vector. This vector indicates the primary direction in the N-dimensional attribute space that leads to inconsistencies in temporal changes, i.e., conflicts.

[0037] S4, based on a probability model, maps the conflict index to anomaly amplitude.

[0038] The inter-attribute conflict index is mapped to a scalar anomaly magnitude using a probabilistic model trained on historical normal data.

[0039] A large amount of known location data collected under normal working conditions was collected, and the above steps were repeated on this historical data to calculate a series of inter-attribute conflict indices. The distribution of historical normal conflict indices is fitted using the generalized Pareto distribution model from extreme value theory to determine the distribution's threshold, scale, and shape parameters. The generalized Pareto distribution is a probabilistic model used to describe the distribution of extreme data exceeding a certain high threshold, and is commonly used for extreme event analysis and anomaly detection. For the currently calculated conflict index... The cumulative distribution function of the trained generalized Pareto distribution model is used for mapping, and the right tail probability of the conflict index or the conditional probability exceeding the preset threshold is calculated. The probability value or the result of the nonlinear function transformation is used as the dimensionless scalar anomaly amplitude s.

[0040] In an alternative embodiment, mapping the inter-attribute conflict index to a scalar anomaly magnitude using a probabilistic model trained on historical normal data includes: Using historical normal location datasets, a series of inter-attribute conflict indices are calculated, and the probability distribution of the inter-attribute conflict indices is fitted using the maximum likelihood estimation method to obtain a gamma distribution model. For the currently calculated inter-attribute conflict index, the corresponding anomaly confidence level is calculated using the cumulative distribution function of the gamma distribution model; The anomaly confidence level is input into a preset nonlinear mapping function to calculate the scalar anomaly magnitude. The nonlinear mapping function outputs zero when the input is zero.

[0041] The probabilistic model is a gamma distribution model, trained offline using historical normal data. The gamma distribution model takes the form of a distribution with shape parameter k and scale parameter k. The defined continuous probability distribution has the following probability density function: ,in It is the gamma function; this model uses the cumulative distribution function for inference calculations. The real-time conflict index to be detected is used. In the online application phase, the input to the gamma distribution model is the scalar inter-attribute conflict index calculated at the current time. The output is a scalar anomaly confidence level. The confidence level is between 0 and 1.

[0042] During the offline phase, a large amount of location data, such as more than 10,000 data points, is collected from the terminal under normal operating conditions. The conflict index between time point attributes is calculated for each time point in this dataset. This yields a set of historical conflict indices. Because the conflict indices... Since the numbers are non-negative real numbers, a gamma distribution is used to fit the probability density. The optimal parameters k and are estimated from this set of historical conflict indices using the maximum likelihood estimation method. For example, we get k=2.1. =0.04, model training complete. In the online application phase, for the conflict index calculated at the current moment... The corresponding P-value, i.e., the anomaly confidence level, is calculated using the cumulative distribution function of the pre-trained gamma distribution model. This value represents the probability that the current conflict index occurs within the normal data distribution; the closer the value is to 1, the more likely it is to be an anomaly. For example, It equals 0.25, calculated as follows It equals 0.98. A nonlinear mapping function is input to generate a scalar anomaly magnitude s. A preferred mapping function is: A and It is an adjustable hyperparameter, for example, A=200. =2.5. This function outputs 0 when the input is 0, and can amplify confidence levels close to 1 into abnormal amplitudes, for example, s≈189.3.

[0043] S5 integrates anomalous amplitude, direction, and noise to construct an anomalous covariance matrix.

[0044] Based on the scalar anomaly amplitude, the dominant anomaly direction vector, and the standard deviation of measurement noise corresponding to each positioning attribute, an anomaly covariance matrix with physical dimensions is constructed in the positioning attribute space.

[0045] The standard deviation of the measurement noise for each positioning attribute Construct as a diagonal matrix The diagonal elements are , , , Through the formula Constructing the anomaly covariance matrix The matrix is ​​an N-row, N-column symmetric positive semi-definite matrix, and its dimensions are consistent with the dimensions of the covariance of the location attributes.

[0046] In an optional embodiment, based on the scalar anomaly amplitude, the dominant anomaly direction vector, and the standard deviation of measurement noise corresponding to each location attribute, an anomaly covariance matrix with physical dimensions in the location attribute space is constructed, including: Obtain the dominant anomaly direction vector and its transpose; The outer product of the dominant anomaly direction vector and its transpose is performed to obtain an N×N direction matrix. Scalar multiplication of the scalar anomaly magnitude and the direction matrix yields a dimensionless covariance matrix; A diagonal standard deviation matrix is ​​constructed using the measurement noise standard deviation corresponding to each positioning attribute. The diagonal standard deviation matrix, the dimensionless covariance matrix, and the diagonal standard deviation matrix are then multiplied sequentially to obtain an anomaly covariance matrix with physical dimensions.

[0047] The N×1 dominant anomaly direction vector obtained from the above calculation and the transpose of 1×N Calculate the outer product of the two to generate an N×N rank-1 direction matrix. This matrix defines the distribution direction of anomalies in the attribute space. The scalar anomaly magnitude s obtained in the previous stage, for example 189.3, is compared with the direction matrix. Performing scalar multiplication yields a dimensionless outlier covariance matrix. Obtain the standard deviation of measurement noise for each attribute. The standard deviation value is usually obtained from the sensor specifications or offline calibration, for example, =15ns, =8ns, =0.02 radians. Construct an N×N diagonal standard deviation matrix. The diagonal elements are Using sandwich multiplication By scaling the dimensionless covariance matrix to the correct physical scale, an anomaly covariance matrix with physical dimensions is obtained. , The diagonal elements are the squares of the units of each attribute, such as nanosecond squared or radian squared, while the off-diagonal elements represent the correlation between attributes caused by anomalies.

[0048] S6, obtain the positioning solution function and solve its current Jacobian matrix.

[0049] Obtain the positioning solution function used to solve the raw positioning data to calculate the measurement position, and calculate the Jacobian matrix of the positioning solution function at the current working point.

[0050] Define a localization function x = h(p) to transform an N-dimensional attribute space into an m-dimensional location space, where p is an N×1 localization attribute vector and x is an m×1 location coordinate vector. This function is derived based on the geometric relationships between the terminal and each base station, such as polygonal and multi-angle measurements. At the current moment, the measured attribute vector... As the working point, compute the Jacobian matrix of the function h at that point. Its dimensions are m×N, and its elements are... The Jacobian matrix represents the partial derivative of the i-th location component with respect to the j-th location attribute. It is a matrix composed of all the first-order partial derivatives of a multivariate vector-valued function and is used to describe the local linear change relationship of the output vector caused by a small change in the input variable at a specific point.

[0051] S7 generates a location space noise adjustment matrix through Jacobian matrix projection.

[0052] Using the Jacobian matrix, the anomaly covariance matrix with physical dimensions in the attribute space is projected onto the location space to generate the measurement noise adjustment matrix in the location space.

[0053] In an optional embodiment, the anomaly covariance matrix with physical dimensions in the attribute space is projected onto the location space using the Jacobian matrix to generate a measurement noise adjustment matrix in the location space, including: The positioning solution function is determined to be a function that maps N positioning attributes to m-dimensional position coordinates, where m is the dimension of the position space; The first-order partial derivatives of the positioning solution function with respect to each input positioning attribute are obtained, and the partial derivatives are arranged into an m×N matrix to obtain the Jacobian matrix. The measurement noise adjustment matrix is ​​obtained by performing matrix multiplication operations sequentially on the Jacobian matrix, the anomaly covariance matrix with physical dimensions, and the transpose of the Jacobian matrix.

[0054] The localization solution function x=h(p) is used to convert an N-dimensional attribute space to an m-dimensional location space, where p is an N×1 localization attribute vector and x is an m×1 location coordinate vector. For example, in two-dimensional localization, m=2. The function h is typically a nonlinear function derived from the geometric relationships of multilateral or multi-angle measurements. To linearly propagate the covariance, it is necessary to calculate the attribute value of this function at the current operating point, i.e., the currently measured value. Jacobian matrix at the location Jacobian matrix It is an m×N matrix with elements This refers to the partial derivative of the i-th component of the location coordinate with respect to the j-th component of the positioning attribute. For example, in a TOA-based two-dimensional three-base station positioning system, N=3, m=2, It's a 2×3 matrix. Obtain the Jacobian matrix. Then, using the linear error propagation law, the abnormal covariance matrix in the attribute space is... Projected into the location space, the calculation formula is: The operations with dimensions of m×N, N×N, and N×m result in an m×m matrix. . This indicates an anomaly at the attribute level, affecting the calculated position coordinates (e.g.) , The contribution of uncertainty.

[0055] S8 merges the adjustment matrix and covariance matrix to update the measurement noise covariance.

[0056] The measurement noise adjustment matrix is ​​added to the measurement noise covariance matrix preset by the Kalman filter to obtain the updated measurement noise covariance matrix.

[0057] Obtain the preset measurement noise covariance matrix in the Kalman filter. The matrix is The location space matrix is ​​typically set as a diagonal matrix, with the diagonal elements representing estimates of the variance of the location measurement noise. Perform matrix addition. ,in The position space measurement noise adjustment matrix generated in the previous step is used to obtain the updated measurement noise covariance matrix. .

[0058] In an optional embodiment, the measurement noise adjustment matrix is ​​added to the measurement noise covariance matrix preset by the Kalman filter to obtain an updated measurement noise covariance matrix, including: Obtain the preset measurement noise covariance matrix at the current time from the Kalman filter; The measurement noise adjustment matrix is ​​added to the preset measurement noise covariance matrix to obtain the updated measurement noise covariance matrix. In the measurement update step of the Kalman filter, the Kalman gain is calculated and the state estimate is updated using the updated measurement noise covariance matrix.

[0059] A preset m×m measurement noise covariance matrix is ​​defined. It represents the statistical characteristics of measurement noise under normal conditions. For example, for a two-dimensional positioning system, Possibly set as ,in, Represents a diagonal matrix. and It is the expected standard deviation of the position measurement in the x and y directions.

[0060] S9 utilizes updated covariance matrix filtering to reduce the impact of anomalies on localization.

[0061] A Kalman filter with an updated measurement noise covariance matrix is ​​used to update the state of the measurement position, thereby reducing the impact of abnormal data on the positioning results.

[0062] In the update step of the extended Kalman filter, the updated measurement noise covariance matrix is ​​used. Perform the following calculations: Calculate the new information covariance matrix. , ,in, It is the observation matrix of the filter. This represents the position uncertainty predicted in the previous step. Calculate the Kalman gain. , Update state estimates , Among them, new information This is the difference between the actual measured value and the predicted measured value calculated based on the predicted state. This represents the predicted position from the previous step. Update the state covariance matrix using either the Joseph form or the standard form. This completes the current iteration.

[0063] Ablation experiments were conducted to verify the rationality of the cross-attribute temporal anomaly detection and filtering scheme.

[0064] The experimental conditions were set as a simulated two-dimensional positioning scenario in an urban environment, with the target moving along a preset trajectory for a total duration of 300 seconds and a data sampling frequency of 10Hz. The positioning system fused three attributes: Time of Arrival (TOA), Time Difference of Arrival (TDOA), and Angle of Arrival (AOA). Between 150 and 152 seconds, a step bias fault with an amplitude of 100 meters was injected into the TOA measurement to simulate a non-Gaussian anomaly. The experiment compared four schemes: a standard extended Kalman filter baseline scheme, a scheme based on a fixed threshold for minimum eigenvalue elimination, a scheme using only the anomaly amplitude for isotropic noise covariance adjustment, and the complete scheme proposed in this invention, combining anomaly amplitude and direction vector.

[0065] The root mean square error (RMSE) of positioning for the four schemes was recorded throughout the entire operating cycle. The RMSE of positioning for the standard extended Kalman filter baseline scheme was 16.82 meters; the RMSE of positioning for the scheme based on the minimum eigenvalue fixed threshold elimination was 6.45 meters; the RMSE of positioning for the scheme that only uses abnormal amplitude for isotropic adjustment was 4.13 meters; and the RMSE of positioning for the complete scheme proposed in this invention was 1.58 meters.

[0066] The complete scheme proposed in this invention offers improved performance compared to other schemes. The standard extended Kalman filter baseline scheme cannot handle non-Gaussian anomalies, leading to severe contamination of state estimation by anomalous data and resulting in the largest error. The fixed threshold rejection scheme improves performance by identifying and discarding anomalous data points, but data rejection leads to information loss, affecting the continuity and accuracy of filtering. The scheme using only the anomalous amplitude further improves performance by adjusting the noise covariance matrix to reduce the weight of anomalous measurements, avoiding complete information loss. However, since the adjustment is isotropic, it fails to accurately identify the source of the anomaly, thus also suppressing normal attribute information to some extent. The complete scheme of this invention not only represents the severity of the anomaly but also locates the direction of the anomaly in the multidimensional attribute space by dominating the anomaly direction vector. This achieves accurate and directional adjustment of the measurement noise covariance matrix, suppressing the negative impact of TOA attribute anomalies while retaining useful information in TDOA and AOA attributes, thus achieving the lowest localization error and the strongest system stability.

[0067] The above are all preferred embodiments of the present invention and are not intended to limit the scope of protection of the present invention. Therefore, all equivalent changes made in accordance with the structure, shape and principle of the present invention should be covered within the scope of protection of the present invention.

Claims

1. A method for cleaning abnormal location data based on multi-base station positioning, characterized in that, Includes the following steps: S1, acquire the original positioning data of the terminal to be positioned collected by multiple base stations within a preset time window and the standard deviation of measurement noise corresponding to each positioning attribute; S2, for continuous data points within the time window, calculate the temporal gradient sequence of each localization attribute and construct a cross-attribute temporal gradient correlation matrix; S3, based on the cross-attribute temporal gradient correlation matrix, solve for the inter-attribute conflict index and extract the dominant anomaly direction vector; S4 maps the inter-attribute conflict index to a scalar anomaly magnitude using a probability model trained on historical normal data. S5. Based on the scalar anomaly amplitude, the dominant anomaly direction vector, and the standard deviation of measurement noise corresponding to each positioning attribute, an anomaly covariance matrix with physical dimensions is constructed in the positioning attribute space. S6, obtain the positioning solution function used to solve the original positioning data to calculate the measurement position, and calculate the Jacobian matrix of the positioning solution function at the current working point; S7. Using the Jacobian matrix, the anomaly covariance matrix with physical dimensions in the attribute space is projected onto the location space to generate the measurement noise adjustment matrix in the location space. S8. Add the measurement noise adjustment matrix to the measurement noise covariance matrix preset by the Kalman filter to obtain the updated measurement noise covariance matrix. S9 employs a Kalman filter with an updated measurement noise covariance matrix to update the state of the measurement position, thereby reducing the impact of abnormal data on the positioning results.

2. The abnormal location data cleaning method based on multi-base station positioning according to claim 1, characterized in that, For continuous data points within a time window, calculate the temporal gradient sequence for each localization attribute and construct a cross-attribute temporal gradient correlation matrix, including: For each positioning attribute at the current moment, the backward difference method is used to calculate the time gradient, and the time gradient sequence of each attribute is obtained. For the original temporal gradient sequence of each attribute, the gradient value at each time step is divided by the standard deviation of the measurement noise corresponding to that attribute to obtain the normalized gradient. Set a sliding time window of length W, and extract the normalized gradient sequence fragments of the current time step and the previous W-1 consecutive time steps. Within the sliding time window, calculate the Pearson correlation coefficient between normalized gradient sequence segments corresponding to any two localization attributes; The correlation coefficients of all attribute pairs are combined into an N×N symmetric matrix, where N is the number of categories of location attributes. This matrix is ​​the cross-attribute temporal gradient correlation matrix.

3. The abnormal location data cleaning method based on multi-base station positioning according to claim 1, characterized in that, Based on the cross-attribute temporal gradient correlation matrix, the inter-attribute conflict index is obtained, and the dominant anomaly direction vector is extracted, including: Eigenvalue decomposition is performed on the cross-attribute temporal gradient correlation matrix to obtain N eigenvalues ​​and their corresponding eigenvectors or feature spaces. The feature value with the smallest numerical value is selected as the inter-attribute conflict index; Extract any normalized feature vector in the feature space corresponding to the minimum feature value, and use it as the dominant anomaly direction vector.

4. The abnormal location data cleaning method based on multi-base station positioning according to claim 1, characterized in that, The inter-attribute conflict index is mapped to a scalar anomaly magnitude using a probabilistic model trained on historical normal data, including: Using historical normal location datasets, a series of inter-attribute conflict indices are calculated, and the probability distribution of the inter-attribute conflict indices is fitted using the maximum likelihood estimation method to obtain a gamma distribution model. For the currently calculated inter-attribute conflict index, the corresponding anomaly confidence level is calculated using the cumulative distribution function of the gamma distribution model; The anomaly confidence level is input into a preset nonlinear mapping function to calculate the scalar anomaly magnitude. The nonlinear mapping function outputs zero when the input is zero.

5. The abnormal location data cleaning method based on multi-base station positioning according to claim 3, characterized in that, Based on the scalar anomaly amplitude, the dominant anomaly direction vector, and the standard deviation of measurement noise corresponding to each location attribute, an anomaly covariance matrix with physical dimensions is constructed in the location attribute space, including: Obtain the dominant anomaly direction vector and its transpose; The outer product of the dominant anomaly direction vector and its transpose is performed to obtain an N×N direction matrix. Scalar multiplication of the scalar anomaly magnitude and the direction matrix yields a dimensionless covariance matrix; A diagonal standard deviation matrix is ​​constructed using the measurement noise standard deviation corresponding to each positioning attribute. The diagonal standard deviation matrix, the dimensionless covariance matrix, and the diagonal standard deviation matrix are then multiplied sequentially to obtain an anomaly covariance matrix with physical dimensions.

6. The abnormal location data cleaning method based on multi-base station positioning according to claim 1, characterized in that, Using the Jacobian matrix, the anomaly covariance matrix with physical dimensions in the attribute space is projected onto the location space to generate the measurement noise adjustment matrix in the location space, including: The positioning solution function is determined to be a function that maps N positioning attributes to m-dimensional position coordinates, where m is the dimension of the position space; The first-order partial derivatives of the positioning solution function with respect to each input positioning attribute are obtained, and the partial derivatives are arranged into an m×N matrix to obtain the Jacobian matrix. The measurement noise adjustment matrix is ​​obtained by performing matrix multiplication operations sequentially on the Jacobian matrix, the anomaly covariance matrix with physical dimensions, and the transpose of the Jacobian matrix.

7. The abnormal location data cleaning method based on multi-base station positioning according to claim 1, characterized in that, The measurement noise adjustment matrix is ​​added to the preset measurement noise covariance matrix of the Kalman filter to obtain the updated measurement noise covariance matrix, which includes: Obtain the preset measurement noise covariance matrix at the current time from the Kalman filter; The measurement noise adjustment matrix is ​​added to the preset measurement noise covariance matrix to obtain the updated measurement noise covariance matrix. In the measurement update step of the Kalman filter, the Kalman gain is calculated and the state estimate is updated using the updated measurement noise covariance matrix.

8. The method for cleaning abnormal location data based on multi-base station positioning according to claim 1, characterized in that, The raw positioning data includes the time difference of arrival, angle of arrival, and reference signal received power.

9. The method for cleaning abnormal location data based on multi-base station positioning according to claim 3, characterized in that, Eigenvalue decomposition of the cross-attribute temporal gradient correlation matrix is ​​performed using the Jacobi iteration method or the QR algorithm.

10. The method for cleaning abnormal location data based on multi-base station positioning according to claim 1, characterized in that, The Jacobian matrix is ​​a matrix composed of all the first-order partial derivatives of a multivariate vector-valued function. It is used to describe the local linear change relationship of the output vector caused by a small change in the input variable at a specific point in the localization solution function.

Images