Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information

Available online at www.sciencedirect.com

ScienceDirect Advances in Space Research xxx (2019) xxx–xxx www.elsevier.com/locate/asr

Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information Fabian Schiemenz a,⇑, Jens Utzmann a, Hakan Kayal b a b

Airbus Defence and Space GmbH, 88039 Friedrichshafen, Germany Julius-Maximilians-Universita¨t Wu¨rzburg, 97074 Wu¨rzburg, Germany

Received 21 March 2019; received in revised form 1 July 2019; accepted 16 August 2019

Abstract Errors in neutral atmospheric density are the dominant contributor to unrealistic orbital state-vector covariances in low Earth orbits (LEO). Density uncertainty is caused by model-uncertainty at spatial scales below and within the model resolution, as well as inputuncertainty of the environmental parameters supplied to the semi-empirical atmospheric model. The paper at hand provides multiple contributions. First, analytic equations are derived to estimate the relative density error due to an input parameter uncertainty in any of the environmental parameters supplied to the model. Second, it is shown on the example of uncertain geomagnetic activity information, how to compute the required inputs to facilitate the accurate estimation of the relative density error. A clamped cubic splining approach for the conversion from geomagnetic amplitude (ap) to the kp index is postulated to perform this uncertainty propagation, as other algorithms were found unsuitable for this task. Results of numerical simulations with three popular semi-empirical models are provided to validate the set of derived equations. It is found that geomagnetic input uncertainty is especially important to consider in case of low global geomagnetic activity. The findings seamlessly integrate with prior work by the authors to perform density-uncertainty considering orbit estimation. Ó 2019 COSPAR. Published by Elsevier Ltd. All rights reserved.

Keywords: Geomagnetic activity uncertainty; Atmospheric density uncertainty; Error propagation; NRLMSISE-00; DTM-2012; DTM-2013

1. Introduction Atmospheric drag is the main contributor to orbital uncertainty in Low Earth Orbits (altitude range 300 km 6 z 6 1200 km) (Emmert et al., 2014). Knowledge of the acceleration due to drag is limited mainly by two factors: knowledge of the ballistic coefficient and knowledge of atmospheric density. Both parameters vary with time, however the uncertainty in the atmospheric density is generally greater than the uncertainty in the ballistic coefficient (Hejduk, 2017). When performing orbit estimation, the atmospheric density estimate is obtained from semi-empirical atmo-

spheric models such as NRLMSISE-00 (Picone et al., 2002), JB2008 (Bowman et al., 2008), DTM-2012 (Bruinsma, 2013) or DTM-2013 (Bruinsma, 2015). If data of calibration satellites is available, also atmospheric calibration techniques, such as the High Accuracy Satellite Drag Model (HASDM) can be used to debias the underlying atmospheric model,1 which removes its static nature and leads to a significant improvement in the accuracy of the density estimate (Storz et al., 2005). Density uncertainty arises at spatial scales below and within the model algorithm resolution, which is in the order of
4000 km for NRLMSISE-00, DTM-2012 and DTM2013. Submodel scale spatial density uncertainty

⇑ Corresponding author.

E-mail address: [email protected] (F. Schiemenz).

1

Currently HASDM is only deployed against JB2008.

https://doi.org/10.1016/j.asr.2019.08.023 0273-1177/Ó 2019 COSPAR. Published by Elsevier Ltd. All rights reserved.

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

2

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

(‘‘subgrid-scale density uncertainty”) has been analyzed for example by Bruinsma and Forbes (2008) and Anderson et al. (2009). Grid-scale density uncertainty is also known as average model uncertainty. It is commonly discussed in the respective model reports (e.g. Bruinsma (2015) or Picone et al. (2002)). Despite model-uncertainty at grid- and subgrid-scales, also errors in the inputs to the atmospheric models cause an additional uncertainty-component, which we call ‘‘input-uncertainty”. This type of uncertainty has only received little treatment in the literature. Semi-empirical models typically make use of a spherical harmonics function to model periodic and non-periodic variations in the species partial number densities, exospheric temperature and, in the case of NRLMSISE-00, also in the baseline altitude temperature (Bruinsma, 2015). This function is commonly denoted Gð LÞ, where L is the vector of the environmental parameters that is supplied to the model and contains the required information about latitude, longitude, local solar time, solar activity and geomagnetic activity. For more information on this model-specific function, see for example (Bruinsma, 2015). Of the environmental parameters that are required to compute a density estimate with a semi-empirical atmospheric model, the geomagnetic activity and solar flux information are worthy of an input-sensitivity analysis. Any uncertainty in these parameters will propagate to an uncertainty in atmospheric density, which subsequently propagates to an uncertainty in the orbital mean motion and mean anomaly. In previous research the authors demonstrated the approximate propagation of solar flux errors to relative atmospheric density errors (Schiemenz et al., 2019a). Emmert et al. (2017) have shown how to make use of relative density errors to obtain estimates of the resulting mean motion and mean anomaly uncertainty. Recent work by the authors also demonstrates how to obtain a density uncertainty covariance matrix in the Geocentric Celestial Reference Frame (GCRF) that can be used for density uncertainty-considering orbit determination (Schiemenz et al., 2019b). The paper at hand is divided into two parts: first we present improved equations for the accurate estimation of the relative density error and its variance due to an error/a variance in the spherical harmonics function of the environmental parameters (Gð LÞ). The error in Gð LÞ (denoted dGð LÞ) and its variance (VarðdGð LÞÞ) can be caused by any of the environmental parameters that populate L, for example the solar flux information or the geomagnetic activity information. As we have already demonstrated the propagation of solar flux errors/variances to errors/variances in Gð LÞ in Schiemenz et al. (2019a), the second part of this paper is devoted to the derivation of dGð LÞ and VarðdGð LÞÞ due to geomagnetic index uncertainty. The combination of both contributions allows to perform the estimation of relative density uncertainty for the case of geomagnetic activity uncertainty. The gained

insight can then be used, for example, to form a physicsbased process noise covariance matrix for the purpose of orbit propagation or determination, as shown in Schiemenz et al. (2019b). Geomagnetic activity parameters are used in semiempirical density models as an indicator of either the global state of the magnetosphere or specific parts of it. The most popular indices for orbit estimation are Ap; Kp; Km and DST, however many more exist. An excellent description of the various geomagnetic activity parameters and their derivation can be found in Mandea and Korte (2010, Section 8). All of these are subject to various types of uncertainty, such as station uncertainty, physical uncertainty and data processing uncertainty, as has been elaborated by Xu (2008). Consequently the geomagnetic parameters obtained from space-weather databases contain an inherent uncertainty. This is especially true if real-time orbit estimation is performed with quick-look or forecast indices. Some authors use capitalized first letters to denote daily average parameters (e.g. Ap) and lower-case letters to denote measurements or n-hourly indices (e.g. ap). We stick to this convention, however not without noting the indices used by the models considered in this paper, which are: daily Ap for NRLMSISE-00 and interpolated 3 h-delayed kp as well as 24 h average Kp for DTM-2012 and DTM-2013. To the authors’ knowledge, the work of Xu (2008) is the only published study which both discusses magnetic index uncertainty and gives quantitative indications of the resulting magnitude. Motivated by the work of Xu, we assume the magnetic index errors to follow a Gaussian distribution with a standard deviation of 0.5 Ap units in geomagnetic amplitude, well noting that future studies may be helpful to further improve the characterization of the uncertainty in the magnetic indices. The remainder of the paper is structured as follows: Section 2 is devoted to the derivation of the equation describing the relative exospheric density error due to an absolute error in the spherical harmonics function (dGð LÞ), including the corresponding variance computation. Subsequently, Section 3 demonstrates how obtain the required estimates of dGð LÞ and VarðdGð LÞÞ for the case of magnetic index uncertainty and the semi-empirical models NRLMSISE-00, DTM-2012 and DTM-2013. We validate our findings with the help of numerical Monte-Carlo simulations in Section 4. A conclusion finalizes the paper. 2. Relative exospheric mass density error In a mixture of gases, mass density is defined as the product of the mean mass and the overall number density of the species contained in the gas: q ¼ nM Mean mass is defined according to Eq. (2): Pnspecies i¼1 ni M i M¼ P nspecies i¼1 ni

ð1Þ

ð2Þ

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

where M i denotes the mass of the species and ni denotes the number density of the species. In the following we use the term density to refer to mass density. The relative exospheric density error is defined as the absolute error (‘‘truth - estimate”), divided by the density estimate, as given in Eq. (3). A hat-sign is used to indicate estimated quantities, i.e. the output obtained from a semi-empirical model. q ¼

^ dq q qq ¼ ¼ 1 ^ ^ ^ q q q

ð3Þ

Similarly for the relative error in the number density and the mean mass: ^  M n q=M q þ 1  1 n ¼  1 ¼ 1¼ ^ ^ M n ^ =M q ^ ^ M M ¼ q þ  1 M M ^ M M 1 M ¼  1 ) ¼ ^ M M þ 1 M

ð5Þ

Inserting Eq. (5) into Eq. (4) yields: 

n ¼ M qþ1 þ M1þ1  1 () q ¼ n þ M þ n M

 number density and mean mass with respect to the temperature at the baseline altitude  number density and mean mass with respect to the number density at the baseline altitude In Schiemenz et al. (2019a) the temperature is assumed to be fixed at the exospheric temperature, which is one of the reasons why a tuning-function is needed. All of the models considered in this paper set the reference altitude to 120 km and assume a Bates temperature profile, as introduced in Bates (1959)):   T 00 f T ðfÞ ¼ T ex  ðT ex  T 0 Þ exp  T ex  T 0 ¼ T ex  ðT ex  T 0 Þ expðrfÞ

ð4Þ

ð6Þ

Eq. (6) allows to compute the relative density error based on the underlying relative errors in the mean mass and number density. Since the masses of the individual species that are considered in a semi-empirical density model are known, we only need a description of the number density per species and the sensitivities of this number density with respect to three quantities in order to evaluate Eq. (6). The quantities of interest that determine the density computation are:  exospheric temperature (T ex )  temperature at the baseline altitude (T 0 )  number density at the baseline altitude (n0 ) Each one of these depends on the spherical harmonics function of the environmental parameters Gð LÞ. For the LEO altitude range of interest, the equation governing the number density estimate at geopotential height f (which corresponds to altitude z) is given by Schiemenz et al. (2019a):   Z     T ð fÞ g0 M i f 1 0 ln nðfÞi ¼ ln ni0  ð1 þ ai Þ ln  0 df T0 k 0 T ðf Þ ð7Þ where  ai : thermal diffusion coefficient of the species  g0 : gravitational acceleration at the baseline altitude The six required derivatives are therefore:  number density and mean mass with respect to exospheric temperature

