A Multi-Fault Detection and Exclusion Method for GPS / BDS / INS Tightly Coupled Navigation Based on Inter-Satellite Differential

Abstract

The present invention discloses a multi-fault detection and rejection method for GPS / BDS / INS tightly-coupled navigation based on inter-satellite difference, comprising the following steps: Step 1, initialize to obtain the inter-system bias between GPS and BDS; Step 2, establish a satellite / inertial tightly-coupled model based on inter-satellite difference; Step 3, perform the first fault detection and rejection on the pseudorange observables including the reference satellite based on the innovation; Step 4, after filtering, perform the second fault detection and rejection based on the residual value and output the final positioning result. The method of the present invention can handle the faults of multi-satellite observables including the reference satellite and maximize the utilization of correct satellite observables.

Description

Technical Field

[0001] The present invention belongs to the technical field of satellite positioning and navigation, and particularly relates to a multi-fault detection and rejection method for GPS / BDS / INS tightly coupled navigation based on inter-satellite differential orientation. Background Technique

[0002] Satellite / inertial integrated navigation is the main positioning means for current outdoor navigation. Among them, the tightly coupled navigation technology is an integrated navigation technology at the satellite observation quantity end. Compared with the loosely coupled one, it has higher positioning accuracy and can work when the number of satellites is less than 4. With the establishment and improvement of satellite navigation systems in various countries, the number of on-orbit navigation satellites is increasing continuously, and more observation quantities can be applied to tightly coupled navigation, which will greatly improve the positioning accuracy of integrated navigation. However, it also increases the probability of simultaneous failures of multiple satellite observation quantities. Failures are often caused by the influence of multi-path effects, artificial interference, spoofing, random software and hardware failures of satellite navigation, etc. Especially in urban environments, high-rise buildings stand, and occlusion and multi-path effects are particularly obvious. These failures are usually unpredictable and cannot be modeled. Therefore, real-time multi-observation quantity fault detection and rejection are of great significance for inertial / satellite tightly coupled navigation systems, which can ensure the accuracy and reliability of navigation and positioning.

[0003] Currently, the satellite information fault detection algorithms applied to integrated navigation systems include Receiver Autonomous Integrity Monitoring (RAIM) and Receiver Autonomous Integrity Extrapolation (AIME). These methods are highly sensitive to single faults but cannot effectively detect multiple faults. Multi-Hypothesis Solution Separation (MHSS) separates different faults into different subsets through combinatorial operations on all satellites, and constructs a test statistic based on the difference between the subset solution and the full-set solution. It is an effective multi-fault detection method and is applied to Advanced Receiver Integrity Monitoring (ARAIM). However, when the number of observation quantities or faults to be monitored is large, it is prone to the situation of computational explosion. At the same time, the above methods cannot effectively distinguish between reference observation quantity faults and non-reference observation quantity faults when facing inter-satellite differential satellite / inertial tightly coupled navigation observation quantity faults. Summary of the Invention

[0004] In order to solve the problems in the prior art, the present invention provides a multi-fault detection and rejection method for GPS / BDS / INS tightly coupled navigation based on inter-satellite differential orientation. This method can handle multi-satellite observation quantity faults including reference satellites and maximize the utilization of correct satellite observation quantities.

[0005] To achieve the above object, the technical solution adopted by the present invention is as follows:

[0006] A multi-fault detection and rejection method for GPS / BDS / INS tightly coupled navigation based on inter-satellite differential orientation, comprising the following steps:

[0007] Step 1, initialize and obtain the inter-system bias between GPS and BDS;

[0008] Step 2, establish a satellite / inertial tight integration model based on inter-satellite differencing;

[0009] Establishment of the measurement equation in the satellite / inertial tight integration model based on inter-satellite differencing in Step 2:

[0010] Step 21, correct the satellite pseudorange observations obtained at the current epoch for satellite clock error, ionospheric error, and tropospheric error, and subtract the inter-system bias from the pseudorange observations in the BDS system. The final obtained quantity is the pseudorange observation value used for subsequent solution;

