Verification and Accuracy Analysis of Single-Frequency Occultation Processing Based on the BeiDou Navigation System

Abstract

GNSS single-frequency occultation processing technology has the advantage of simple instrumentation, but it is not clear about the accuracy of the Beidou-based single-frequency occultation processing. This paper verifies the single-frequency occultation processing algorithm of the BeiDou navigation system (BDS) and analyzes its accuracy based on occultation observation data from the FY3E satellite. The research aimed to verify the single-frequency ionospheric relative total electron content (relTEC), analyze the accuracy of the reconstructed second frequency ${\mathrm{B}}_3^{\ast }$’s excess phase Doppler, and analyze the accuracy of the refractive index products. Results: (1) As for relTEC and excess phase Doppler, the correlation coefficient between single-frequency occultation processing and dual-frequency occultation processing is greater than 0.95. (2) The relative average deviations of the excess phase Doppler of ${\mathrm{B}}_3^{\ast }$are mostly less than 0.2%, and the relative standard deviations are mostly around 0.5%. (3) The bias index and root mean square index of single/dual-frequency inversion have good consistency compared with ERA5 data. All the results show that the single- and dual-frequency inversion refractive index products have comparable accuracies, and the accuracy of the standard deviation of single-frequency inversion refractive index products over 25 km being slightly lower than that of dual-frequency inversion refractive index products.

1. Introduction

Global navigation satellite system (GNSS) radio occultation was originally developed to measure the planetary atmosphere when exploring the solar system [1,2]. Then, in the mid-1960s, Fishbach put forward the GNSS radio occultation theory for remote sensing of the Earth’s atmospheric parameters [3,4], which provides accurate atmospheric parameter profiles worldwide. The theory is useful for numerical weather forecasting, climate monitoring, and atmospheric research. Because this theory requires putting multiple signal transmitting devices and signal receiving devices in space, until the late 1980s the United States Global Positioning System (GPS) tended to be perfect, providing the necessary satellite constellation for GNSS radio occultation [5,6].

GNOS, a global navigation satellite radio occultation detector carried by China’s Fengyun 3 series meteorological satellites, can simultaneously receive GPS navigation signals and BeiDou navigation system (BDS) navigation signals for occultation detection [7]. It is the world’s first radio occultation detector compatible with BDS and GPS systems [8]. The BDS occultation observation data in this paper comes from the dual-frequency occultation observations of BDS B1I signals and BDS B3I signals carried by the second-generation radio occultation detector (GNOS II) on the FY-3E satellite (FY-3E) in the Fengyun series of missions. The operating frequency of BDS B1I signals is 1561.098 MHz. The signal adopts multifrequency transmission technology and time-division frequency multiple-access technology, which have strong anti-interference ability and can eliminate errors in navigation and positioning worldwide. The BDS B1I signal has higher accuracy and reliability. The operating frequency of BDS B3I signals is 1268.52 MHz. In this paper, the carrier phase (${\mathrm{B}}_3^{\ast }$) of a quasi-B3I signal is reconstructed, and the results are compared with the actual BDS B3I signal in dual-frequency measurement.

In GNSS atmospheric occultation data processing, the ionospheric error will greatly affect the accuracy of atmospheric occultation retrieval products [9,10]. Dual-frequency occultation processing means that when the radio-wave signal transmitted by the LEO satellite is received by the high-precision dual-frequency radio occultation receiver carried by the LEO satellite, the radio wave signal will have a certain amount of bending and delay effects through the lower atmosphere and the ionosphere due to atmospheric refraction, scattering, and other reasons, and the bending and delay effects will cause the frequency of the receiver signal to change [11,12]. In general, the ionosphere-free combination using dual-frequency observations can effectively eliminate the influence of the first-order term of the ionosphere, thus improving the product accuracy [13,14]. After eliminating the ionospheric delay through the ionospheric combination between the dual-frequency BDS-based carrier-phase observations, the Doppler frequency shift of the BDS signal is calculated using the high-precision velocity, position, and clock-difference information of the receivers of the BDS satellite (reference satellite) and the low-orbit satellite (occultation satellite), so as to invert the curved angular profile of the atmosphere [15,16,17,18].