3

ð8Þ

T0

0 can be considered as an inverse scale height and r ¼ T ex T 0 is also known as ‘‘shape factor” (Emmert, 2015). The quantity T 00 represents the temperature gradient at the baseline altitude (commonly in units of [K/km]). The solution of the integral of the inverse temperature profile with respect to geopotential height, as required in Eq. (7), is derived in Appendix A and allows to obtain a closed-form description of the number density versus altitude. The solution (Eq. (A.8)) reads:   T0 Z f ln T ðf Þ 1 f 0  ð9Þ 0 df ¼ T ex rT ex 0 T ðf Þ

Introducing Eq. (9) to Eq. (7) then yields the desired description of ni ðfÞ, as for example published in Chamberlain and Hunten (1990) or Picone et al. (2013):   T0 g Mi ni ¼ ni0 T ðfÞ1þai þci exp  0 f ð10Þ kT ex where ci ¼

g0 M i kT ex r

ð11Þ

with k denoting the Boltzmann constant. Eq. (10) is valid for all models assuming hydrostatic equilibrium and a Bates temperature profile (e.g. NRLMSISE-00, DTM2012 and DTM-2013). In all other cases (e.g. JB2008), its usage to derive relative density error estimates introduces undesired errors, however is still more accurate than the assumption of a constant thermospheric temperature, as applied in Schiemenz et al. (2019a). In the following we make use of Eq. (10) to derive the sensitivities of the density with respect to n0 ; T ex and T 0 . 2.1. Relative exospheric number density error The relative number density error is given to first order by

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

4

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

Pspecies i

n ffi

@ni @ni0

dni0

^ n

Pnspecies þ

i

@ni @T ex

^ n

dT ex

Pnspecies þ

i

@ni @T 0

^ n

dT 0 ð12Þ

Eq. (12) expresses the relative number density error as the sum of the individual error contributors, namely the exospheric temperature error, baseline temperature error and baseline number density error. Each of the terms needs to be formulated in order to derive an expression for n . 2.1.1. Sensitivity of number density with respect to number density at the baseline altitude The models considered evaluate the number density at the baseline altitude according to Eq. (13): X

nspecies

n0 ¼

ni00 expðGi ð LÞÞ

ð13Þ

2.1.2. Sensitivity of number density with respect to exospheric temperature The sensitivity of the number density with respect to exospheric temperature requires the computation of the derivative   g M @ @ T0 g0 M i 1þai þkT0ex ri ni ¼ ni T ð fÞ exp  f @T ex @T ex 0 kT ex gð xÞ

This can be written as dxd af ð xÞ expðcð xÞÞ. The corresponding rule for the differentiation is (see Appendix B for the derivation): d gð xÞ af ð xÞ expðcð xÞÞ dx

d f ð xÞ d d gð xÞ dx cð xÞ þ lnð f ð xÞÞ gð xÞ þ gð xÞ ¼ af ð xÞ expðcð xÞÞ f ð xÞ dx dx ð18Þ

i¼1

where i denotes the ith species (helium, hydrogen, atomic oxygen, atomic nitrogen, molecular oxygen and molecular nitrogen for the DTM-models and additionally argon, as well as ionized oxygen for NRLMSISE-00) and ni00 the global and temporal average number density of the species at the baseline altitude, prior to modification via the exponential of the spherical harmonics function. ni00 is a modeldependent coefficient. An error in the environmental parameters results in an error in n0 according to:        ^ þ dL  exp Gi L ^ dni0 ¼ ni00 exp Gi L    ^ ðexpðdGi Þ  1Þ ¼ ni00 exp Gi L |fflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflffl}

Using x ¼ T ex ; a ¼ ni0 ; f ðT ex Þ ¼

gðT ex Þ ¼ 1 þ ai þ

ð14Þ

M i g0 ðT ex  T 0 Þ ¼ 1 þ ai þ c i ; kT ex T 00

we therefore obtain:

ð15Þ

@ @T ex

   1þai þci T0 M i g0 f ni ¼ ni0 exp  T ðfÞ kT ex |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} ni

"

The absolute error is then given by the introduction of Eq. (15) and summation over the species:

M i g0 f kT 2ex

¼ ni

nspecies

dn0 ¼

ni ðexpðdGi Þ  1Þ

;

M i g0 f kT ex

cðT ex Þ ¼ 

From Eq. (10) it follows that

X

T 00 f T ex T 0

and

:¼ni0

@ni ni ¼ @ni0 ni0

T0 T0  ¼ T ðfÞ T ex  ðT ex  T 0 Þ exp 

ð16Þ

i¼1

Inserting Eq. (16) with the model-estimated number density into the first summand of Eq. (12), we obtain: Pspecies ^ ni ðexpðdGi Þ  1Þ i n0 :¼ ^ n Pspecies ^ ni expðdGi Þ 1 ð17Þ ¼ i ^ n Eq. (17) describes the sensitivity of the number density at the altitude of interest with respect to the number density at the baseline altitude, divided by the absolute number density estimate.

h

M i g0 kT 2ex

þ MkTi g0 0TT2 0 ln 0 ex





T0 T ðfÞ

 ð 1 þ ai þ c i Þ

T 0 ½1exp ðrfÞð1þrfÞ T ðfÞ2 T0 T ðfÞ

#

ð19Þ

i    i þci f þ TT 00 ln TTð0fÞ  1þa ð1  expðrfÞ½1 þ rfÞ T ðfÞ 0

From the definition of the Bates temperature profile (Eq. (8)) it follows that expðrfÞ ¼

T ex  T ðfÞ T ex  T 0

Hence, Eq. (19) can be written as: "    @ M i g0 T0 T0 ni ¼ n i f þ 0 ln @T ex T0 T ð fÞ kT 2ex  

1 þ ai þ c i T ex  T ðfÞ 1  ð1 þ rfÞ T ex  T 0 T ð fÞ

ð20Þ

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

Defining the second summand of Eq. (12) as N T ex , we obtain the contribution of an exospheric temperature error to the overall particle number density by N T ex

¼

i

^ n

@T ex

¼ dT ex kTg02



ex

 

i

dT ex

f þ TT 00 ln 0

Pnspecies i

 

T0 T ðf Þ T

ni M i

^ n

Pnspecies

a¼

i

ni ai

n



 T ðfÞ



ð21Þ

dT ex 1 Texex T ð1þrfÞ 0

T ðf Þ

ni ð1þai þci Þ ^ n

;

ð22Þ

and the (number density weighted) mean ratio of temperature and species scale heights as Pnspecies ni c i c¼ i ; ð23Þ n yields the desired description of N T ex upon evaluation with the estimated quantities: "

N T ex

!! T^ 0 T^ 0 f þ 0 ln T^ 0 T^ ðfÞ  # T^ ex  T^ ðfÞ 1þ^ a þ ^c ^fÞ dT ex 1  ð1 þ r T^ ðfÞ T^ ex  T^ 0

^ gM ¼ 0 k T^ 2ex

ð24Þ

2.1.3. Sensitivity of number density with respect to baseline temperature Similarly as for the exospheric temperature, we make use of Eq. (18) with the definitions x ¼ T 0; a ¼ ni 0 ; f ðT 0 Þ ¼

ni

 kTM0iTg0ex 0

Making use of the definition of mean mass (Eq. (2)) and defining the (number density weighted) mean thermal diffusion coefficient as Pnspecies

  1þai þci T0 M i g0 f ni ¼ ni0 exp  kT ex T ð fÞ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} 

@ @T 0

"

Pnspecies @ni

T0 T0  ¼ T ðfÞ T ex  ðT ex  T 0 Þ exp 

T 00 f T ex T 0

5

 ln

T0 T ðf Þ



þ ð 1 þ ai þ ci Þ

T ðfÞT 0 expðrfÞ½1þrf T ðfÞ2

which can be simplified to   @ M ig T0 ni ¼ ni0  0 0 ln @T 0 kT 0 T ex T ð fÞ T ex T ðfÞ 1 þ ai þ ci T ðfÞ  T 0 T ex T 0 ½1 þ rf þ T0 T ð fÞ

#

T0 T ðfÞ

# ð25Þ

With the help of Eq. (25) we can now evaluate the third term of Eq. (12), which we denote N T 0 : Pnspecies @ni dT 0 N T 0 ¼ i ^n @T 0   Pnspecies M n i i i ð26Þ ¼ dT 0 kTgex0T 0 ln TTð0fÞ ^ n 0 Pnspecies T ex T ðfÞ T ðfÞT 0 T ex T ð1þrfÞ ni ð1þai þci Þ i 0 þ TdTðf0Þ ^ T0 n The evaluation of Eq. (26) with estimated quantities from an atmospheric model and Eqs. (2), as well as (22) and (23) with the estimated number density, results in: "   ^ T^ ðfÞ g0 M 1 þ ^a þ ^c N T 0 ¼ ln  k T^ ex T^ 00 T^ 0 T^ ðfÞ 

T^ ex  T^ ðfÞ T^ ðfÞ ^ fÞ   ð1 þ r ð27Þ dT 0 T^ 0 T^ ex  T^ 0 Inserting Eqs. (17), (24) and (27) into Eq. (12) yields a first-order estimate of the relative error in the total number density, caused by an absolute error in the function of the environmental parameters (dGi ð LÞ): h ^  i   ^ ^ T^ ex T^ ðfÞ aþ^c ^ 1  ð 1 þ r f Þ n ffi kgT0^M2 f þ TT^ 00 ln T^Tð0fÞ  1þ^ ^ ^ ^ T T T ðf Þ ex

0

ex

0

dT ex ðdGT ex ð LÞÞ h ^ ^   i T^ ðfÞ aþ^c T^ ex T^ ðfÞ ^ þ kT^g0 MT^ 0 ln TT^ðfÞ  1þ^ ð 1 þ r f Þ  ^ ^ ^ ^ T0 T ex T 0 T ðf Þ ex 0 0   dT 0 dGT ð LÞ Pnspecies 0 ^ ni ðexpðdGi ÞÞ þ i 1 ^ n ð28Þ

;

and

