US 10,768,312 B1 · Integrity analysis method based on kinematic-to-kinematic relative positioning scenario

Abstract

The present invention provides an integrity analysis method based on a kinematic-to-kinematic relative positioning scenario, including the following steps: a) establishing a kinematic-to-kinematic relative positioning model, and inputting navigation data; b) calculating a float solution of an integer ambiguity; c) detecting and correcting cycle slips based on a total electron content rate; d) calculating a probability of correct fix and a probability of incorrect fix for the integer ambiguity; e) determining a fault to be detected and a satellite fault probability; 0 calculating a standard deviation δ_(v|CF) and a position domain deviation b_m; and g) calculating an integrity risk value of a carrier phase. The present invention provides an integer ambiguity calculation algorithm for a kinematic-to-kinematic positioning system in the case of a long baseline, to calculate carrier phase integrity.

An integrity analysis method based on a kinematic-to-kinematic relative positioning scenario, wherein, the analysis method comprises the following steps:

a) establishing a kinematic-to-kinematic relative positioning model for a navigation system, and inputting navigation data;
b) calculating a float solution of an integer ambiguity;
c) detecting and correcting cycle slips based on a total electron content rate, comprising:
step 1:
performing linear combination of wide-lane phase minus narrow-lane pseudorange and detecting an ionospheric total electron content rate for the foregoing solution data respectively:
[
Δ
⁢
N
1
P
⁡
(
k
)
-
Δ
⁢
N
2
P
⁡
(
k
)
]
=
N
WL
P
⁡
(
k
-
1
)
-
N
WL
P
(
⁢
k
)
=
[
⁢
Φ
1
P
⁡
(
k
-
1
)
-
Φ
2
P
⁡
(
k
-
1
)
-
f
1
⁢
P
1
P
⁡
(
k
-
1
)
+
f
2
⁢
P
2
P
⁡
(
k
-
1
)
λ
WL
⁡
(
f
1
+
f
2
)
]
-
[
⁢
Φ
1
P
⁡
(
k
)
-
Φ
2
P
⁡
(
k
)
-
⁢
f
1
⁢
P
1
P
⁡
(
k
)
+
f
2
⁢
P
2
P
⁡
(
k
)
λ
WL
⁡
(
f
1
+
f
2
)
]
;
⁢
[
λ
1
⁢
Δ
⁢
N
1
P
⁡
(
k
)
-
λ
2
⁢
Δ
⁢
N
2
P
⁡
(
k
)
]
=
4
⁢
0
.
3
×
1
⁢
0
16
⁢
(
γ
-
1
)
⁢
Δ
⁢
⁢
t
·
TECR
Φ
⁡
(
k
)
f
1
2
-
λ
1
⁡
[
Φ
1
P
⁡
(
k
)
-
Φ
1
P
⁡
(
k
-
1
)
]
+
λ
2
⁡
[
Φ
2
P
⁡
(
k
)
-
Φ
2
P
⁡
(
k
-
1
)
]
;
where f1, f2 represent two frequencies of a received double-frequency satellite signal;
λ1, λ2 represent wavelengths at corresponding frequencies;
γ
=
f
1
2
f
2
2
;
λ
WL
=
c
(
f
1
-
f
2
)
;
TECRΦ(k) represents an ionospheric total electron content rate at an epoch K;
Φ1P(k) and Φ1P(k−1) represent carrier phase measurement values between two adjacent epochs at the frequency f1;
P1P(k) and P1P(k−1) represent pseudorange measurement values between two adjacent epochs at the frequency f1;
ΔN1P(k) and ΔN2P(k) represent cycle slips that are detected at the frequencies f1 and f2 using a cycle as a unit respectively;
c represents a light velocity;
NWLP(k) and NWLP(k−1) represent integer ambiguities between two adjacent epochs;
step 2:
calculating influence of a pseudo random error on cycle slip detection;
SPIEmax=max{|δP1P(k)|,|δ2P(k)|},δP1P(k),δP2P(k)∈(−∞,+∞);
where δP1P(k) and δP2P(k) represent pseudorange measurements at f1 and f2 at the epoch K respectively, and SPIEmax represents a maximum value of a pseudo random error;
step 3:
calculating cycle slip search areas;
SSΔN1P(k):ΔN1P(k)−δΔN1P(k)≤ΔN1P(k)cand≤ΔN1P(k)+δΔN1P(k);
SSΔN2P(k):ΔN2P(k)−δΔN2P(k)≤ΔN2P(k)cand≤ΔN2P(k)+δΔN2P(k);
where ΔN1P(k)cand represents a correct cycle slip search area at L1, and ΔN2P(k)cand represents a correct cycle slip search area at L2;
step 4:
calculating a wide-lane cycle slip variation by using a cycle slip search rule 1;
and calculating a total electron content rate threshold by using a cycle slip search rule 2;
TECRΦ(k)=H{λ1[Φ1P(k)−Φ1P(k−1)]+λ2[Φ2P(k−1)−Φ2P(k)]+[λ1ΔN1P(k)cand−λ2ΔN2P(k)cand]};
wherein
H
=
f
1
2
4
⁢
0
.
3
×
1
⁢
0
16
⁢
(
γ
-
1
)
⁢
Δ
⁢
t
;
d) calculating a probability of correct fix PCF and a probability of incorrect fix PIF for the integer ambiguity;
calculating a probability of a candidate ambiguity value Z:
P
⁡
(
a
⋓
=
z
)
=
∏
i
=
1
m
⁢
[
Φ
⁡
(
1
-
2
⁢
I
i
T
⁡
(
a
-
z
)
2
⁢
σ
i
/
I
)
+
Φ
⁡
(
1
+
2
⁢
I
i
T
⁡
(
a
-
z
)
2
⁢
σ
i
/
I
)
-
1
]
;
where ă represents an integer ambiguity estimated using a bootstrap method;
a represents a correctly fixed integer ambiguity;
σi/I represents a standard deviation of ambiguities, and the standard deviation depends on a set I={1, 2, . . . , i−1};
Ii represents the ith column vector of a unit-lower-triangular matrix L−T generated from decomposition of LDLT of a fuzzy covariance matrix;
m represents the number of fixed integer ambiguities, and
Φ
⁡
(
x
)
=
∫
-
∞
x
⁢
(
1
/
2
⁢
π
)
⁢
exp
⁡
(
-
1
2
⁢
v
2
)
⁢
dv
;
a−z represents a difference obtained by subtracting an integer ambiguity calculated using a correlation method from an integer ambiguity calculated using a bootstrap method;
when the correctly fixed integer ambiguity a is equal to the candidate ambiguity value Z, the probability of correct fix PCF is calculated by:
P
⁡
(
a
⋓
=
a
)
=
∏
i
=
1
m
⁢
[
2
⁢
Φ
⁡
(
1
2
⁢
σ
i
/
I
)
-
1
]
;
the probability of incorrect fix PIF is calculated by using an exclusion relation between a correct fix event and an incorrect fix event;
e) determining a fault to be detected and a satellite fault probability;
f) calculating a standard deviation δv|CF and a position domain deviation bm;
wherein the following models are established:
y=Gx+ε and
{
x
^
=
Sy
S
=
(
G
T
⁢
WG
)
-
1
⁢
GW
;
and an ambiguity float solution-relative position vector matrix Ŝ and a related estimated error covariance matrix is obtained according to the foregoing equations:
S
^
=
[
x
^
N
^
]
;
P
=
[
P
x
^
P
x
^
⁢
N
^
P
N
^
⁢
x
^
P
N
^
]
;
where G represents a coefficient matrix;
ϵ represents a total satellite error;
{circumflex over (x)} represents an estimated value in a solution state;
W represents a relation matrix between a pseudorange and an error;
{circumflex over (N)} represents an estimated value of an integer ambiguity;
S represents a state vector;
P represents a covariance matrix of position and float ambiguity state estimated error;
step 1:
calculating the kth integer ambiguity element Nz(k);
Nz(k)=round([01×3Zk]S(k-1));
where Zk represents the kth row of a transition matrix Z;
Z1:k=[01×(k-1) 1 01×(n-k)];
S(k-1) represents an updated state vector S based on fixed ambiguities 1→K;
step 2:
calculating the standard deviation δv|CF;
continuously updating position vector estimation and remaining ambiguity real values S;
K=P(k-1)[01×3Zk]T([01×3Zk]TP(k-1)[01×3Zk]T)−1;
S(k)=S(k-1)+K(NZ(k)−[01×3Zk]TS(k-1));
P(k)=(I−K[01×3Zk])P(k-1);
where K represents a gain;
P represents a covariance matrix of position and float ambiguity state estimated error;
P(k) represents a covariance matrix of an updated state vector after the kth ambiguity is fixed;
and it is assumed that δv|CF=√{square root over (P(3,3)(m))} after m iterations;
step 3:
calculating the position domain deviation bm;
K=P[0m×3Z1:m]T([0m×3Z1:m]P[0m×3Z1:m]T)−1;
S(m)=Ŝ+K(NZ(1:m)−[0m×3Z1:m]TŜ);
where Z1:m represents a transition matrix Z of a row vector 1 to m;
Z1:m=[Im×m 0m×(n-m)];
NZ(1:m) represents a fixed solution for fix of 1 to m integer ambiguities;
the position domain deviation bm caused by a relative position error of the mth incorrect candidate pair (a−z)m is:
bm=SCF(m)−SIF(m)=K(a−z)m;
where am represents a correctly fixed ambiguity vector;
zm represents an incorrectly fixed ambiguity vectors;
SCF(m) represents a correct state obtained after fix of m ambiguities;
SIF(m) represents an incorrect state obtained after fix of m ambiguities;
SCF(m)=Ŝ+K(am−[0m×3Z1:m]Ŝ);
SIF(m)=Ŝ+K(zm−[0m×3Z1:m]Ŝ);
g) calculating an integrity risk value IH0 of a carrier phase for the GNSS;
I
H
⁢
0
=
1
-
(
1
-
P
⁢
{

x
^
v
-
x
v

>
VAL
|
C
⁢
F
}
)
⁢
PCF
-
∑
n
=
1
k
⁢
(
1
-
P
⁢
{

x
^
v
-
x
v

>
VAL
|
IF
n
}
)
⁢
PIF
n
;
where PVALCF=P{|{circumflex over (x)}v−xv|>VAL|CF} represents a probability that a vertical error is greater than an alarm threshold under the condition of correct fix;
PCF represents a probability of correct fix;
PIFn represents a probability of incorrect fix;
PVALIF=P{|{circumflex over (x)}v−xv|>VAL|IFn} represents a probability that a vertical error is greater than an alarm threshold under the condition of incorrect fix;
and
h) applying the integrity risk value IH0 in the navigation system to improve position precision of the navigation system on the kinematic-to-kinematic relative positioning scenario by comparing the integrity risk value IH0 to a zero-fault integrity value IH0req, wherein a geometric mechanism of the navigation system is considered to be available when the integrity risk value IH0 is less than the zero-fault integrity value.

635 more elements

Integrity analysis method based on kinematic-to-kinematic relative positioning scenario

Life of this patent

Integrity analysis method based on kinematic-to-kinematic relative positioning scenario

Life of this patent

Abstract

Claims11