Galileo is the European Global Navigation Satellite System (GNSS) that is currently at the deployment stage. Two early test satellites, called Galileo In-Orbit Validation Element-A (GIOVE-A) and GIOVE-B, were launched in 2005 and in 2008, respectively. These test satellites were used to validate critical technology components, as well as to secure the Galileo frequencies (Montenbruck et al.2015a). GIOVE-A and GIOVE-B were decommissioned in mid-2012 after the successful launch of their successors, i.e., the first pair of full-featured In-Orbit Validation (IOV) Galileo spacecraft (Steigenberger and Montenbruck2017). In total, four IOV spacecraft were launched between 2011 and 2012.

Groundtracks of Galileo satellite positions for the arc length of 26 h with a 15-min sampling interval for the nominal orbit (top) and incorrect orbit (bottom). Color bar corresponds to satellite height above the Earth’s surface. A clear difference between satellite’s velocity in perigee (dark blue) and apogee (red) can be observed for E14 (dots are plotted every 15 min)

Full size image

In 2014, the first pair of Full Operational Capability (FOC) Galileo satellites was launched on a Soyuz ST rocket, however, into wrong, highly eccentric orbits. Instead of a nominal altitude of 23,225 km above the Earth’s surface with an inclination of\(56^{\circ }\) and a revolution period of 14h05m, the satellites were orbiting at heights between 13.713 and 25,900 km with an inclination of 50.1\(^{\circ }\) and a revolution period of 11h42m. A series of satellite maneuvers allowed for a correction of the eccentricity at the cost of substantial fuel consumption. Currently the satellites orbit with a revolution period of 12h56m at heights between 17,178 and 26,019 km (see Fig. 1), which reduces the radiation exposure from the Van Allen radiation belts. Although the perigee was substantially increased, the orbits are still highly eccentric with\(e= 0.1585\), whereas the orbital inclination angle could not be corrected at all. These two FOC satellites cannot be included in the almanac of the broadcast navigation message, because the eccentricity and the semi-major axis exceed the limit of the deviation w.r.t. the nominal orbit, which amounts to\(\Delta e_{\max }=0.03125\) and\(\Delta a_{\max }= 87\) km (Steigenberger and Montenbruck2017).

In 2015 another six and in May 2016 another two FOC Galileo satellites were launched, all of which reached their target orbits. As a result, today there are 14 active Galileo satellites: 12 in three nominal orbital planes: A, B, and C; and 2 in the fourth ‘abnormal’ or ‘extended’ plane.

1.1SLR tracking of Galileo spacecraft

As opposed to the current constellation of GPS satellites, all Galileo spacecraft are nominally equipped with laser retroreflector arrays (LRA) designed for satellite laser ranging (SLR). Galileo LRAs consist of 84 and 60 corner cube retroreflectors for IOV and FOC, respectively. The fused silica corner cubes are uncoated on the rear reflecting side for both IOV and FOC satellites, whereas the front side of IOV cubes is coated by the anti-reflection indium tin oxide^Footnote1 for 532 nm (Dell’Agnello2010). The uncoated corner cubes are preferable, as they increase the return rate signal strength so that the target is easier to acquire by the laser stations (Wilkinson and Appleby2011). The test GIOVE-A and GIOVE-B satellites were equipped with aluminum-coated retroreflectors; thus, their return rates were typically lower than currently for the IOV and FOC satellites.

Tracking high-orbiting GNSS satellites poses a challenge for laser stations: on the one hand due to a weak returning signal in the presence of background noise especially during the daytime and on the other hand due to an increasing number of targets scheduled for SLR tracking. However, SLR stations with a proper automation and scheduling can track more than 140 different satellites (Kirchner and Koidl2015).

The International Laser Ranging Service (ILRS) (Pearlman et al.2002) initiated three GNSS intensive tracking campaigns between 2014 and 2015 in order to improve the scheduling, to increase the data volume especially in the daytime, and to find the best possible tracking strategy which is a trade-off between the maximum number of registered satellite pass segments and the maximum number of tracked satellites (ILRS2014). The GNSS tracking campaigns were coordinated by the ILRS Study Group LARGE (LAser Ranging to GNSS s/c Experiment) which was appointed tasked in 2014 to promote SLR tracking of GNSS satellites.

The SLR observations to GNSS satellites allows validating the GNSS microwave orbits (Pavlis1995; Urschl et al.2007; Montenbruck et al.2015b), combining SLR and GNSS techniques for the improvement of GNSS orbits’ quality (Hackel et al.2015) or using co-locations in space for the estimation of the satellite microwave antenna offsets, LRA offsets, and the scale transfer from the SLR to GNSS solutions (Thaller et al.2011,2014,2015).

1.2IGS MGEX campaign

In 2012, the International GNSS Service (IGS) (Dow et al.2009) established the Multi-GNSS experiment (MGEX) (Montenbruck et al.2013) in order to prepare operational services providing precise orbit and clock products for new and upcoming GNSS, such as the European Galileo, Chinese BeiDou, and regional systems, such as Japanese Quasi-Zenith Satellite System (QZSS) and Indian NAVigation with Indian Constellation (NAVIC) system. The multi-GNSS solutions are a challenge for MGEX analysis centers due to the increasing number of satellites, new signals and observations types, and orbit modeling difficulties caused by different attitude steering modes as well as different satellites’ sensitivities to perturbing forces. Currently seven MGEX Analysis Centers provide orbit and clock products of at least one new GNSS: (1) Center for Orbit Determination in Europe (CODE), Switzerland (Prange et al.2016), (2) Centre National d’Études Spatiales (CNES), France (Loyer et al.2012), (3) Deutsches GeoForschungsZentrum (GFZ), Germany (Uhlemann et al.2016), (4) Technische Universität München (TUM), Germany (Steigenberger et al.2013), (5) Wuhan University (WU), China (Guo et al.2015), (6) Japan Aero-space Exploration Agency (JAXA), Japan (Kasho2014), and (7) European Space Agency (ESA), Europe (Springer et al.2012). Comparisons of orbit and clocks provided by different MGEX Analysis Centers can be found in Steigenberger et al. (2015); Guo et al. (2017), and Prange et al. (2017).

1.3The objective and the structure of this paper