In some flight missions, single-frequency receivers may also be used for GNSS atmospheric occultation detection. Since sometimes even single-frequency receiver occultation detection cannot eliminate the ionospheric error through combination, the half-sum combination of the single-frequency carrier phase and pseudorange has been developed to eliminate ionospheric error [19,20,21]. For example, a GPS radio occultation receiver named Turbooge carried on the Danish Ørsted satellite launched in February 1999 only received the single-frequency occultation signal normally [22], so the scientific researchers used the double-difference single-frequency occultation processing method for GNSS atmospheric occultation processing. However, the double-difference single-frequency occultation processing results used in the follow-up calculations cannot be compared with the double-frequency occultation processing results because the Turbooge receiver did not receive the dual-frequency signal. The single-frequency occultation processing method has the advantage of simple instrumentation, but it is not clear about the accuracy of the Beidou-based single-frequency occultation processing. At present, the latest GNSS radio occultation technology uses single-difference or non-difference technology for data processing [23,24], which requires algorithm verification and accuracy analysis of the single-frequency occultation processing method based on BDS observations under the current differential technology.

This paper introduces a method of using ${\mathrm{B}}_1$ single-frequency carrier phase and pseudorange observations to calculate the relative ionospheric total electron content (TEC). It then reconstructs the algorithm of ${\mathrm{B}}_3$ carrier-phase observations, which is described in detail in Section 2. The sources of the data are described in Section 3. Then, Section 4 reconstructs the second frequency ${\mathrm{B}}_3^{\ast }$ of occultation data observed by the FY3E satellite, compares the consistency of the ${\mathrm{B}}_3^{\ast }$ with the actual ${\mathrm{B}}_3$, and ends up with inverting the atmospheric refractive index product. In Section 4.1, the relTEC change rate results are compared, and in Section 4.2, the correctness of the single-frequency occultation processing method is verified by comparing the excess phase Doppler. In Section 4.3, this paper compares the BDS single-frequency occultation refractive index product with the ERA5 reanalysis field data to evaluate the accuracy of the single-frequency occultation product.

2. Single-Frequency Occultation Data Processing Method Based on BDS

The inversion algorithm based on BDS single-frequency radio occultation data means that when there is only observation of frequency band ${\mathrm{B}}_1$, the pseudorange and carrier phase are used to calculate the ionospheric delay, and further reconstruct the inversion of occultation products (atmospheric refractive index) in the frequency band ${\mathrm{B}}_3^{\ast }$. In this paper, the BDS observation data used to retrieve atmospheric occultation products is the ${\mathrm{B}}_1$ carrier-phase measurement at 50 Hz sampling rate and B1I pseudorange measurement (${\mathrm{C}}_1$) at the 1–50 Hz sampling rate. Due to the lack of data for the second frequency, the error caused by the ionospheric effect cannot be eliminated by difference combination. However, when only single-frequency measurement is used, the linear combination of ${\mathrm{B}}_1$ and ${\mathrm{C}}_1$ data can be used to construct ${\mathrm{B}}_3^{\ast }$, but the result will be affected by the noise of the observation value ${\mathrm{C}}_1$ itself [25,26,27].

This paper adopts a special filtering method to solve this problem. During data processing, with the help of a low-pass filter, the ${\mathrm{B}}_1-{\mathrm{C}}_1$ data are filtered. At the same time, the regularization method [28] that minimizes the second derivative is used to interpolate to 50 Hz and repair the data gap so as to reduce the noise impact from the observations. The filter can be described as matrix $\mathrm{F}$, given by the following formula:

$${\mathrm{F}=(\tilde{\mathrm{I}}+\mathsf{\gamma }\mathrm{S}}^{\mathrm{T}}{\mathrm{S})}^{-1}$$

(1)

