© S. Bruinsma et al., Published by EDP Sciences 2017

1. Introduction

Semi-empirical thermosphere models are used in the computation of the atmospheric drag force in satellite orbit determination and prediction, as well as in atmospheric studies. They predict pointwise temperature and (partial) density as a function of location (altitude, latitude, longitude, local solar time), solar and geomagnetic activities, and season. The Committee on Space Research (COSPAR) selected three thermosphere models in 2012, and they are described in the COSPAR International Reference Atmosphere (CIRA-2012) report. These models, NRLMSISE-00 (Picone et al. 2002), JB2008 (Bowman et al. 2007), and DTM2009 (Bruinsma et al. 2012), are also the ISO (International Organization for Standardization) models for the neutral upper atmosphere. The Drag Temperature Model (DTM) was updated recently, and the current version is DTM2013 (Bruinsma 2015). It was developed in the framework of the Advanced Thermosphere Modelling and Orbit Prediction project (ATMOP; http://www.atmop.eu), which was a European Union 7th Framework project. DTM2013 and the pre-ATMOP benchmark DTM2009 were evaluated by comparing to total density data in order to demonstrate and quantify the significant improvements made, especially in the 250–500 km altitude range (Bruinsma et al. 2012; Bruinsma 2015).

The satellite GOCE (Gravity field and steady-state Ocean Circulation Explorer; ESA 1999), was the first Earth Explorer mission of the European Space Agency (ESA). It was launched in March 2009 in a 96.5° inclination, quasi dusk-dawn orbit. The Science Mission phase, maintained at constant mean altitudes between 275 and 230 km thanks to drag compensation using ion propulsion, lasted from 1 November 2009 through 20 October 2013. GOCE neutral densities, with a resolution of 80 km along the orbit and a precision of a few percent (Bruinsma et al. 2014), are made available by ESA (https://earth.esa.int/web/guest/missions/goce/goce-thermospheric-data). Whereas the Bruinsma et al. (2014) paper compared GOCE densities (data version 1_0) to models mainly for validation purposes of the density data up to May 2012, in this study the validated and complete GOCE density dataset (data version 1_4) is exploited to evaluate the models thoroughly.

The GOCE dataset is the only high-resolution, high-accuracy dataset at low altitude. As such, it is the best dataset for thermosphere model evaluation at low altitude. The Xenon for the ion propulsion was exhausted on 21 October 2013, and GOCE re-entered the atmosphere three weeks later on 11 November 2013. The Global Positioning System (GPS) receiver operated almost to the end, and the accelerometers measured the drag acceleration until they saturated on 8 November 2013. Density was derived from these measurements, which provides additional information for model testing at altitudes between 230 and 170 km. This density dataset was compared to the models separately from the ESA density dataset (described in Sect. 2.1.1), because it was computed with different software as part of this study (see Sect. 2.1.2), and because the GOCE satellite was no longer in Science mode.

The objective of this study is to quantify the performance of the CIRA models at low altitudes of 270 km and less, i.e. when re-entry of objects in near-circular orbits becomes imminent, using GOCE density data. Performance is evaluated on long time scales of the order of years, which is important in case of mission design and lifetime (i.e. required thrust and fuel), to typical operational time scales of a month, as well as time scales of the order of days, which is important in case of mission operations and re-entry predictions. Model evaluations have been published in the past, but they may no longer be pertinent because the CIRA models have since been updated and because the quality of the density data is not comparable to GOCE (e.g. Marcos et al. 1983; Hedin 1987; Marcos 1990). Model evaluations are also given in the NRLMSISE-00, DTM2009, and DTM2013 papers, but not in such detail at low altitude. Other studies focused on correcting model error (model calibration), which required quantifying the error, i.e. model evaluation (e.g. Doornbos et al. 2008).

The next section briefly reviews the models and the density data used in the study. The evaluation results are presented in Section 3. A summary of findings and the conclusions of the analyses are given in Section 4.

2. Description of models and data

2.1. Relevant information on NRLMSISE-00, JB2008, and DTM2013

The thermosphere varies on global scales with annual, semiannual, solar rotation (about 27 days), solar active region (a few months), and solar cycle (11 years on average) periodicities. These variations are represented in current models of the thermosphere. The CIRA models were constructed by fitting to their respective underlying density databases as well as possible in the least-squares sense. They are climatology (or “specification”) models of the thermosphere with a low spatial and temporal resolution of the order of thousands of kilometers and hours, respectively. Density variations with smaller spatial scales or shorter time scales are sources of geophysical noise (e.g. gravity waves, upward propagating tides; Bruinsma & Forbes 2008), which, in the absence of all other sources of error, represent the ultimately achievable model accuracy. The solar and geomagnetic indices limit the temporal resolution of these models to 24 and 1–3 h, respectively. The main data sources and solar and geomagnetic indices of each model are listed in Table 1; more detailed information for each model is given in Picone et al. (2002), Bowman et al. (2007), and Bruinsma (2015). The dates and cadences for some of the NRLMSISE-00 datasets could not be retrieved, but the purpose of this table is to show the differences in the databases. DTM2013 is the only model that has partly assimilated the GOCE density dataset, and consequently it is expected to compare best. NRLMSISE-00 has not assimilated any recent total density data. Figure 1 shows predictions of the CIRA models at 250 km altitude for low (mean F10.7 = 75 sfu; 14 December 2009) and medium (mean F10.7 = 144 sfu; 14 December 2011) solar activity in the Science Mission phase of GOCE. The geomagnetic activity kp was less than 2 on both dates. It is a qualitative example of model predictions from which typical differences can be gleaned that are due to the underlying database of each model. DTM2013 and JB2008 predictions are rather similar in shape, with differences in amplitude of about 10–15%. The amplitudes predicted with NRLMSISE-00 lead to about the same differences, but its structure is rather different. This is due to the very large amount of recent satellite-drag inferred densities assimilated in DTM2013 and JB2008 as opposed to incoherent scatter radar data (i.e. a temperature measurement) combined with a small amount of satellite-drag inferred densities assimilated in NRLMSISE-00, and a more (NRLMSISE-00) or less (JB2008) sophisticated modeling of the diurnal variation, e.g. with or without cross coupling with solar and geomagnetic activities.

Fig. 1. Model predictions, in g/cm3, at 250 km altitude close to winter solstice for low (left frames) and medium solar activity.

DTM2013 and NRLMSISE-00 can be used for the entire space age (ap, F10.7 and F30 measurements start before 1957), and they predict temperature, density, and composition. JB2008 is constructed using a combination of solar- and geomagnetic proxies and indices that have been available since 1998 (i.e. it cannot be used before that year), and it predicts temperature and density. An additional weakness of JB2008 is that atmospheric composition, which is required together with temperature to calculate the aerodynamic coefficient of a satellite, is not standard output. Figure 2 shows the 81-day mean F10.7 (NRLMSISE-00 and JB2008), F30 (DTM2013), and S81c (JB2008; integrated 26–34 nm) proxies.

Fig. 2. 81-day mean of the solar activity proxies F10.7, F30, and S81c (left axis), and planetary geomagnetic index ap (right axis).

The way the solar and geomagnetic indices should be used (NRLMSISE-00 and DTM2013) or are used (JB2008) in the models is not entirely clear. Specifically, the models can be run with indices that are constant or interpolated over the day in case of solar proxies, or constant or interpolated over the 3-hour intervals in case of planetary geomagnetic indices ap. There is no consensus or standard method, and this may lead to differences in model predictions of a few percent. In this study, mean (daily) solar proxies were interpolated linearly over the current day (the current day delayed by 24 h) for DTM2013 and NRLMSISE-00. For DTM2013, the ap, delayed by 3 h and the mean of the last 24 h, were interpolated linearly in the 3-hour intervals and then converted to kp. NRLMSISE-00 was run in precise mode, with seven (interpolated) values of ap as specified in the subroutine header. JB2008 code was run with the proxy files (SOLFSMY and DTCFILE) made available by Space Environment Technologies. The user only has to specify the date, and the JB2008 subroutine then returns a density.

2.2. The total density data

The neutral density data that were used in this study are inferred from GOCE data, with additional high-resolution density data from Challenging Mini-Satellite Payload (CHAMP) for specifically testing the performance as a function of local solar time.

2.2.1. Total density dataset of GOCE in Science Mode

Densities were inferred from thruster data for the Science Mission (Doornbos et al. 2014), which lasted from 1 November 2009 through 20 October 2013. GOCE neutral densities have a 0.1 Hz sampling, which leads to approximately 80 km resolution along track, with a precision of a few percent at worst and better than 1% after 2010 (Bruinsma et al. 2014). Due to the quasi dusk-dawn orbit, the local solar time coverage of GOCE densities is limited to 6–8 am & pm (the orbit precessed from 6 to 8 at the end of the mission). The density data is freely available on the dedicated ESA GOCE website (https://earth.esa.int/web/guest/missions/goce/goce-thermospheric-data). In this study, the GOCE data were scaled by a factor of 1.25, as was done in the ESA Final Report (Doornbos et al. 2014), effectively assuming that the High Accuracy Satellite Drag Model (HASDM) of the US Air Force Space Command (Storz et al. 2005) is unbiased. This is not necessarily true because the scale of a thermosphere model is closely related to the aerodynamic coefficients that were used in the derivation of atmospheric density. As a consequence, DTM2013 and JB2008 (which can be considered the climatology version of the HASDM model) are expected to be unbiased with respect to the scaled GOCE densities. Applying a different scale factor than 1.25 will lead to different results of the evaluation, because values of observed and modeled density are analyzed in the form of mean of density ratios (i.e. “model bias”), and not only as standard deviation and correlation. The daily-mean densities and the spherical (above 6378 km) altitude are shown in Figure 3 (the slight scatter in the altitude is due to averaging over incomplete days). The increasing density in 2012–2013 is due to the staggered lowering of the orbit altitude as well as solar activity.

Fig. 3. Daily-mean GOCE densities and mean spherical altitude above 6378 km (red; right axis).

2.2.2. Total density dataset of GOCE in the re-entry phase

The last three weeks of the GOCE mission are not part of the Science Mission, and the GOCE accelerometer measurements in a lower resolution mode were provided by ESA as a special dataset. Densities were derived from the special dataset using the methodology described in Bruinsma et al. (2004) from 22 October to 8 November 2013, and they were also scaled to HASDM for consistency. Their relative precision was estimated at a few percent, i.e. comparable to the thruster-inferred densities. This high precision was possible thanks to the ever-increasing drag acceleration experienced as GOCE decayed, which offset the lower accelerometer resolution in the re-entry phase. Precise kinematic orbit positions were used in a least-squares dynamic orbit adjustment procedure in which the accelerometer bias was estimated every 12 h. The scale factors have been fixed to the values provided in Visser & van den IJssel (2015). Actual scale factors in the re-entry phase are likely different by a few percent, but scaling the densities to HASDM minimizes the already small effect of the error. The satellite macro model for GOCE was simply a rectangular box with a frontal area of 0.7 m² and sides of 10.77 and 5.90 m² in the cross track and radial directions, respectively, and a mass of 1000 kg. The drag coefficient was computed using the Sentman model (Sentman 1961) for flat plates. The macro model was correctly oriented in inertial space thanks to the level-1b attitude quaternions (GO_CONS_EGG_IAQ_2C; Gruber et al. 2012).

The precise kinematic GOCE orbits were obtained by processing the dual-frequency GPS observations collected by the on-board 12-channel Lagrange GPS receiver (Intelisano et al. 2008). The undifferenced carrier phase data were used in the ionosphere-free linear combination for a least-squares adjustment of the epochwise kinematic positions and receiver clock corrections, as well as the carrier phase ambiguities (Svehla & Rothacher 2005). The pseudo-range observations were used only for the initial receiver clock synchronization and to obtain a very first a priori orbit from a kinematic code point positioning. Using the latest version of the Bernese GNSS Software v5.3 (Dach et al. 2015), the GPS data were processed at the full sampling of 1 Hz, together with the attitude data delivered by the on-board star cameras according to Bock et al. (2014). The GPS orbits, the 5 s GPS clock corrections, and the Earth rotation parameters of the Center for Orbit Determination in Europe (CODE; Bock et al. 2009; Dach et al. 2009) were used, to guarantee the highest consistency between the GOCE orbit determination and the GPS orbit determination, which is also performed with the Bernese GNSS Software.

For the performance of the GOCE POD during the Science Mission phase, the reader is referred to Bock et al. (2014). In particular, the independent orbit validation by means of Satellite Laser Ranging (SLR) yielded a mean value of 1 mm and a root mean square (RMS) of 24.2 mm for the SLR residuals of the kinematic orbits over the entire Science Mission phase. During the last three weeks, GOCE was tracked by SLR on three days only, due to the challenge imposed by the low orbital altitude to provide reliable orbit predictions to the tracking network (Jäggi et al. 2011). Table 2 shows the daily statistics of the SLR validation for the kinematic GOCE orbits of these three days. It has to be noted, however, that the number of SLR observations is very limited and the statistics, therefore, not very robust. Figure 4 shows the daily RMS values of the ionosphere-free carrier phase residuals for the last three weeks of GOCE, which is comparable to results of the Science Mission. The quality of the GOCE orbits in the last three weeks is only very slightly lower than during the Science Mission; it allows accurate accelerometer bias estimation to the end of the mission. As a result, the densities are only minimally affected by calibration error, and negligibly by orbit position error.

Fig. 4. RMS of the GOCE GPS ionosphere-free phase residuals.

2.2.3. Total density inferred from accelerometer measurements on CHAMP

Neutral densities were derived from accelerometer measurements on the Challenging Mini-Satellite Payload (CHAMP; Reigber et al. 1996) in the altitude range 450–250 km using the methodology described in Bruinsma et al. (2004). The CHAMP dataset used in this study covers the period 20 May 2001 through 2 September 2010, i.e. up to 17 days before its re-entry in the atmosphere. Due to its precessing orbit plane, 24-hour local solar time sampling was achieved approximately every four months. The accelerometer provided high-resolution measurements from which densities were inferred with about 80 km along track resolution and a precision of a few percent (Bruinsma et al. 2004) nearly from pole-to-pole (87° inclination).

3. Results

The models are evaluated by computing the mean of the density ratios, defined as observed-to-calculated ratios (“O/C”; unity for an unbiased model), the standard deviation and RMS of the density ratios minus one (simply “standard deviation of the density ratios” and “RMS of the density ratios” in the following). The numbers reflect the relative precision of the models; the absolute precision, i.e. in kg/m³, is rather difficult to interpret due to its (orders of magnitude) difference as a function of altitude or solar cycle. Model bias, i.e. the mean of the density ratios minus one, is most damaging in orbit extrapolation because it causes position errors that increase in time. The standard deviation represents a combination of the ability of the model to reproduce the observed variations and the geophysical and instrumental noise in the observations, inclusive of bias in case of RMS. The correlation coefficients R are also computed. Contrary to RMS, they are insensitive to model bias and R ² represents the fraction of observed variance captured by the model. In the next sections, the relative precisions of the models are evaluated on three typical time scales, ranging from years (i.e. solar cycle or part of), to months (i.e. solar rotation), to days.

3.1. Performance on long time scales

The model performance on long time scales is evaluated by computing the density ratios averaged over the entire GOCE Science Mission of nearly 4 years, as well as per calendar year. The results are listed in Table 3. The bias of the models over the entire Science Mission is small, less than 5% for all three models. This remains true only for DTM2013 when the averaging is done over one year intervals (its bias is smallest, 2% or less), whereas the NRLMSISE-00 bias is 11% for 2013. The number of observations varies per year due to data gaps, which is why the results for the entire Science Mission are not the average of the annual results. The standard deviations of the density ratios are smallest for DTM2013, notably in 2010 during the low solar activity conditions and the onset of solar cycle 24; it is less than 10% in 2012 and 2013. JB2008 performs worst when solar activity was low, with largest bias and standard deviation in 2010, and smallest standard deviation in 2013. NRLMSISE-00 on the other hand performs worst when solar activity was highest, with largest bias and standard deviation in 2013. The standard deviation of the density ratios is expected to decrease as solar activity increases (Bruinsma & Forbes 2008). The correlation coefficient of the models is 0.97 or larger for the entire Science Mission. DTM2013 has rather similar correlations per year, whereas NRLMSISE-00 presents the largest variations in annual correlations. JB2008 has the smallest correlation in 2010, confirming its worst performance when solar activity was low. The three models have the highest correlation in 2011 for unexplained reasons.

3.2. Performance on monthly time scales

The monthly evaluation is done because it is an operational time scale (mission planning) as well as a physical one, close to the solar rotation period of approximately 27 days. The performance of the models is expected to decrease when averaging is done over shorter time spans. This is confirmed in Figure 5, which displays the monthly mean density ratios (Fig. 5a), standard deviations and correlations (bottom). The DTM2013 (JB2008; NRLMSISE-00) bias shown in Figure 5a, which was 2% (6.3%; 11.5%) or less for the yearly averages, increases to 8.4% (12.3%; 20.1%) or less. It does not present a clear variation that hints at model error. The same holds for JB2008 except for 2010, which displays a semiannual variation in the density ratios. NRLMSISE-00 displays an increasing bias in 2013, which is probably due to the significantly weaker F10.7 proxy for that year compared to F30 and S81c (Fig. 2). The arithmetic mean and standard deviation of the time series of the monthly values is given for each model, both of which are best for DTM2013.

Fig. 5. Monthly mean density ratios (a), standard deviation (std) of the density ratios (b), and correlations (c). The last (solid) symbol is with the re-entry dataset. JB2008 computed with proxies from June 2016.

The monthly mean standard deviations of the density ratios (Fig. 5b) all present a negative slope in 2010–2011, after which they stabilize for DTM2013 and JB2008 but increase again for NRLMSISE-00 in 2013. DTM2013 and JB2008 achieve very low standard deviations of less than 8% for several months. For the low solar activity conditions in the first few months of 2010, DTM2013 and NRLMSISE-00 have much smaller standard deviations than JB2008, and this may be due to the different variation in activity between S81c, F10.7, and F30 that is visible in Figure 2.

The monthly correlation coefficients (Fig. 5c) are highest with DTM2013 for most of the months. NRLMSISE-00 has the least high correlations, and an annual minimum around January is visible. JB2008 correlations become higher when solar activity increased and orbit altitude decreased. It has the smallest correlation of the three models in the Science Mission for the month January 2010, when it drops to just below 0.86. Periods of diminished correlation are visible for all three models. This will be shown in more detail in the next section.

The last month (#47), plotted in solid symbols, is actually the 3-week re-entry period dataset described in Section 2.1.2. The models achieve the highest correlation (0.98) in that period. JB2008 has the smallest bias and standard deviation with the re-entry dataset, followed by DTM2013 and NRLMSISE-00. The models’ performances compared with the Science Mission and re-entry datasets are the same, i.e. they are the same in the altitude range 275–170 km; the best results are obtained with DTM2013 and JB2008.

3.3. Performance on daily time scales

The model performance on daily time scales is most important for (imminent) re-entry predictions, as well as spacecraft conjunction predictions. The performance of the models decreases again compared to the monthly averages. This is shown in Figure 6, which displays the daily-mean density ratios with DTM2013 (Fig. 6a), NRLMSISE-00 and JB2008 (Figs. 6b and 6c). The re-entry data is plotted in black for the three models. The mean and standard deviation of the time series of daily-mean density ratios are also given for each model. The maximum model biases, and the standard deviations of the time series are much larger than the corresponding numbers of the monthly averages (Fig. 5a). This is mainly due to weaker model performance during periods of enhanced geomagnetic activity or storms, which typically have a duration of hours to a few days; these events are the main cause of the spikes in Figure 6. Inspection of time series of density ratios reveals that the largest biases (spikes) for DTM2013 (JB2008; NRLMSISE-00) on daily time scales are 16–24% (20–28%; 30–38%), and the standard deviation of its time series increases by 59% (57%; 46%) compared to monthly averages. However, the standard deviations are still rather small for DTM2013 and JB2008 (5.9% and 6.4%, respectively), and larger for NRLMSISE-00 (9.4%) due to the drift in the density ratios. The density ratios of the three models computed with the re-entry dataset present large variations. It is largest with DTM2013, which has a clear step about in the middle of the period (31 October), and smallest with NRLMSISE-00 (but biased by 15% on average). The performance of the three models is (visually as well as in terms of mean and standard deviation of the time series, which are given for each model in Fig. 6) in agreement with results for the entire Science Mission.

Fig. 6. Daily-mean density ratios with DTM2013 (a), NRLMSISE-00 (b), and JB2008 (c) computed with proxies from June 2016. The re-entry dataset is plotted in black for all three models.

Due to the small number of geomagnetic storms over the GOCE mission, we cannot glean robust and accurate statistical results on model performance with respect to geomagnetic activity level. DTM2013 has a tendency to overestimate GOCE densities for enhanced geomagnetic and storm conditions by 10–20% (not shown), whereas NRLMSISE-00 by about the same amount underestimates them; JB2008 does not appear to have a bias that is proportional to geomagnetic activity. The reduced model performance with geomagnetic activity can however be demonstrated by sorting the standard deviations of the density ratios and correlations over kp (the daily mean). This is displayed in Figure 7 (a zoom on days for which kp > 2), which shows that for enhanced and storm conditions, the standard deviations of the density ratios increase while correlations decrease. Averages for low activity and enhanced to storm activity were calculated and provided in the figure frames.

Fig. 7. Daily-mean standard deviations (std) of the density ratios (a) and correlations (b) sorted over kp (black; right axis). JB2008 computed with proxies from June 2016.

The daily-mean RMS of the density ratios of the three models and the 81-day mean F30 flux are displayed in Figure 8. The time series are shown for NRLMSISE-00 (Fig. 8a), JB2008 (Fig. 8b), and DTM2013 (Fig. 8c); the average daily RMS of the density ratios of DTM2013 (JB2008 and NRLMSISE-00) is 0.10 (0.12 and 0.14). The relative variability of the thermosphere is highest under low solar activity conditions (Bruinsma & Forbes 2008), which is confirmed by the high RMS in 2009–2010 (for which model biases are small, i.e. RMS is comparable to standard deviation). In 2011, when the solar activity increased sharply, DTM2013 and JB2008 become significantly more accurate, their RMS dropping below 0.10. It is less pronounced for NRLMSISE-00 due to its bias (Fig. 6b); in 2013, the bias becomes larger than 20% and the RMS is largest here too. Solar cycle 24 is weak with a small maximum reached in April 2014. Model performance (relative) is expected to be even better for average and strong cycles, for which the standard deviation (and RMS in case bias is of the order of a few percent) of the density ratios is expected to drop again by 0.01–0.03 at GOCE altitude.

Fig. 8. Daily-mean RMS of the density ratios of NRLMSISE-00 (a), JB2008 (b), and DTM2013 (c), and the 81-day mean F30 (black; right axis). JB2008 computed with proxies from June 2016.

The daily correlations, shown for DTM2013 in Figure 9, reveal several 2–4 month intervals that stand out. The most important ones are centered on December–January (in the black ellipses) and present a clear minimum. The same intervals are detected in the NRLMSISE-00 and JB2008 daily correlations (not shown). The 81-day mean F30 and the daily minus the mean F30 (i.e. how the proxies are actually used in the models) are also plotted in Figure 9 in order to visually verify a possible relation with solar activity. None is evident except for the peak in activity at the end of 2011. A similar inspection with regard to geomagnetic activity (kp; not shown) also does not indicate a relation with the intervals. The minimum reached at the end of 2011 is less profound with JB2008 than with DTM2013 and NRLMSISE-00, which means that at least part of the effect is due to the solar activity proxies used. A spacecraft-related issue (thruster, attitude, and/or accelerometers) affecting the densities can be ruled out because it would have had direct consequences on the orbit altitude for example. That leaves another geophysical source of perturbations, which is waves generated in the lower atmosphere that propagate to the thermosphere. Verifying this hypothesis requires a Whole Atmosphere Model, which is beyond the scope of this study. The effect of waves and tides in the GOCE data has been analyzed (e.g. Gasperini et al. 2015; Häusler et al. 2015), although not in the context of a modification of the mean state of the thermosphere.

Fig. 9. Daily-mean correlations with DTM2013, and the 81-day mean F30 and daily minus the mean F30 (purple; right axis).

3.4. Evaluation over 24-hour local solar time with CHAMP

GOCE was in a quasi Sun synchronous dusk-dawn orbit and, as a consequence, the densities covered only 4 h in local solar time from 6 to 8 am & pm. Therefore, the model evaluations presented thus far may not be representative of performance for all local solar times and consequently be overly optimistic, or pessimistic. As stated earlier, GOCE delivered the only high-resolution and high-precision dataset at low altitude and a dataset with full local time coverage and comparable quality is only available at higher altitudes. The CHAMP densities, in the altitude range 450–250 km (but only the last month below 275 km), were used to verify the model error as a function of local solar time. Figure 10 presents the mean density ratios and the correlations using the nine-year CHAMP dataset binned in local solar time. Except for considerable but quite stable biases of JB2008 and NRLMSISE-00, a typical modeling error (e.g. a diurnal variation) is not evident in the density ratios. The correlations reveal minima for the GOCE orbit configuration and are smallest in the dawn sector. Therefore, we conclude that model evaluation with GOCE densities provides performance results that are representative or even slightly pessimistic, when applied to all local solar times.

Fig. 10. CHAMP density ratios (a) and correlations (b) binned in local solar time. JB2008 computed with proxies from June 2016.

4. Summary and conclusions

GOCE thruster (Science Mission phase) and accelerometer-inferred (3-week re-entry phase) densities were used to evaluate the performance of the CIRA models DTM2013, JB2008 (using the June 2016 version of solar activity proxies) and NRLMSISE-00 in the 275–170 km altitude range. Three typical time scales – annual and longer, monthly and daily – were used in the evaluation in order to address different applications of the models. The range of possible values for bias, standard deviation, RMS, and correlation becomes larger when the evaluation interval decreases; this study provides these statistics based on data in the altitude range of 275–170 km (and probably applicable to a slightly wider range) that heretofore were not documented. DTM2013 performed best of the three models on all time scales. This was expected because it is the only model that has assimilated (part of) the GOCE density dataset as well as the entire CHAMP dataset. NRLMSISE-00 performs worst, with considerable bias in 2009 and 2013, and systematically higher standard deviations (lower correlations) than the other two models. Nonetheless, its performance is exceptional taking the absence of recent and high accuracy density data in its construction into account. JB2008 performance is presently, thanks to the release 4_2g solar activity proxies, very close to that of DTM2013. However, it still presents some weakness under the lowest solar activity conditions in 2009 and 2010.

On the longest (4 years) time scale, all three models have small biases of less than 5% and correlations of 0.97. When the evaluation is performed on annual time scales, the range of calculated biases increases by a few percent, and the correlations remain high but decrease to 0.93–0.96. The standard deviation of the density ratios becomes smaller than 10% with DTM2013 for the years 2012 and 2013.

On monthly and daily time scales, the model biases cover an ever-wider range, roughly doubling from about 10% to 20%, respectively; standard deviation of the density ratios reaches about 8% (6%) for the best months (days) and about 15% (30%) for the worst with DTM2013 and JB2008; the minimum correlations drop below 0.90 and 0.75, respectively. The sometimes significantly diminished model performance on daily time scales is due to geomagnetic enhanced and storm conditions, which have typical durations of hours to a few days.

Model comparisons using the entire Science Mission dataset or the re-entry dataset, with densities from 230 to 170 km, showed statistically the same performance according to our metric. While overall most accurate, the DTM2013 performance is trailing JB2008 in the re-entry phase.

The CHAMP densities were used here to establish the model errors as a function of local solar time, which is not possible with only GOCE due to its quasi Sun-synchronous dusk-dawn orbit with 4 h of local time coverage (6–8 am & pm). The correlations were smallest for the GOCE orbit configuration, i.e. the model errors are not underestimated when comparing to the GOCE densities.

The CIRA model accuracy at low altitude from 275 to 170 km is presently quantified.

Acknowledgments

This study received funding from ESA under Contract 4000115173/15/F/MOS. The authors wish to thank the reviewers for their suggestions to improve the paper. The editor thanks Eelco Doornbos and an anonymous referee for their assistance in evaluating this paper.

References

All Tables

All Figures

	Fig. 1. Model predictions, in g/cm3, at 250 km altitude close to winter solstice for low (left frames) and medium solar activity.
In the text

	Fig. 2. 81-day mean of the solar activity proxies F10.7, F30, and S81c (left axis), and planetary geomagnetic index ap (right axis).
In the text

	Fig. 3. Daily-mean GOCE densities and mean spherical altitude above 6378 km (red; right axis).
In the text

	Fig. 4. RMS of the GOCE GPS ionosphere-free phase residuals.
In the text

	Fig. 5. Monthly mean density ratios (a), standard deviation (std) of the density ratios (b), and correlations (c). The last (solid) symbol is with the re-entry dataset. JB2008 computed with proxies from June 2016.
In the text

	Fig. 6. Daily-mean density ratios with DTM2013 (a), NRLMSISE-00 (b), and JB2008 (c) computed with proxies from June 2016. The re-entry dataset is plotted in black for all three models.
In the text

	Fig. 7. Daily-mean standard deviations (std) of the density ratios (a) and correlations (b) sorted over kp (black; right axis). JB2008 computed with proxies from June 2016.
In the text

	Fig. 8. Daily-mean RMS of the density ratios of NRLMSISE-00 (a), JB2008 (b), and DTM2013 (c), and the 81-day mean F30 (black; right axis). JB2008 computed with proxies from June 2016.
In the text

	Fig. 9. Daily-mean correlations with DTM2013, and the 81-day mean F30 and daily minus the mean F30 (purple; right axis).
In the text

	Fig. 10. CHAMP density ratios (a) and correlations (b) binned in local solar time. JB2008 computed with proxies from June 2016.
In the text