A method and system for anti-trajectory-pulling deception positioning oriented to urban canyon environment

Abstract

The application discloses a kind of anti trajectory traction deception positioning technical field for urban canyon environment's anti trajectory traction deception positioning method and system, method includes the following steps: positioning system obtains SLAM data and GNSS data, and obtains pose estimation value and GNSS observation value after pre-processing;Obtain prior state estimation and prior state covariance matrix;Based on GNSS observation value, prior state estimation and prior state covariance matrix, state machine determination is carried out and positioning result is obtained.The application carries out signal and physical layer pre-inspection when state machine determination is determined, and whether it is suffered from trajectory traction deception based on Mahalanobis distance is judged;If suffer from deception, then cut off GNSS signal, and update GNSS observation value, when pre-inspection passes through and trajectory traction deception disappears, based on damping confidence weight update GNSS observation value, reduce the rate of missed detection, avoid system long time isolation lock, improve the stability of positioning system.

Description

Technical Field

[0001] This invention relates to the field of anti-trajectory traction deception positioning technology, specifically to an anti-trajectory traction deception positioning method and system for outdoor unmanned equipment in urban canyon environments. Background Technology

[0002] With the large-scale deployment of autonomous driving and mobile robots, outdoor unmanned equipment (such as logistics vehicles and security patrol vehicles) heavily relies on multi-source combined navigation schemes using GNSS and local laser SLAM. GNSS is the abbreviation for Global Navigation Satellite System, and SLAM is the abbreviation for Simultaneous Localization and Mapping, used to solve the problem of robot localization and map building when moving in unknown environments. However, civilian GNSS signals have an open structure and weak power, making them extremely vulnerable to generative spoofing interference. Among them, slowly varying traction spoofing injects positional deviations slowly through extremely low induction rates, exhibiting strong concealment and harmfulness. Existing defense systems against slowly varying traction spoofing have the following theoretical and engineering shortcomings: Firstly, the deceptive drifting during the incubation and quarantine periods creates a dual antagonistic effect: Slowly varying deception causes the single-step residual to always fall within the confidence threshold of conventional filtering, leading to the failure of robust filtering based on single-epoch residual verification. The system needs to rely on higher-dimensional verification mechanisms (such as comparing the shape consistency between the GNSS trajectory and the SLAM-calculated trajectory) to detect the deception. However, in this process, the system accumulates double drift: during the deception latency period, the global pose is pulled off course due to the fusion of false signals; during the isolation period, after GNSS is cut off, the system relies solely on local SLAM calculations, accumulating relative drift.

[0003] Secondly, single-step state transition and operational degradation during the recovery period: When the actual GNSS signal is re-filtered, the accumulated spatial residuals will cause a drastic single-step state jump within the conventional EKF framework, disrupting the planning layer. and The continuity constraint caused a surge in tracking errors, forcing the system to slow down or suspend operations. For position constraints, For speed constraints.

[0004] Third, the inherent flaws of existing smooth recovery schemes: High-performance global optimization (such as factor graph optimization): Although it can achieve global smoothness, its huge overhead of state maintenance and iterative solution makes it difficult to deploy directly on embedded mobile platforms with limited computing power.

[0005] Stepped weight recovery: This scheme is based on stepped weight recovery of system-level state machine scheduling, but it will still cause discrete jitter in some areas and increase the coupling of the system.

[0006] Covariance adaptive adjustment based on statistical domain: This scheme falls into the "smoothing-convergence" paradox. If the R matrix is ​​too large, it will lead to slow convergence and decreased sensitivity of secondary detection; if the R matrix is ​​too small, it will be impossible to avoid jumps.

[0007] Single-step state limiting: This scheme directly truncates the state increment after the EKF measurement update. EKF is the abbreviation for Extended Kalman Filter. After truncation, the covariance and state become mismatched, destroying the filter consistency and causing subsequent spoofing detection to fail.

[0008] This shows that the existing anti-fraud system lacks a low-overhead security recovery mechanism that balances statistical consistency and operational continuity. Summary of the Invention

[0009] The purpose of this invention is to provide an anti-trajectory traction deception positioning method and system for urban canyon environments. It solves the technical problems of existing robust filtering algorithms when dealing with GNSS gradual traction deception of outdoor unmanned equipment, such as the easy contamination of the prior state in the early stage of deception, the easy divergence of pure odometer calculation during deception, and the kinematic step caused by the sudden change in positioning during the signal recovery period after deception, which destroys the continuity of the motion trajectory.

[0010] The present invention achieves the above objectives through the following technical solutions: An anti-trajectory traction deception positioning method for urban canyon environments, applied to outdoor unmanned equipment, includes the following steps: Outdoor unmanned equipment receives GNSS signals, acquires SLAM / oddsight pose data and GNSS data, and performs data preprocessing to obtain pose estimation values ​​and GNSS observation values. Based on the pose estimation value, obtain the prior state estimate and the prior state covariance matrix; State machine determination is performed based on GNSS observations, prior state estimates, and prior state covariance matrix. Positioning results are obtained based on prior state estimates, Kalman gain matrix, and equivalent observations. During the state machine determination process, multi-epoch cumulative residual characteristics are obtained based on pose estimation values ​​and GNSS observation values, and the Mahalanobis distance is used to determine whether trajectory traction deception has occurred. If the determination is negative, the GNSS signal is trusted and the GNSS observation value is used as the equivalent observation value. If the determination is positive, the GNSS signal is cut off and the equivalent observation value is updated based on the prior state estimation.