where $\tilde{\mathrm{I}}$ is an incomplete identity matrix, which means the only diagonal element corresponding to the position where both ${\mathrm{B}}_1$ and ${\mathrm{C}}_1$ data are available is 1. Since in this paper the interpolation needs to be 50 Hz, the size of $\tilde{\mathrm{I}}$ is about 500 × 500. Although there are many elements in the matrix $\tilde{\mathrm{I}}$, only about 10–500 elements on the diagonal are 1, and the other elements are all 0. The specific value of $\tilde{\mathrm{I}}$ is determined by the sampling rate of the carrier phase and the pseudorange. Matrix S is a second-order derivative finite-difference operator [29], ${\mathrm{S}}^{\mathrm{T}}$ is its transposition, and $\mathsf{\gamma }$ is a regularization parameter that determines the smoothing strength. Here, $\mathsf{\gamma }$ = 10⁶, corresponding to a low-pass filter with a cut-off frequency of ~0.05 Hz, the cut-off frequency is selected to effectively eliminate the uncertain high-frequency noise in ${\mathrm{C}}_1$ data from the signal BDS B1I. Therefore, the result of the filtered ${\mathrm{B}}_1$ minus ${ \mathrm{C}}_1$ is obtained by the following formula:

$$\overline{{\mathrm{B}}_1{-\mathrm{C}}_1}{=\mathrm{F}(\mathrm{B}}_1{-\mathrm{C}}_1)$$

(2)

where ${(\mathrm{B}}_1{-\mathrm{C}}_1)$ is a vector with about 500 elements, and the missing data point of ${\mathrm{B}}_1$ or ${\mathrm{C}}_1$ is 0. The expression $\mathrm{k}(\overline{{\mathrm{B}}_1{-\mathrm{C}}_1})$ is the relative total electron content of the ionosphere, relTEC, in which the coefficient k is the conversion coefficient between the ionospheric delay and the total electron content. The carrier-phase measurement of ${\mathrm{B}}_3^{\ast }$ is based on the carrier-phase measurement of ${\mathrm{B}}_1$ and the result of filtered ${\mathrm{B}}_1{-\mathrm{C}}_1$. The formula is as follows:

$${\mathrm{B}}_{3 }^{\ast }{=\mathrm{B}}_1-0.5(1-\frac{{\mathrm{f}}_1^2}{{\mathrm{f}}_3^2}\left)\right(\overline{{\mathrm{B}}_1{-\mathrm{C}}_1})$$

(3)

where ${\mathrm{f}}_1=1561.098 \mathrm{MHz}$ is the observed signal frequency, and ${\mathrm{f}}_3$ is the frequency of the observed signal to be constructed. In order to verify and analyze the reconstructed ${\mathrm{B}}_3^{\ast }$ and single-frequency occultation processing method, this paper proposes constructing the B3I signal with a frequency of 1268.52 MHz. The carrier-phase measurement of the reconstructed ${\mathrm{B}}_3^{\ast }$ without ionospheric influence can be directly calculated through the carrier-phase measurement ${\mathrm{B}}_1$ and pseudorange measurement. Follow-up calculations based on ${\mathrm{B}}_1$ and ${\mathrm{B}}_3^{\ast }$ (the single-frequency occultation processing method) to invert of occultation products and based on ${\mathrm{B}}_1$ and ${\mathrm{B}}_3$ (the dual-frequency occultation processing method) are the same [30,31,32,33].

3. Data Sources

The occultation observation data are selected from the ROEX-format occultation observation data received by the GNOS radio occultation detector on Fengyun 3E satellite on 4 January 2022. According to the instrument design, the GNOS instrument’s occultation data sampling rate is 50 Hz for closed-loop and 100 Hz for open-loop, and the accuracy of the occultation data is comparable to that of similar international occultation products [34].

As for precision orbit determination, the positioning data also come from the BDS positioning data of GNOS. The GNSS precise ephemeris selects the GBM fast product released by the GFZ analysis center of IGS, and uses the simplified dynamic precise orbit determination to obtain the FY3E precise orbit and clock deviation. Precision orbit determination of low-orbit satellites is a repeated iterative process. Firstly, the initial orbit of the low-orbit satellite is obtained by pseudorange and carrier phase measurement, which could calculate the initial orbit parameters. Then, the orbital parameters are further corrected by the least squares LSQ. A more precision orbit is then calculated with new orbital parameters. Through multiple iterations and residual judgment, the precise orbit and orbital acceleration parameters are finally obtained. In addition, the time arc of each precision orbit determination is 30 h, so for the adjacent two days, there is a 6 h orbit overlap, which is used to evaluate the accuracy.