Eq. (28) only depends on the outputs of the atmospheric model and model-dependent coefficients (e.g. the thermal diffusion coefficients). It can be evaluated with a single model call, once dGi ð LÞ is available and dT ex and dT 0 are computed as shown in Appendix C.

cðT 0 Þ ¼ 0

2.2. Relative exospheric mean mass error

Therefore the derivative of the number density profile with respect to the baseline temperature is given by:

To complete the estimation of the relative density error, we also need to derive the sensitivities of the mean mass

gðT 0 Þ ¼ 1 þ ai þ

M i g0 ðT ex  T 0 Þ ¼ 1 þ ai þ c i ; kT ex T 00

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

6

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

with respect to T ex ; T 0 and n0 . The relative error in the mean mass is given to first order by Eq. (29): Pnspecies @M @M dni0 @M dT ex @T dT 0 i @ni0 þ @T ex þ 0 ð29Þ M ffi ^ ^ ^ M M M

!

X

nspecies @ @T ex

@ @T ex

M i ni

n

Pnspecies

¼

i

@ n @T ex

M i ni

i

M¼

!

X

nspecies

i n2

M i @T@ex ni

n

X

nspecies

 Mn

@ @T ex

ni

i

2.2.1. Sensitivity of mean mass with respect to number density at the baseline altitude Using the definition of mean mass (Eq. (2)), we can compute the partial derivative with respect to the species number density at the baseline altitude: ! ! nspecies nspecies X X @ M j nj n M j nj @n@ n Pnspecies @ni i0 0 M j nj j j j @M @ ¼ ¼ 2 @ni0 @ni0 n n Pnspecies @ nspecies X M j @n nj j i0 @ ¼  Mn nj n @ni j

0

The derivative with respect to the i th species is only nonzero if j ¼ i, otherwise it is zero. We can therefore drop the sum and set the index equal to i:

The derivative of the number density with respect to exospheric temperature is already available (Eq. (20)) and needs to be inserted. Hence: Pnspecies

The derivative @M ¼ @ni0

M i nnii

@ @ni0

ni is given in Eq. (15). Hence:

ex

i

n

The contribution of an error in the number density at the baseline altitude therefore becomes:

h i Pnspecies M i nni0i M ni  n ni dni0 Pnspecies M i ni  M ni dni i n i 0 0 ni0 ni0 ¼ M n0 ¼ ^ ^ M nM

¼

i

Pnspecies i

M i n M n i0

i0

ni0 ðexpðdGi Þ1Þ

^ nM

M i ni expðdGi ÞM i ni Mni expðdGi ÞþMni ^ nM

which, using estimated quantities and the definition of mean mass, leads to Pnspecies Pnspecies ^ M i^ ni expðdGi Þ ni expðdGi Þ i M n0 :¼  i ð30Þ ^ ^ n ^ nM 2.2.2. Sensitivity of mean mass with respect to exospheric temperature To obtain the contribution of an error in the exospheric temperature to a relative error in the mean mass, we pro@M ceed similarly by computing the derivative @T : ex

T0 T ðfÞ

1þai þci T ðfÞ





T

T ðfÞ 0

i

1 Texex T ð1þrfÞ

0

T

þ

T ðfÞ

1 Texex T ð1þrfÞ 0

T ð fÞ

  Pnspecies Pnspecies M i ni ai M i ni ci i i  Ma þ Mc   n n The contribution of the absolute exospheric temperature error to the relative mean mass error is given by: @M @T ex

dT ex ^ M

which, upon evaluation with the estimated quantities, results in:  Pnspecies    M 2i ^ ni g0 T^ 0 T^ 0 i ^ P M T ex ¼ f þ T^ 0 ln T^ ðfÞ kT^ 2 M nspecies Mi^ ni ex 0 i   T^ T^ ðfÞ ð1þ^ rfÞ T ex T^ 0

1 ^ex

þ

ð31Þ

T^ ðfÞ

 

Pnspecies Pnspecies Mi^ ni ai Mi^ ni ci ^a þ ^c  i ^nM^  i ^nM^ dT ex

Introducing the definition of the absolute error in the species number density from Eq. (14), results in Pnspecies h ni ni i M n0 ¼

 

T0 ln T0 0

fþ

which can be simplified to Pnspecies 2     M i ni g0 T0 T0 2 @ i M ¼ f þ ln  M @T ex T 00 T ðf Þ n kT 2ex  

M T ex :¼

M ni  n ni0

0



M i g0 kT 2ex

¼ n nspecies h    i X Mg  T ex T ðfÞ i þci  Mn ni kTi 2 0 f þ TT 00 ln TTð0fÞ  1þa 1  ð 1 þ rf Þ T ð fÞ T ex T 0

@

@M M i @ni0 ni M @  ¼ ni @ni0 n @ni0 n

i

@ M @T ex

h

M i ni

2.2.3. Sensitivity of mean mass with respect to baseline temperature Finally, we need to compute the derivative of the mean mass with respect to the baseline temperature: Pnspecies species M i @T@ 0 ni M nX @ @ i  M¼ ni @T 0 n i @T 0 n The derivative of the number density with respect to the baseline temperature is given in Eq. (25), which leads to:

  T T ðfÞ T ðfÞT 0 ex ½1þrf Pnspecies T ex T 0 M i g0 T0 1þai þci @ M @T 0

M i ni0 

i

¼

ln

T ðfÞ

þ

T ðf Þ

T0

n

X

nspecies



kT 0 T ex 0

M n

i



  T ex T ðfÞ M i g0 1þai þci T ðfÞT 0 T ex T 0 ½1þrf T0 ni0  kT 0 T ex ln T ðfÞ þ T ðfÞ T0 0

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

Simplification yields  Pnspecies 2  @ g0 T ðfÞ M i ni i 2 M M¼ 0 ln @T 0 T0 kT 0 T ex n Pnspecies T ðf ðfÞ  TTexexT ð1 þ rfÞ M i ni ai T 0 T0 Ma þ Mc  i  T ðfÞ n ! Pnspecies M i ni c i  i ð32Þ n

"

!! T^ 0 T^ 0 f þ 0 ln T^ 0 T^ ðfÞ   T^ ex  T^ ðfÞ 1 þ ^a þ ^c ^ fÞ 1  ð1 þ r T^ ðfÞ T^ ex  T^ 0 !!Pn  species T^ 0 T^ 0 g0 M 2i ^ni ^ i þ f þ 0 ln M Pnspecies M i ^ni T^ 0 T^ ðfÞ k T^ 2ex i

^ gM q ffi 0 ^ k T 2ex

^

The contribution of an absolute error in the baseline temperature to a relative error in the mean mass is to first order given by Eq. (33): M T 0 ¼

@M @T 0

dT 0 ^ M

ð33Þ

which, upon evaluation with estimated quantities, yields "   Pnspecies 2 T^ ðfÞ g0 Mi ^ ni ^ i M T 0 ¼ ln  M P nspecies M i^ ni T^ 0 k T^ 00 T^ ex i 

T^ ðf T^ 0

^

^

^fÞ  TT^ex TT^ðfÞ ð1 þ r ex

0

T^ ðfÞ !# Pnspecies Pnspecies M i^ M i^ ni ai ni c i i i  ^ a þ ^c   dT 0 ð34Þ ^ ^ ^ ^ nM nM Inserting Eqs. (30), (31) and (34) into Eq. (29, allows to obtain a single expression for the relative error in the mean mass caused by an absolute error in the spherical harmonics function of the set of environmental parameters: " ! !! Pnspecies 2 ^ T^ 0 T^ 0 g0 M n i ^ M ffi f þ 0 ln Pinspecies i  M M i^ ni T^ 0 T^ ðfÞ k T^ 2ex i ^

^

!# Pnspecies ^ ^ M a M c n n i i i i i i  ^ a þ ^c  i  i ^ ^ ^ ^ nM nM " ! T^ ðfÞ g0 ln dT ex ðdGðLÞÞ þ T^ 0 k T^ 00 T^ ex Pnspecies

Pnspecies 2 M ^ ni ^  Pinspecies i  M ^ M n i i i



T^ ðf T^ 0



T^ ex T^ ðfÞ ð1 T^ ex T^ 0

B B B B T^ B B B @

 T^ ðfÞ

Pnspecies þ

i

^ aþ^c

7 7 7 f 7 7dT 0 ^ ex T^ ðfÞ T 7 T^ 0  ^ ^ ð1þ^ rfÞ T ex T 0 7 Pnspecies Pnspecies 5 Mi^ ni ai Mi^ ni ci i

^ ^ nM



i

^ ^ nM

M i ^ni expðdGi Þ  1 þ n M ^ ^nM

ð36Þ

2.3. Variance of relative exospheric mass density error

^fÞ þr

T^ ðfÞ

!# Pnspecies ^ ^ M a M c n n i i i i i i  ^ a þ ^c  i  i ^ ^ ^ ^ nM nM Pnspecies M i^ ni expðdGi Þ dT 0 ðdGðLÞÞ þ i ^ ^ nM Pnspecies ^ ni expðdGi Þ  i ^ n

^

^ fÞ 1  TT^exexTT^ðfÞ ð1 þ r 0 þ ^ T ð fÞ Pnspecies Pnspecies  

M i ^ni ai M i ^ni ci i i ^ ^  aþc  dT ex ^ ^ ^nM ^nM "   ^ T^ ðfÞ g0 M 1 þ ^a þ ^c þ ln  0 k T^ ex T^ 0 T^ 0 T^ ðfÞ     T^ ex  T^ ðfÞ T^ ðfÞ T^ ðfÞ g ^ fÞ   ð1 þ r þ 0 0 ln T^ 0 T^ 0 T^ ex  T^ 0 k T^ 0 T^ ex  nspecies 2  um M ^ni ^  Pni species i  M M i ^ni i 3 0

Note that the product of the relative errors in Eq. (36) is by far smaller than the individual contributions of the number density and the mean mass errors (typically
j106 j vs.
j103 j). It should however not be neglected for the computation of the relative density error, as it is important in cases where n and M are of similar magnitude but opposite sign.

^fÞ 1  TT^ex TT^ðfÞ ð1 þ r ex 0 þ T^ ðfÞ

!

7

Pnspecies

ð35Þ

Eqs. (28) and (35) can now be inserted into Eq. (6) to compute the overall relative density error, which results in:

We devote this section to the approximation of the variance of Eq. (36). In the following, we assume that dGi follows a zero-mean Gaussian distribution. VarðdT ex Þ and VarðdT 0 Þ are linearly related to VarðdGi Þ, as shown in Appendix C, and are therefore also zero-mean Gaussian. Since variances are always positive, we can safely drop the product of the relative errors in the number density and mean mass. Masking the pre-factors with scalars a and b, we can write Eq. (36) as: Pnspecies M i ^ni expðdGi Þ q ffi adT ex þ bdT 0 þ i 1 ^ ^nM The relative density error variance is therefore given by

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

8

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

X

!

nspecies

  Var q ¼ a2 VarðdT ex Þþb2 VarðdT 0 Þþ

Var

Mi^ ni expðdGi Þ

i

ð^nM^ Þ ! nspecies X ð37Þ ni expðdGi Þ þ 2abCovðdT ex ;dT 0 Þþ2 ^naM^ Cov dT ex ; M i^ X

!

nspecies

þ 2 ^nbM^ Cov dT 0 ;

2

i

M i^ ni expðdGi Þ

i

VarðdT ex Þ and VarðdT 0 Þ are given by Eqs. (C.4) and (C.5). The remaining components that are to be computed, are therefore:    

CovðdT ex ; dT 0 Þ Pnspecies  Var M i^ ni expðdGi Þ i   Pn ni expðdGi Þ Cov dT ex ; i species M i ^   Pnspecies Cov dT 0 ; i M i^ ni expðdGi Þ

Therefore Eq. (42) becomes    1 Cov expðdGi Þ; exp dGj ffi e2ðVarðdGi ÞþVarðdGj ÞÞþCovðdGi ;dGj Þ ð43Þ Eq. (40) can now be evaluated using Eqs. (41) and (43), which yields: X

!

nspecies

Var

M i^ ni expðdGi Þ ¼

i

ð38Þ

ð39Þ

þ2

2.3.2. Variance of the error in the particle density at the baseline altitude The variance of the sum over all considered species is defined as: ! n nspecies species X X Var M i^ M 2i ^ ni expðdGi Þ ¼ n2i VarðexpðdGi ÞÞ i

XX

nspecies nspecies

þ2

i

   M i ni M j nj Cov expðdGi Þ; exp dGj

ð40Þ

j–i

Since we assume dGi to follow a zero-mean normal distribution, expðdGi Þ is lognormal with VarðexpðdGi ÞÞ ¼ eVarðdGi Þ ðVarðdGi Þ  1Þ

ð41Þ

Using the definition of covariance, we obtain:       Cov expðdGi Þ; exp dGj ¼ E expðdGi Þexp dGj    E½expðdGi ÞE exp dGj    ð42Þ ¼ E exp dGi þ dGj    E½expðdGi ÞE exp dGj

X

nspecies

M 2i ^ n2i eVarðdGi Þ ðVarðdGi Þ  1Þ

i

n X n h 1 i X 1 1 M i ni M j nj e2ðVarðdGi ÞþVarðdGj ÞÞþCovðdGi ;dGj Þ  e2VarðdGi Þ e2VarðdGj Þ i

j–i

ð44Þ

2.3.3. Covariance of exospheric temperature error and species number density errors at the baseline altitude   Pn Next, we compute Cov dT ex ; i species M i ^ni expðdGi Þ . From the zero-mean assumption of dT ex it follows that the covariance is equal to the expectation of the product: X

nspecies

i

1

1

As dT ex and dT 0 are assumed to follow zero-mean Gaussian distributions, the covariance is equal to the expected value of the product of the absolute errors in the baseline and exospheric temperatures. Using Eqs. (C.2) and (C.3), we obtain:  CovðdT ex ; dT 0 Þ ffi E½dT ex dT 0  ¼ T 00 T ex0 E dGT 0 dGT ex

   dGi þ dGj
N ldGi þ ldGj ;r2dGi þ r2dGj þ 2Cov dGi ;dGj Þ      ffi N 0:0; VarðdGi Þ þ Var dGj þ 2Cov dGi ;dGj

 e2VarðdGi Þ e2VarðdGj Þ

2.3.1. Covariance of exospheric temperature error and baseline temperature error The definition of covariance implies that CovðdT ex ; dT 0 Þ ¼ E½dT ex dT 0   E½dT ex E½dT 0 

Each expectation may be computed from the first moment of the lognormal distribution. In the first term, we have the addition of two correlated zero-mean Gaussians. In this case it holds that:

Cov dT ex ;

! M i^ ni expðdGi Þ ¼ E

"n species X

i

# M i^ ni dT ex expðdGi Þ

i

X

nspecies

¼

ð45Þ

M i^ ni T ex0 E½dGT expðdGi Þ

i

In Eq. (45) we need to formulate the expected value of the product of a normal and lognormal distribution. This type of distribution is sometimes denoted ‘‘Normal-Lognormal Mixture” (NLNM); however, it has only received little treatment in the literature. One of the very few papers on this type of distribution is (Yang, 2008). In his work Yang states the first two moments of the NLNM created from a standard normal distribution with unit variance and a generic normal distribution without giving a derivation. Both are connected via an arbitrary correlation coefficient q. We seek a description of the first moment of the NLNM which is created from two correlated normal distributions with non-unit variance. The computation of the required first moment is given in Appendix D and extends the work of Yang for the case of non-standard-normal Gaussians. Using the final result, Eq. (D.9), we obtain:

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx 1r 2 dGi

E½dGT expðdGi Þ ¼ rdGT rdGi qe2

Introducing the definition of the correlation coefficient, it follows that: E½dGT expðdGi Þ ¼ rdGT rdGi CovðdGT ; dGi Þ 1  pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffipffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi e2VarðdGi Þ VarðdGT Þ VarðdGi Þ ¼ CovðdGT ; dGi Þe2VarðdGi Þ 1

ð46Þ

Introducing Eq. (46) into Eq. (45) allows computation of the covariance as the sum of the individual covariances of the species: ! nspecies X Cov dT ex ; M i^ ni expðdGi Þ i

X

9

case of NRLMSISE-00, also in the baseline temperature. The vector of the environmental parameters, L, contains latitude, longitude, local solar time, solar flux and geomagnetic activity parameters. An error in each of these parameters will cause a corresponding error and uncertainty/variance2 in Gð LÞ. In this section we consider uncertain geomagnetic activity information, i.e. uncertain magnetic amplitudes. In case of NRLMSISE-00, the uncertainty in the daily Ap index can be propagated to VarðdGi Þ without intermediate conversion steps. For the DTM models we first need to propagate the uncertainty in the 3 h-delayed ap data to the 3 h-delayed kp index, as this is the geomagnetic activity information required by the DTM models. For a treatment of uncertain solar flux data the interested reader is referred to (Schiemenz et al., 2019a).

nspecies

¼

M i^ ni T ex0 CovðdGT ; dGi Þe2VarðdGi Þ 1

ð47Þ

i

Eq. (47) is also valid for the last term to be computed, i.e.   Pn Cov dT 0 ; i species M i ^ ni expðdGi Þ , as dT ex and dT 0 only differ by the scalar coefficients. Therefore it holds that: ! nspecies X M i^ ni expðdGi Þ Cov dT 0 ; i

X

nspecies

¼

M i^ ni T 00 CovðdG0 ; dGi Þe2VarðdGi Þ 1

ð48Þ

i

3.1. NRLMSISE-00: error/variance propagation from Ap to G(L) In NRLMSISE-00, the relationship between the change in Gi ð LÞ and a change in the daily magnetic amplitude, is given by: dGi ¼ ðc1 þ c2 p40 þ c3 p20 þ cc Þ dAð dApÞ where

h i ÞÞ1 Að ApÞ ¼ ð Ap  4Þ þ ðc5  1Þ  Ap  4 þ expðc4 ðcAp4 4

This completes the equations to evaluate the variance of the relative density error. Inserting Eqs. (39), (44) and (48) into Eq. (37) leads to:

¼ Ap þ ðc5  1ÞAp  4c5  c5c1 þ c5c1 expðc4 Ap þ 4c4 Þ 4 4 ¼ c5 Ap þ expðc4 ApÞ

  Var q ¼ a2 VarðdT ex Þ þ b2 VarðdT 0 Þ Pnspecies 2 2 VarðdG Þ i ðVar ðdGi Þ1Þ Mi ^ ni e i þ 2 ð^nM^ Þ nspecies 1  Pnspecies X 1Var dG 1Var dG Var ðdGi ÞþVar ðdGj ÞÞþCovðdGi ;dGj Þ M i ni M j nj e2ð e2 ð i Þ e2 ð j Þ

c5  1 expð4c4 Þ 4c5  c5c1 4 c4 |fflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflffl} :¼ d

¼ c5 Ap þ d expðc4 ApÞ  4c5  c5c1 4 ð51Þ

i

j–i

þ 2

2 ð^nM^ Þ

  þ 2abT 00 T ex0 Cov dGT 0 ; dGT ex nspecies X 1 M i ^ni T ex0 CovðdGT ex ;dGi Þe2VarðdGi Þ þ 2 ^naM^

ð49Þ

i nspecies

þ 2 ^nbM^

X

ð50Þ

M i ^ni T 00 CovðdG0 ; dGi Þe2VarðdGi Þ 1

i

Eq. (49) yields an approximate description of the relative density error variance. It only depends on outputs of a single atmospheric model call and model-dependent coefficients. The variances and covariances of the errors in the function of the environmental parameters are modeldependent and subject of the next section.

and c1 ; c2 ; c3 ; c4 and c5 are model-dependent coefficients. p20 and p40 represent Legendre-polynomials of degree two and four and cc additional cross-correlation terms in latitude and longitude. To complete Eq. (50), we need to derive an expression ^ þ dAp, we obtain: for dAð dApÞ. Defining Ap ¼ Ap   ^ þ d expðc4 ApÞ  d exp c4 Ap ^ dAð dApÞ ¼ c5 Ap  c5 Ap     ^ expðc4 dApÞ  d exp c4 Ap ^ ¼ c5 dAp þ d exp c4 Ap |fflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflffl} :¼f

¼ c5 dAp þ f expðc4 dApÞ  f ð52Þ

3. Absolute error, variance and covariance in G i ðLÞ due to uncertain magnetic index information

Inserting Eq. (52) into Eq. (50), results in

As explained in Section 1, the model-specific function Gð LÞ is used to model periodic and non-periodic variations in the number density, exospheric temperature and, in the