[0011] As a further optimization of the invention, when trusting the GNSS signal, the steady-state physical drift prior is continuously updated. After the GNSS signal is cut off, the system enters an isolated state and updates the equivalent observations to virtual observations based on the prior state estimation at the time of signal cutoff and the steady-state physical drift prior. When the trajectory traction deception ends, the system enters a recovery state and updates the equivalent observations to the fused observations based on the damped confidence weight fusion of the GNSS observations and virtual observations.

[0012] As a further optimization of the invention, the steady-state physical drift prior is iteratively updated through the following steps: In a normal fusion state without deception interference, the difference between the GNSS pose increment and the SLAM pose increment during the synchronization period is dynamically calculated, and the difference is divided by the time interval to obtain the instantaneous physical drift velocity. The instantaneous physical drift velocity is limited by the chassis kinematic physical boundary, and values ​​exceeding the limit are eliminated; The instantaneous physical drift velocity after clipping is input into a low-pass filter for iterative updates, and the steady-state physical drift prior is extracted.

[0013] As a further optimization of the invention, the multi-epoch cumulative residual feature is the squared Mahalanobis distance between the GNSS observation value and the pose estimation value within a fixed-length sliding window. When the squared Mahalanobis distance exceeds a preset safety threshold, it is determined to be due to GNSS deception interference. The formula for calculating the squared Mahalanobis distance is:

[0014] in, The squared Mahalanobis distance index for the current period. and These are the components of the innovation accumulation vector in the horizontal and vertical directions, respectively. The noise variance is measured as a benchmark. This represents the actual number of samples within the current sliding window.

[0015] As a further optimization of the invention, in the isolated state, GNSS observations are updated to virtual observations based on the following formula:

[0016] in, These are virtual observations. For prior state estimation, For the steady-state physical drift prior when the isolated state is triggered, For a single cycle.

[0017] As a further optimization of the invention, in the recovered state, GNSS observations are iteratively updated to fused observations based on the following formula:

[0018]

[0019] in, The first and Confidence weights for the period The damping coefficient is... To merge observations, These are GNSS observations.

[0020] As a further optimization of the invention, in the process of obtaining the positioning results: Constructing an adaptive expansion matrix based on an exponential adaptive robust model:

[0021] in, The preset reference measurement noise covariance matrix, It is an adaptive expansion matrix; When trusting GNSS signals, =1; When the GNSS signal is cut off, ; in, The squared Mahalanobis distance index for the current period. As a safety threshold, This is the weight decay rate factor. This is the lower bound of the weight.

[0022] As a further optimization of the invention, the preprocessing includes dynamic alignment of the global GNSS coordinate system and the local SLAM coordinate system, which includes: Coarse alignment stage: Acquire synchronization frame sequences, calculate the geometric centroid difference between the GNSS observation point set and the SLAM inferred point set, and extract the initial translation vector; Fine alignment stage: When the movement trajectory of the outdoor unmanned equipment reaches the preset kinematic excitation threshold, the collected GNSS observation point set and SLAM calculated point set are decentrifuged, the covariance matrix after decentrifugation is decomposed into singular values, the rotation matrix and translation vector are solved, and the coordinate system is aligned.

[0023] An anti-trajectory traction spoofing positioning system for urban canyon environments, comprising: The data acquisition module is used to receive GNSS signals, SLAM / odometer pose data and GNSS data and perform data preprocessing to obtain pose estimation values ​​and GNSS observation values; The state prediction module is used to obtain the prior state estimate and the prior state covariance matrix based on the pose calculation value. The positioning update module performs state machine determination based on GNSS observations, prior state estimates, and prior state covariance matrix, and obtains positioning results based on prior state estimates, Kalman gain matrix, and equivalent observations. The location update module includes: The deception detection unit is used to obtain multi-epoch cumulative residual characteristics based on pose estimation values ​​and GNSS observation values, and to determine whether the trajectory traction deception has been suffered based on Mahalanobis distance; The normal state processing unit is used to trust GNSS signals when not being deceived, continuously update the steady-state physical drift prior, use GNSS observations as equivalent observations, and obtain positioning results; The isolated state processing unit is used to cut off the GNSS signal when it is deceived, and update the equivalent observations based on the prior state estimate at the time of signal cutoff and the prior steady-state physical drift, and obtain the positioning results; The recovery state processing unit is used to fuse GNSS observations and virtual observations based on damped confidence weights when the pre-detection passes and the trajectory traction deception is determined to have ended, update the equivalent observations, and obtain the positioning results.

[0024] A storage medium having computer program instructions stored thereon, which implement the above method when executed by a processor.

[0025] The beneficial effects of this invention are as follows: This invention uses a SLAM front-end as an electromagnetically independent dead reckoning benchmark, combined with a dual-track detection mechanism of sliding window accumulated Mahalanobis distance and single-frame innovation squared. Multi-epoch sliding window accumulation can effectively identify sub-meter level slowly changing position drift that is extremely difficult to detect with traditional single-frame algorithms. At the same time, it can quickly break the alarm deadlock by judging the single-frame innovation at the moment of real signal return, reduce the false negative rate, and prevent the system from falling into a long-term isolation deadlock. This invention proposes a dynamic robustness mechanism based on asymmetric exponential response. During the deception period, the continuous differentiability of the exponential function is used to replace the traditional hard threshold switching, eliminating the sawtooth trajectory jitter caused by abrupt changes in filter gain. During the signal recovery period, a dual mechanism of locking the two-dimensional physical difference and the measurement covariance benchmark is introduced to slowly align the virtual calculation with the real observation. At the same time, the Kalman filter is forced to absorb a smooth trajectory without oscillation at the mathematical level, fundamentally eliminating the cliff-like position discontinuity and sawtooth jitter during mode switching, and ensuring the absolute stability of the robot control system. This invention designs a two-stage online spatiotemporal calibration strategy that progresses from coarse to fine. During the robot motion initialization stage, singular value decomposition and Helmert transform are used to automatically solve and compensate for rotational and translational deviations between GNSS and SLAM systems. This eliminates the reliance on precision tooling fixtures and manual calibration, significantly reducing the hardware integration threshold and maintenance complexity. Attached Figure Description