The position accuracy of the final FY3E precision orbit is at the centimeter level, and the speed accuracy is on the order of millimeters per second. This accuracy can usually be used for subsequent occultation data processing. It is worth noting that the single- and dual-frequency occultation processing use the same FY3E orbit data, so we believe that the influencing factors of the occultation products in the study are only reflected in the difference between the single- and the dual-frequency occultation processing.

4. Verification and Accuracy Analysis of Single-Frequency Occultation Processing

4.1. relTEC Correctness Verification

In Section 2, the total relative electron content of the single-frequency ionosphere, relTEC, is measured by the difference between the pseudorange data of frequency ${\mathrm{B}}_1$ observed from occultation and the carrier-phase observation data. The relTEC in the dual-frequency observations in this paper is the difference between the ${\mathrm{B}}_1$ obtained from occultation observation data and the ${\mathrm{B}}_3$ carrier-phase observation data and the correctness of relTEC is verified based on the relTEC change rate observed by dual-frequency observations. Because the excess phase Doppler solution (in Section 4.2) in the occultation inversion reflected in the TEC is the relTEC change rate, the relTEC change rate at different heights will better reflect the influence of single- and dual-frequency occultation processing on occultation products. Taking two periods of occultation time as an example, Figure 1 shows the calculation results of the relTEC change rate of the BDS single- and dual-frequency at 0–140 km altitudes of occultation tangency points, where the horizontal axis is the relTEC change rate in the unit mM/s, and the vertical axis is the tangent height in the unit km. Blue and pink dots represent the change rate of single-frequency and dual-frequency relTEC, respectively. The results show that the relTEC change rate between single- and dual-frequency occultation processing is in the same order of magnitude, and the two have good trend consistency. Compared with the pink dots (the dual-frequency relTEC change rate), the blue point (the single-frequency relTEC change rate) shows a slight divergence, which means that the relTEC change rate calculated by the measurement difference of the pseudorange and carrier phase of only one frequency still has residual ionospheric error, and it is also affected by the pseudorange noise. Even though the introduction of pseudorange observations also brings some noise, and there is a slight difference between the two at this stage, the overall trend of relTEC change rate is the same, and the effect is ideal. The relTEC change rate calculated by the differential solution of the pseudorange and carrier phase of the single-frequency and calculated by the differential solution of the carrier phase of the dual-frequency are basically consistent at different altitudes within 0–120 km.

To better quantify the consistency of single-frequency and dual-frequency occultation processing results, first-order linear fitting should be done on all single-frequency relTEC change rates calculated by an occultation event, and then the fitted function could be compared with the function y = x. Figure 2 is the comparison diagram of the relTEC change rate function of between the BDS single-frequency observations and the BDS dual-frequency observations under two occultation events with different data volumes. In the figure, the horizontal axis is the relTEC change rate calculated by the dual-frequency carrier-phase observation combination, and the vertical axis is the relTEC change rate calculated by the single-frequency pseudorange-carrier phase combination. The blue dots represent the relTEC change rate of each single-frequency observation. The fitting function of the single-frequency relTEC change rate is represented by a red line, and the function y = x, representing the consistency of the two, is represented by a blue line. As the figure shows, the more data there is for an occultation event, the more consistent the results of the single-frequency occultation processing and the dual-frequency occultation processing. Moreover, the concentration is higher when the data volume is greater. The two occultation events show that the blue dots are distributed near the curve of the function y = x, and the red fitting curve is very close to the curve position and curve trend of the function y = x. The effect is highly ideal.

The analysis in Figure 2 shows the consistency of the intensity of the linear relationship between the single-frequency and dual-frequency occultation processing. It shows that the two keep pace to a large extent in terms of their increasing or decreasing trends, even if they do not maintain a proportional relationship; that is, the characteristics of the position sequence in the occultation data processing are not considered. Therefore, the consistency of the relTEC results calculated by single-frequency and dual-frequency observations will be analyzed from the perspective of order.