This paper evaluates systematic effects and the quality of Galileo orbits using SLR data with a main focus on Galileo satellites launched into incorrect orbits. Section2 provides general information on the SLR and GNSS data collecting and processing strategies. Section3 describes the current space segment of the Galileo system with differences between satellites and the orbit characteristics. Sections 4 and6 provide an analysis of SLR residuals to individual Galileo satellites, different orbital planes, different types of Galileo spacecraft, and different empirical orbit models. Section5 discusses the issues related to multi-photon and single-photon SLR detector types. The results are summarized and concluded in the last section.

2.1CODE MGEX solutions

CODE, as one of the MGEX analysis centers, was providing three system solutions based on GPS, GLONASS, Galileo in 2012. BeiDou (MEO and IGSO) was additionally included in late 2013. Since 2014 CODE has been providing five satellite system solutions based on GPS, GLONASS, Galileo, BeiDou, and QZSS (Prange et al.2017) on an operational basis. CODE orbit solutions are based on double-difference microwave GNSS phase observations using the ionosphere-free linear combination. For Galileo, the linear combination is formed on a basis of E1 and E5a signals. The orbits are always referred to the middle day of 3-day arcs, because the 3-day solutions are typically associated with much smaller misclosures at the orbit boundaries, as well as better estimates of Earth rotation parameters when new or incomplete GNSS systems are involved (Lutz et al.2016). The observation sampling is 180 s with the elevation cutoff angle equal to 3\(^{\circ }\).

In this paper, we analyze the SLR observation residuals with respect to microwave Galileo orbits generated by CODE for the period between January 2014 and July 2016. In January 2015, CODE changed the orbit model from the classical Empirical CODE Orbit Model (ECOM1) (Beutler et al.1994; Springer et al.1999) to the extended ECOM2 (Arnold et al.2015). In order to keep a full consistency throughout the entire analyzed time series, for 2015 and 2016 the operational CODE MGEX solutions are used, whereas for 2014 the reprocessed orbits using ECOM2 are employed, instead of operational MGEX products based on ECOM1.

The orbit definition consists of 6 Keplerian parameters for the initial epoch of a satellite arc and up to 9 empirical orbit parameters, decomposed into three orthogonal directions: theD axis pointing from the satellite toward the Sun, theY axis lying along the solar panels, and theB axis completing the right-handed coordinate orthogonal frame (see Fig. 2 for an explanation).

Sun–Earth–satellite reference frame showing the elevation angle of the Sun over orbital plane\( \beta \), argument of satellite latitude with respect to the Sun\(\Delta u\) and the elongation angle\(\varepsilon \)

Full size image

The classical ECOM1 has been used for generating high-precision GPS orbits by most of the IGS Analysis Centers. ECOM1 includes the estimation of 5 empirical parameters:\(D_0\),\(Y_0\),\(B_0\),\(B_C\), and\(B_S\) (Beutler et al.1994). Recent studies revealed, however, that for new GNSS satellites, such as GLONASS-M (Arnold et al.2015) or Galileo (Montenbruck et al.2015b; Prange et al.2017) the classical ECOM1 introduces spurious signals into orbit and clock products. Therefore, the new model, i.e., ECOM2, was developed for a better absorption of the solar radiation pressure acting upon the elongated satellite bodies (Arnold et al.2015). ECOM2 was originally derived for GLONASS-M as it remarkably reduces spurious periodic signals in the combined GPS \(+\) GLONASS products.

When using ECOM2, a set of 9 empirical parameters is estimated:\(D_0\),\(Y_0\),\(B_0\),\(D_{C2}\),\(D_{S2}\),\(D_{C4}\),\(D_{S4}\),\(B_C\), and\(B_S\) (Arnold et al.2015). All parameters are estimated without any constraints. No a priori solar radiation pressure model, no albedo model, and no antenna thrust are used for Galileo orbits generated at CODE. The penultimate section of this paper (Sect. 6) provides the comparison between Galileo orbits based on ECOM1 and ECOM2.

It is assumed that all Galileo satellites are subject to the yaw steering attitude control strategy, despite that for very low\(| \beta |\) angles of less than 2\(^{\circ }\) a patented dynamic yaw steering is in fact employed (Ebert and Oesterlin2008). In the dynamic yaw steering mode, the yaw angle is replaced by a smoothed version to keep the maximum slew rate at noon and midnight within the defined threshold, which also avoids discontinuities in the yaw rate and its higher-order derivatives (Montenbruck et al.2015a). Neglect of the dynamic yaw steering modeling may lead to some orbit errors, but only for very low\(| \beta |\) angles.

Besides the satellite orbits, CODE solutions include the following parameters: GNSS station coordinates, 2 h tropospheric zenith path delays, 24 h tropospheric horizontal gradients, and 24 h Earth rotation parameters (X andY pole coordinates, their drifts, and the excess of length of day). The no-net-rotation minimum constraint is imposed on a subset of stable GNSS stations. The same software, i.e., the Bernese GNSS Software 5.3 (Dach et al.2015) is used along with the same a priori models for microwave as well as for SLR solutions. For more details on CODE MGEX solutions, see Prange et al. (2017).

2.2SLR solutions

The SLR observation residuals are computed as differences between laser ranges and the microwave-based positions of GNSS satellites. SLR provides the accuracy assessment mainly for the radial component, because the maximum nadir angle of SLR observations for the Galileo satellites in nominal orbits is 12.3\(^{\circ }\), whereas for the Galileo satellites in highly eccentric orbits the maximum nadir angle varies between 11.3\(^{\circ }\) and 15.5\(^{\circ }\) in apogee and perigee, respectively.

The SLR station coordinates are fixed to the a priori reference frame SLRF2008.^Footnote2 GNSS-derived Earth rotation parameters from CODE MGEX solutions are used for the transformation between the Earth-fixed and inertial reference frames. The LRA offset values provided by ESA have been applied in order to refer the SLR observations to the spacecraft’s center-of-mass.^Footnote3 The station displacement models, including solid Earth tides, ocean tidal loading, and the mean pole definition are consistent with the International Earth Rotation and Reference System Service (IERS) Conventions 2010 (Petit and Luzum2011) and thus also with the microwave-based GNSS solutions.