[0026] Figure 1 This is a flowchart of the anti-trajectory traction deception positioning method of the present invention; Figure 2 This is a schematic diagram comparing the effects of online spatiotemporal calibration of heterogeneous sensors based on simulation data according to the present invention before and after calibration. Figure 3 This is a simulation diagram of the deception detection method according to the present invention; Figure 4 This is a simulation diagram of the asymmetric robust weight adjustment of the present invention; Figure 5 This is a 2D trajectory simulation diagram under the GNSS deception attack of the present invention; Figure 6 This is a simulation diagram of the absolute error suppression and stability of the present invention. Detailed Implementation

[0027] The present application will now be described in further detail with reference to the accompanying drawings. It should be noted that the following specific embodiments are only used to further illustrate the present application and should not be construed as limiting the scope of protection of the present application. Those skilled in the art can make some non-essential improvements and adjustments to the present application based on the above application content. Example

[0028] like Figure 1-4 As shown, this embodiment relates to an anti-trajectory traction deception positioning method for urban canyon environments. This method is applied to outdoor unmanned equipment and includes the following steps: Step S100: The outdoor unmanned equipment receives GNSS signals through the positioning system, acquires SLAM data and GNSS data, preprocesses the SLAM data and GNSS data, and obtains the pose estimation value. and GNSS observations This enables multidimensional feature extraction and robust filtering. Step S200, this step corresponds to Figure 1 The time update process involves calculating values ​​based on pose. Perform SLAM state prediction and prior state covariance matrix prediction to obtain prior state estimates. and the prior state covariance matrix ; Step S300, this step corresponds to Figure 1 Measurement updates in the process, based on GNSS observations Prior state estimation and prior state covariance matrix Perform state machine determination and equivalent measurement to obtain equivalent observations. Based on prior state estimation Kalman gain matrix and equivalent observations Get the location results.

[0029] The initial state of the aforementioned state machine is the normal state. During the state machine's decision-making process, a signal and physical layer pre-check is performed. After the pre-check passes, the pose value is calculated. and GNSS observations The multi-epoch cumulative residual features are obtained, and the Mahalanobis distance is used to determine whether trajectory traction deception has occurred. The determination process is as follows: If the determination is negative, the GNSS signal is trusted, the state machine remains in the normal state, and the steady-state physical drift prior is continuously updated. GNSS observations As equivalent observations And repeat the time update and measurement update; If the determination is yes, the GNSS signal is cut off, and the state machine enters the isolated state, based on the prior state estimate at the time of signal cutoff. and steady-state physical drift prior Obtain virtual observations Virtual observations As equivalent observations And repeat the time update and measurement update; When the pre-detection passes and the trajectory traction deception is determined to have ended, the state machine enters the recovery state, based on damped confidence weights. Fusion of GNSS observations With virtual observations To obtain fused observations , will merge observations As equivalent observations And repeat time updates and measurement updates.

[0030] In the recovery state, Gradually increase, when When the recovery state ends, the state machine re-enters the normal state and repeats the time update and measurement update.

[0031] To address the issue of unknown initial installation errors between the GNSS receiver and the wheeled odometer (or SLAM system) in low-cost mobile robots, step S100 of this embodiment employs a two-stage coordinate system alignment strategy, from coarse to fine, to overcome this problem. For example... Figure 2 As shown, step S100 is mainly executed in the positioning data processing flow, and the specific flow is as follows: Step S110, State Initialization: The positioning system first initializes the calibration state, setting it to an unaligned initial state. Simultaneously, it initializes the GNSS data buffer and SLAM / odot pose data buffer, used for storing historical trajectory data, and initializes the coordinate transformation parameters (rotation matrix and translation vector) to the default identity matrix and zero vector. The SLAM / odot pose data buffer stores SLAM / odot pose data. This SLAM / odot pose data can be either SLAM data or odot pose data.

[0032] When a frame of raw GNSS data is received from a GNSS receiver, its geodetic coordinates are first converted to Cartesian coordinates in the local tangent plane coordinate system. Since the alignment with the SLAM / mileage pose data coordinate system has not yet been completed, this data is temporarily used as GNSS observation data to be aligned.

[0033] Step S120, Static coarse alignment stage: During the initial robot startup, the positioning system enters data acquisition mode, collecting GNSS and SLAM data, and performing static coarse alignment of the GNSS and SLAM data. The specific process is as follows: Step S121, Coarse Alignment Trigger Determination: The positioning system continuously collects GNSS observation points and corresponding SLAM / odometer pose points, forming valid data in the buffer queue consisting of a set of GNSS observation points and a set of SLAM / odometer pose points. GNSS observation points are local coordinate points obtained after coordinate transformation of the original GNSS data, and SLAM / odometer pose data outputs local pose points at the same time. When the number of valid data frames in the buffer queue reaches a preset coarse calibration threshold, coarse alignment calculation is automatically triggered.