In statistics, Spearman’s correlation coefficient quantifies the correlation between independent variable X and dependent variable Y. It relies fully on the data sample ranking value. If Y tends to increase when X increases, Spearman’s correlation coefficient is positive, while if Y tends to decrease when X increases, Spearman’s correlation coefficient is negative. When X and Y are more and more close to a complete monotonic correlation, Spearman’s correlation coefficient will increase in absolute value, and when X and Y are completely monotonically correlated, the absolute value of Spearman’s correlation coefficient is 1. The fully monotonically increasing relationship means that for any two pairs of data (${\mathrm{X}}_{\mathrm{i}}$, ${\mathrm{Y}}_{\mathrm{i}}$) and (${\mathrm{X}}_{\mathrm{j}}$, ${\mathrm{Y}}_{\mathrm{j}}$), the results of ${\mathrm{X}}_{\mathrm{i}}{-\mathrm{Y}}_{\mathrm{i}}$ and ${\mathrm{Y}}_{\mathrm{i}}{-\mathrm{Y}}_{\mathrm{j}}$ always share the same plus or minus sign. Spearman’s correlation coefficient is calculated as

$$\mathsf{\rho }=\frac{\frac{1}{\mathrm{n}}{{\displaystyle \sum }}_{\mathrm{i}=1}^{\mathrm{n}}\left(\left(\mathrm{R}\left({\mathrm{X}}_{\mathrm{i}}\right)-\overline{\mathrm{R}\left(\mathrm{x}\right)}\right)\cdot \left(\mathrm{R}\left({\mathrm{Y}}_{\mathrm{i}}\right)-\overline{\mathrm{R}\left(\mathrm{y}\right)}\right)\right)}{\sqrt{\left(\frac{1}{\mathrm{n}}{{\displaystyle \sum }}_{\mathrm{i}=1}^{\mathrm{n}}{\left(\mathrm{R}\left({\mathrm{X}}_{\mathrm{i}}\right)-\overline{\mathrm{R}\left(\mathrm{x}\right)}\right)}^2\right)\cdot \left(\frac{1}{\mathrm{n}}{{\displaystyle \sum }}_{\mathrm{i}=1}^{\mathrm{n}}{\left(\mathrm{R}\left({\mathrm{Y}}_{\mathrm{i}}\right)-\overline{\mathrm{R}\left(\mathrm{y}\right)}\right)}^2\right)}}$$

(4)

where $\mathrm{R}\left({\mathrm{X}}_{\mathrm{i}}\right)$ represents the ith relTEC change rate sequence calculated by the occultation data at the BDS single-frequency received by FY3E; $\mathrm{R}\left({\mathrm{Y}}_{\mathrm{i}}\right)$ represents the ith relTEC change rate sequence calculated by the occultation signal at the dual-frequency of BDS received by FY3E; and $\overline{\mathrm{R}\left(\mathrm{x}\right)}$ and $\overline{\mathrm{R}\left(\mathrm{y}\right)}$ represent the median of single-frequency and dual-frequency relTEC change rates, respectively.

Figure 3 shows the distribution histogram of Spearman’s correlation coefficient between single- and dual-frequency occultation processing of an occultation event on 4 January 2022. The dual-frequency carrier phase combination is used to calculate the relTEC change rate (as an independent variable X), and the single-frequency pseudorange-carrier phase combination is used to calculate the relTEC change rate (as a dependent variable Y), with an interval of 0.05. It is shown that the correlation coefficient between the rate of change of the dual-frequency solution and the rate of change of the single-frequency solution is greater than 0.95 in most (nearly 90%) of the occultation events of the BDS B1I signal. Moreover, the proportion of data with Spearman’s correlation coefficient greater than 0.8 is 96.8%, showing a strong correlation. The Spearman’s correlation coefficient of all data is greater than 0.7. Therefore, the relTEC change rate results of single- and dual-frequency observation have strong consistency related to order.

The analysis of the relTEC change rate form above three angles—direct plot the relTEC change rate with the altitude of 0–120 km (Figure 1), y = x consistency comparison (Figure 2), Spearman correlation coefficient histogram (Figure 3)—shows that in the occultation observations of BDS B1I signals, the relTEC calculated by the combination of single-frequency pseudorange carrier phase shows good consistency with the relTEC calculated by the dual-frequency carrier phase, and the result of the single-frequency occultation processing is ideal.

4.2. Reconstruction Excess Phase Doppler Correctness Verification