[0011] Step 22, select the satellite with the highest satellite elevation angle as the reference satellite L in the inter-satellite differencing, and calculate the satellite pseudorange inter-satellite difference value of the i-th satellite and the inertial navigation predicted pseudorange inter-satellite difference value

[0012]

[0013] where: and are the pseudorange values of the i-th satellite and the reference satellite L, and are the predicted pseudoranges of the i-th satellite and the reference satellite L solved by inertial navigation, and their values are the geometric distances between the corresponding satellites and the receiver. Assuming the total number of satellites is N, the measurement quantity Z k and the measurement matrix H k / k-1 are respectively:

[0014]

[0015] where: are the inertial navigation predicted pseudorange inter-satellite difference values of the 1st, i-th, and (N - 1)-th satellites at epoch k, are the pseudorange inter-satellite difference values of the 1st, i-th, and (N - 1)-th satellites at epoch k; correspond to the 1st, i-th, and (N - 1)-th satellites respectively, and are defined as the difference in the unit observation vectors between each satellite and the reference satellite, that is and are the unit observation vectors of the i-th satellite and the reference satellite L respectively, is the rotation matrix from the navigation system n to the earth system e;

[0016] Furthermore, Step 3, perform the first fault detection and rejection on the pseudorange observations including the reference satellite based on the innovation quantity, including the following steps:

[0017] Step 31, calculate the innovation quantity Υ k , the innovation quantity Υ k is defined as the difference between the measurement and the predicted measurement, that is:

[0018]

[0019] where: X k / k-1 = F k / k-1 X k-1 is the one-step predicted value of the filter state at the k-th epoch, F k / k-1 is the state transition matrix, X k-1 is the state quantity of the previous epoch. Since in closed-loop filtering, previous errors will be continuously estimated and corrected, so X k-1 is 0, then X k / k-1 = 0, that is, in closed-loop tight integration, the innovation quantity Υ k is equal to the measurement Z k ;

[0020] Step 32, calculate the first-step fault detection quantity F Υ,k , F Υ,k consists of N - 1 elements, that is T is the transpose symbol, and are the fault detection quantities corresponding to the 1st, the i-th, and the (N - 1)-th satellites respectively, and are defined as:

[0021]

[0022] where: is the innovation quantity of the i-th satellite, equal to the i-th element in Υ k , is the variance value of the innovation quantity of the i-th satellite, equal to the diagonal element located in the i-th row and i-th column of the innovation quantity covariance matrix ∑ Υ,k , ∑ Υ,k The calculation method is:

[0023]

[0024] where: P k / k-1 is the one-step prediction covariance matrix of the state quantity, R k is the covariance matrix of the measurement noise;

[0025] Step 33, analyze the distribution law of F Υ,k , screen and record the values in F Υ,k that conform to the standard normal distribution in the set A0, that is And count the number of element values that conform to the standard normal distribution, denoted as m0, i.e., m0 = length(A0);

[0026] Step 34, screen record F Υ,k The values that conform to the non-zero mean normal distribution in it are recorded in set a1, i.e., μ is the mean of the normal distribution, μ≠0. Count the number of element values that conform to the non-zero mean normal distribution, denoted as m1, i.e., m1 = length(A1); Since the value of μ is unknown, only when there are no less than three elements that all satisfy the non-zero mean distribution can judgment and recording be carried out, i.e., m1≥3 or m1 = 0;

[0027] Step 35, when m0≥2, that is, there are more than two innovation measurement quantities that satisfy the standard normal distribution. At this time, it is considered that the reference star pseudorange observable is fault-free, and go to Step 36; when m0<2, consider the reference star pseudorange observable and go to Step 37;

[0028] Step 36, the non-reference star pseudorange observables corresponding to all test statistics that do not satisfy the standard normal distribution are detected as faulty and excluded, and go to Step 310;