[0034] Step S122, Coarse Alignment Calculation: Considering that the robot may be stationary or in a slightly moving state at the initial stage of startup, making it difficult to accurately analyze the rotational components, this stage mainly focuses on eliminating translational deviations. The positioning system calculates the geometric centroids of the GNSS observation point set and the SLAM / odot pose point set in the buffer queue, respectively. By comparing the positional differences between the two centroids, the initial translation vector is calculated.

[0035] Step S123, First State Transition: Update the calculated initial translation vector to the positioning system parameters and switch the positioning system state to coarse alignment mode. At this point, the positioning system has eliminated the main coordinate origin deviation.

[0036] Step S130, Dynamic Refinement Stage: After entering the fine alignment mode, the positioning system begins monitoring the robot's motion state to calculate rotational deviation. This rotational deviation is then eliminated through dynamic fine alignment, completing the calibration. The specific process of dynamic fine alignment is as follows: Step S131, Coarse Alignment Trigger Determination: The positioning system calculates the robot's current cumulative movement distance in real time. To prevent numerical instability caused by forced calculation during stationary rotation or minor jitter, this invention sets a motion excitation threshold, which is preferably 2.0 meters in this embodiment. When the displacement of the odometry trajectory exceeds this threshold, it indicates that the robot has generated a sufficiently long effective trajectory segment, triggering fine alignment calculation based on the Helmert transform. The Helmert transform is used to solve for the rotation matrix and translation vector between the GNSS local coordinate system and the SLAM / odometry pose local coordinate system.

[0037] Step S132: Perform dynamic fine-pair calculation based on SVD decomposition: Constructing SLAM / Oregon pose set and GNSS observation point set First, the two point sets are decentroided to obtain the decentroided SLAM / odot pose point set. and decentroided GNSS observation set Construct the covariance matrix H:

[0038] Perform singular value decomposition (SVD) on matrix H to obtain .

[0039] The optimal rotation matrix R is calculated as follows: In the formula, in the formula, For the centroid-free SLAM / milestone pose set, the first One point, For the centroid-decentrifuged GNSS observation point set, the first... One point, for The transpose of . It is a left singular vector matrix. It is a right singular vector matrix. It is a diagonal matrix, representing and Scaling factor between corresponding directions, and All are orthogonal matrices. for The transpose of .

[0040] If the calculated determinant |R| < 0, then correct the sign of the last column of V. Let the GNSS observation point set be... SLAM / otometry pose set is The centers of mass of the two are respectively and After obtaining the rotation matrix R using SVD, the translation vector is finally updated. .

[0041] Step S133, Secondary state transition: The calculated high-precision rotation matrix R and translation vector t are solidified as positioning system calibration parameters, the data buffer is cleared to release memory, and the state is switched to calibration completion mode, marking the end of the entire calibration process.

[0042] Step S140, Real-time coordinate unification: After entering calibration completion mode, for each newly input GNSS raw coordinate frame, the positioning system performs a transformation using the latest calibration parameters to obtain aligned GNSS observations. SLAM / Odometry Pose Estimation Value Obtained recursively from SLAM / odometry pose data. The centroid-decentered GNSS observation set is then used. As GNSS observations .

[0043] The calibration parameter transformation is performed based on the following formula:

[0044] In the formula, The GNSS local coordinates before alignment. This is the GNSS observation point after rotation and translation. Through the above steps, the absolute positioning data is mapped to the robot's local control coordinate system in real time. This step ensures the accuracy of the GNSS observation values. With pose estimation value Under the same spatiotemporal reference, it provides a precise geometric basis for the subsequent construction of an independent kinematic dead reckoning reference and the calculation of Mahalanobis distance features.

[0045] This embodiment details the multidimensional feature extraction and robust filtering process in step S100. To address the highly deceptive nature of slow-moving, induced deception, the positioning system does not employ a traditional fixed threshold elimination strategy. Instead, it constructs a Mahalanobis distance-sensitive detection model based on dead reckoning and dynamically adjusts the fusion weights using an asymmetric exponential response mechanism to counter complex deception interference.

[0046] Step S200 specifically includes the following steps: Step S210, SLAM state prediction: Based on the obtained local pose increment, further prediction is performed through the state transition equation to obtain the prior state estimate. .

[0047] The local pose increment is calculated based on the pose estimation values ​​of two adjacent frames, and this local pose increment is used to recursively extrapolate the posterior state estimate of the previous time step to obtain the prior state estimate of the current time step, thus realizing SLAM / odot pose estimation. Specifically:

[0048]

[0049] in, For the first The pose estimate at time t is the first time. The pose estimate at time t. This represents the local pose increment between two adjacent frames. This is the posterior state estimate from the previous time step. This is the prior state estimate for the current moment.

[0050] Step S220: Prediction of the prior state covariance matrix: After obtaining the prior state estimate at the current moment, the prior state covariance matrix at the current moment is calculated based on the posterior state covariance matrix of the previous moment and the noise matrix of the SLAM / odometry pose estimation process:

[0051] in, Let be the prior state covariance matrix at the current moment. Let be the posterior state covariance matrix of the previous time step. This is the noise matrix for the SLAM / odometry pose estimation process.

[0052] Step S300 includes the following steps: Step S310, Filter initialization and relative coordinate construction: During the positioning system startup or filter reset phase, the origin of the coordinate system is first determined, and then a unified local coordinate system is established. The specific process is as follows: Step S311, Origin determination: The positioning system uses the first frame of valid GNSS absolute coordinates received as the reference offset.