In Section 2, after the ionospheric relTEC step is calculated by combining the pseudorange and carrier phase of frequency ${\mathrm{B}}_1$, the quasi-carrier-phase observation ${\mathrm{B}}_3^{\ast }$ of frequency band ${\mathrm{B}}_3$ is then reconstructed with the combination of frequency ${\mathrm{B}}_1$. In the double difference method for occultation processing, excess phase Doppler is a very important parameter [35]. The excess phase Doppler $\Delta \mathrm{f}$ of frequency ${\mathrm{B}}_1$ is calculated by the excess phase delay:

$$\Delta \mathrm{f}=-\frac{1}{\mathsf{\lambda }}·\frac{\mathrm{d}\tilde{\mathrm{L}}}{\mathrm{dt}}$$

(5)

where $\tilde{\mathrm{L}}$ indicates the excess phase delay of the single-frequency carrier phase and the quasi-carrier-phase (or the excess phase delay of the dual-frequency double difference), and $\mathsf{\lambda }$ indicates the wavelength of frequency ${\mathrm{B}}_1$.

It is helpful to study the excess phase Doppler in the occultation process. Figure 4 shows the comparison of the excess phase Doppler calculated from the reconstructed quasi-carrier-phase observation ${\mathrm{B}}_3^{\ast }$ and from actual dual-frequency observation, in which the horizontal axis is the excess phase Doppler, in Hz, and the vertical axis is the tangent point height, in km. Pink and blue dots represent the excess phase Doppler of the reconstructed ${\mathrm{B}}_3^{\ast }$ and the excess phase Doppler of the actual observation ${\mathrm{B}}_3$. Due to the coincidence of many data dots in Figure 4, it also shows the result of 100× magnification at some positions. At the beginning of occultation, the tangent point of the ray passes through the ionosphere and the upper stratosphere, where the medium is relatively thin, the amplitude of the signal is relatively constant, and the excess phase Doppler is close to 0. As the ray drops in the atmosphere, the vertical gradient with stronger refractive index will cause the ray to become more curved, resulting in excess phase Doppler delay changes. The figure also shows that the overall trend of the excess phase Doppler of reconstructed ${\mathrm{B}}_3^{\ast }$ and that of the actual observation ${\mathrm{B}}_3$ is the same and highly consistent. No matter how many observations are conducted in an occultation event, the excess phase reconstructed by a single frequency is consistent with that of the real observation ${\mathrm{B}}_3$.

Next, the trend change between the reconstructed excess phase and the actual excess phase are compared, a first-order linear fitting of all the reconstructed excess phase Doppler of an occultation event is performed, and then the fitted function is compared with the function y = x. Figure 5 is the comparison of the excess phase Doppler function between BDS single-frequency observation reconstruction and BDS dual-frequency observation under two occultation events. The enlarged details are shown in the small box. The horizontal axis is the excess phase Doppler of real observation ${\mathrm{B}}_3$, and the vertical axis is the excess phase Doppler of the reconstructed ${\mathrm{B}}_3^{\ast }$ (the blue dots). At the same time, the fitting function of the single-frequency reconstructed ${\mathrm{B}}_3^{\ast }$’s excess phase Doppler is represented by a blue curve, and the function y = x, which represents the 100% consistency of the two, is represented by a red curve. As can be seen in Figure 5, when there is more data of an occultation event (Figure 5A), the fitting function of the single-frequency reconstructed ${\mathrm{B}}_3^{\ast }$ excess phase Doppler is closer to the red curve. The two occultation events show that the blue data dots are distributed near y = x, and the red fitting curve is almost identical to y = x. The result of ${\mathrm{B}}_{3 }^{\ast }$ excess phase Doppler reconstruction is very promising.

Similar to the consistency analysis of relTEC change rate in Figure 3, Spearman’s correlation coefficient is also used as a statistical indicator in the consistency analysis of reconstructed excess phases Doppler. Figure 6 shows the statistical results of Spearman’s correlation coefficient histogram of the excess phase Doppler between the reconstructed ${\mathrm{B}}_3^{\ast }$ and the actual ${\mathrm{B}}_3$ observed by FY3E satellite on 4 January 2022, with the histogram interval of 0.05. The figure shows that the Spearman correlation coefficients of the data are all above 0.5, and most of the correlation coefficients are concentrated near 1.0, showing strong correlations. The proportion of Spearman correlation coefficients greater than 0.8 (that is, those showing extremely strong correlations) exceeds 90%, while those with correlation coefficients greater than 0.95 make up nearly 78%. Therefore, the correlation between the reconstructed excess phase Doppler and the actual excess phase Doppler is significant and highly consistent.

