Characterizing BDS signal-in-space performance from integrity perspective

2.1 Datasets and data cleaning

This study is carried out by employing the historical broadcast and precise ephemerides of BDS. The performance evaluation for BDS-2 is conducted over a 4.5-year period from 2016 to mid-2020 (June 21). And the analysis for BDS-3 only involves the recent 1.5-year data starting in 2019 due to the limited availability of BDS-3 precise products.

International GNSS Service (IGS) Multi-GNSS Experiment (MGEX) provides daily combined broadcast ephemerides (BCEs) by merging observed data from global GNSS tracking stations. This product, called BRDM, shows satisfactory continuity over the analysis period and thus is adopted as the primary BCE data source in this study. Besides, we utilize another two BCE products for data augmentation: the product from the Test and Assessment Research Center (TARC) of CSNO (available at: http://www.csno-tarc.cn/en/support/downloads) and the product offered by Innovation Academy for Precision Measurement Science and Technology (APM), Chinese Academy of Sciences (CAS) (available at: ftp://159.226.162.31).

A sanity check is an important preprocessing step to remove data logging errors from the original records, which was implemented through a majority voting approach (Heng et al., 2012). However, this approach may be biased because the error patterns are highly sensitive to the hardware and software used in the stations (Montenbruck et al., 2018). In this work, the BRDM product is mostly employed for the analysis, which can overcome this shortcoming by adopting the following data cleaning strategy (O. Montenbruck, personal communication, July 10, 2019): first, it performs the least significant bit recovery for converting back to broadcast bits (Heng et al., 2012) and then removes the incorrect numbers due to rounding and truncation; second, the incomplete records are isolated; finally, different records at the same epoch are all retained in the final file if their correctness cannot be unambiguously judged.

To perform a reliable SISRE evaluation, we also need to identify which ephemeris is actually available to the navigation users at a specific epoch. For this purpose, the BCE products provide the transmission time of message (𝑡_tm) to indicate the time when an ephemeris was first received by the tracking stations. With this parameter, the ephemeris selection process was implemented as follows (Montenbruck et al., 2018): for a given epoch 𝑡, one should select the record that minimizes (𝑡_tm −𝑡) with the condition that (𝑡_tm −𝑡) ≥ 0. However, this strategy is not well suited for integrity evaluation because it does not take incorrect records into account. In this work, we propose a new strategy by applying the following criteria (CSNO, 2019a):

1. Normally, the orbit parameters are updated at the start of every BDT hour, and the reference time of ephemeris (𝑡_oe) is an integral point. The value of 𝑡_oe shall change if any of the orbit parameters change. Under abnormal conditions, the parameters may be updated at nonintegral points, and the value of 𝑡_oe will accordingly change. These rules also apply to the clock correction parameters.

2. The health flag of a satellite can be updated at any time, and this will not influence the values of 𝑡_oe and 𝑡_oc (the reference time of clock parameters).

Based on the criteria above, one can judge whether there are incorrect records at a given epoch. If it proves that wrong records do exist, we then use a posterior approach to identify the correct message: multiple sets of orbit and clock parameters are separately derived using different records, and the record that is most consistent with the precise product is regarded as the correct one. This method is initially acceptable because it is very unlikely that a faulted ephemeris (i.e., the satellite broadcast inaccurate information) coincidently becomes consistent with the precise one after being disturbed by random recording errors. However, please note, this approach may lead to a selection bias that might hide true faults. A reliable approach should not depend on precise ephemerides, which will be further investigated in the future.

When selecting the ephemeris, we also take the applicability interval into consideration by applying the following rule: the ephemeris for use should be the one that minimizes (𝑡_oe −𝑡) subject to 0≤(𝑡_oe −𝑡) ≤1 hour. An epoch should be labeled as “no BCE” if there is no correct record satisfying this condition. To reveal the possible causes of the missing BCEs and to prove the rationality of the rule above, we provide some further analyses as follows.

Normally, the broadcast ephemeris is updated at the start of every BDT hour, while it is missing for a few hours in the dataset. The absences may be caused by (a) the satellite encountered a signal outage (i.e., stopped signal transmission), (b) the tracking stations failed to record that message or the BCE products failed to restore it, or (c) the navigation message was not normally updated. To preliminarily identify the cause of each missing message, we employ the raw observation files from the MGEX stations for judging the pseudorange availability when the ephemeris was absent. The result suggests that most of the missing BCEs are likely caused by satellite signal outages because the pseudoranges were also unavailable over those periods. However, this is not true for a few absences occurring before mid-2017. At the current stage, we cannot unambiguously find the reason behind these absences, and we will conduct a more comprehensive analysis to account for every absence in the future.

In this study, the application interval of BDS ephemerides is set to one hour. This is justified because, as mentioned above, the missing BCEs are mostly caused by signal outages. Besides, it is abnormal that the ephemeris is not updated at the start of a BDT hour (CSNO, 2019a). Therefore, it is suggested that the expired ephemeris should not be used in safety-critical applications. This will potentially reduce the integrity risk at the cost of a slight increase in continuity risk.

The BDS precise ephemerides (PCEs) are routinely generated by several MGEX Analysis Centers (ACs) (Montenbruck et al., 2017). For the period of interest, two ACs continuously provide the final PCEs of all the BDS-2 satellites: WHU (Wuhan University) and GFZ (Geo-ForschungsZentrum) (IGS MGEX, 2020a). In this study, WHU products are mostly employed for the BDS-2 SIS performance evaluation, and GFZ products are regarded as a backup when the former has short breaks. As for BDS-3, only WHU can continuously provide the final PCEs over the analysis period. In this work, the sample interval of the precise products is 15 minutes. For an epoch, the PCE product for a satellite is deemed as invalid if the associated precise orbit or clock information is absent or invalid. The epochs when the two products are simultaneously invalid are labeled as “no PCE.” Appendix A gives some further analyses on the missing PCEs, including their causes and the potential effects on the evaluation.

Figures 1 and 2 show the data status for individual satellites of BDS-2 and BDS-3. Due to the significant volume of missing data, the integrity performance of the satellites C03/C17, C16, C18, and C31-C37 is not evaluated in this study. The absence of precise products may be caused by the fact that the MGEX observations are likely to be limited for new satellites. And the continuous absence of BCEs for some satellites (e.g., C03/C17) may be due to one of the following reasons: (a) the MGEX stations may not track the newly launched satellites and those experimental satellites (e.g., C17, C18, C31) that are outside the nominal constellation, or (b) these satellites may not broadcast open-service navigation messages. The latter case exerts little effect on the integrity analysis because the navigation solution will not be affected by a satellite that does not broadcast navigation messages. Besides, it makes sense that the satellites outside the nominal constellation shall not be used in safety-critical applications. In this work, we also discard the observations of C13/C15 before 2017, considering that the stability of the satellite itself and of the precise products may be poor in the first several months after launch.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 1

Summary of observations for BDS-2 satellites, where green indicates good observations, blue and yellow indicate missing WHU and GFZ precise ephemerides respectively, black means the absence of PCEs, and red indicates no BCE. Note: the satellite C018 switched its PRN from C17 to C03 on 2018/09/29, and C017 switched its PRN from C15 to C13 on 2016/10/11

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 2

Summary of observations for BDS-3 satellites, where green indicates good observations, black means missing PCE, and red indicates the absence of BCEs. Note: the satellites C18 and C31 are outside the nominal BDS-3 constellation

For the rest of the satellites, Figure 3 and Figure 4 present a direct view on their data unavailability percentage. The results suggest that despite the efforts devoted to data augmentation, a relatively large percentage of BCE and PCE products are absent. This will negatively impact the evaluation of the BDS integrity performance. Therefore, it is worth noting that the following analyses can only coarsely capture the characteristics of BDS SISREs, and we will conduct more rigorous evaluations by accounting for every absence in future work.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 3

The percentages of missing BCE and PCE products for BDS-2 satellites

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 4

The percentages of missing BCE and PCE products for BDS-3 satellites

The precise products must be much more accurate than the BCEs to serve as the reference for SISRE evaluation. For BDS-2, prior studies have shown the feasibility of utilizing precise ephemerides in SISRE evaluation (Montenbruck et al., 2018; Yang, Xu, Nie, Gao, & Guan, 2019). As for BDS-3, the precise products of some satellites (including C19∼C22 and C29∼C30) achieve a 7-cm level radial orbit accuracy, while the products of other satellites may have poor quality due to limited tracking data (J. Guo, personal communication, August 1, 2019). Therefore, only these six satellites will be involved in the preliminary BDS-3 SISRE evaluation.

2.2 Time, coordinate, and corrections

Because the BCEs and PCEs are referred to the different time and coordinate systems, adequate care must be taken to ensure the consistency in the comparison. Aside from time and coordinate system alignment, several additional corrections are required in the data preprocessing procedures.

The BDS BCEs are referred to the Beidou time (BDT) scale, while the PCEs are given in the GPS time (GPST). In this study, GPST is selected as the basic time scale, and a 14-s BDT-GPST time offset is applied to the broadcast ephemerides. The synchronization error between GPST and BDT is in the order of 50 ns (CSNO, 2019a; Wang & Du, 2018), which can be ignored in the comparison between the broadcast orbits and the precise ones.

The broadcast and precise orbits are respectively specified in CGCS2000 and ITRF2008. However, the two frames are commonly considered to agree at few centimeters level (Montenbruck, Steigenberger, et al., 2015). Since this difference is well below the uncertainty of the BCE, it can be ignored in the following analysis. It is worth mentioning that BDS and GPS adopt a different geocentric gravitational constant and Earth’s rotation rate, and incorrect substitution could lead to orbital errors up to several meters (Yang et al., 2019).

The satellite phase center offset (PCO) also significantly impacts the SISRE evaluation, which describes the location of the satellite’s antenna phase center (APC) relative to its center of mass (CoM). The PCEs generally provide the CoM orbit positions and the APC clock offsets, while the orbit and clock parameters derived from the BCEs are both referred to the APC. As an exception to this rule, the reference point of the broadcast orbits for BDS-2 was CoM before Day-Of-Year (DOY) 17, 2017 and was changed to the APC of B3I signal thereafter (B. Wang, Chen, Wang, & Zhou, 2018). Therefore, PCO corrections should be applied in the comparison between the broadcast and precise orbits after that day.

The PCO values of BDS-2 and BDS-3 satellites are generally provided in the satellite-fixed coordinate frame: the origin is the CoM of the satellite, the z-axis points toward the geocenter, the y-axis corresponds to the cross product of the z-axis with the satellite-Sun vector, and the x-axis completes the right-handed system (Montenbruck, Schmid, et al., 2015). Different PCOs are adopted by different MGEX ACs. For both BDS-2 and BDS-3, WHU used the manufacturer-provided PCOs (CSNO, 2019c) before GPS week 2072 and has been adopting the IGS official values included in igs14.atx (IGS MGEX, 2020b) since then. In the meantime, GFZ has been using the European Space Agency (ESA)-estimated BDS-2 PCO values (Dilssner, Springer, Schönemann, & Enderle, 2014). The broadcast orbit errors in Earth-centered Earth-fixed (ECEF) frame, [ Δ𝑥 Δ𝑦 Δ𝑧 ]^T, are

1

where denotes the attitude matrix describing the transformation from the satellite-fixed frame to ECEF; 𝑑𝑥, 𝑑𝑦, and 𝑑𝑧 are the PCOs; and subscripts “b” and “p” respectively represent the BCEs and PCEs. In this work, is calculated by assuming that the satellites adopt a nominal yaw-steering attitude mode, although the satellites sometimes switch to the orbit-normal attitude mode (Montenbruck, Schmid, et al., 2015). This is acceptable because adopting incorrect attitude modes can only cause a projected range error of about 15 cm at most for any point on the Earth (Montenbruck et al., 2018).

The BDS broadcast clock offset of satellite i is referred to the APC of the B3I signal (CSNO, 2019a), while the precise clock offset is referred to the equivalent APC of the dual-frequency ionosphere-free (IF) combination. Therefore, time group delay (TGD) correction should be applied to the broadcast clock offset for comparison with the precise product. In the BCEs, the TGD corrections for B1I/B3I and B2I/B3I are given in and , respectively. Please note, the errors of the broadcast TGDs contribute to the broadcast clock errors, and their effects on the user ranging errors are accounted in the above process. Therefore, we do not additionally calculate the TGD errors using differential code biases (DCBs), and this is justified for the evaluation of the aforementioned dual-frequency users. However, in future work, we will evaluate the BDS TGD accuracy as Martini, Sgammini, and Boyero (2020b) do to support the users using single frequency or other dual-frequency combinations. Table 1 shows the associated TGD corrections for different precise products. Note that, since the relativistic effect is not corrected in the PCEs, the correction should also NOT be applied to the broadcast clock offsets. Accordingly, the difference between the broadcast and precise clock offsets is given by the following:

2

Besides, there is a time-varying and satellite-independent timescale bias 𝜇 (relative to GPST) in the reference clock of the PCEs. Therefore, at epoch k, can be modeled as

3

where Δ𝑇^𝑖(𝑘) denotes the broadcast clock error for satellite i. To obtain Δ𝑇^𝑖(𝑘), the bias term 𝜇(𝑘) needs to be estimated. We adopt the robust iterative-weighted-average method (Wu et al., 2017) to realize the estimation because of its good tolerance to the potential clock outliers.

2.3 Ephemeris status: Age-of-data, URA index, and health flag

The SIS accuracy is usually affected by the elapsed time since the generation of the orbit and clock information. The broadcast ephemeris employs age-of-data-ephemeris (AODE) and age-of-data-clock (AODC) to indicate how recent this navigation message is. Tables 2 and 3 summarize the percentages of different AODE and AODC values for BDS satellites over the analysis period. The results suggest a noticeable difference between the values of AODC and AODE. This is caused by the fact that the generation of clock information is separate from the orbit determination process (Wu et al., 2017). The data also reveal that the BDS-2 MEOs have the longest AODs due to the limited ground tracking and uploading stations. As for BDS-3, the AODs are mostly shorter than an hour, which highlights the benefits of the ISLs.

URA Index (URAI) is broadcast as an indicator of URA, which provides a conservative RMS estimate of the SIS-induced pseudorange errors. A part of the relationship between URAI and URA is briefly shown as follows, and the full version is given in the interface control document (ICD) (CSNO, 2019a): URAI = 0, 0.00 m < URA≤2.40 m; URAI = 1, 2.40 m < URA≤3.40 m; URAI = 2, 3.40 m < URA ≤4.85 m; URAI > 2, URA > 4.85 m. Table 4 shows that the URAI value is nearly always equal to two, i.e., the broadcast URA is usually bounded by 4.85 meters.

For safety purposes, the broadcast ephemeris indicates the health status of a satellite using the health flag. When this flag is set to “1,” the associated satellite is deemed as unhealthy and encounters an outage. Table 5 shows the percentages of unhealthy ephemerides for different types of satellites. The results indicate that the BDS-2 MEOs are most frequently subject to failures.

In this work, we perform the SIS performance evaluation with all healthy ephemerides, without partitioning the data by AODC, AODE, and URAI. Therefore, we can obtain a unique overbounding URA for each satellite, which can cover all nominal situations. However, it should be noted that this may lead to optimistic estimates for some AOD values because the SIS accuracy might vary as a function of AODs. Appendix B provides a quantitative analysis on this relationship, and the results suggest that the variation in accuracy with AODs is considerable for a few BDS-2 satellites (mainly the MEOs) but not very obvious for the others. The analysis also reveals that the partition by AODs may fail because the samples at some AODs are not sufficient to be statistically meaningful.

3.1 Definitions of SIS faults

SIS anomalies (i.e., faults) are mainly caused by incorrect navigation data, orbital maneuvers, abnormal orbital disturbances, clock jumps, clock drifts, and other unusual events (Gunning et al., 2018). SIS faults can potentially lead to large navigation errors and constitute the major threats to navigation safety in aviation applications.

In this study, we adopt two different definitions of SIS faults and term them as Definition a and Definition b, respectively. The former is a new definition that takes the fault shape into consideration, and the latter is a snapshot approach recommended by Walter et al. (2019). Both definitions depend on a predefined threshold, and the rationale of the threshold determination is given as follows. In prior studies, the threshold is usually given as , denotes the broadcast URA (Walter et al., 2019). However, we advocate selecting a set of thresholds for the following reasons: first, we are not sure whether the broadcast URA can properly describe the true error distributions; second, this threshold is an important parameter that significantly influences the estimated URA and fault probability; third, this approach can offer a more comprehensive characterization of the SISREs than solely adopting one threshold.

The two definitions are illustrated in detail as follows. Under Definition a, an SIS fault is defined as an event occurring when there is an abnormal value or trend in the SIS UREs. The abnormal value refers to the IURE whose magnitude exceeds the threshold. And an abnormal trend is said to exist when (a) the IURE data at several consecutive epochs form a curve with a ramp, sinusoidal, or any nonrandom shape, and (b) this curve contains at least one abnormal value. Under Definition b, an SIS fault is defined to occur whenever the IURE is greater than the threshold (Walter et al., 2019).

Figure 6 shows the difference between the two definitions by simulating a ramp fault. This figure shows that the fault duration is different under these two definitions. At the beginning of this fault, the URE magnitudes are below the threshold, but there is clearly an abnormal trend. These samples are regarded as faulted under Definition a, while they are considered unfaulted under Definition b. Therefore, the fault duration is extended under Definition a, and thus, the fault probability is increased, which is the trade-off for potentially reducing the URA.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 6

Graphical illustration of the SIS fault definitions

Prior studies usually adopt Definition b because (a) it avoids any ambiguity about whether a fault exists and maintains consistency with GPS performance standards and (b) the below-threshold segment imposes little threat to the users. However, the new definition is also reasonable and meaningful, although it is somewhat ambiguous. This is because, first, it is safe (i.e., conservative) to regard a sample as faulted in the determination of ISM parameters. This is actually a trade-off between URA and fault probability. Second, as an alternative, this definition can broaden the understanding of the distinction between “faulted” and “nominal” samples, which will be further illustrated in Section 4.2.

3.2 Anomaly detection

This section presents how to detect the anomalies in the observed data. The worst-case URE (MPE) metric is employed in this task because the most vulnerable user must be protected in integrity-related applications (Walter & Blanch, 2015). It is straightforward to implement anomaly detection under Definition b (Walter et al., 2019). However, it is hard to unambiguously perform this task under Definition a due to the difficulty in accurately identifying all the samples that belong to a fault. In this paper, we present a simple and empirical scheme to roughly estimate the start and end time of a fault as follows. This scheme is preliminary and should be viewed as an example. And due to the manual operation in the final step, it cannot automate in an ISM generator (ISMG). In future work, it will be improved to support autonomous and reliable anomaly detection.

1. Detect the faults based on Definition b, and record the faulted epochs.

2. If a fault lasts over 45 minutes (changeable) and its first or last data differs from the threshold by less than 10 meters (changeable), go to step (3); otherwise, this fault is regarded as an outlier or a step fault.

3. Plot the curve of this fault with the surrounding data, and manually estimate its start and end epochs based on the curve.

3.3 Determination of fault probabilities

Walter et al. (2019) have presented the precise definitions of prior fault probabilities and have used them for the assessments of GPS and GLONSS. These definitions are also adopted in this study, and they are briefly illustrated as follows.

According to the ARAIM baseline algorithm, we consider two nonoverlapping classes of SIS fault probability that account for satellite faults and constellation faults, respectively. The probability of satellite fault, 𝑃_sat, is defined as the probability that an individual satellite is in a fault state at any given time. And the probability of constellation fault, 𝑃_const, is the probability that multiple satellites in a constellation are simultaneously faulted due to a common cause.

Fault rate 𝑅 is defined as the probability that a healthy satellite or constellation becomes faulted within a period of time. Let 𝑁_𝐹 be the number of fault occurrences within a total time interval 𝑇_𝑡, evaluating 𝑅 is equivalent to deriving the expected value of 𝑅 given 𝑁_𝐹 and 𝑇_𝑡, i.e., 𝐸(𝑅𝑁_𝐹). Rigorous derivations have been established by Walter et al. (2019) with key assumptions laid out and all the steps carefully justified. They examined a priori probabilities of 𝑅 that are proportional to 𝑅^𝑚, where 𝑚 is an exponent. Based on this, the final estimate of is (Walter et al., 2019)

7

This equation applies 𝑚 = −1∕2 in the final term, and the justifications for this choice are provided by Walter et al. (2019). In brief, the choice of 𝑚 = −1∕2 ensures that the evaluation results are consistent across different parameterizations (i.e., fault rate, 𝑅, versus mean time between faults, 1/ 𝑅). In this paper, we will directly employ Equation (7) for fault rate evaluation, although this is not the unique approach to do this.

Fault probability 𝑃 is related to fault rate 𝑅 through the following equation:

8

MTTN denotes the expected average fault duration. A fault is said to be over when it is corrected or the user is notified by the broadcast ephemeris. MTTN is an important indicator that shows the control segment’s ability to monitor and repair SIS faults. Since CSNO has not yet provided any commitment on this parameter, in this work, we preliminarily estimate the MTTN values based on the historical observations.

4.1 Unfaulted SIS UREs and URA definition

Section 4 focuses on evaluating the nominal SIS performance with unfaulted SIS UREs. If an IURE is not identified as a fault, then it is under nominal conditions. Therefore, the evaluation will be directly affected by the fault definition and the threshold. This evaluation aims at properly estimating the URA values for each monitored satellite. URA provides a conservative estimate (i.e., an integrity overbound) of the uncertainty on the SIS URE (Walter, Blanch, & Gunning, 2019). Although this parameter is available in the navigation message, it must be validated or adjusted (if necessary) by air navigation service providers (ANSP) before being provided to aviation users. For this purpose, Section 4.2 presents the methodology to estimate a safety-assured URA based on the observed data.

4.2 Overbounding SIS UREs

In ARAIM algorithms, range errors are usually described by Gaussian models. However, the actual error distributions are usually not ideally Gaussian due to complex disturbances (Walter et al., 2018). Fortunately, various overbounding techniques have been proposed to derive a conservative Gaussian distribution with the observed one. The definition of overbounding can be mathematically described as follows (Larson, 2018):

9

where denotes the probability of the event inside the parentheses, 𝑋_ob is the random variable of the overbounding distribution, 𝑋 is the random variable of the actual distribution, and is the error-free value of 𝑋. Equation (9) is an expression in the form of a cumulative distribution function (CDF). In GNSS integrity domain, single-CDF overbounding (DeCleene, 2000), paired-CDF overbounding (Rife, Pullen, Enge, & Pervan, 2006; Rife, Walter, & Blanch, 2004), two-step overbounding (Blanch, Walter, & Enge, 2018), and three-condition overbounding (Martini et al., 2020a) are four Gaussian overbounding approaches. A peculiar single-CDF approach, called folded CDF (FCDF), is adopted in this work. In this approach, is equal to the median of 𝑥, and thus, (Larson, 2018).

It is acceptable to use a single-CDF approach for the analysis, although the other three approaches show advantages in (a) avoiding the requirements on symmetry and unimodality and (b) ensuring the conservativeness after convolution. This is because, first, the overbounding results might be excessively conservative when applying the latter three approaches (Walter et al., 2018). Second, Appendix C shows that for most of the satellites, the orbit and clock errors approximately conform to symmetric and unimodal distributions. Moreover, prior studies on GPS have examined that the range errors do not convolve in a hazardous way before the position error exceeds the protection level (Walter et al., 2018). In future work, we will analyze the effect of error convolutions in BDS and will also consider the three advanced approaches for more reliable evaluations.

Figure 7 provides a graphical depiction of the FCDF overbounding approach. In this figure, the actual heavy-tailed distribution is simulated by mixing 40 points uniformly distributed in (-6, -4] ∩ [4, 6) with 10⁵ points from the Gaussian distribution ℕ(0, 1²), and the overbounding distribution is ℕ(0, 1.5²). This figure indicates that merely overbounding the core distribution may lead to an overly optimistic estimate, and it also suggests that the beginning segment of a ramp fault is a possible cause of heavy tails.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 7

Graphical depiction of the FCDF overbounding concept

As mentioned in Section 2.4, the user-grid URE metric is used to evaluate the nominal SIS performance. To conduct the evaluation for a satellite, we first select a mesh of users that are evenly distributed in its coverage footprint. For BDS, the orbit altitude is 35,786 kilometers for the GEOs and IGSOs and 21,528 kilometers for the MEOs (CSNO, 2019a). Therefore, if the user is near the Earth’s surface and its mask angle is zero, the minimum latitude 𝜃_min is 8.73^◦ for the GEO/ IGSO satellites and 13.18^◦ for the MEO satellites. Table 6 gives the locations of the selected users with respect to the O-XYZ frame shown in Figure 5. These users are fixed in this frame to keep their LOS vectors toward the analyzed satellite invariant with time. As a reminder, the definitions of latitude and longitude here are different from those in the geographic coordinate systems (e.g., WGS84).

Finally, the following steps show how to determine the overbounding URA (𝜎_URA) for a satellite.

1. Calculate the IUREs for every selected user at each epoch when there are both healthy BCE and valid PCE.

2. Perform anomaly detection by comparing the MPEs with the threshold; then, if epoch i is labeled as “faulted,” change all the IUREs of this epoch to ∞×sign(MPE[i]).

3. Form the IURE time series for each user so that we can evaluate the IURE distributions for each user separately.

4. For each user, (a) estimate the mean of the corresponding IURE time series, (b) remove this value to form the zero-mean IURE distribution, and (c) generate an FCDF curve based on this IURE distribution.

5. Find a Gaussian distribution that can strictly overbound the FCDF curves of ALL users for ∀ 𝐼𝑈𝑅𝐸 ∈ [−𝑇_𝑙,𝑇_𝑟], where 𝑇_𝑙 = min{IURE | IURE < 0 and IURE is finite} and 𝑇_𝑟 = max{IURE | IURE > 0and IURE is finite}.

Figure 8 provides the overbounding results of C09 to further illustrate the method above. In this figure, each colored thin line represents the IURE distribution of a user, and the heavy blue line denotes the overbounding Gaussian distribution. Under the 10-meter threshold, the number of faulted samples is 86 under Definition a and 41 under Definition b. As a trade-off, the nominal distributions that are to be bounded for evaluating the URA are obviously less heavy-tailed under the former definition as compared to the latter one. Specifically, |𝑇_𝑙| and 𝑇_𝑟 are about 8 meters under Definition a and 10 meters under Definition b. This results in a smaller overbounding URA under the former definition at the cost of an increase in fault probability. The trade-off shown here proves the significances of the new definition (Definition a): it broadens the understanding of SIS faults and serves as a supplement to the original definition (Definition b).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 8

Bounding the user-grid UREs for C09 under different fault definitions, where the green shadow indicates the “nominal” distribution to be bounded

5.1 Example analysis

Based on the proposed SISRE evaluation scheme, this subsection shows some results as the examples to further illustrate the associated methods and to identify some special BDS SIS features. This part is a complement to the results in Section 5.2 where we quantitatively evaluate the integrity performance of every monitored satellite.

First, Figure 9 compares the error behaviors among different types of satellites. In this figure, the error distributions are formed without removing any sample. And to allow for a clearer view, the x-axis is limited to ±30 meters and the y-axis ranges from 10^–5 to 1. Note that, the errors of C21 are all smaller than 8 meters. This figure proves the superior SIS accuracy of BDS-3 as compared to BDS-2. And it also shows that these error distributions (except MPE) appear to have a Gaussian core while their tail data seem to have a larger variance than the majority. This is clear evidence of the occurrence of heavy tails, and this phenomenon is also reported in GPS and other constellations (Walter et al., 2018).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 9

Radial (R), along-track (a), cross-track (c), clock, and the maximum projected error distributions in the form of FCDF for various satellites

Figure 10 shows the profiles of several SIS fault types observed in the dataset. The shape of a fault is related to its underlying causes. And the causes can be identified by analyzing the surrounding broadcast parameters. For example, the ramp fault in this figure was caused by clock instability, and the sinusoidal fault was due to incorrect orbit parameters (specifically, the rate of right ascension and the sine harmonic correction term to inclination angle 𝐶_𝑖𝑠 were incorrect).

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 10

The profiles of various types of SIS faults observed in the dataset

Next, we use an example of an SIS fault occurring in an MEO to illustrate how a fault is corrected in orbit. Figure 11 presents the orbit and clock errors over the corresponding period, and Table 7 gives the broadcast navigation parameters. As seen, this is a ramp clock fault whose magnitude increases from a small value to 60 meters. And the table indicates that this fault was caused by incorrect clock parameters. This fault had lasted 18 hours before the user was notified by a new navigation message (the fourth row in Table 7), where the health flag was set to one. Note that, the clock correction parameters were not updated in this message, and thus, the AODC was unchanged. The reason behind the long fault duration is the limited distribution of BDS-2′s ground tracking and uploading stations (Wu et al., 2017). Finally, Table 7 shows that the satellite recovered from the fault after the new clock parameters were uploaded to the satellite and broadcast to the users.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 11

Orbit and clock errors when a fault occurs in C11 (MEO)

Special attention should be paid to the anomalous SIS performance evaluation for BDS-2 MEOs. Since BDS-2 only serves the Asia-Pacific region (55^◦S∼55^◦N, 70^◦E∼150^◦E), it is reasonable to remove the faulted data that do not affect any user in this region from the evaluation. And in safety-critical applications, it is strongly suggested that the ARAIM software and/or the receiver disable the BDS-2 satellites when the user is outside this region. Figure 12 summarizes all the faulted data of C11 (45 faults in total) and shows the corresponding subsatellite points. There is a green curve in this figure, and the area inside this curve is called the “dark zone.” When an MEO flies over this zone, it is completely invisible to all the users in the BDS-2′s service region (the user mask angle is 15^◦). This figure indicates that many faulted data correspond to this zone and they should be removed from the determination of fault probability, fault rate, and MTTN.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 12

Subsatellite points of C11 when it is under faulted conditions, where black indicates clock faults, red denotes orbit faults, and the green curve shows the dark zone

5.2 Preliminary estimates of URA and fault probability for BDS satellites

This section presents detailed BDS SISRE characteristics and compares the SIS performance among different satellites. First, we provide the estimated values of URA and fault probability for each satellite as follows. In the evaluation, Definition b is adopted as the primary definition while Definition a is only applied to some satellites where the two definitions can lead to significantly different estimates of URA. Pease note, the following results can only be viewed as preliminary due to (a) the relatively frequent absence of PCEs and BCEs and (b) the limited data volume.

Figures 13–15 present the estimated values of URA (𝜎_URA) and fault probability for various types of BDS-2 satellites. We can draw the following conclusions from the results. First, as the detection threshold is raised, 𝜎_URA tends to increase while the fault probability will decrease. Their values will probably level off if the threshold is large enough. Besides, for some satellites (e.g., C09-C11), 𝜎_URA differs significantly under different fault definitions. This implies that the tail error distributions of these satellites are greatly impacted by ramp faults. Similarly, Figure 16 shows the results for the monitored BDS-3 satellites. It is obvious that the conclusions above are also valid for BDS-3. Please note, these figures demonstrate the variation in the results with fault definitions and thresholds, but we do not know which set of results can best describe the long-term planned behaviors. This remains an open problem for future investigations.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 13

The estimated values of URA and fault probability (P) for BDS-2 GEO satellites under Definition b

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 14

The estimated values of URA and fault probability (P) for BDS-2 IGSO satellites, where the black dashed lines indicate Definition a and the red solid lines denote Definition b

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 15

The estimated values of URA and fault probability (P) for BDS-2 MEO satellites, where the black dashed lines indicate Definition a and the red solid lines denote Definition b

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 16

The estimated values of URA and fault probability (P) for BDS-3 satellites, where the black dashed lines indicate Definition a and the red solid lines denote Definition b

Then, we further compare the SIS performance among different satellites as follows. Hereafter, the threshold is set to 10 meters for BDS-2 satellites and 4 meters for BDS-3 satellites. Please note, although Table 4 offers the broadcast URAs, we do not choose the thresholds based on them because we are not confident that (a) they can accurately describe the actual SISRE behaviors (the results above show that they might be overly conservative), and (b) the Constellation Service Provider (CSP) will make the commitments based on them (because the CSP may make commitments in written specifications, e.g., Open Service Performance Standard). Figures 9 and A8-A9 (Appendix C) demonstrate that for most BDS-2 satellites, the UREs are mostly below 10 meters and almost all the UREs above this number correspond to some fault events (e.g., ramp faults). This conclusion is also true for BDS-3 satellites when the threshold is 4 meters. Therefore, it is initially acceptable to use these two thresholds for the preliminary performance evaluation and comparison.

Figure 17 presents the estimated URAs (𝜎_URA) for all the monitored satellites, and Figure 18 shows the corresponding fault probabilities. The results suggest that most of the BDS-2 satellites hold the 𝜎_URA of 2.4∼3.2 meters, except C11 whose 𝜎_URA is greater than 3.4 meters. Their fault probabilities also vary a lot with a range from 1.5×10^–4 (C01 and C08) to 4.5×10^–3 (C11). Please note, these results are derived under the assumption that there is no fault during any of the data outages. Given the large variability in performance over satellites, it is strongly suggested to set different values for individual satellites in the ISM rather than merge them together.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 17

URA values of various satellites, where triangles indicate Definition a and circles denote Definition b

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 18

Fault probabilities for various satellites, where triangles indicate Definition a and circles denote Definition b

Moreover, the results suggest the superior SIS accuracy of the BDS-3 satellites to the BDS-2 satellites. For most of the BDS-3 satellites, the 𝜎_URA is about 1 meter when Definition a is considered. This mainly benefits from new spaceborne clocks and the equipped ISL payloads. The former helps improve the precision of short-term clock prediction, and the latter mitigates the challenges caused by the limited ground stations by enabling satellite-to-satellite communication and ranging.

Although Figures 16 and 18 also show the estimated fault probabilities of the BDS-3 satellites, we are not confident in these numbers because the fault probabilities of some satellites (e.g., C20, C29, and 30) are significantly influenced by some prolonged (e.g., over 10 hours) faults. And the long fault duration may be attributed to the maturation phase that BDS-3 has been undergoing. It is reasonable to expect that the fault duration will not be such long after this phase. Please also note that the current results may not provide a fully reliable picture of the future BDS-3 SIS performance. This is because, first, there might be some satellite faults during the data outages. Second, the data volume is not sufficiently large for a robust evaluation. Moreover, the integrity performance may be significantly improved after the constellation becomes mature. Therefore, it is reasonable that BDS-3 can make commitments into the SARPs with better performance than the results obtained in this work. Toward making BDS-3 available to ARAIM users, we will continuously conduct performance monitoring for BDS-3 in future work.

Figure 19 shows the observed fault rates of each monitored satellite. Note that, the fault rate is exactly the same under the two fault definitions because they only impact the fault duration. The preliminary results suggest that the fault rate varies from 1.3×10^–4 to 10^–3 per hour among different BDS-2 satellites and ranges from 4×10^–5 (0 fault in 1.5 years) to 3.5×10^–4 (4 faults in 1.5 years) per hour for BDS-3 satellites. Besides, there are no concurrent faults and thereby no constellation fault observed in the dataset of BDS-2 and BDS-3. Based on Equation (7) (Walter et al., 2019), the constellation fault rate can therefore be preliminarily set to 1.3×10^–5 per hour for BDS-2 and 4×10^–5 per hour for BDS-3. We cannot provide the estimate of constellation fault probability because the MTTN of constellation faults is unknown. This result is obtained by assuming that there is no constellation fault hidden in the missing data. We will examine this assumption and update the result (if necessary) in future work.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 19

Estimated fault rates for BDS-2 and BDS-3 satellites

Due to the limited data volume of BDS-3, we will only focus on BDS-2 in the following analyses. MTTN is an important performance indicator of a navigation system, which shows how quickly it can respond to an SIS fault. Figure 21 presents the preliminary estimates of MTTN for each BDS-2 satellite. The MTTN values here are simply estimated by calculating the average fault duration. As can be seen, for most of the satellites, the MTTN appears to be below 2.5 hours under Definition b and below 3 hours under Definition a. However, there is also an obvious variation in this number among different satellites: the MTTN varies from 30 minutes (C01) to about 4 hours (C11). On the whole, the GEOs have the lowest MTTN while the MEOs have the longest one. This figure also suggests that the MTTN is impacted by the fault definition, and this is because the fault duration is extended under Definition a.

However, the above MTTN estimation method may not provide a reliable estimate with sufficient accuracy for the following reasons. First, the number of observed faults is not large enough to accurately capture the SIS fault characteristics. Second, the calculation relies on the assumptions that (a) the duration of a fault is at least 15 minutes and (b) a fault occurring in one time interval does not affect the faults occurring in any other intervals, but both of them may not be true. Besides, the individual fault durations for the MEOs would be significantly affected by when the fault initiates compared to the visibility of the ground monitoring stations. Therefore, the results in Figure 20 are preliminary and should be further examined in the future. Integrity Failure Modes and Effects Analysis (IFMEA) is an effective approach to help assess the MTTN, which intends to identify satellite failure modes, determine their probability of occurrence, and analyze their effects on the users (Ober & Harriman, 2006). Therefore, we provide some initial assessments on the fault characteristics of BDS-2 as follows. In future work, more efforts will be taken to complete the IFMEA for BDS and to further assess the MTTN with the analysis results.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 20

Estimated MTTN for BDS-2 satellites, where triangles indicate Definition a and circles denote Definition b

Figure 21 summarizes the sources of the observed BDS-2 SIS faults, and Figure 22 presents the relative frequency histogram of their maximum magnitudes. The results indicate that over 70% of the SIS faults are attributed to the clock component. Clock anomalies even occupy more than 88% for MEOs, while this proportion is only about 50% for GEOs and IGSOs on the whole. As for the fault magnitude, Figure 22 suggests that the magnitude of most faults is below 100 meters, although large (over 1 kilometer) faults also account for a significant proportion.

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 21

The fault sources of BDS-2 satellites, where the bottom indicates clock and the top denotes orbit

• Download figure
• Open in new tab
• Download powerpoint

FIGURE 22

The relative frequency histogram of the maximum magnitudes of BDS-2 SIS faults