[0053] Step S312, Relative coordinate calculation: All subsequent input GNSS observations Subtract this reference offset from all values ​​to convert them into relative coordinates with the first frame as the origin. This step avoids the precision loss caused by large floating-point calculations, while ensuring that the local coordinate systems of GNSS data and SLAM / odometer pose data are on the same order of magnitude.

[0054] Step S313, Threshold setting: Set the target weight exponential decay starting threshold for the Mahalanobis distance squared index. Weight decay rate factor Weight lower bound and the recovery damping coefficient In this embodiment, The preferred value is 3.0. The preferred value is 2.0. The preferred value is 0.001. The preferred value is 0.015. Since the positioning system obtains dead reckoning baseline based on the SLAM front end, the prediction residual is extremely small. Therefore, the extremely low threshold can ensure that detection is triggered in the early stages of deception.

[0055] Step S320, Steady-state physical drift prior iterative update: Step S321, Instantaneous physical drift velocity acquisition: based on GNSS observations Obtain the GNSS pose increment, dynamically calculate the difference between the GNSS pose increment and the local pose increment within the synchronization time period, and divide the difference by the time interval to obtain the instantaneous physical drift velocity.

[0056] Step S322: Steady-state physical drift prior extraction. The instantaneous physical drift velocity is limited by a preset drift velocity boundary. This boundary is determined based on the mobile robot chassis kinematic constraints, the normal GNSS drift range, and the SLAM short-time error range. Instantaneous physical drift velocities exceeding the chassis kinematic physical boundary limits are removed. The filtered instantaneous physical drift velocities are then input into a low-pass filter for iterative updates to extract the robot's steady-state physical drift prior. The positioning system can also be used in unmanned vehicles and other unmanned outdoor equipment.

[0057] It should be noted that the convergence time constant of the low-pass filter is several times the cumulative time of the subsequent fixed-length sliding window, ensuring that it is not contaminated by spoofing signals when learning the steady-state physical drift velocity. The fixed-length sliding window detects the spoofing signal first, triggers the anomaly judgment, and freezes the update of the steady-state physical drift prior.

[0058] Step S330: Two-way detection based on sliding window Mahalanobis distance and single-frame innovation: The positioning system performs the following dual-detection process during each measurement update: Step S331, Signal and Physical Layer Pre-check: First, check the trace and altitude change of the observation quality covariance matrix output by the GNSS receiver. If the trace of the matrix is ​​greater than the quality threshold or the altitude change is greater than the altitude threshold, the positioning system marks it as a positioning anomaly, indicating that the positioning system has been spoofed. In this embodiment, the quality threshold is preferably 0.5, and the altitude threshold is preferably 2.0m.

[0059] Step S332: Detect spoofing based on sliding window Mahalanobis distance: Set a fixed-length sliding window, use the local odometry to estimate the state as a benchmark, and calculate the current GNSS observation value. The formula for calculating the new information vector is as follows:

[0060] in, For the first Frame information vector, For the first The independent local pose estimation reference at time step. The independent local pose estimation reference is obtained recursively from the SLAM / odot pose increment, i.e. , .in, For the first k SLAM / odot pose estimation value at time 10:00. This represents the local pose increment between two adjacent frames.

[0061] Write the innovation vectors from multiple consecutive frames into a fixed-length sliding window, and calculate the cumulative innovation vector within the fixed-length sliding window. The Mahalanobis distance squared index is calculated based on the aforementioned cumulative vector of innovation:

[0062] in, This refers to the Mahalanobis distance squared index for the current period, which is the cumulative residual characteristic across multiple epochs. and These are the components of the innovation accumulation vector in the horizontal and vertical directions, respectively. The noise variance is used as the baseline measurement. The innovation accumulation vector is the sum of the cross-source residual vectors across multiple frames within the sliding window, and is further normalized based on this accumulation residual vector. , which serves as the multi-epoch cumulative residual characteristic and the squared Mahalanobis distance index for the current period. In this embodiment, the actual number of samples within the current sliding window is [number]. The preferred value is 10.

[0063] index Reflects current GNSS observations The degree of deviation from the dead reckoning baseline in the statistical space. When trajectory spoofing occurs, although the position deviation of a single point may be small, its statistical distribution will deviate significantly from the prediction domain of SLAM, resulting in an abnormally large Mahalanobis distance.

[0064] In addition, in some other embodiments, the signal and physical layer pre-detection may also detect only the trace of the covariance matrix. When the trace of the matrix is ​​greater than a quality threshold, the positioning system marks it as a positioning anomaly, indicating that the positioning system has been deceived.

[0065] Step S340, State Machine Determination: ① Normal State Determination: If the state machine is in a normal state before encountering a slow-change traction deception, or in the initial stage of encountering a slow-change traction deception, and the Mahalanobis distance does not exceed the safety threshold, then the state machine is determined to be in a normal state. The conditions for determining the normal state are as follows: ,in, As a safety threshold, this embodiment .

[0066] ② Isolation State Determination: After encountering a slow-changing traction deception for a period of time, if the Mahalanobis distance exceeds a preset safety threshold, the state machine enters an isolation state. The determination condition for the isolation state is as follows: Similarly, the =3.0.

[0067] When a positioning system encounters a gradually changing traction spoofing, the deviation in a single frame is small, but the accumulated Mahalanobis distance will rapidly increase and exceed the preset safety threshold. The positioning system can quickly detect this gradually changing traction spoofing, improving positioning reliability and security.