In addition to the above-mentioned analysis, the relative average deviation and relative standard deviation of the reconstructed excess phase Doppler are also calculated. Figure 7 shows the histogram of the relative average deviation and relative standard deviation distribution of the excess phase Doppler between the reconstructed ${\mathrm{B}}_3^{\ast }$ and the actual ${\mathrm{B}}_3$. The horizontal axis of Figure 7A is the percentage of the deviation between the excess phase Doppler of the reconstructed ${\mathrm{B}}_3^{\ast }$ and that of ${\mathrm{B}}_3$. The horizontal axis of Figure 7B is the relative standard deviation percentage of the excess phase Doppler of the reconstructed ${\mathrm{B}}_3^{\ast }$ and the actual ${\mathrm{B}}_3$. The figure shows that the relative average deviation of the single-frequency measurement of the reconstructed ${\mathrm{B}}_3^{\ast }$’s excess phase Doppler is less than 1%, the absolute value of most of the relative average deviation is less than 0.2%, and the deviation amplitude is low. The relative standard deviation is less than 1%, and most deviations are around 0.5%. Compared with the excess phase Doppler measured by dual-frequency, the data quality of the excess phase Doppler reconstructed by single frequency is nearly ideal.

To sum up, the reconstructing excess phase ${\mathrm{B}}_3^{\ast }$ is verified by considering the excess phase Doppler. By the analysis of four aspects: direct plot the excess phase Doppler (Figure 4), y = x consistency comparison (Figure 5), Spearman correlation coefficient histogram (Figure 6), and plot the relative average deviation and relative standard deviation (Figure 7). All results show that the reconstructed excess phase Doppler and the dual-frequency excess phase Doppler not only have good consistent frequencies but also have strong monotonic correlations. The standard deviation and average deviation of the reconstructed excess phase Doppler are very low, indicating that the quality of the reconstructed excess phase Doppler is good.

4.3. Accuracy Analysis of Refractive Index Products

The vertical profile of the refractive index is important in occultation data processing. In the ionosphere, the refractive index is directly related to the electron density [36,37]. In a neutral atmosphere, the refractive index is a function of pressure, temperature, and water vapor density [38,39]. In this section, the accuracy of refractive index products measured by single-frequency occultation of BDS signal is analyzed. The reference atmospheric refractive index product is the ERA5 reanalysis field data provided by the European Centre for Medium-Range Weather Forecasts. Specifically, the pressure, temperature, and humidity data in ERA5 reanalysis field data are converted into a refractive index, and then as a function of geopotential height, the consistency between the refractive index of single-frequency occultation inversion and dual-frequency inversion, as well as the deviation between the refractive index of ERA5 reanalysis field data and the refractive index of dual-frequency inversion, are compared through the average deviation index and the root mean square error index, so as to evaluate the accuracy of single-frequency inversion products.

The single-/dual-frequency refractive index inversion results of the BDS B1I and B3I signal from FY3E in the first two weeks of January 2022 are selected. The ERA5 reanalysis field data product is used as the reference value for accuracy statistics. The results are shown in Figure 8. Figure 8A compares the single-frequency and dual-frequency refractive index accuracy in the first week of 2022, and Figure 8B compares the single-frequency and dual-frequency refractive index accuracy in the second week of 2022. The solid line represents the average deviation bias between the inversion result and ERA5 data, the blue solid line represents the percentage deviation of the inversion result of the single-frequency reconstruction, and the black solid line represents the percentage deviation of the dual-frequency inversion result. It can be seen that between the geopotential heights of 10 km and 40 km, the average deviation index of the refractive index produced by single-frequency and dual-frequency inversion is less than 1%, and both are distributed around 0, with a very high degree of consistency. The average deviation index of the refractive index near the Earth’s surface is 2–5% for both single-frequency and dual-frequency inversions, and it quickly converges to 0. The dotted line represents the distribution of root mean square standard deviation between the inversion results and ERA5 data, the red dotted line represents the percentage standard deviation of the inversion results of single-frequency reconstruction, and the blue dotted line represents the percentage standard deviation of the dual-frequency inversion results. It can be seen from the figure that the single-frequency and dual-frequency refractive index inversion results are consistent with the standard deviation root mean square of ERA5 data. From the surface to 25 km, the root mean square error index is highly consistent. The root mean square error index of the single-frequency and dual-frequency inversion results of the part with the geopotential height greater than 25 km has some slight deviation, but the deviation does not exceed 1%, and the standard deviation from ERA5 data is generally less than 6%.