Since we assume the errors in GðLÞ to follow zero-mean Gaussians, the variance is the only moment describing the distribution and hence can be treated as a synonym for the overall uncertainty. 2

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

10

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

dGi ¼ ðc1 þ c2 p40 þ c3 p20 þ cc Þð c5 dAp þ f expðc4 dApÞ  f Þ |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} c

   Cov dGi ;dGj ffi E dGi dGj

 ¼ ci c5i cj c5j E dAp2 |fflfflfflffl{zfflfflfflffl}

ð53Þ

¼ cc5 dAp þ cf expðc4 dApÞ  cf

Var:ofGuassian

   þ ci c5i cj f j E dApexp c4j dAp |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

Eq. (53) relates an absolute error in Ap to an absolute error in the spherical harmonics function. The index i also applies to the set of coefficients and denotes the species under consideration, which can be either a chemical species (i.e. the gas constituents) or a ‘‘temperature species” (i.e. the exospheric temperature and the baseline temperature). The model-coefficients for exospheric temperature in NRLMSISE-00 are: c1 ¼ 4:93933  103 ; c2 ¼ 2:50802  103 ; c3 ¼ 5:72562  103 ;c4 ¼ 8:47001  102 ;c5 ¼ 1:70147  101 . Eq. (53) can now be used to derive the mean and variance of the resulting distribution. Since the relationship between dAp and dGi is not linear, the third and higher moments are non-zero, however we will continue with the ^ is expected to be very assumption of Gaussianity, as Ap close to Ap, which means that we are considering dAp to be confined to a small linear neighborhood. This implies ^ From that rdAp is expected to be small in comparison to Ap. Eq. (53) it follows that

Inserting the corresponding moments of the lognormal distribution and the NLNM (Eq. (D.9)) we obtain:

VarðdGi Þ ¼ c2 c25 Varð dApÞ þ c2 f 2 Varðexpðc4 dApÞÞ

3.2. DTM-2012/2013: error/variance propagation from ap to kp

þ 2c2 c5 fCovð dAp; expðc4 dApÞÞ

ð54Þ

Since dAp is assumed to be zero-mean Gaussian, expðc4 dApÞ follows a lognormal distribution. The covariance therefore matches the first moment of the resulting NLNM with a correlation coefficient of 1.0. Using Eq. (D.9) results in  1 2 E dApec4 dAp ¼ c4 Varð dApÞe2c4 Varð dApÞ ð55Þ After application of the second moment of the lognormal distribution, we obtain the following expression for VarðdGi Þ:   2 2 VarðdGi Þ ¼ c2 c25 Var dAp þ c2 f 2 eðc4 Varð dApÞÞ eððc4 Varð dApÞÞ1Þ 1 2 2c2 c c fVarð dApÞeð2c4 Varð dApÞÞ ð56Þ 5 4

Similarly, we can evaluate the first moment to  12  E½dGi  ¼ cf e2c4 Varð dApÞ  1

ð57Þ

which is non-zero. Since, however, c4 is small (8:47001  102 for the exospheric temperature) and taken to the square, the resulting mean is generally very close to zero, such that the assumption of dGi being zero-mean Gaussian can be justified in the close neighborhood around ^ Ap. Using  Eqs. (53),  we can now also compute an expression for Cov dGi ; dGj , as required by Eq. (49). Continuing with the assumption of dGi being zero-mean, we can express the covariance as the expected value of the product of the random variables:

NormalLognormalMixture

þ ci f i cj c5j E½ dApexpðc4i dApÞ

ð58Þ

    ci f i cj f j E½expðc4i dApÞ ci f i cj f j E exp c4j dAp |fflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflffl} Lognormal

   þ ci f i cj f j E expðc4i dApÞexp c4j dAp þci f i cj f j |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl} Lognormal

  Cov dGi ; dGj ffi

h 1c2 Var ð dAp Þ Varð dApÞ ci cj c5i c5j  ci cj c5i f j c4j e2 4j i 1c2 Varð dAp Þ ci cj c5j f i c4i e2 4i 2 2 ð59Þ 1 1 þci cj f i f j e2ðc4i c4j Þ Varð dApÞ  e2ðc4i Þ Varð dApÞ

2 1 e2ðc4j Þ Varð dApÞ þ 1

In case of the DTM-models we treat the error in the 3 hdelayed magnetic amplitude as a zero-mean Gaussian random variable. As the relationship between ap and kp is not linear, the resulting error process in kp cannot be Gaussian, meaning that also higher moments have non-zero contribution. If the following conditions however are met, it is possible to estimate a near-Gaussian error process by only considering the resulting first and second moments, as the propagation can be considered linear in the neighborhood of ^ap:  rdap must be small in comparison to ldap  If a logarithmic mapping is considered, it must be ensured that ^ap > 0, as the logarithm is only defined for positive arguments First of all, we require a continuous function to map ap to kp .3 The exact mapping of the linear magnetic amplitude to the quasi-logarithmic kp=Kp index is only available for discrete values of ap=Ap via conversion tables introduced by Bartels ((Bartels, 1949; Bartels and Veldkamp, 1953)). Wertz defines the following mapping, which is quoted to be accurate to within
10 % (Wertz, 2012): kp ¼ 1:75 ln ap  1:6 This mapping is also not suitable, as its conversion results are too inaccurate for error propagation purposes in case of small values of ap. Also it should be noted that 3

The same function is also valid for the daily magnetic indices.

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

errors in kp are more pronounced for small values of ap, as the knot size of the discrete ap to kp conversion table increases with increasing magnetic activity. Hence, we seek for a very accurate conversion from ap to kp, especially for small values of ap. Vallado postulates a cubic splining approach for its accuracy and consistency (Vallado, 2013). His algorithm however is also not suitable, as it results in undesired bumps at ap ¼ 2 and ap ¼ 300, which are the consequence of the undefined additional two degrees of freedom in the outer knots when performing the cubic splining (c.f. Fig. 1). We therefore postulate a clamped cubic splining approach which makes sure that the transitions at ap ¼ 2 and ap ¼ 300 remain smooth. The respective coefficients of the third-order polynomials are given in Appendix E. Independently of the boundary conditions, a cubic splining approach generally results in a table of third-order coefficients that allow to obtain the matching cubic polynomial for the current value of ap: kp ¼ a3 ap3 þ a2 ap2 þ a1 ap þ a0

ð61Þ

An absolute error in ap is defined as dap ¼ ap  ^ ap. Using Eq. (61), we therefore find the corresponding error in kp via: dkp ¼ a3 dap3 þ ð3a3 ^ ap þ a2 Þdap2   þ 3a3 ^ ap2 þ 2a2 ^ ap þ a1 dap

ð62Þ

Eq. (62) relates an error in ap to an error in kp. The associated mean and variance in dkp then allow to formulate a Gaussian approximation of the resulting distribution. Denoting c3 ¼ a3 ; c2 ¼ 3a3 ^ and ap þ a2 c1 ¼ 3a3 ^ ap2 þ 2a2 ^ ap þ a1 , the variance in dkp is defined as:     VarðdkpÞ ¼ c23 Var dap3 þ c22 Var dap2   þ c12 VarðdapÞ þ 2c3 c2 Cov dap2 ; dap3   þ 2c3 c1 Cov dap3 ; dap   þ 2c2 c2 Cov dap2 ; dap ð63Þ Recall that we assume dap to follow a zero-mean Gaussian distribution. The variances and covariances thus equate to:     2 3 Var dap3 ¼ E dap6  E dap3 ¼ 15VarðdapÞ ð64Þ     2 Var dap2 ¼ 2VarðdapÞ Cov dap3 ; dap     ¼ E dap4  E dap3 E½dap ¼ 3Var dap2 ð65Þ      3 2 5 3 2 Cov dap ; dap ¼ E dap  E dap E dap  ¼ E dap5 ¼ 0:0 ð66Þ     2 3 2 Cov dap ; dap ¼ E dap  E dap E½dap ¼ 0:0 ð67Þ Hence, the variance in dkp is given by: VarðdkpÞ ¼ 15a23 VarðdapÞ3    ap2 þ 2a2 ^ap þ a1 þ 2½3a3 ^ap þ a2 2 VarðdapÞ2 þ 6a3 3a3 ^  2 þ 3a3 ^ap2 þ 2a2 ^ap þ a1 VarðdapÞ ð68Þ

11

Fig. 1. Differences of the ap to kp cubic splining algorithm given in Vallado (2013) and a clamped cubic spline. The Vallado-algorithm shows undesired strong changes in the slope at the outer two ap knots.

Due to the linearity of expectation, the mean is computed from Eq. (62) as  ð69Þ E½dkp ¼ ð3a3 ^ap þ a2 ÞE dap2 ¼ ð3a3 ^ap þ a2 Þr2dap which is non-zero. A little consideration with the coefficients of Table E.2 however shows that for an assumed standard deviation of rdap ¼ 0:5 the mean in kp remains below j0:05j and can therefore be neglected in allmost all cases. Summarizing:   dap : N 0:0;r2dap !  3 dkp : N ð3a3 ^ap þ a2 Þr2dap ;15a23 VarðdapÞ   2 2 þ 6a3 ½3a3 ^ap2 þ 2a2 ^ap þ a1  þ 2½3a3 ^ap þ a2  VarðdapÞ 2

þð3a3 ^ap2 þ 2a2 ^ap þ a1 Þ VarðdapÞ

 ð70Þ

where the coefficients a3 ; a2 and a1 are to be taken from the row of Table E.2 that belongs to ^ap. 3.3. DTM-2012/2013: error/variance propagation from kp to G(L) Now that the approximate Gaussian in the kp error has been determined, the next step is to relate the error and its uncertainty in kp to an error and its variance in the spherical harmonics of the environmental parameters. In the following we analyze dGðdkpÞ for DTM-2012 and DTM-2013. Both models contain the same description of Gð L), however use difference parametrizations. 15 terms model the effects of the magnetic index inputs to the models. Only six of the terms however change if the 3 h-delayed kp index input changes. The dependency reads:

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

12

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

dGðdkpÞ ¼ ðc1 þ c2 p20 þ c5 p40 Þdkp þ ðc3 þ c4 p20 Þ      kp2  ^kp2 þ c6 kp4  ^kp4

4

ð71Þ