[0068] ③ Recovery State Determination: When the state machine enters the isolated state, the positioning system continuously monitors the latest GNSS observations. Independent local pose estimation benchmark single-frame information square When the value falls back to When the trace of the covariance matrix output by the synchronously detected GNSS receiver is less than the set quality threshold within the confidence threshold, it indicates that the real signal and dead reckoning information have reached physical consistency, and the state machine transitions from the isolated state to the recovered state. The formula for calculating the square of the single-frame innovation is as follows:

[0069] In the formula, and These represent the residuals of the actual GNSS observation position and the dead reckoning position in the east and north directions, respectively, for the current k frames. The noise variance is used as a benchmark measurement.

[0070] Step S350, Normal State Location Update: When the state machine is in a normal state, the following steps are executed: Step S351: Trust the GNSS signal and continuously update the steady-state physical drift prior. .

[0071] Step S352: Perform normal state equivalent measurements and convert GNSS observations into... As an equivalent measurement result, i.e. ,in, These are equivalent observations.

[0072] Step S353: Obtain normal localization results: Constructing an adaptive expansion matrix based on an exponential adaptive robust model:

[0073] In the formula, The preset reference measurement noise covariance matrix, This is the adaptive expansion matrix for the current period. The confidence weights are the equivalent observations for the current period. As the lower bound of the weights, in the normal state, The equivalent observation is the GNSS observation, i.e. , ,and .

[0074] Since both the state estimate and the equivalent observation in this embodiment represent position quantities in a unified local coordinate system, the observation matrix can be taken as the identity matrix, and the Kalman update formula is:

[0075] in, For the first The periodic Kalman gain matrix, Let be the prior state covariance matrix. For prior state estimation, This serves as a posterior state estimate and is used as the positioning result for the current cycle. Let be the posterior state covariance matrix. It is an identity matrix.

[0076] Step S360, Isolation-state location update: When the state machine is in an isolated state, the following steps are executed: Step S361: The positioning system instantly cuts off GNSS trust, freezes and locks the steady-state physical drift prior, and records the current steady-state physical drift prior as the locked steady-state physical drift prior. .

[0077] Step S362: In the isolated state, based on the independent pure dead reckoning benchmark and the locked steady-state physical drift prior. The virtual observations of a single period are reconstructed and used as equivalent measurement results. The independent pure dead reckoning baseline is the prior state estimate. The equivalent observation expression for the isolated state is as follows:

[0078] in, These are virtual observations. For prior state estimation, For steady-state physical drift priors, This is a single-cycle time interval.

[0079] In the isolated state, the virtual observations and the reference measurement noise covariance matrix are synchronously input into the filter for state estimation updates, maintaining the statistical consistency of the filter covariance matrix.

[0080] Step S363: Construct the benchmark measurement noise covariance matrix based on the exponential adaptive robust model: In isolated mode, the positioning system calculates GNSS observations based on the squared Mahalanobis distance index for the current period. The confidence weights. Let the squared Mahalanobis distance index for the current period be... Then GNSS confidence weight Calculate using the following formula:

[0081] In the formula, This is the starting threshold for the exponential decay of the confidence weights, i.e., the aforementioned safety threshold. This is the weight decay rate factor. This is the lower bound of the weight. In this embodiment, in the isolated state, The preferred value is 3.0. The preferred value is 2.0. The preferred value is 0.001.

[0082] Constructing an adaptive expansion matrix based on confidence weights:

[0083] In the formula, The preset reference measurement noise covariance matrix, Given the adaptive expansion matrix for the current period, in the isolated state, with the squared Mahalanobis distance index... The increase of confidence weight Rapidly decrease, adaptive expansion matrix Increase, making the deceived GNSS observations The impact on filtering updates is reduced. In the isolated state, the original GNSS observations are not directly used as equivalent observations in the state update. The adaptive expansion matrix is ​​mainly used to characterize and suppress the influence of deceived GNSS observations. The actual state update is performed using virtual observations as equivalent observations.

[0084] Robust state update and dynamic compensation in isolated states: The isolated state virtual observations constructed in step S362 Equivalent observations for the current period and equivalent observations and adaptive expansion matrix Substituting into the Kalman update formula, we update the posterior state estimate. At this point, the original GNSS observations... The influence of spoofed GNSS is weakened, and the positioning results are mainly maintained by SLAM / odometry pose estimation and locking steady-state physical drift prior, thereby reducing the impact of spoofed GNSS on the positioning results.

[0085] When deception is detected Rapid decay approaching 0, leading to As it approaches infinity, the Kalman filter gain matrix becomes... The values ​​tend to zero, thus reducing the impact of spoofed GNSS observations on the filter update. In the isolated state, the original GNSS observations are not directly used as equivalent observations in the state update. Instead, a unified Kalman update is performed using virtual observations constructed from prior state estimates and locked steady-state physical drift priors to maintain the continuity of the positioning results. The calculation formula for the isolated state positioning results is the same as that for the normal state.

[0086] Step S370, Recovery State Location Update: When the state machine is in the recovery state, the following steps are executed: Step S371: The positioning system instantly clears the sliding window, breaks the alarm deadlock, and enters the smooth recovery phase, smoothly recovering the GNSS signal. At this time, physical space two-dimensional geometric interpolation is performed for smooth transition, and the confidence weight is adjusted by iteratively updating the formula for smooth transition. The iterative update formula is as follows:

[0087] in, The damping coefficient is... For the first Confidence weights for the period For the first Confidence weights for the period.

[0088] Step S372: In the recovered state, apply the confidence weight to the GNSS observations. With virtual observations Linear weighted fusion is performed to generate fused observations, which are then used as equivalent observations. The expression for the recoverable-state equivalent observation is as follows:

[0089] In the formula, To merge observations. Confidence weights. During iterative updates, the value gradually approaches 1, thus making the equivalent observation value gradually approach the GNSS observation value. .

[0090] In the recovery state, the equivalent observations obtained by linear weighted fusion are input into the state machine to prevent the Kalman filter from oscillating and diverging due to changes in confidence weights during the transition period.

[0091] Step S373: Obtain the recovered state localization result: Forced locking of the adaptive expansion matrix.

[0092] Robust state update and dynamic compensation in isolated states: equivalent observations obtained by linearly weighted fusion of recovery states and adaptive expansion matrix Substitute the values ​​into the Kalman update formula to update the posterior state estimate, and repeat the alternating time and measurement updates. The formula for calculating the recovered state localization result is the same as that for the normal state.

[0093] Please see Figure 2 and Figure 3 , Figure 3 The variation trend of the innovative Mahalanobis distance based on dead reckoning was simulated under a trajectory-driven deception attack scenario, based on the multidimensional deception features, and the corresponding dynamic response curve of the asymmetric robust weight. Figure 3 The asymmetric robust weights in the equation are the confidence weights. To verify the dynamic robustness of this invention under complex interference environments, this example establishes a mobile robot simulation platform based on the Gazebo physics engine. Gazebo is an open-source physics engine. In the simulation, a slow, deceptive traction interference lasting 20 seconds was set from 46 to 66 seconds. The experimental comparison results are shown in the attached figures in the manual. ① Regarding deception detection feature response: like Figure 3 As shown in the upper part, in the early stage of deception injection, the novel Mahalanobis distance feature constructed by the system increases with the accumulation of deviation and exceeds the set detection threshold. This result shows that the feature can effectively respond to covert traction interference and trigger an early warning.

[0094] ② Regarding the asymmetric weighting adjustment mechanism: like Figure 3 As shown in the lower half, when the Mahalanobis distance exceeds the threshold, the measurement confidence weight of the system instantly drops from 1.0 to a minimum value. After the deception interval ends, the weight does not undergo a transient jump, but instead exhibits a smooth damped rise curve, gradually recovering to 1.0. This result demonstrates the execution logic of the asymmetric response mechanism (fast truncation, smooth recovery) of this invention.

[0095] ③ Regarding resistance to damage: like Figure 4 As shown in the lower half, during the attacked gray time period, the absolute positioning error of the standard EKF diverges rapidly; both the comparison algorithm MRKF and the method of this invention cut off the contamination signal. Due to the activation of kinematic compensation, the absolute error of the method of this invention increases slowly during this blind walking phase (the peak value is controlled at about 1.1 meters), maintaining effective tracking of the real trajectory on the two-dimensional trajectory.

[0096] ④ Regarding the smoothing effect of signal recovery: At the moment the deception ends, the absolute error curve of the compared MRKF algorithm shows a sharp, vertical drop, resembling a cliff. In contrast, the error curve of the method in this invention exhibits a continuous and smooth downward trend without any vertical drops, eventually converging gently to a normal level. This method does not experience abrupt positional changes but instead generates a smooth fusion transition curve, effectively eliminating trajectory jumps and underlying oscillations during the recovery period. Example

[0097] This embodiment relates to an anti-trajectory traction spoofing positioning system for urban canyon environments, including: The data acquisition module is used to acquire SLAM and GNSS data, preprocess the SLAM and GNSS data, and obtain pose estimation values. and GNSS observations ; The state prediction module is used to calculate values ​​based on pose. Obtain prior state estimates and prior state covariance matrix ; The positioning update module is used to update GNSS observations. Prior state estimation and prior state covariance matrix Perform state machine determination and equivalent measurement and obtain positioning results; The location update module includes: The pre-inspection unit is used to perform signal and physical layer pre-inspection. The deception detection unit is used to infer values ​​based on pose after the pre-detection is passed. and GNSS observations Obtain multi-epoch cumulative residual features and determine whether trajectory traction deception has occurred based on Mahalanobis distance; The normal state processing unit is used to trust GNSS signals when not being deceived and to continuously update the steady-state physical drift prior. And obtain the location results; An isolated state processing unit is used to cut off GNSS signals when spoofed, based on prior state estimation at the time of signal cutoff. and steady-state physical drift prior Update GNSS observations and obtain positioning results; The recovery state processing unit is used to process GNSS observations based on damped confidence weights when the pre-detection passes and the trajectory traction deception is determined to have ended. With virtual observations The data is then fused, GNSS observations are updated, and positioning results are obtained. Example

[0098] This embodiment relates to a storage medium storing computer program instructions, which, when executed by a processor, implement the aforementioned anti-trajectory traction deception positioning method for urban canyon environments.

[0099] The embodiments described above are merely examples of several implementations of the present invention, and while the descriptions are relatively specific and detailed, they should not be construed as limiting the scope of the present invention. It should be noted that those skilled in the art can make various modifications and improvements without departing from the concept of the present invention, and these modifications and improvements all fall within the scope of protection of the present invention.

Claims