Daily statistics of SLR observations to Galileo satellites:a number of SLR observations to all Galileo satellites collected by SLR stations as a function of time;b number of SLR observations to specific Galileo satellites as a function of time;c number of SLR observations to specific Galileo satellites as a function of SLR stations

Full size image

We employ a two-stage residual screening for SLR data. In the first step, only the largest residuals exceeding the 1-m level are excluded. Then, themean values and standard deviations (sigma) are calculated. Finally, those SLR residuals which exceed the threshold of\(mean\pm 3\)sigma are excluded and newmean andsigma values are calculated again. More details on SLR validation of GNSS orbits at CODE can be found in Sośnica et al. (2015).

One has to keep in mind that the SLR observations, despite their high precision, may not be fully free of biases (Appleby et al.2016). Therefore, the ILRS Analysis Standing Committee established a pilot project in order to investigate the SLR range biases over the long time spans for LAGEOS and Etalon satellites. Initial results indicate differences between biases to mid-orbiting LAGEOS and high-orbiting Etalons of the order of 20–30 mm for SLR stations equipped with the multi-photon detectors and virtually no biases for SLR stations equipped with single-photon detectors (Appleby and Rodríguez2016). Similar differences in SLR offsets between multi-photon and single-photon stations were also found when analyzing GLONASS observations (Sośnica et al.2015).

Figure 3 shows statistics concerning the SLR observations collected to different Galileo satellites as a function of time and laser station. An outstanding performance, in terms of number of collected data, can be observed for Yarragadee (7090), Mt. Stromlo (7825), Changchun (7237), Graz (7839), Wettzell (8834), Herstmonceux (7840), Zimmerwald (7810), and Monument Peak (7110). The largest number of data was collected to the IOV satellites: E11, E12, E19, especially between August and October 2015, i.e., during the third ILRS GNSS intensive tracking campaign.

The Galileo satellite constellation is planned to be a Walker constellation with 27 satellites in circular orbits distributed over 3 orbital planes separated by 120\(^{\circ }\) with an additional inactive spare satellite in each plane, all of which orbit at the inclination angle\(i=56^{\circ }\) (see Table 1). The satellites’ altitude of about 23,225 km above the Earth’s surface corresponds to a semi-major axis of about 29,600 km and revolution period of 14h05m. Every satellite completes 17 orbital revolutions in 10 sidereal days. The draconitic year of Galileo satellites, i.e., a period of two consecutive passes of the Sun through the orbital plane in the same direction, equals 355.6 days, which is by about 3 days longer than that for GPS satellites. The maximum\(\beta \) angle for Galileo is plane- and time-dependent, but never exceeds 80\(^{\circ }\) (see Table 1). The area of the IOV solar panels and of theX andZ sides of the satellite body is 10.8, 1.3, and 3.0 m\(^2\), respectively. Assuming a mean exposure to the Sun and a mass of 700 kg, the mean area-to-mass ratio equals to 0.019 m\(^2\) kg\(^{-1}\), which is about 21% more than in case of GPS Block IIA satellites. Moreover, large differences between the area of theZ side (pointing toward the Earth) and theX side of the satellite may introduce some systematic variations and thus can be problematic for precise orbit modeling (Arnold et al.2015; Montenbruck et al.2015b).

Currently, plane B of the Galileo constellation is populated by 2 IOV: E11 and E12, and 2 FOC satellites: E26, E22, launched in October 2011 and March 2015, respectively. Plane C contains also 2 IOV: E19 and E20, and 2 FOC: E08 and E09, launched in October 2012 and December 2015, respectively. However, E20 has been transmitting only on one frequency E1 since September 2014 due to a serious power outage. Therefore, the CODE MGEX orbits are not generated for this satellite as the application of the ionosphere-free linear combination is not possible and the SLR validation is also not performed after 2014 due to missing microwave-based orbits. Plane A is occupied by 4 FOC satellites: E24 and E30 launched in September 2015, and E01, E02, launched in May 2016, which are still (August 2016) in the commissioning phase.

Median values and interquartile ranges of SLR residuals for individual Galileo satellites, different orbital planes, and different satellite types

Full size image

Satellites E18 and E14 form another orbital plane, whose orientation changes with respect to the other Galileo orbital planes due to a different revolution period of the ascending node, i.e., 26 years instead of 37 years. The shorter nodal revolution period is caused by the smaller semi-major axis, larger eccentricity, and most importantly, the lower inclination angle. The satellites complete 13 orbit revolutions in 7 days. The masses of E18 and E14 are smaller by about 46 kg with respect to other FOC satellites, due to the fuel loss expended on orbit maneuvers when reducing the orbital eccentricity after the launch. The smaller mass results in increase in the area-to-mass ratio by 7% and thus also in a larger sensitivity to non-gravitational orbit perturbations. The draconitic year of E18 and E14 is about 351.6 days, which is very similar to that of GLONASS orbits. Most of the deficiencies in GNSS orbit modeling appear as spurious variations of the draconitic year period or its harmonics in the GNSS-derived parameters (Rodriguez-Solano et al.2014; Lutz et al.2016). Draconitic years for nominal and abnormal Galileo orbits differ only slightly (by 4 days), which will not be very helpful in separating the draconitic from geophysical signals in the GNSS-derived parameters.

The maximum\(| \beta |\) angle is relatively low for the incorrect orbital planes, equaling to 48\(^{\circ }\). The orbit modeling for eclipsing satellites and satellites at low\(| \beta |\) is typically most challenging, as the Sun illuminates the three surfaces of the satellite body (\(+Z, +X\), and\(-Z\)),^Footnote4 the illuminated cross-sectional area varies most, and the satellite has to slew quickly in order to keep the proper yaw orientation. Finally, one has to keep in mind that neither ECOM1 nor ECOM2 were designed for the solar radiation pressure modeling of the non-circular GNSS orbits or attitude modes deviating from strict yaw attitude, which happens when dynamic yaw steering is applied during deep eclipses.

This section evaluates the quality of Galileo orbits as a function of orbital planes, satellite types, and\(| \beta |\) angles (see Fig. 2). Here, all SLR stations are taken into account without a distinction between the stations equipped with different detector types.

4.1Mean offsets and the RMS of SLR residuals