The reason for the slight deviation of the single- and dual-frequency inversion processing results may be that it is necessary to assume the symmetry of the large weather balloon and use the Abel integral transform [40,41] to determine the refractive index from the bending angle in the inversion of the refractive index. Because the Abel transform suppresses high-frequency noise to a certain extent [42], the fractional error of the refractive index is often several times smaller than the corresponding fractional error of bending angle, so in the troposphere, the residual ionospheric noise (not filtered out below 15 km) is not expected to have a significant impact on the error in the retrieved refractive index profile. However, the root mean square error of partial single-frequency and dual-frequency inversion of refractive index above 25 km may come from the smoothed residual ionospheric noise and the prior climatologic situation near the transition layer, which have a certain impact. However, on the whole, the refractive index trend of single/dual-frequency inversion is basically the same, with good consistency.

To sum up, the standard deviation of the single-frequency inversion of refractive index is slightly lower than that of dual-frequency inversion over 25 km, but the system deviation of the single-frequency reconstruction and that of the dual-frequency observation inversion of the refractive index are basically the same.

Through the correctness verification of the single-frequency relTEC change rate, reconstruction excess phase Doppler and the accuracy analysis of the refractive index products, the single-frequency occultation inversion algorithm based on the BDS B1I signal in this paper is shown to be correct and feasible.

5. Conclusions and Discussion

This paper introduces a method of using the ${\mathrm{B}}_1$ single-frequency carrier phase and pseudorange observations to calculate the relative ionospheric TEC, then reconstructs the algorithm of carrier-phase observation ${\mathrm{B}}_3$. The quasi-carrier-phase observation ${\mathrm{B}}_3^{\ast }$ is reconstructed from the pseudorange and carrier-phase observation of the frequency ${\mathrm{B}}_1$. The refractive index product is inverted, and the accuracy of the relTEC, the reconstructed excess phase ${\mathrm{B}}_3^{\ast }$, and the refractive index product in the single-frequency inversion algorithm is evaluated. On the whole, the reconstructed excess phase Doppler and dual-frequency excess phase Doppler not only have good frequency consistency but also have a strong monotonic correlation. The relative average deviation of the reconstructed excess phase Doppler is less than 1%, and most deviations are less than 0.2%. The relative standard deviation is less than 1%, and most standard deviations are less than 0.4%, indicating that the reconstructed excess phase Doppler quality is good. The refractive index of single-frequency and dual-frequency inversion shows good consistency, and the deviation from ERA5 results is less than 1%. The single- and dual-frequency inversion products have comparable accuracies, and the inversion results are satisfactory. The accuracies of dual- and single-frequency inversion products are on the same order of magnitude, with overall comparable accuracy of products, except for the slight decrease in the standard deviation accuracy of single-frequency inversion products over 25 km.

The use of dual-frequency receivers for single-difference or non-difference technology for occultation data processing is maturing. The BDS B1I signal is a high-precision and robust signal that covers the whole world. This paper studies the correctness of the single-frequency occultation inversion algorithm and product availability by taking the occultation observation of the BDS B1I signal by the FY3E satellite as an example. On the one hand, single-frequency occultation processing has the advantage of simple instruments. On the other hand, in practical applications, due to cost considerations or other accidents, the occultation receiver can only rely on single-frequency occultation measurement for GNSS atmospheric occultation detection. At this time, the ionospheric error cannot be eliminated by dual-frequency combination, which requires algorithm verification and accuracy analysis of the single-frequency processing method under current differential technology based on BDS occultation observations.

Overall, the development of a single-frequency receiver for occultation data processing using reconstruction methods can reduce satellite-related costs. With the further advancement of global navigation satellite signals, low-cost and fast atmospheric occultation detection data processing will provide more support for meteorological business applications.