[0029] Step 37, the reference star pseudorange observable is faulty and the reference star needs to be replaced: when m1≥3, there are no less than 3 correct non-reference star pseudorange observation values, and go to Step 38; when m1<3, it is impossible to judge the non-faulty non-reference star pseudorange observation values, and go to Step 39;

[0030] Step 38, remove the original reference star, arbitrarily select a satellite in set A1 as the reference star, and the remaining satellites in A1 as non-reference stars, and calculate the measurement Z k and the measurement matrix H k / k-1 and the covariance matrix R of the measurement noise k , and go to Step 310;

[0031] Step 39, remove the original reference star, arbitrarily select one of the remaining satellites as the reference star, record the number of times K of replacing the reference star. When K = N - 2, go to Step 311, otherwise go to Step 22;

[0032] Step 310, perform extended Kalman filtering, output the positioning result, and the first-step fault detection scheme based on the innovation ends;

[0033] Step 311, output the inertial navigation recursive result as the positioning result, and the fault detection scheme based on the innovation ends.

[0034] In the said Step 1, in the initialization stage, the prior inter-system deviation between GPS and BDS is obtained before the system officially starts working. The method for solving the inter-system deviation at epoch k is:

[0035] δ tτ,k = δt uG,k -δt uB,k (1)

[0036] Where: δt uG,k and δt uB,k are the GPS receiver clock error and BDS receiver clock error obtained by pseudorange single-point positioning solution at epoch k respectively, and δt τ,k is the inter-system bias at the k-th epoch. The mean value of the inter-system bias obtained in the initialization stage is the prior inter-system bias δt τ :

[0037]

[0038] Where: j represents an epoch from epoch 1 to k, and δt τ,j is the inter-system bias at the j-th epoch.

[0039] Furthermore, in step 4, after filtering, perform a second fault detection and rejection based on the residual value and output the final positioning result;

[0040] Step 41, calculate the measurement residual s k :

[0041] s k = Z k - H k X k / k (10)

[0042] Where: X k / k is the posterior estimation of the state quantity, Z k and H k are the measurement quantity and the measurement equation after the first fault detection,

[0043] Step 42, construct a test statistic λ based on the residual k :

[0044]

[0045] Where: is the transpose matrix of s k ; is the n-th row element of s k ; is the element in the n-th row and n-th column of the residual covariance matrix, κ is the dimension of the measurement quantity and the measurement matrix after the first fault detection and rejection. If no fault occurs, E() is the mathematical expectation function, and χ 2 (κ) is the central chi-square distribution with k degrees of freedom. λ k follows the central chi-square distribution with k degrees of freedom, that is, λk ~X 2 (k); If a fault occurs, λ k follows a non - central chi - square distribution, that is is the non - centrality parameter;

[0046] Step 43, perform residual - based fault detection:

[0047]

[0048] where: F d is the detection threshold calculated in real - time, calculated according to the preset false - alarm rate P f = 1 / 1000:

[0049]

[0050] where: P(λ k <F d ) represents the probability of λ k <F d , and x is an arbitrary variable;

[0051] When the second - step fault detection scheme based on residuals fails to pass, it is considered that the first - step fault detection scheme based on innovations is invalid. At this time, the positioning will no longer rely on satellites, but directly output the inertial navigation recurrence result; if the second - step test statistic based on residuals passes, the first - step fault detection statistic scheme based on innovations is effective, and the integrated navigation filtering result is output.

[0052] Compared with the prior art, the present invention has the following beneficial effects:

[0053] Through inter - satellite differencing, the GPS / BDS / INS tight integration in the present invention no longer needs to estimate the receiver clock error, reduces the state variables to be estimated, avoids the problem of inaccurate receiver clock error modeling, and improves the positioning accuracy;

[0054] In the present invention, the two - step fault detection method is carried out before and after filtering respectively, which reduces the possibility of missed detection to the greatest extent and ensures the reliability of positioning;