Figure4 and Table 2 provide the statistics for individual Galileo satellites with the mean and median SLR offsets, mean RMS and interquartile ranges, maximum\(| \beta |\) angle for orbital planes (between 2014.0 and 2016.5), as well as the dependency of SLR residuals on the elongation of the Sun, which characterizes the quality of satellite orbits.

The mean offset of SLR residuals is smallest for E18 and E14, i.e., both in incorrect orbital planes and for E24 and E30, i.e., recently activated FOC satellites. For these satellites, the mean offset is between\(-22\) and\(-26\) mm, whereas for remaining Galileo satellites the mean offset is in the range between\(-36\) and\(-49\) mm. The largest offset and the largest RMS are for E20, i.e., an IOV satellite with the microwave GNSS solutions only in 2014, as it has been transmitting the signal only on one frequency since 2014. In 2014, the ground network of IGS MGEX stations tracking the Galileo signal was inferior as compared to the situation from 2015 and 2016, which is related to the poorer quality of E20 orbits. The RMS of residuals is above 40 mm for most satellites, except for E30, E18, and E12. E30 is a new satellite with a relatively small number of SLR observations; E18 is a satellite launched into incorrect orbit, whereas E12 is an IOV satellite with a relatively large number of SLR observations at high\(| \beta |\) angles. The RMS values of SLR residuals are smaller when considering data only from selected high-performing stations as given in Table 3 (see discussion in Sect. 5).

When comparing different orbital planes, one would expect that for the planes with lower\(\beta _{max }\), the RMS of SLR residuals should be larger due to solar radiation pressure modeling issues. However, this is not the case, as the orbital plane A has the smallest RMS and the smallest mean offset of 28 and\(-25\) mm, respectively, despite having the lowest\(\beta _{max }\). The orbital plane A has recently been populated; thus, the quality of microwave Galileo orbits is high due to a good distribution of the GNSS ground observing network. For planes C and B, the mean offset is\(-44\) mm with the RMS of about 42 mm.

From the comparison between IOV, FOC, and FOC in incorrect orbits, the RMS is 41, 45, and 40 mm and the mean offset is\(-45, -35, \,\mathrm{and}\, -22\) mm, respectively. The satellites in incorrect orbits have reduced masses by 46 kg, due to the fuel loss that was needed for maneuvers which corrected the orbital eccentricity. This raises the question of whether the different offsets in the radial direction, which coincides with theZ axis of the satellite, can be explained by the displacement of the satellite center-of-mass. ESA has published time-dependent center-of-mass values for all Galileo satellites.^Footnote5 According to this information, the difference of the center-of-mass for Galileo E18 between November 2014 and February 2016 is 3.0, 0.0, and 0.0 mm for theX,Y, andZ components, whereas the difference between E18 and other FOC amounts to about 57.4, 4.2, and 0.8 mm for theX,Y, andZ components, respectively. This means that the fuel consumption changes the satellite center-of-mass mostly along theX axis, whereas theZ axis pointing to the Earth and lying along the radial direction, i.e., the direction of the largest sensitivity to SLR observations, changes only insignificantly. In no case can the reduced satellite mass account for the different SLR offset of the order of 13–22 mm.

SLR observation residuals to Galileo E11 as a function of time and the Sun elevation angle\(| \beta |\) above the orbital plane

Full size image

We consider it likely that the SLR offsets shown in Fig. 4 might well be reduced through a proper modeling of the Earth radiation pressure and the transmitter antenna thrust. For example, for GPS and GLONASS, the albedo modeling reduces the offset by about 10 mm, whereas the antenna thrust reduces the offset further by 5–10 mm (Rodriguez-Solano et al.2012; Fritsche et al.2014). Steigenberger et al. (2017) found that the missing antenna thrust accounts for a radial offset of 15 and 26 mm for IOV and FOC, respectively, and 17–32 mm for FOC in incorrect orbits depending whether the satellites are in perigee or apogee. Using the a priori solar radiation pressure model may also reduce the SLR offset (Steigenberger and Montenbruck2017).

4.2Orbit modeling issues

Figure5 shows the time series of SLR residuals to E11 and the absolute value of the Sun elevation angle\(| \beta |\) above the orbital plane. The Sun crosses through the orbital plane every 177.8 days, whereas the passage in the same direction takes place every 355.6 days, i.e., once per draconitic year. The maximum beta angle changes over time due to a slow precession of the orbital plane, whose full revolution lasts on average 37 years (see Table 1).

Despite using the improved model ECOM2, which substantially reduces the systematic\( \beta \)-dependent effects, the spread of SLR residuals is still much larger for the lowest\(| \beta |\) angles. This is on the one hand due to the orbit modeling deficiencies and on the other hand due to the dynamic yaw steering of Galileo satellites for\(| \beta |<2^{\circ }\) (Ebert and Oesterlin2008).

4.3Dependency on satellite–Sun elongation angle

Now, we examine the SLR residuals as a function of the satellite–Sun elongation angle\( \varepsilon \) and the satellite latitude\(\Delta u\) w.r.t. the Sun that are defined as (see Fig. 2):

Dependency of SLR residuals on the satellite elongation angle w.r.t. the Sun when using ECOM2

Full size image

Figure6 shows that all satellite types, IOV, FOC incorrect, and FOC, exhibit a residual dependency on the elongation angle. For high\(| \beta |\) (red and yellow dots in Fig. 6), the elongation angle is close to 90\(^{\circ }\) and the spread of SLR residuals becomes much smaller than for low\(| \beta |\) having typically a large dispersion (blue dots in Fig. 6). When the\(| \beta |\) angle is high, the Sun illuminates the\(+X\) surface of the satellite box for most of the time (\(+X_\mathrm{IGS}\) in Fig. 9 from Montenbruck et al.2015a). For low\(| \beta |\) angles, the Sun illuminates\(+X, +Z\), and\(-Z\) sides, which results in large temporal variations of the surface area exposed to the Sun for elongated satellite bodies, and thus causes difficulties in a proper solar radiation pressure modeling. The twice-per-revolution components in the D direction of the new ECOM2 were added to absorb in particular the effect of the temporal exposure variability (Arnold et al.2015).