The cx -coefficients are model- and species-dependent, meaning that the appropriate values that match the species under study need to be loaded from the matching model coefficients file. p20 and p40 are Legendre polynomials and represent kp-/latitude-/longitude cross-couplings. Defining dkp ¼ kp  ^kp and resolving the differences in the powers results in:   dGi ðkpÞ ¼ c1 þ c2 p20 þ c5 p40 þ 2^kpðc3 þ c4 p20 Þ þ 4^kp3 c6 dkp   þ c3 þ c4 p20 þ 6c6 ^kp2 dkp2 þ 4c6 ^kpdkp3 þ c6 dkp4 ð72Þ Based on  Eq. (72),  we can now compute E½dGi ; VarðdGi Þ and Cov dGi ; dGj , which are needed by Eq. (49). To simplify notation, define b1 ¼ c1 þ c2 p20 þ c5 p40 þ 2^kpðc3 þ c4 p20 Þ þ 4^kp3 c6 ; b2 ¼ c3 þ c4 p20 þ 6c6 ^kp2 ; b3 ¼ 4c6 ^kp and b4 ¼ c6 . Then:       VarðdGi Þ ¼ b24 Var dkp4 þ b23 Var dkp3 þ b22 Var dkp2 þ b21 VarðdkpÞ     þ 2b4 b3 Cov dkp4 ;dkp3 þ 2b4 b2 Cov dkp4 ; dkp2     þ 2b4 b1 Cov dkp4 ;dkp þ 2b3 b2 Cov dkp3 ; dkp2     þ 2b3 b1 Cov dkp3 ;dkp þ 2b2 b1 Cov dkp3 ; dkp ð73Þ For the variance terms it holds that     2 Var dkp4 ¼ E dkp8  E dkp4 4

2

¼ 105VarðdkpÞ  9VarðdkpÞ 

 3

¼ 96VarðdkpÞ

4

ð74Þ

3

2

Var dkp ¼ 15VarðdkpÞ  ð0:0Þ ¼ 15VarðdkpÞ   2 Var dkp2 ¼ 2VarðdkpÞ

3

And for the covariances we obtain:   Cov dkp4 ; dkp3 ¼ 0:0   3 3 Cov dkp4 ; dkp2 ¼ 15VarðdkpÞ  3VarðdkpÞ 

 1

¼ 12VarðdkpÞ

3

Cov dkp4 ; dkp ¼ 0:0   Cov dkp3 ; dkp2 ¼ 0:0   2 2 Cov dkp3 ; dkp1 ¼ 3VarðdkpÞ  0:0 ¼ 3VarðdkpÞ   Cov dkp2 ; dkp ¼ 0:0

ð75Þ ð76Þ

ð77Þ

ð78Þ ð79Þ ð80Þ ð81Þ ð82Þ

Combining the variances and covariances with Eq. (73) results in:

VarðdGi Þ ¼ 96b24 VarðdkpÞ   3 þ 15b23 þ 24b4 b2 VarðdkpÞ   þ 2b22 þ 6b3 b1 VarðdkpÞ2 þ b21 VarðdkpÞ

ð83Þ

Similarly we obtain the mean via:   E½dT ex  ¼ b4 E dkp4 b2 E dkp2 2

¼ 3b4 VarðdkpÞ þ b2 VarðdkpÞ

ð84Þ

Combining Eqs. (84) and (83) now allows to formulate the approximate Gaussian in the spherical harmonics function error:    2 dkp : N 0:0; r2dkp ! dGi :
N 3b4 VarðdkpÞ 4

þ b2 VarðdkpÞ; 96b24 VarðdkpÞ   3 þ 15b23 þ 24b4 b2 VarðdkpÞ     þ 2b22 þ 6b3 b1 VarðdkpÞ2 þ b21 Var dkp4 Þ

ð85Þ

The model-coefficients cx for DTM-2012 and DTM-2013 in the case of exospheric temperature are: c1 ¼ 0:128036  101 ðDTM  2012Þ; 0:124511  101 ðDTM  2013Þ c2 ¼ 0:154477  101 ðDTM  2012Þ; 0:163993  101 ðDTM  2013Þ c3 ¼ 0:792957  103 ðDTM  2012Þ; 0:184118  102 ðDTM  2013Þ c4 ¼ 0:0 ðDTM  2012Þ; 0:0 ðDTM  2013Þ c5 ¼ 0:0 ðDTM  2012Þ; 0:0 ðDTM  2013Þ c6 ¼ 0:457116  104 ðDTM  2012Þ; 0:234802  104 ðDTM  2013Þ

The model coefficients of the other species can be inferred form the respective model-coefficients file. The mean in dGi is non-zero. To test its magnitude/relevance, we perform a study of the model-coefficients for a representative species. For the case of exospheric temperature, the magnitude of the resulting E½dT ex  is given in Table 1 for a couple of representative scenarios: Table 1 indicates that the mean decreases with increasing ^ap and decreasing rdap . Only for large rdap and small ^ap it exceeds 0:1K. Hence, the Gaussian approximation of dT ex is also approximately zeromean.   We need yet to determine Cov dGi ; dGj . As reasoned above, we assume that the mean in dGi can be neglected. Therefore the covariance is equal to the expectation of the product, which results in:     Cov dGi ;dGj ffi E b4i dkp4 þ b3i dkp3 þ b2i dkp2 þ b1i dkp    b4j dkp4 þ b3j dkp3 þ b2j dkp2 þ b1j dkp  which is a polynomial in dkp of order eight. Since we assume dkp to be zero-mean, the resulting expectations are given by the central moments of the normal distribution, which are non-zero in case of even powers of dkp and zero otherwise. Hence:

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

13

Table 1 E½dT ex  for different ^ ap and rdap combinations. ^ ap ¼ 2:0 ^ ap ¼ 15:0 ^ ap ¼ 300:0

rdap ¼ 0:5

rdap ¼ 1:0

0.01762 (DTM-2012)/ 0.03371 (DTM-2013) 0.01073 (DTM-2012)/ 0.01023 (DTM-2013) 0.00006 (DTM-2012)/ 0.00004 (DTM-2013)

0.14341 (DTM-2012)/ 0.31673 (DTM-2013) 0.08867 (DTM-2012)/ 0.08449 (DTM-2013) 0.00057 (DTM-2012)/ 0.00033 (DTM-2013)

   Cov dGi ; dGj ffi b4i b4j E dkp8    þ b4i b2j þ b3i b3j þ b2i b4j E dkp6    þ b3i b1j þ b2i b2j þ b1i b3j E dkp4  þ b1i b1j E dkp2 which equates to   Cov dGi ; dGj ffi 105b4i b4j VarðdkpÞ4   3 ð86Þ þ15 b4i b2j þ b3i b3j þ b2i b4j VarðdkpÞ   2 þ3 b3i b1j þ b2i b2j þ b1i b3j VarðdkpÞ þ b1i b1j VarðdkpÞ 4. Validation The validation of the analytic error and uncertainty propagation equations is the subject of this section. For the error propagation (Eq. (36)), we establish a geomagnetic activity/altitude grid, which covers the full ap scale and the altitude range 300 6 z 6 1200 km, which corresponds to the LEO altitude range. The variance  estimation  (Eq. (49)) is evaluated using VarðdGi Þ and Cov dGi ; dGj of the respective atmospheric model and is validated with the help of a Monte-Carlo (MC) analysis. 4.1. Validation of relative density error 

Fig. 2 depicts the ratios ^qq over the magnetic amplituide/ altitude-grid, where ^q is computed using Eq. (36) and q from a model-based propagation. A value of 1.0 indicates a perfect match of estimation and truth. The grid is comparable to those presented in Schiemenz et al. (2019a), with the exception that geomagnetic amplitudes and not solar flux information is analyzed. As a value of 1.0 is obtained basically over the entire grid, no additional tuning is needed. For NRLMSISE-00, each grid-point resides at a value of 1.0, which indicates a perfect estimation of the relative density error via Eq. (36). The same also holds true for the DTM-models, except for ^ ap ¼ 2. For this very low magnetic amplitude, the components of ^q , i.e. ^n and ^M , are of equal magnitude but opposite signs. In the computation of Fig. 2 we explicitly took the product of the errors in Eq. (36) into account. Nevertheless, small errors in the estimation of ^n and ^M can cause errors in ^q in the order of multiple 10%, depite the absolute error in the density being small. For practical applications the outliers do not lead to

any issues, as the orbital motion depends on the double integral of the relative density error (Emmert et al., 2017) and the altitude/magnetic index combination is commonly not at one of the ill-conditioned grid-points for longer periods of time. 4.2. Validation of variance of relative density error Each point of the grid in Fig. 2 corresponds to a unique altitude/magnetic amplitude-combination for which the corresponding variance is estimated via Eq. (49). In the following we present for each model the variance validation based on the representative scenarios:  Altitude: 300 km (low LEO), magnetic amplitude: 9 (low)  Altitude: 500 km (medium LEO), magnetic amplitude: 300 (very high)  Altitude: 700 km (high LEO), magnetic amplitude: 2 (very low)  Altitude: 1100 km (very high LEO), magnetic amplitude: 50 (medium) For each model and case, 200 samples are computed via time-consuming individual numerical orbit propagations, for which we evaluate the true density using Ap (NRLMSISE-00)/ap (DTM-2012 and DTM-2013) and ^ ap. dAp and dap are treathe variance estimation using Ap/^ ted as zero-mean Gaussian random variables. The resulting variance in the relative density error is then compared to the actual spread of the Monte-Carlo iterations. In the case of NRLMSISE-00, the subplots depict the MC-iterations of the baseline temperature (first column), the exospheric temperature (second column) and finally the relative density error (third column). In case of DTM-2012 and DTM-2013 the baseline temperature is modeled constant. Hence, the first column depicts the variance propagation from the 3 h-delayed magnetic amplitude to the quasi-logarithmic kp index. Judging from the results depicted in Figs. 3–5, the estimated variance of the relative density error and all intermediate stages well matches the corresponding variances of the Monte-Carlo iterations, which validates Eq. (49) and the underlying assumptions. Nevertheless, it has to be noted that for very low magnetic amplitudes a bias may be introduced in the ^ap to ^kp mapping, which propagates to the relative density errors. While the zero-mean assumption in this case is still valid for the variance computations

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

14

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx



Fig. 2. Plots of the ratio ^qq for all considered models with respect to altitude and magnetic amplitude. The estimates match the true values. For the DTM models outliers are possible in case of very low magnetic amplitudes.