[0055] In the first - step fault detection method of the present invention, a new fault detection strategy is proposed according to the distribution of the detection quantity, which can detect faults of multiple observation quantities including the reference satellite observation quantity, maximize the utilization of correct satellite observation quantities, and thus ensure the accuracy and stability of positioning. Brief Description of the Drawings

[0056] Figure 1 is the algorithm structure diagram of the present invention;

[0057] Figure 2It is a comparison diagram of the positioning error between the fault detection and rejection scheme of the present invention and the scheme without fault detection under step error;

[0058] Figure 3 It is a comparison diagram of the positioning error between the fault detection and rejection scheme of the present invention and the scheme without fault detection under ramp error. Detailed implementation manners

[0059] The present invention will be further described below in conjunction with embodiments.

[0060] Embodiment 1

[0061] As Figure 1 shown, a multi-fault detection and rejection method based on GPS / BDS / INS tight integration navigation oriented to inter-satellite difference includes the following steps:

[0062] Step 1, initialize to obtain the inter-system bias between GPS and BDS;

[0063] Step 2, establish a satellite / inertial tight integration model based on inter-satellite difference;

[0064] Establishment of the measurement equation in the satellite / inertial tight integration model based on inter-satellite difference in step 2:

[0065] Step 21, correct the satellite clock error, ionospheric error, and tropospheric error for the satellite pseudorange observation values obtained at the current epoch, and subtract the inter-system bias from the pseudorange observation values in the BDS system. The finally obtained quantity is the pseudorange observation value used for subsequent calculation;

[0066] Step 22, select the satellite with the highest satellite elevation angle as the reference satellite L in the inter-satellite difference, and calculate the satellite pseudorange inter-satellite difference value and the inertial navigation predicted pseudorange inter-satellite difference value

[0067]

[0068] where: and are the pseudorange values of the i-th satellite and the reference satellite L respectively, and are the predicted pseudoranges of the i-th satellite and the reference satellite L calculated by inertial navigation respectively. Their values are the geometric distances between the corresponding satellites and the receiver. Assuming the total number of satellites is N, the measurement quantity Z k and the measurement matrix H k / k-1 are respectively:

[0069]

[0070] where: The inertial navigation predicted pseudo-range inter-satellite differential values of the 1st, the i-th, and the (N-1)-th satellites at epoch k, respectively, The pseudo-range inter-satellite differential values of the 1st, the i-th, and the (N-1)-th satellites at epoch k, respectively; Corresponding to the 1st, the i-th, and the (N-1)-th satellites respectively, which is defined as the difference of the unit observation vectors between each satellite and the reference satellite, that is and are the unit observation vectors of the i-th satellite and the reference satellite L respectively, is the rotation matrix from the navigation system n to the earth system e;

[0071] Step 3, perform the first fault detection and rejection on the pseudo-range observations including the reference satellite based on the innovation, including the following steps:

[0072] Step 31, calculate the innovation γ k , the innovation γ k is defined as the difference between the measured value and the predicted measured value, that is:

[0073]

[0074] where: X k / k-1 = F k / k-1 X k-1 is the one-step prediction value of the filter state at the k-th epoch, F k / k-1 is the state transition matrix, X k-1 is the state quantity of the previous epoch. Since in the closed-loop filtering, the previous errors will be continuously estimated and corrected, so X k-1 is 0, then X k / k-1 = 0, that is, in the closed-loop tight integration, the innovation γ k is equal to the measured value Z k ;

[0075] Step 32, calculate the first fault detection quantity F Υ,k , F Υ,k consists of N-1 elements, that is T is the transpose symbol, and are the fault detection quantities corresponding to the 1st, the i-th, and the (N-1)-th satellites respectively, defined as:

[0076]

[0077] where: is the innovation of the i-th satellite, equal to the i-th element in γ k , is the variance value of the innovation of the i-th satellite, equal to the innovation covariance matrix ∑Υ,k The diagonal elements in the i-th row and i-th column of Υ,k The calculation method is:

[0078]

[0079] Where: P k / k-1 is the one-step prediction covariance matrix of the state quantity, R k is the covariance matrix of the measurement noise;

[0080] Step 33, Analyze F Υ,k Distribution law, filter records F Υ,k The values ​​that conform to the standard normal distribution are recorded in the set A0, that is, And count the number of element values ​​that conform to the standard normal distribution, recorded as m0, that is, m0 = length (A0);

[0081] Step 34, filter record F Υ,k The values ​​in that conform to the normal distribution with non-zero mean are recorded in set A1, that is, μ is the mean of the normal distribution, μ≠0, and the number of element values ​​that meet the non-zero mean normal distribution is counted, recorded as m1, that is, m1=length(A1); since the value of μ is unknown, it can only be judged and recorded when there are no less than three elements that meet the non-zero mean distribution, that is, m1≥3 or m1=0;

[0082] Step 35, when m0≥2, that is, there are more than two new information detection quantities that satisfy the standard normal distribution, at this time, it is considered that the reference star pseudorange observation quantity is fault-free, and go to step 36; when m0<2, it is considered that the reference star pseudorange observation quantity is fault-free, and go to step 37;

[0083] Step 36, all non-reference star pseudorange observations corresponding to the test statistics that do not satisfy the standard normal distribution are detected as faulty and excluded, and the process goes to step 310;

[0084] Step 37, the reference star pseudorange observation value fails, and the reference star needs to be replaced: when m1≥3, there are no less than 3 non-reference star pseudorange observation values ​​that are correct, go to step 38; when m1<3, it is impossible to judge the non-faulty non-reference star pseudorange observation value, go to step 39;

[0085] Step 38: Eliminate the original reference star, select any satellite in set A1 as the reference star, and use the remaining satellites in A1 as non-reference stars to calculate the measurement Z k , measurement matrix H k / k-1 and the covariance matrix R of the measurement noise k , and go to step 310;

[0086] Step 39: Exclude the original reference star, select any one of the remaining satellites as the reference star, record the number of times K of replacing the reference star. When K = N - 2, go to Step 311; otherwise, go to Step 22.

[0087] Step 310: Perform extended Kalman filtering and output the positioning result. The first - step fault detection scheme based on innovation ends.

[0088] Step 311: Output the inertial navigation recursion result as the positioning result. The fault detection scheme based on innovation ends.

[0089] In the above - mentioned Step 1, in the initialization stage, the prior inter - system deviation between GPS and BDS is obtained before the system officially starts working. The method for solving the inter - system deviation at epoch k is as follows:

[0090] δ tτ,k = δt uG,k - δt uB,k (1)

[0091] Where: δt uG,k and δt uB,k are the GPS receiver clock error and BDS receiver clock error obtained by pseudo - range single - point positioning calculation at epoch k respectively, and δt τ,k is the inter - system deviation at the k - th epoch. Calculate the mean value of the inter - system deviation obtained in the initialization stage to get the prior inter - system deviation δt τ :

[0092]

[0093] Where: j represents a certain epoch from epoch 1 to k, and δt τ,j is the inter - system deviation at the j - th epoch.

[0094] As Figure 2 shown, in Step 4, after filtering, perform the second - step fault detection and exclusion based on the residual value and output the final positioning result;

[0095] Step 41: Calculate the measurement residual s k :

[0096] s k = Z k - H k X k / k (10)

[0097] Where: X k / k is the posteriori estimation of the state quantity, Z k and H k are the measurement quantity and measurement equation after the first - step fault detection,

[0098] Step 42: Construct a test statistic λ based on the residualk :

[0099]

[0100] Wherein: is the transpose matrix of s k , is the n-th row element of s k , is the element in the n-th row and n-th column of the residual covariance matrix, κ is the dimension of the measurement vector and the measurement matrix after the first fault detection and rejection. If no fault occurs, E() is the mathematical expectation function, χ 2 (κ) is the central chi-square distribution with degree of freedom k, λ k follows the central chi-square distribution with degree of freedom k, that is, λ k ~X 2 (k); if a fault occurs, λ k follows the non-central chi-square distribution, that is is the non-centrality parameter;