The mean slope of the elongation dependency is\(-0.50, -0.43\), and\(-0.80\) mm/\(^{\circ }\) for IOV, incorrect, and FOC satellites, respectively (see Fig. 6). The smallest slope results for incorrect orbital planes, despite that the empirical orbit models were elaborated for circular and not for highly eccentric satellite orbits with variable angular and linear velocities. For the incorrect orbital planes,\(\beta _{max }\) is relatively low as compared to IOV and FOC, and does not exceed 48\(^{\circ }\), which should cause a higher sensitivity of these planes to orbit modeling deficiencies. Despite lower\(| \beta |\), smaller masses, and larger area-to-mass ratios, the orbit quality of Galileo satellites in incorrect planes is remarkably good with the smallest slope of the elongation dependency when compared to other Galileo satellites. High orbit quality of Galileo satellites in incorrect orbits is not only essential for precise geodetic and geophysical applications, but also is a prerequisite for the test of the gravitational redshift employing stable Galileo orbits and clocks (Delva et al.2015).

Prange et al. (2017) found that the interquartile ranges of orbit misclosures in the radial direction are 15 and 11 mm for IOV and FOC in incorrect orbits when using ECOM2, and 38 and 34 mm when using ECOM1, respectively. Smaller orbit misclosure errors in the radial direction for the incorrect orbits confirm again the high quality of Galileo orbits in the incorrect planes. On the other hand, the quality of the along-track component is inferior for the incorrect orbits with a median discontinuity of 26 mm when compared to 3 mm discontinuity of IOV (Prange et al.2017) and it may even increase when the 1-day satellite arcs are generated instead of 3-day arcs. However, SLR is not sensitive to the along-track orbit errors; thus, SLR residuals do not demonstrate any orbit misclosure issues in along-track.

The SLR stations typically employ 3 different detector types:

• Micro-Channel Plate (MCP),

• Photo-Multiplier Tube (PMT),

• Compensated Single-Photon Avalanche Diode (C-SPAD).

MCP and PMT detectors have the capability of registering many photons and thus belong to the type of high-energy detectors. SPAD detectors register single to a few photons with the fast rise time of the avalanche to give good epoch timing and high stability. The improved version of the SPAD detectors is the C-SPAD with a better temperature stabilization and the compensation for time walk. However, in single-photon detectors for each firing of the laser only a single event can be recorded by the C-SPAD, which means that the detector may instead register background noise or daylight as an observation. Noise is reduced by close gating, small fields of view and narrow spectral filters tuned to the laser wavelength.

Satellite signature effect for multi-photon detectors:a the nadir angle is 0 and the LRA centroid coincides with the mean observation,b the LRA is inclined and the detection timing is defined at some threshold level at the leading edge of the return pulse; as a result, the mean registered range is shorter than the actual one

Full size image

Most of the NASA SLR stations are equipped with high-energy MCP detectors, e.g., Monument Peak (7110, USA), Greenbelt (7105, USA), McDonald (7080, USA), Yarragadee (7090, Australia), Tahiti (7124, French Polynesia), as well as Matera (7941, Italy) and Wettzell (8834, Germany). The PMT detectors are typically employed by the Russian stations, such as: Arkhyz (1886), Svetloe (1888), Altay (1879), Zelechukskaya (1889), Baikonur (1887), Komsomolsk (1868, all in Russia), Brasilia (7407, Brazil), as well as by Hartebeesthoek (7501, South Africa). SPAD or C-SPAD detectors are used by the majority of the European SLR stations: Herstmonceux (7840, UK), Graz (7839, Austria), Zimmerwald (7810, Switzerland), Grasse (7845, France), as well as by recently upgraded Chinese stations: Changchun (7237), Beijing (7249), Shanghai (7821, all in China), San Juan (7403, Argentina), and by Mt Stromlo (7825, Australia). Some SLR stations are equipped with multiple detectors, e.g., PMT and C-SPAD in Potsdam (7841, Germany).

MCP and PMT detectors may be affected by the satellite signature effect when the energy is not properly controlled. The satellite signature effect is defined as a spread of optical pulse signals due to reflection from multiple reflectors, which causes variations of mean SLR residuals as a function of satellite’s nadir angle (Otsubo et al.2001). The effective array size, which is the measure of the spread of optical pulse signals due to the reflection from multiple reflectors, is higher for high-energy detection systems, because the detection timing is defined at some threshold level at the leading edge of the return pulse. When the flat LRA is inclined, and not perpendicular with respect to the observer, the first response comes from the nearest edge of the LRA, and thus, the ranges registered are shorter when the threshold level does not allow registering all incoming photons (see Fig. 7).

RMS of SLR residuals and the number of observations to Galileo collected by different SLR stations

Full size image

In the single-photon systems, operating in the low-energy mode, the probability of a photon reflection is the same for every corner cube, and thus, the average reflection point, after several hundreds of reflections, coincides with the center of LRA which in turn corresponds to the centroid of the residual histogram. As each detected photon may originate from any of the retroreflectors, the spatial distribution of the whole array is mapped over many detections (Otsubo et al.2001). Thus, the SLR stations operating in the single-photon mode should not be affected by the issues related to different incidence angles of a laser beam for flat LRAs. However, the single-photon stations require post-processing rejection of high-rate data and responsive energy control by signal reduction filtering when tracking at different elevation angles, because the optical path of a laser pulse through the troposphere varies substantially between the zenith and low-elevation observations, thus leading to intrinsically variable return signal strengths. In addition, the signal variability during a pass may be caused by the telescope and laser pointing, either due to imperfect tracking or variability due to atmospheric turbulence.

Dependency of SLR residuals on the incidence angle of the laser beam with respect to flat LRAs onboard Galileo satellites

Full size image