via Eq. (49), the resulting mean in the relative density error should be considered when propagating the relative density uncertainty to the orbital uncertainties. We want to emphasize that the Monte-Carlo derived uncertainties, albeit free of any of the assumptions made in Sections 2 and 3, are suitable for validation purposes, however not for practical operations, as they require the time-consuming evaluation of hundreds of individual numerical orbit propagations in order to build the resulting probability distribution. The   presented analytic equations for q (Eq. (36)) and Var q (Eq. (49)), however, can be evaluated without any numerical error propagations or

additional model calls, as the required information can be extracted from the same model-call, which is made to obtain the atmospheric density estimate. They are therefore compatible with operational scenarios, such as density uncertainty considering orbit estimation. 5. Summary and conclusion Equations have been derived to compute accurate estimates of the relative density error resulting from an input-uncertainty in the set of environmental parameters, which needs to be supplied to popular semi-empirical

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

15

Fig. 3. Comparison of NRLMSISE-00 derived absolute baseline temperature errors (first column), absolute exospheric temperature errors (second column) and relative exospheric density errors (third column). In each case the blue line depicts the estimated 3r deviations, whereas the red line gives the Monte-Carlo (MC) 3r deviations. The 200 MC-iterations are plotted in green color. The thin black line depicts the mean of the MC iterations. In all figures the estimated uncertainty agrees with the MC-derived uncertainty, which validates Eq. (49). Also the zero-mean assumptions hold in any case. The resulting relative density error is largest for small magnetic amplitudes (5% in case of (c)).

atmospheric models, such as NRLMSISE-00, DTM-2012 and DTM-2013. Subsequently, we showed how to derive the corresponding variance in the relative density error, which is the key parameter needed to perform orbit estimation including density uncertainty consideration, as has been worked out in Schiemenz et al. (2019b). For the error and variance estimations we assumed a Bates temperature profile of the exosphere. To validate the derived equations, we presented a mapping from an error/uncertainty in the geomagnetic activity information to the spherical harmonics function used inside the respective semi-empirical model. In case of the DTM-models this mapping makes use of a clamped cubic

splining-approach in order to perform the transition from geomagnetic amplitude to the kp index, since it turned out that presently published algorithms are unsuitable for this purpose. It has been shown, that both, the equations describing the relative density error and its variance, are able to accurately predict the actual error/uncertainty in the relative density as obtained from a model-based Monte-Carlo computation. We also found, that geomagnetic activity uncertainties are especially important to consider in case of low global geomagnetic activity. With the help of the uncertainty estimation equations presented in this work, it is now possible to analytically account for the input-uncertainty of the environmental

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

16

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

Fig. 4. Comparison of DTM-2012 derived absolute exospheric temperature errors (second column) and relative exospheric density errors (third column). The DTM-models consider T 0 to be constant. The first column therefore shows the results of the ^ap to ^kp propagation. In each case the blue line depicts the

3r deviations, whereas the red line gives the Monte-Carlo 3r deviations. While the variance is generally correctly estimated, (c) (very low ^ ap) depicts a bias which is introduced in the propagation from the ^ap error to the ^kp error. This bias does not degrade the subsequent variance computations (which assume zero bias), but propagates to the relative density error and therefore also eventually to the orbital state-vector, which should be kept in mind when performing orbit estimation simulations. In all cases the standard deviation in the 3 h delayed magnetic amplitude was set to rdap ¼ 0:5.

parameters that need to be supplied to atmospheric models. The need of a tuning-function, as it was necessary in Schiemenz et al. (2019a), has been removed. Nevertheless, it has to be noted that the presented propagation from geomagnetic activity to relative density errors is, most likely, only directly applicable to the measurement uncertainty of the geomagnetic activity information. The more interesting case, i.e. geomagnetic forecast uncertainty, is expected to be too large to justify the assumption of a near-linear uncertainty propagation as applied in this work. For this case, non-linear error/uncertainty

propagations need to be considered to capture all moments of the resulting distribution. The validation of such equations could then also be performed using Monte-Carlo analysis. In a future work we plan to combine our previous efforts on solar flux input-uncertainty and density uncertainty considering orbit estimation with the work at hand, in order to perform orbit determination that is accounting for the density uncertainty due to geomagnetic and solar flux input-uncertainty. It is hoped that the incorporation of this physics-based GCRF density-uncertainty covari-

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

17

Fig. 5. Comparison of DTM-2013 derived absolute exospheric temperature errors (second column) and relative exospheric density errors (third column). The left plot shows the ^ ap to ^kp propagation. The results are similar to those obtained with DTM-2012 as atmospheric model in the sense that the variance is correctly estimated in all cases, however for very low magnetic amplitudes a bias is introduced in the ^ap to ^kp propagation, which subsequently propagates to the relative density error.

ance matrix will further improve covariance realism and avoid filter smugness.

Appendix A. Inverse Bates-profile integral To obtain a description of the number density vs. altitude, it is required to solve the integral of the inverse temperature profile with respect to geopotential height. For a Bates temperature profile, this integral can be solved analytically, as shown in the following.

Inserting the Bates temperature profile (Eq. (8)) results in: Z f Z f 1 1 0 0 ðA:1Þ 0 df ¼ 0 df T f T ð Þ  ð T  T Þ exp ð rf Þ ex ex 0 0 0 Substitute u ¼ rf ) df0 ¼  r1 du0 . Then the integral becomes Z Z f 1 1 rf 1 0 du0 ¼  df 0 r 0 T ex  ðT ex  T 0 Þ expðu0 Þ 0 T ðf Þ ðA:2Þ

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

18

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

dv Next, substitute expðuÞ ¼ v ) du ¼ v, which yields: Z expðrfÞ Z f 1 1 1 0 dv0 0 df ¼  0 Þv0 r T ð  ð T  T Þv T f ð Þ ex ex 0 0 1

Introducing Eq. (B.1) results in: d af ð xÞgðxÞ dx

ðA:3Þ Next, replace the integrand with its partial fraction expansion: Z Z f 1 1 expðrfÞ 1 T 0  T ex 0 þ dv0 0 df ¼  0 0 r T v T f ð Þ T ex v ðT ex  T 0 Þ  T 2ex ex 0 1

¼

Z

1 T ex r

expðrfÞ

1

1 0 T ex  T 0 dv  v0 rT ex

ðA:4Þ

Z

expðrfÞ 1

1 dv0  T ex  v0 ðT ex  T 0 Þ

ðB:3Þ Hence, the required equation of the derivative is d gð xÞ af ð xÞ expðcð xÞÞ dx

d f ð xÞ d d gð xÞ cð xÞ þ lnð f ð xÞÞ gð xÞ þ gð xÞ dx ¼ af ð xÞ expðcð xÞÞ f ð xÞ dx dx ðB:4Þ

ðA:5Þ

¼ ðT ex  T 0 Þ Finally, substitute w ¼ T ex  vðT ex  T 0 Þ ) dw dv to obtain: Z

expðcð xÞÞ h   d gð xÞ d dxf ð xÞ g ð x Þ ln ð f ð x Þ Þ þ g ð x Þ ¼ a f ð xÞ expðcð xÞÞ dx f ð xÞ i þf ð xÞgðxÞ dxd expðcð xÞÞ

Appendix C. Dependence of exospheric and baseline temperature errors/variances on an absolute error in the spherical harmonics function

f

1 1 0 ½rf  lnð1Þ 0 df ¼  T T f ð Þ ex r 0 Z T ex ðT ex T 0 Þ expðrfÞ T ex  T 0 1 dw0 ðA:6Þ þ w0 rT ex ðT ex  T 0 Þ T 0 Z T ð fÞ f 1 1 f 1 þ dw0 ¼ þ ½lnðT ðfÞÞ  lnðT 0 Þ ðA:7Þ T ex rT ex T 0 w0 T ex rT ex

Therefore the solution to the integral of the inverse Bates temperature profile with respect to geopotential height is:   Z f ln TTð0fÞ 1 f 0  ðA:8Þ 0 df ¼ T ex rT ex 0 T ðf Þ

NRLMSISE-00, DTM-2012 and DTM-2013 express exospheric temperature according to Eq. (B.1) T ex ¼ T ex0 þ T ex0 Gð LÞ

where T ex0 is a model-coefficient parameter and represents the global and temporal average exospheric temperature if all dependencies are neglected. An error in Gð LÞ is therefore linearly related to an error in exospheric temperature: dT ex ¼ T ex0 dGT ex ð LÞ

gðxÞ

expðcðxÞÞ with respect

The derivative in Eq. 2.1.2 belongs to the class expðcð xÞÞ. Hence we need to derive the appropri-

gð xÞ d af ð xÞ dx

ate rule for the differentiation. First, define y ð xÞ ¼ f ð xÞgðxÞ and take the natural logarithm on both sides: lnð y ð xÞÞ ¼ gð xÞ lnð f ð xÞÞ ()

d dx

()

d dxy ð xÞ

()

d f ð xÞgð xÞ dx

lnð y ð xÞÞ ¼ dxd gð xÞ lnð f ð xÞÞ þ gð xÞ

y ð xÞ

gð xÞ d af ð xÞ dx

¼

d gð xÞ lnð f ð xÞÞ dx

¼ f ð xÞgðxÞ



þ gð xÞ

d dxf ð xÞ

f ð xÞ

ðB:1Þ

d dxf ð xÞ

f ð xÞ

d gð xÞ lnð f ð xÞÞ dx

d f ð xÞ

þ gð xÞ dxf ðxÞ



expðcð xÞÞ can now be evaluated by a simple application of the product rule:  d d  gð xÞ gð xÞ af ð xÞ expðcð xÞÞ ¼ a f ð xÞ expðcð xÞÞ dx dx gð xÞ d expðcð xÞÞ ðB:2Þ þaf ð xÞ dx

ðC:2Þ

In NRLMSISE-00, the same dependency is also used to compute the temperature at the baseline altitude: dT 0 ¼ T 00 dG0 ð LÞ

Appendix B. Derivative of af ðxÞ to x

ðC:1Þ

ðC:3Þ

Comparing Eq. (C.3) with Eq. (C.2), we note that they are of identical form. The only difference is the scaling factor (T 00 vs. T ex0 ). The relationships of the variances read: VarðdT ex Þ ¼ T 2ex0 VarðdGi ð LÞÞ

ðC:4Þ