[0101] Step 43, perform residual-based fault detection:

[0102]

[0103] Wherein: F d is the detection threshold calculated in real time, and is calculated according to the preset false alarm rate P f = 1 / 1000:

[0104]

[0105] Wherein: P(λ k < F d ) represents the probability of λ k < F d , and x is any variable;

[0106] When the second-step residual-based fault detection scheme fails, it is considered that the first-step innovation-based fault detection scheme is invalid. At this time, the positioning will no longer rely on satellites, but directly output the inertial navigation recurrence result; if the second-step test statistic based on residuals passes, the first-step fault detection statistic scheme based on innovation is effective, and the integrated navigation filtering result is output.

[0107] For the fault detection part: The test statistics are constructed respectively through the innovation and the residual. First, according to the distribution law of the test statistic based on the innovation, the correct non-reference satellite observables are identified, and whether the reference satellite has a fault is detected. If there is a fault, a new reference satellite is selected from the correct non-reference satellites for replacement and the process starts over; if there is no reference satellite fault, the faulty non-reference satellite observables are removed and filtering is performed. After the filtering is completed, the second step of fault detection is carried out. The used satellite observables are detected according to the test statistic based on the residual. If it passes, it proves that the solution detected based on the innovation in the first step is effective, and the filtered positioning result is output; if it does not pass, it means that the solution detected based on the innovation in the first step is invalid, and the inertial navigation recursive result is directly output.

[0108] Embodiment 2

[0109] Simulation description

[0110] The following factors are considered in the fault simulation: fault type (step fault, ramp fault), the number of satellite observables with faults, whether the reference star has a fault, the magnitude of the fault, and the fault duration. In this Embodiment 2, different types of step fault and ramp fault simulations are set respectively, and the specific simulation situations are shown in Table 1 and Table 2 respectively.

[0111] Step fault

[0112] The detailed situation of the step fault type is shown in Table 1. The main aspects of each fault situation are the fault magnitude, the number of faults, and whether there is a reference star observable fault. Among them, the G23 satellite is the reference star in this experiment, that is, cases 2 and 4 contain reference star faults.

[0113] Table 1. Simulation situation of step fault type

[0114]

[0115] Reference Figure 2 , step faults of different situations are added to the observed pseudorange. The positioning error obtained by the present invention and the positioning error of the scheme without fault detection are as Figure 2 shown. The method described in the present invention can effectively identify and eliminate step faults of multiple satellite observables including the reference star, can identify and eliminate step faults as low as 5m, maintain the positioning error within the normal standard, and ensure the integrity and reliability of positioning.

[0116] Ramp fault

[0117] The detailed situation of the ramp fault type is shown in Table 2. The main aspects of each fault situation are the fault slope, the fault duration, the number of faults, and whether there is a reference star observable fault. Among them, the G23 satellite is the reference star in this experiment, that is, cases 2 and 4 contain reference star faults.

[0118] Table 2. Simulation situation of ramp fault type

[0119]

[0120] Reference Figure 3 , different forms of step faults are added to the observed pseudorange. The positioning error obtained by the present invention and the positioning error without the fault detection scheme are shown in the appendix Figure 3 as shown. When fault detection and rejection are not used, the positioning error caused by the ramp fault is huge and the original positioning accuracy cannot be restored in time after the fault disappears. The method described in the present invention can effectively identify and reject the ramp faults of multi-satellite observables including reference stars, significantly improving the positioning accuracy and effectively ensuring the integrity and reliability of positioning.

[0121] The above are only the preferred embodiments of the present invention. It should be pointed out that for those of ordinary skill in the art, without departing from the principle of the present invention, several improvements and refinements can be made, and these improvements and refinements should also be regarded as the protection scope of the present invention.

Claims