Figure8 shows the RMS of SLR observation residuals to all Galileo satellites with a distinction between the stations equipped with different detector types. For most of the SLR stations, the RMS values are between 29 and 39 mm. The smallest RMS is obtained for stations, Potsdam, Grasse, and Herstmonceux, and equals to 29, 30, and 32 mm, respectively. All of these stations are equipped with C-SPAD or SPAD detectors or operate in the low-energy return rate regime (Potsdam). PMT stations are typically characterized by larger RMS than MCP and C-SPAD stations. This can be related to the fact that the a priori station coordinates in SLRF2008 are based only on a relatively short time series for newly established stations, e.g., Arkhyz, Altay, Brasilia, and that some PMT stations are characterized by a high RMS of single-shot observations, e.g., Komsomolsk. The mean RMS of SLR residuals for PMT stations is 54.0 mm for all PMT stations and 51.2 mm when excluding Komsomolsk, respectively (see Table 3). For MCP stations, the mean RMS of SLR residuals is 38.4 mm for all stations and 38.2 mm when excluding the Tahiti station which is characterized by the highest RMS value and employs interval timers, whereas most of remaining stations upgraded already to event timers. The RMS values for C-SPAD detectors do typically not exceed 40 mm, except for three Chinese stations: Changchun, San Juan, and Beijing. The mean RMS values for all C-SPAD stations is 41.6 mm for all stations and 35.1 mm when excluding three Chinese stations (see Table 3).

The mean offset equals to\(-43.4, -43.9\), and\(-29.5\) mm for selected MCP, PMT, and C-SPAD stations, respectively (see Table 3). The difference of 14 mm between the offsets for multi-photon and single-photon stations is mostly due to the satellite signature effect. The remaining offset is most likely due to orbit modeling deficiencies, in particular the missing transmitter antenna thrust and Earth radiation pressure models.

5.1Dependence on the incident angle

SLR observation residuals for IOV, FOC, and incorrect Galileo orbits as well as for MCP, C-SPAD, and PMT detectors

Full size image

Figure9 shows four examples of the dependence of the SLR observation residuals on the incident angle of a laser beam. Wettzell is a good representative of all high-energy stations equipped with MCP or PMT: The offset in the nadir direction is about\(-20\) mm and the slope w.r.t. nadir angle is about\(-2.0\) mm/\(^{\circ }\). The slope is a consequence of the satellite signature effect due to the first reflection of laser pulses from the nearest edge of the flat retroreflector. Herstmonceux is one of the SLR stations operating strictly in the single-photon regime. For Herstmonceux and for most of the C-SPAD stations, the nadir offset equals to\(-20\) mm, whereas the slope is slightly smaller than that for high-energy stations and is about\(-1.5\) mm/\(^{\circ }\). The slope for C-SPAD detectors is quite large when compared, e.g., to the GLONASS-M satellites equipped with uncoated corner cubes, for which the slope did not exceed 0.1 mm/\(^{\circ }\). The larger slope may suggest the existence of intensity-dependent bias or some orbit modeling issues as the quality of Galileo orbits is still slightly worse than the GLONASS-M orbits, i.e., the RMS of SLR residuals is 32 and 36 mm for GLONASS-M with uncoated and coated retroreflectors, respectively (Sośnica et al.2015).

Figure9 shows examples of two single-photon stations with an anomalous pattern as compared to other single-photon stations. Changchun has a large negative offset in the nadir direction of\(-93\) mm and the slope of\(+2.6\) mm/\(^{\circ }\), whereas Grasse, equipped with an uncompensated SPAD, has an offset of\(-66\) mm and the slope of\(+2.3\) mm/\(^{\circ }\). The slope is about three times larger than previously observed for C-SPAD stations tracking GLONASS-M with coated corner cubes (Sośnica et al.2015). However, the Galileo satellites have no coating on the rear side of the corner cubes and only the front side of IOV is coated with an anti-reflection layer. Therefore, the anomalously large positive slope rather suggests that there may be some issues with proper background noise filtering or intensity-dependent bias at these SPAD stations. Moreover, Dong et al. (2016) reported recently discovered issues with the stability of temperature around the Changchun station.

5.2IOV, FOC, and Galileo satellites in incorrect orbits

The SLR residuals are a function of the detector type as well as the retroreflector type. The LRA installed on IOV satellites comprises 80 corner cubes of a size of\(33.0\times 23.3\) mm, whereas FOC satellites carry LRA with 60 corner cubes of a size of\(28.2\times 19.1\) mm. LRAs were not designed to be carried on satellites in eccentric orbits whose angular velocity is not constant during a single pass. However, the velocity aberration problems as seen for lunar laser ranging (Otsubo et al.2011) should not be directly observable, as the SLR is mostly sensitive to the radial orbit direction. Nevertheless, all these issues may introduce different characteristics for different satellite types and detector types.

Figure10 shows the SLR residuals for IOV, FOC, and satellites in incorrect orbits collected by selected stations employing MCP, C-SPAD, and PMT detectors. The selected MCP stations include: Yarragadee, Wettzell, Matera, Monument Peak, Greenbelt, and McDonald. The C-SPAD stations include Herstmonceux operating strictly in the single-photon mode, Graz and Zimmerwald collecting one to a few photons and Potsdam operating in the low-energy regime. All of these C-SPAD stations operate in at least 100 Hz mode (Zimmerwald) or kHz mode and employ event timers. PMT stations include Hartebeesthoek and all the stations from the Russian network except for Komsomolsk.

The RMS of SLR residuals is largest for PMT stations and amounts between 48 and 89 mm for IOV and FOC, respectively (see Table 4). The number of collected SLR observations is also remarkably smaller for PMT stations as compared to MCP and C-SPAD, despite that this group includes the largest number of SLR stations.

Interestingly, the smallest offset amounting to only\(-6.7\) mm is obtained for C-SPAD detectors when tracking Galileo in incorrect orbits (see Table 4). The RMS is at the level of 23 mm, which is to some extent similar to the quality of SLR observations to historical GPS-35/36 satellites (offset:\(-12\) mm, RMS: 23 mm) (Sośnica et al.2015). However, only few of these observations were collected for\(| \beta |\) angles lower than 2\(^{\circ }\), whereupon the solar radiation pressure modeling is more demanding and the dynamic yaw steering is not entirely handled in a proper way. The smallest RMS of SLR residuals to incorrect Galileo orbits is 22.3, 23.0, and 25.1 mm for Graz, Herstmonceux, and Zimmerwald, respectively. The RMS for Potsdam tracking incorrect Galileo is only 8.6 mm; however, the number of observations amounts to only 33; thus, the RMS calculation may not be fully reliable. The slope dependency between the SLR residuals and nadir angle is rather detector-specific than satellite type-specific (see Table 4).