and VarðdT 0 Þ ¼ T 200 VarðdGi ð LÞÞ

ðC:5Þ

where the index i is used to denote the ‘‘temperature species” and determines the relevant model-coefficients to be used. For NRLMSISE-00, T 00 ¼ 364:7105 K and T ex0 ¼ 1027:318 K. The DTM-models assume a T 00 that is fixed at T 00 ¼ 380:0 K, independently of the value of Gð LÞ. Hence the variance of the baseline temperature is only non-zero for NRLMSISE-00. The parametrization of T ex0 for DTM-2012 and DTM-2013 is T ex0 ¼ 1029:77 K. Appendix D. First moment of generic Normal-Lognormal mixture Define the distributions X and Y as:

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

19

Table E.2 Ap to Kp coefficients obtained by performing a clamped cubic splining of the discrete Bartels Ap to Kp grid. 0 6 ap < 2 2 6 ap < 3 3 6 ap < 4 4 6 ap < 5 5 6 ap < 6 6 6 ap < 7 7 6 ap < 9 9 6 ap < 12 12 6 ap < 15 15 6 ap < 18 18 6 ap < 22 22 6 ap < 27 27 6 ap < 32 32 6 ap < 39 39 6 ap < 48 48 6 ap < 56 56 6 ap < 67 67 6 ap < 80 80 6 ap < 94 94 6 ap < 111 111 6 ap < 132 132 6 ap < 154 154 6 ap < 179 179 6 ap < 207 207 6 ap < 236 236 6 ap < 300 300 6 ap < 400+



X Y




N

 0:0 ; 0:0

a3

a2

a1

a0

0.012821285570668300 0.035903617898679900 0.010281138261386200 0.005220935146865060 0.010602602326074200 0.037189474157432600 0.014790252966883600 0.000193356837978427 0.000254028476054232 0.000822757066238519 0.000334053475314261 0.000188285245951570 0.000229028798183291 8:02535376886814  105 5:75256037028785  105 7:96728884985953  105 3:12985956717728  105 5:00757865332135  107 2:55002325897578  106 8:47046305237753  107 2:11242016197263  106 1:05179943442250  106 2:39014661504994  107 7:18791807451520  107 1:55381471501983  106 6:34152268878335  107 1:27594489473349  107

0.007690762191996690 0.300040183008086000 0.115622622432509000 0.070402258466506300 0.166950803627583000 0.693306573075540000 0.398267696535099000 0.006289768196174410 0.009816103109001310 0.038639246294172500 0.023828522949677600 0.014207819811740100 0.019594617763183700 0.010096486480525600 0.007437318204186690 0.012319264672825500 0.006323944667796300 0:0000677253931617976 0.000559549087636273 0.000398424529471943 0.000819794023764658 0.000433236936407817 0.000163119175910688 0.000420759503283973 0.000990529147170736 0.000558551477429165 0.000127020605087351

0.10000000000000000 0.4846988416321790 0.76228957468960600 0.01819005109354530 1.20495536156399000 3.9565888986547500 3.68443098861973000 0.04341380603826460 0.23668426170037300 0.4901459793472340 0.63427386704206800 0.42261839800744200 0.4900474165154990 0.46006791928319800 0.35636035650598000 0.5919556215906070 0.45206410148421600 0.02382220740002320 0.0155236881579348 0.07452583185023730 0.12129784571672900 0.0441022410260379 0.04773660027101190 0.09385421887082980 0.1982825317732950 0.16730049563228200 0.0383711291226729

0.0000000000000000 0.3897992277547860 0.857189188566999 0.1349435095610820 1.842998674556330 8.4800898458811500 9.348956557759290 1.5740949899851000 0.8010131673366630 4.4351643725747000 2.311354705761110 0.759214599507187 7.4547777311992800 2.679785850653490 1.331587534549650 13.841468114995700 5.646900049067620 3.9171689188126900 4.9663928003582400 2.1448411734355200 0.4142766603753310 7.6918804770570600 2.9774866238085100 0.2258020473527070 20.383237841797300 8.375960314108050 12.191202161387400

r2X

qrX rY

qrX rY

r2Y

! ðD:1Þ

and the resulting mixture as:

h rY 1 i 2 2 E½u ¼ E½ XZ  ¼ E XerX qX e2ð1q ÞrY h rY i 1 2 2 ¼ e2ð1q ÞrY E XerX qX For the exponential function it holds that

u ¼ XeY ¼ XZ

ðD:2Þ

To derive the expected value of u we need the following tools:

rY d rrY qX q rrY qX rX d rrY qX eX ¼ Xe X () XerX qX ¼ eX drY rX q drY

 Conditional distribution  function of bivariate correlated 

From Eqs. (D.5) and (D.6) it follows that rX 12ð1q2 Þr2Y d h rrY qX i E ½ u ¼ e E eX drY q

 First

X

Gaussians:Y jX ¼ N lY þ rrXY qð X  lX Þ; ð1  q2 Þr2Y moment

of

the

ðD:5Þ

lognormal

distribution:

r2 Z

E½eZ  ¼ elZ þ 2  Law of total expectation: E½ XY  ¼ E½ XE½Y jX 

ðD:7Þ

rY r

qX is normal with zero mean and variance 2 rY q r2X ¼ r2Y q2 . Applying the first moment of the logrX

normal distribution a second time results in

Applying the conditional distribution function results in     2 rY 2 Y jX ¼ N qX ; 1  q rY ðD:3Þ rX The law of total expectation then yields:   E½ XZ  ¼ E½ XE½ZjX  ¼ E XE eY jX   ¼ E XE eY jX

ðD:6Þ

E ½ u ¼

rX 12ð1q2 Þr2Y d 1r2 q2 e e2 Y drY q

Hence:  1 2 E XeY ¼ rX rY qe2rY

ðD:8Þ

ðD:9Þ

Eq. (D.9) describes the first moment of the generic NLNM. Appendix E. Ap to Kp polynomial coefficients ðD:4Þ

The first moment of the lognormal distribution results in:

The Ap to Kp third-order polynomial coefficients of the proposed cubic splining-approach are given in Table E.2.

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023

20

F. Schiemenz et al. / Advances in Space Research xxx (2019) xxx–xxx

The coefficients are to be used with Eq. (61) in order to perform the conversion from magnetic amplitude to the kpindex. The resulting curve closely follows Vallado’s algorithm (Vallado, 2013), however is free of abrupt slope changes at the outer ap-knots. References Anderson, R.L., Born, G.H., Forbes, J.M., 2009. Sensitivity of orbit predictions to density variability. J. Spacecraft Rockets 46 (6), 1214– 1230. Bartels, J., 1949. The standardized index, Ks, and the planetary index, Kp. IATME Bull. 12b 97 (2010), 2021. Bartels, J., Veldkamp, J., 1953. International data on magnetic disturbances, second quarter, 1953. J. Geophys. Res. 58 (4), 543–545. Bates, D.R., 1959. Some problems concerning the terrestrial atmosphere above about the 100 km level. Proc. R. Soc. Lond. A 253 (1275), 451– 462. Bowman, B., Tobiska, W.K., Marcos, F., Huang, C., Lin, C., Burke, W., 2008. A new empirical thermospheric density model JB2008 using new solar and geomagnetic indices. In: AIAA/AAS Astrodynamics Specialist Conference and Exhibit, p. 6438. Bruinsma, S., 2013. The semi-empirical thermosphere model DTM2012. In: 6th European Conference on Space Debris, vol. 723.. Bruinsma, S., 2015. The DTM-2013 thermosphere model. J. Space Weather Space Clim. 5, A1. Bruinsma, S.L., Forbes, J.M., 2008. Medium-to large-scale density variability as observed by CHAMP. Space Weather 6 (8). Chamberlain, T.P., Hunten, D.M., 1990. Theory of Planetary Atmospheres: An Introduction to Their Physics and Chemistry, vol. 36. Academic Press. Emmert, J., 2015. Thermospheric mass density: a review. Adv. Space Res. 56 (5), 773–824.

Emmert, J., Byers, J., Warren, H., Segerman, A., 2014. Propagation of forecast errors from the Sun to LEO trajectories: how does drag uncertainty affect conjunction frequency? In: Advanced Maui Optical and Space Surveillance Technologies Conference.. Emmert, J., Warren, H., Segerman, A., Byers, J., Picone, J., 2017. Propagation of atmospheric density errors to satellite orbits. Adv. Space Res. 59 (1), 147–165. Hejduk, M., 2017. Collision avoidance short course part I: theory. In: Advanced Maui Optical and Space Surveillance Technologies Conference; 19–22 Sep. 2017; Wailea, HI; United States.. Mandea, M., Korte, M., 2010. Geomagnetic Observations and Models, vol. 5. Springer. Picone, J., Meier, R., Emmert, J., 2013. Theoretical tools for studies of low-frequency thermospheric variability. J. Geophys. Res. Space Phys. 118 (9), 5853–5873. Picone, J.M. et al., 2002. NRLMSISE-00 empirical model of the atmosphere: Statistical comparisons and scientific issues. J. Geophys. Res. Space Phys. 107 (A12). Schiemenz, F., Utzmann, J., Kayal, H., 2019a. Propagating EUV solar flux uncertainty to atmospheric density uncertainty. Adv. Space Res. 63 (12), 3936–3952. Schiemenz, F., Utzmann, J., Kayal, H., 2019b. Least squares orbit estimation including atmospheric density uncertainty consideration. Adv. Space Res. 63 (12), 3916–3935. Storz, M.F., Bowman, B.R., Branson, M.J.I., Casali, S.J., Tobiska, W.K., 2005. High accuracy satellite drag model (HASDM). Adv. Space Res. 36 (12), 2497–2505. Vallado, D., 2013. Fundamentals of Astrodynamics and Applications, fourth ed. Microcosm Press. Wertz, J.R., 2012. Spacecraft Attitude Determination and Control. Springer Science & Business Media. Xu, W., 2008. Uncertainty in magnetic activity indices. Sci. China Ser. E: Technol. Sci. 51 (10), 1659–1664. Yang, M., 2008. Normal log-normal mixture, leptokurtosis and skewness. Appl. Econ. Lett. 15 (9), 737–742.

Please cite this article as: F. Schiemenz, J. Utzmann and H. Kayal, Accurate estimation of relative atmospheric density error on the example of uncertain geomagnetic activity information, Advances in Space Research, https://doi.org/10.1016/j.asr.2019.08.023