1. A multi-fault detection and rejection method for GPS / BDS / INS tightly-coupled navigation based on inter-satellite difference, characterized in that It includes the following steps: Step 1, initialize and obtain the inter-system bias between GPS and BDS; Step 2, establish a satellite / inertial tight integration model based on inter-satellite differences; Step 3, perform the first fault detection and rejection on the pseudorange observations including the reference star based on the innovation; Step 4, after filtering, perform the second fault detection and rejection based on the residual value and output the final positioning result; In the satellite / inertial tight integration model based on inter-satellite differences in Step 2, the establishment of the measurement equation: Step 21, correct the satellite clock error, ionospheric error, and tropospheric error for the satellite pseudorange observations obtained at the current epoch, and subtract the inter-system bias from the pseudorange observations of the BDS system. The finally obtained quantity is the pseudorange observation value used for subsequent solution; Step 22: Select the satellite with the highest satellite elevation angle as the reference satellite L in the inter-satellite difference, and calculate the inter-satellite difference value of the satellite pseudorange of the i-th satellite respectively and the predicted inter-satellite difference value of the inertial navigation pseudorange Wherein: and are the pseudorange values of the i-th satellite and the reference satellite L, and are respectively the predicted pseudoranges of the i-th satellite and the reference satellite L obtained by inertial navigation solution, and their values are the geometric distances between the corresponding satellites and the receiver. Assuming the total number of satellites is N, the measurement Z k and the measurement matrix H k / k-1 are respectively: Wherein: They are respectively the inertial navigation predicted pseudo-range inter-satellite difference values of the 1st, the i-th, and the (N-1)-th satellites at epoch k, They are respectively the pseudo-range inter-satellite difference values of the 1st, the i-th, and the (N-1)-th satellites at epoch k; They respectively correspond to the 1st, the i-th, and the (N-1)-th satellites, and are defined as the difference of the unit observation vectors between each satellite and the reference satellite, that is and They are respectively the unit observation vectors of the i-th satellite and the reference satellite L, is the rotation matrix from the navigation system n to the earth system e; The said Step 3 includes the following steps: Step 31, calculate the innovation quantity Υ k , the innovation quantity Υ k is defined as the difference between the measured quantity and the predicted measured quantity, i.e.: Where: X k / k-1 = F k / k-1 X k-1 is the one-step predicted value of the filter state at the k-th epoch, F k / k-1 is the state transition matrix, X k-1 is the state quantity of the previous epoch. Since in closed-loop filtering, previous errors will be continuously estimated and corrected, so X k-1 is 0, then X k / k-1 = 0, that is, in the closed-loop tight integration, the innovation γ k is equal to the measurement Z k ; Step 32, calculate the first-step fault detection quantity F based on the innovation Υ,k , F Υ,k consists of N - 1 elements, that is T is the transpose symbol, and are the fault detection quantities corresponding to the 1st, the ith, and the (N - 1)th satellites respectively, and are defined as: Wherein: is the innovation of the i-th satellite, equal to the i-th element in γ k , is the variance value of the innovation of the i-th satellite, equal to the diagonal element at the i-th row and i-th column in the innovation covariance matrix ∑ Υ,k , ∑ Υ,k The calculation method is: Where: P k / k-1 is the one-step prediction covariance matrix of the state quantity, and R k is the covariance matrix of the measurement noise; Step 33, analyze F Υ,k for its distribution pattern, and filter and record the values of F Υ,k that conform to the standard normal distribution in set A0, i.e., also count the number of element values that conform to the standard normal distribution, denoted as m0, i.e., m0 = length(A0); Step 34, screen record F Υ,k Values that conform to the non-zero mean normal distribution in are recorded in set A1, that is μ is the mean of the normal distribution, μ≠0. The number of element values that conform to the non-zero mean normal distribution is counted and denoted as m1, that is, m1 = length(A1); Since the value of μ is unknown, only when there are at least three elements that all satisfy the non-zero mean distribution can judgment and recording be carried out, that is, m1≥3 or m1 = 0; Step 35, when m0≥2, that is, there are more than two innovation detection quantities satisfying the standard normal distribution. At this time, it is considered that the reference star pseudorange observation is fault-free, and go to Step 36; when m0<2, consider the reference star pseudorange observation, and go to Step 37; Step 36, the non-reference star pseudorange observations corresponding to all test statistics that do not satisfy the standard normal distribution are detected as faulty and excluded, and go to Step 310; Step 37, the reference star pseudorange observation is faulty and the reference star needs to be replaced: when m1≥3, there are no less than 3 correct non-reference star pseudorange observations, and go to Step 38; when m1<3, it is impossible to judge the non-reference star pseudorange observations without faults, and go to Step 39; Step 38, eliminate the original reference star, arbitrarily select a satellite in set A1 as the reference star, and the remaining satellites in A1 as non-reference stars, and calculate the measurement Z k , the measurement matrix H k / k-1 and the covariance matrix R k of the measurement noise, and go to Step 310; Step 39, eliminate the original reference star, select any one of the remaining satellites as the reference star, record the number of times K of replacing the reference star. When K = N - 2, go to Step 311, otherwise go to Step 22; Step 310, perform extended Kalman filtering, output the positioning result, and the first fault detection scheme based on the innovation ends; Step 311, output the inertial navigation recursion result as the positioning result, and the fault detection scheme based on the innovation ends.