The FOC satellites tracked by C-SPAD stations have larger offsets by 19 mm and RMS values larger by 14 mm than Galileo satellites in incorrect orbits despite the comparable number of observations. Therefore, one can conclude that the quality of Galileo satellite orbits launched into incorrect planes is slightly better than IOV and FOC based on the smaller offset of SLR residuals and a remarkably smaller RMS of residuals for the single-photon detectors.

Figure10 implies as well that the SLR residuals should rather be analyzed separately for different detector types instead of a collective analysis for all SLR stations. The individual treatment of satellite center-of-mass corrections for every detector type was proposed and implemented for LAGEOS, Etalon, and Ajisai (Otsubo and Appleby2003) and recently also for LARES (Otsubo et al.2015) on the basis of simulated detected distributions and target response impulse functions from a range of laser widths, total system noise, and return rate threshold settings. Figure10 suggests that the GNSS satellites equipped with flat LRAs should also have detector-specific corrections with a properly considered satellite signature effect so that the SLR measurements represent the LRAs’ centroid, instead of using one offset value for all SLR stations. The alternative for the current approach may be the detection of the signal incoming from individual LRA corner cubes as proposed by Kucharski et al. (2015) for Ajisai.

In order to evaluate the impact of the a priori empirical orbit model, the operational CODE MGEX products from 2014 generated using ECOM1 are compared to the reprocessed products from 2014 using ECOM2.

6.1Dependency on the Sun elevation angle

Figure11 compares the SLR residuals as a function of the Sun elevation\(| \beta |\) and relative Sun–satellite latitude\(\Delta u\) angles when using ECOM1 and ECOM2 for all Galileo satellites active and tracked by SLR in 2014. For ECOM1, there is a very strong dependency between SLR residuals and\(\Delta u\), especially for\(| \beta |<50^{\circ }\). When\(\Delta u\) is close to 180\(^{\circ }\), i.e., when a satellite is in the opposite direction than the Sun as seen from the Earth, the SLR residuals assume the values of about\(-200\) mm. These large negative SLR residuals do not only occur for eclipsing (\(| \beta |<12.3^{\circ }\)), but for all Galileo satellites when\(\Delta u \approx 180^{\circ }\). On the contrary, when\(\Delta u \approx 0^{\circ }\) or 360\(^{\circ }\), the SLR residuals assume maximum values of about\(+150\) mm. The residuals would almost certainly be even larger if the SLR stations were able to track satellites close to the Sun, i.e., when\(\Delta u\)\(\approx 0^{\circ }\) and\( \beta \approx 0^{\circ }\). This is, however, impossible for optical satellite measurements.

SLR residuals as a function of the Sun elevation\(| \beta |\) and relative Sun–satellite latitude\(\Delta u\) angles when using ECOM1 and ECOM2 for Galileo satellites in 2014

Full size image

Dependency between SLR residuals, local time of an observation acquisition, and the satellite–Sun elongation angle for all Galileo satellites observed in 2014 when using ECOM1 and ECOM2

Full size image

Dependency between SLR residuals, local time of an observation acquisition, and the satellite elongation angle for IOV and FOC in incorrect orbital planes observed in 2014.0–2016.5 by 4 European C-SPAD stations when using ECOM2

Full size image

The pattern of the latitude dependency for Galileo satellites using ECOM1 is similar to that observed for GPS Block IIA satellites and opposite to that of GLONASS from Fritsche et al. (2014) (Fig. 10). This is due to different surface areas for theX andZ sides of the satellite body. The surface area ofX andZ is 2.7 and 2.9; 4.2 and 1.7; and 1.3 and 3.0 m\(^2\); for GPS-IIA, GLONASS, and Galileo, respectively (Flohrer et al.2015). When theZ surface area is larger thanX, as in the case of GPS-IIA and Galileo, the SLR residuals assume negative values for\(\Delta u \approx 180^{\circ }\). On the contrary, for theX area larger thanZ, the SLR residuals assume positive values for\(\Delta u \approx 180^{\circ }\), as in the case of GLONASS (Fig. 8 in Fritsche et al.2014).

The switch from ECOM1 to ECOM2 remarkably reduces the spurious SLR residual dependency on the\(\Delta u\) angles (see Fig. 11, bottom). However, the systematic pattern is not completely eliminated, as the residuals are about\(-65\), and\(-15\) mm for\(\Delta u \approx 180^{\circ }\), and\(\Delta u \approx 0^{\circ }\), respectively. Moreover, large negative SLR residuals occur for\(| \beta |<8^{\circ }\) and\(\Delta u \approx 150^{\circ }\). This can partly be explained by the unmodelled dynamic yaw steering (Ebert and Oesterlin2008). However, dynamic yaw steering should be visible for\(| \beta |<2^{\circ }\), whereas it seems that it is employed already for higher\(| \beta |\) angles. Larger residuals could also be explained by the eclipsing seasons, but the satellites enter the shadow already for\(| \beta |<13^{\circ }\), whereas increased residuals are not visible for 8\(^{\circ }<| \beta |<13^{\circ }\). This issue requires thus further investigations.

An alternative improvement of Galileo orbits, instead of using ECOM2, could be obtained through using an adjustable box-wing model (Rodriguez-Solano et al.2014) or an a priori box-wing model as in ESA (Flohrer et al.2015) and CNES (Loyer et al.2012), or by using an a priori model with a different decomposition of main axes for the empirical accelerations as proposed by Montenbruck et al. (2015a,b). However, all these solutions are based on the estimation of empirical orbit parameters which may absorb some geophysical signals, e.g., the geocenter motion. Therefore, another alternative could be the measurement of orbit perturbations using onboard accelerometers as proposed by Kalarus et al. (2016) for future Galileo spacecraft.

6.2Dependency on the observation time