1. A method for anti-trajectory traction deception positioning in urban canyon environments, characterized in that, For application in outdoor unmanned equipment, the following steps are included: Outdoor unmanned equipment receives GNSS signals, acquires SLAM / oddsight pose data and GNSS data, and performs data preprocessing to obtain pose estimation values ​​and GNSS observation values. Based on the pose estimation value, obtain the prior state estimate and the prior state covariance matrix; State machine determination is performed based on GNSS observations, prior state estimates, and prior state covariance matrix. Positioning results are obtained based on prior state estimates, Kalman gain matrix, and equivalent observations. During the state machine determination process, multi-epoch cumulative residual characteristics are obtained based on pose estimation and GNSS observations, and the Mahalanobis distance is used to determine whether trajectory traction deception has occurred. If the determination is negative, the GNSS signal is trusted and the GNSS observations are used as equivalent observations. If the determination is positive, the GNSS signal is cut off and the equivalent observations are updated based on prior state estimation. When trusting GNSS signals, the steady-state physical drift prior is continuously updated. After the GNSS signal is cut off, the system enters the isolation state. Based on the prior state estimate at the time of signal cutoff and the steady-state physical drift prior, the equivalent observations are updated to virtual observations. When the trajectory traction deception ends, the system enters the recovery state. Based on the damped confidence weight, the GNSS observations and virtual observations are fused, and the equivalent observations are updated to the fused observations. In the recovered state, GNSS observations are iteratively updated to fused observations based on the following formula: in, The first and Confidence weights for the period The damping coefficient is... To merge observations, These are GNSS observations. These are virtual observations; The steady-state physical drift prior is iteratively updated through the following steps: In a normal fusion state without deception interference, the difference between the GNSS pose increment and the SLAM pose increment during the synchronization period is dynamically calculated, and the difference is divided by the time interval to obtain the instantaneous physical drift velocity. The instantaneous physical drift velocity is limited by the chassis kinematic physical boundary, and values ​​exceeding the limit are eliminated; The instantaneous physical drift velocity after clipping is input into a low-pass filter for iterative updates, and the steady-state physical drift prior is extracted.

2. The anti-trajectory traction deception positioning method according to claim 1, characterized in that: The multi-epoch cumulative residual characteristic is the squared Mahalanobis distance between the GNSS observation value and the pose estimation value within a fixed-length sliding window. When the squared Mahalanobis distance exceeds a preset safety threshold, it is determined to be due to GNSS deception interference. The formula for calculating the squared Mahalanobis distance is: in, The squared Mahalanobis distance index for the current period. and These are the components of the innovation accumulation vector in the horizontal and vertical directions, respectively. The noise variance is measured as a benchmark. This represents the actual number of samples within the current sliding window.

3. The anti-trajectory traction deception positioning method according to claim 1, characterized in that: In the isolated state, GNSS observations are updated to virtual observations based on the following formula: in, For prior state estimation, For the steady-state physical drift prior when the isolated state is triggered, For a single cycle.

4. The anti-trajectory traction deception positioning method according to claim 3, characterized in that, During the process of obtaining location results: Constructing an adaptive expansion matrix based on an exponential adaptive robust model: in, The preset reference measurement noise covariance matrix, It is an adaptive expansion matrix; When trusting GNSS signals, =1; When the GNSS signal is cut off, ; in, The squared Mahalanobis distance index for the current period. As a safety threshold, This is the weight decay rate factor. This is the lower bound of the weight.

5. The anti-trajectory traction deception positioning method according to claim 1, characterized in that: The preprocessing includes dynamic alignment between the global GNSS coordinate system and the local SLAM coordinate system, which includes: Coarse alignment stage: Acquire synchronization frame sequences, calculate the geometric centroid difference between the GNSS observation point set and the SLAM inferred point set, and extract the initial translation vector; Fine alignment stage: When the movement trajectory of the outdoor unmanned equipment reaches the preset kinematic excitation threshold, the collected GNSS observation point set and SLAM inferred point set are decentrifuged, the covariance matrix after decentrifugation is decomposed into singular values, the rotation matrix and translation vector are solved, and the coordinate system is aligned.

6. A trajectory traction spoofing positioning system for urban canyon environments, based on the method of any one of claims 1-5, characterized in that, include: The data acquisition module is used to receive GNSS signals, SLAM / odometer pose data and GNSS data and perform data preprocessing to obtain pose estimation values ​​and GNSS observation values; The state prediction module is used to obtain the prior state estimate and the prior state covariance matrix based on the pose calculation value. The positioning update module performs state machine determination based on GNSS observations, prior state estimates, and prior state covariance matrix, and obtains positioning results based on prior state estimates, Kalman gain matrix, and equivalent observations. The location update module includes: The deception detection unit is used to obtain multi-epoch cumulative residual characteristics based on pose estimation values ​​and GNSS observation values, and to determine whether the trajectory traction deception has occurred based on Mahalanobis distance; The normal state processing unit is used to trust GNSS signals when not being deceived, continuously update the steady-state physical drift prior, use GNSS observations as equivalent observations, and obtain positioning results; The isolated state processing unit is used to cut off the GNSS signal when it is deceived, and update the equivalent observations based on the prior state estimate at the time of signal cutoff and the prior steady-state physical drift, and obtain the positioning results; The recovery state processing unit is used to fuse GNSS observations and virtual observations based on damped confidence weights when the pre-detection passes and the trajectory traction deception is determined to have ended, update the equivalent observations, and obtain the positioning results.

7. A storage medium storing computer program instructions thereon, characterized in that, When the computer program instructions are executed by the processor, the method described in any one of claims 1-5 is implemented.

Images