2. The multi-fault detection and rejection method based on GPS / BDS / INS tight integration navigation oriented to inter-satellite differences according to claim 1, characterized in that In the said Step 1, in the initialization stage, the prior inter-system bias between GPS and BDS is obtained before the system officially starts working. The solution method for the inter-system bias at epoch k is: δt τ,k = δt uG,k -δt uB,k (1) where: δt uG,k and δt uB,k are the GPS receiver clock error and BDS receiver clock error obtained by pseudorange single-point positioning solution at epoch k respectively, and δt τ,k is the inter-system bias at the k-th epoch. By averaging the inter-system biases obtained in the initialization stage, the prior inter-system bias δt τ is obtained: where: j represents one of the epochs from epoch 1 to k, and δt τ,j is the inter-system bias at the j-th epoch.

3. The multi-fault detection and rejection method based on GPS / BDS / INS tight integration navigation oriented to inter-satellite differences according to claim 2, characterized in that The said Step 4: After filtering, perform the second fault detection and rejection based on the residual value and output the final positioning result, Step 41, calculate the computational amount measurement residual s k : s k = Z k - H k X k / k (10) Where: X k / k is the posterior estimate of the state quantity, Z k and H k are the measured quantity and the measurement equation after the first fault detection Step 42, construct a residual-based test statistic λ k : Wherein: is the transpose matrix of s k , is the nth row element of s k, is the element at the nth row and nth column of the residual covariance matrix, k is the dimension of the measurement vector and the measurement matrix after the first fault detection and rejection. If no fault occurs, E(·) is the mathematical expectation function, χ 2 (κ) is the central chi-square distribution with degree of freedom κ, λ k follows the central chi-square distribution with degree of freedom κ, that is, λ k ~χ 2 (κ); if a fault occurs, follows the non-central chi-square distribution, that is, λ k ~χ 2 (κ, ζ), where ζ is the non-centrality parameter; Step 43, perform fault detection based on the residual; Where: F d is the detection threshold calculated in real time, and is calculated according to the preset false alarm rate P f = 1 / 1000: Where: P(λ k <F d ) represents the probability of λ k <F d , and x is an arbitrary variable; When the second fault detection scheme based on the residual cannot pass, it is considered that the first fault detection scheme based on the innovation is invalid. At this time, the positioning will no longer rely on satellites, but directly output the inertial navigation recursion result; if the second test statistic based on the residual passes, the first fault detection statistic scheme based on the innovation is valid, and the combined navigation filtering result is output.

Images