Using ECOM1 results in a strong correlation between the SLR residuals and the local time of acquisition of the SLR observations. In case of GLONASS-M, the SLR residuals for daytime tracking are at the level of\(-30\) mm, whereas for nighttime tracking the residuals are at the level of\(+20\) mm (Sośnica et al.2015). The opposite pattern is found for Galileo satellites when using ECOM1 (see Fig. 12, left): The daytime tracking is associated with the positive residuals of about\(+100\) mm, whereas the nighttime tracking is associated with negative residuals of\(-150\) mm. This pattern is directly related to orbit modeling deficiencies as it can substantially be reduced by using ECOM2. Figure 12 (right) illustrates that using ECOM2 reduces the nighttime and daytime dependency by a factor of four: The daytime residuals are reduced to about\(-10\) mm and the nighttime residuals to\(-70\) mm.

However, few SLR residuals with large negative values do still exist, especially for the nighttime hours, as shown in Fig. 12 (right). Only the best SLR stations are able to track high-orbiting targets during the daytime, whereas almost all SLR stations are capable of tracking GNSS during the nighttime. Therefore, Fig. 13 shows only the SLR residuals from the European C-SPAD stations, which are typically characterized by the highest quality of collected observations (see Sect. 5). Figure 13 (left) shows the SLR residuals to IOV satellites, which suggests that the difference between daytime and nighttime tracking exists even for high-performing stations but is smaller when compared to all SLR stations. The dependency of SLR residuals on the local time of observation acquisition for single-photon stations tracking Galileo satellites in incorrect orbits is smallest out of Galileo satellites (see Fig. 13, right), which confirms a high quality of these orbits. However, a longer time series of observations and a larger number of high-quality SLR observations are needed for a reliable analysis focusing on a quality assessment of individual stations and individual Galileo satellites.

In this paper, we compared SLR observations with respect to microwave-based Galileo orbits generated by CODE in the framework of the IGS MGEX campaign for the period 2014.0–2016.5. The resulting RMS of SLR residuals is between 29 and 39 mm for most of the high-performing SLR stations when tracking all Galileo satellites. The mean offset equals to\(-43.4, -43.9\), and\(-29.5\) mm for selected MCP, PMT, and C-SPAD stations, respectively. Most of the kHz stations and stations that upgraded from interval timers to event timers were able to significantly reduce their biases. However, the SLR residuals depend not only on the detector type, but also on the generation of Galileo satellites. The smallest offset, amounting to\(-6.7\) mm, is obtained for the two Galileo satellites in incorrect orbits and when tracked by the European single-photon stations. The mean RMS is then at the level of 23 mm, which is very similar to the quality of SLR observations to GPS-35/36 satellites. The RMS of SLR residuals to the Galileo satellites in incorrect orbital planes equals 22.3, 23.0, and 25.1 mm for Graz, Herstmonceux, and Zimmerwald, respectively. This confirms the very high quality of Galileo orbits, especially for those satellites which were launched into incorrect orbits.

ECOM1 introduces a very strong dependency between SLR residuals and the satellite’s argument of latitude, which results in mean SLR residual variations between\(-200\) mm when the satellite is in the opposite direction in the sky w.r.t. the Sun and\(+150\) mm when the satellite’s argument of latitude is close to zero and the\(| \beta |\) angle is low. A similar dependency is also visible for the nighttime and daytime tracking: The daytime SLR residuals are equal to about\(+100\) mm, whereas the nighttime residuals equal on average to about\(-150\) mm when using ECOM1.

The new ECOM2 results in a superior performance with respect to the previously used ECOM1 as it substantially reduces the dependence between SLR residuals and satellites’ elongation angles with respect to the Sun’s position and the local time of SLR observations. ECOM2 models much better the impact of direct solar radiation pressure acting on a variable surface area of the satellite boxes exposed to the Sun, especially for Galileo satellites whose side surface areas differ remarkably, i.e.,\(X=1.3\) and\(Z=3.0\) m\(^2\). However, some lesser systematic effects still exist in the SLR residuals especially when\(| \beta |<\)8\(^{\circ }\) even when using ECOM2. The mean slope of the elongation dependency w.r.t. SLR residuals is\(-0.50, -0.43\), and\(-0.80\) mm/\(^{\circ }\) for IOV, incorrect, and FOC satellites, respectively, when using ECOM2. The smallest slope, and thus also orbits of superior quality in the radial direction, is obtained for the two satellites in incorrect orbital planes, despite the fact that the empirical orbit models were developed for circular and not for highly eccentric satellite orbits with a variable aberration due to changeable angular and linear velocities.

The correlations between SLR residuals and the beta angles, elongation angles, or the local time of observation acquisition are related to limitations in modeling of the impact of solar radiation pressure and other non-gravitational forces acting upon GNSS satellites. Modeling of Galileo satellite orbits is very challenging, due to a large area-to-mass ratio of these satellites, which causes high sensitivity of satellite orbits to non-gravitational orbit perturbations. With ECOM1, the modeling errors were much larger than in case when using ECOM2, due to a limited number of estimated empirical orbit parameters, which turned out to be insufficient to appropriately absorb the impact of solar radiation pressure. Extending the number of estimated parameters in ECOM2, especially adding twice-per-revolution parameters in theD direction, substantially reduces spurious signals observed in SLR residuals. However, SLR residuals still reveal some imperfectnesses of ECOM2; thus, further investigation related to a proper modeling of forces acting upon Galileo satellites is indispensable in order to achieve Galileo orbit models freed from any systematic errors.

In 2015, the GREAT experiment (Galileo gravitational Redshift test with Eccentric sATellites) has been launched with the goal to test the gravitational redshift using the onboard atomic clocks on Galileo satellites in elliptical orbits (Delva et al.2015). The unintended elliptic orbits induce periodic modulation of the gravitational redshift, while their very good clock stability allows measuring this periodic modulation very accurately. ILRS supports the GREAT experiment by organizing a series of intensive tracking campaigns for Galileo satellites to make significant progress in disentangling the desired redshift signals from systematics due to the clock, orbital modeling, target signatures, and tracking system specific errors.

High quality of orbital ephemerides for the two Galileo satellites in incorrect orbits is thus not only crucial for precise positioning and geodesy, but also constitutes a prerequisite for the validation of general relativity effects using high-quality Galileo clocks. All orbit modeling deficiencies in the radial direction affect the estimation of GNSS clock products. Fortunately, SLR observations are good proxies of the quality of the radial component and thus may help evaluate and improve both the orbit quality in the radial direction and the clocks, whose utmost quality will be essential for the confirmation of relativistic effects.