A (huge) collection of interchangeable empirical models to compute the thermosphere density, temperature and composition

Thermosphere models are usually employed in space activities such as mission planning and analysis, orbit determination, attitude control design, maneuvering planning, fuel estimation and orbit decay studies, among others.
 
The Orbital Dynamics group of INPE coded
during the 80’s several thermospheric computer programs in FORTRAN. Those models were based on Jacchia works (70 and 77) and their analytical versions, since Jacchia assumed numeric integration of the diffusion equation in order to compute density and composition as function of height. Numeric integration normally demands computing time, and therefore analytical means “fast” in this sense. All the models were recently converted to C/C++ by Carrara, because FORTRAN compilers seem to be quickly vanishing nowadays.

Presently there is a set of 6 thermospheric models:
  • Jacchia 70 (in fact it is almost the same as the 64 model)
  • Jacchia-Roberts
  • Jacchia 77
  • LaFontaine-Hughes
  • Jacchia-Lineberry (77 model)
  • MSIS-86 (C code)
All these programs are interchangeable, which means that any model, except MSIS, can be easily changed to any other. Of course there are several models not included in this list, but, with few exceptions, they do not satisfy the requirement of "interchangeability" since they need different input data, like, for instance, the Jacchia-Bowman JB2006 and JB2008 models (Bowman et al. 2008a, 2008b), which considers the effect of solar UV radiation.

The interchangeable models (except MSIS again) present 3 entry points:
  • Static profile (function of the exospheric temperature and height)
  • Static density (same as static profile, but function of solar flux F10.7 and geomagnetic activity besides height).
  • Dynamic model.
In order to make it easy to remember and to exchange between models, they were respectively named as:
  • xmowei
  • xsmade
  • xsdamo
both in FORTRAN and C, where x can be “j” for Jacchia 70, “r” for Robert’s analytical version of Jacchia 70, “i” for Jacchia 77 model (numerical integration), “a” for the analytical version from LaFontaine and Hughes of Jacchia’s 77, and “p” for the Lineberry polynomial analytical version of Jacchia’s 77. “mowei” stands for MOecular WEIght, since the static profile produces mainly the molecular weight as function of height. “smade” means Static Model for Atmospheric DEnsity, and finally “sdamo” is an acronymic for Static and Dynamic Atmospheric MOdel, as seem in Table 1.

Table 1 - Profile, Static and Dynamic functions


Static profile Static density Dynamic model
Jacchia 70 jmowei jsmade jsdamo
Jacchia-Roberts rmowei rsmade rsdamo
Jacchia 77 imowei ismade isdamo
LaFontaine-Hughes amowei asmade asdamo
Jacchia-Lineberry (77 model) pmowei psmade psdamo
MSIS-86 ---
--- msdamo (C only)

FORTRAN versions were compiled and tested with FORTRAN  Powerstation 4 and  Compact Visual Fortran 6 on Windows XP. C versions were implemented and tested with MS Visual C++ Express running on Windows XP.
Parameters and usage of all these functions can be found in comments in the source code, or in Carrara, 1990 (in portuguese). Very short descriptions of these models are presented bellow.

Jacchia’s 70.

Jacchia produced his first model in 1964 (Jacchia, 1964 - subroutine ADEN). An improved version appeared later on CIRA 72.

Jacchia-Roberts

Soon after the first Jacchia’s model, Roberts (Roberts, 1971) produced an analytical version.

Jacchia’s 77

Jacchia improved his model with recently (at that time) satellite observations and orbit measuring. The 77’s integrates the diffusion equation for each individual atmospheric component rather than the mean molecular weight, unlike the preceding model. The version available here was coded by Mattos in 1984, in FORTRAN
.

LaFontaine-Hughes

LaFontaine and Hughes (LaFontaine and Hughes, 1983) presented their analytical version of Jacchia’s 77 model in 83. Again Mattos coded it in FORTRAN.

Jacchia-Lineberry

Lineberry suggested that the static profile could be approximated by a truncated McLaurin serie. The original version made by Mueller (Mueller, 1982) was based on the Jacchia’s 70. In 1984 Mattos produced a Jacchia-Lineberry version based on the 77’s model.
Unfortunately the author didn’t publish it (except by an oral presentation in a conference – Mattos, 1984), but a brief description of the work appeared in a portuguese publication of INPE (Carrara, 1990).

MSIS-86

There are several “MSIS” models, namely the MSIS-86, MSISE-90 and NRLMSIS-00. They can be easily found and downloaded at Goddard NASA’s Model Web . Above 120 km they are essentially the same, since major differences occur at low altitudes. They were coded, naturally, in FORTRAN. Some C/C++ versions can be found on net, but it seems that there is no C version for 86’s model, except by those converted with f2c, which leaves lots of junk in the code. MSIS models were produced by Hedin (CIRA, 1986; Hedin, 1987) and coded in Fortran by D. Bilitza. The version available here is the C version of MSIS-86, converted to C by Carrara (manually, no f2c). Sorry, there is no manual yet available for C version, but the parameters are almost identical to the Fortran version.

Model Parameters

Although some models has it's own parameters and output, an interface
function was created so as to make them interchangeable both in input and output parameters. Tables 2 and 3 presented below show those parameters and a simplified explanation of each. You may refer to specific model documentation for a deeper explanation.

Table 2 - Thermospheric model input parameters

Par Type Meaning static
profile
static
density
dynamic
model
tmp
double
exospheric temperature, Kelvin


alt
double
local geocentric altitude, meters (90000 to 2000000)

sa array, double, 3 components
1- local right ascension, radians (0 to 2.*pi)
2- local declination (geocentric latitude), radians (0 to 2.*pi)
3- local geocentric altitude, meters (90000 to 2000000)



su array, double, 2 components
1- Sun right ascension, radians (0 to 2.*pi).
2-
Sun declination, radians (-pi to pi)


sf array, double, 3 components 1- daily solar flux F10.7, Jensen
2- averaged daily flux, Jensen
3- 3 hourly planetary geomagnetic index Kp, according to model specific definition.
mjd
double modified Julian date referred to 1950.0


dayf double UTC time, seconds


gsti double Greenwich sidereal time, radians (0 to 2.*pi)

       
Table 2 - Thermospheric model output parameters

Par Type Meaning static
profile
static
density
dynamic
model
tp array, double, 2 components
1- exospheric temperature, Kelvin.
2- local temperature, Kelvin

ad
array, double, 6 components
1- log10 of Helium number-density
2-
log10 of molecular Oxigen number-density
3- log10 of molecular Nitrogen number-density
4- log10 of Argon number-density
5- log10 of atomic Oxigen number-density
6- log10 of atomic Hydrogen number-density
wmol
double thermospheric mean molecular weight, kg/kgmol.
dens
double
thermospheric density, kg/m3
 
MSIS dynamic model has its own parameters. To make it compatible to Jacchia models, some changes in its parameters had to be made. The sf parameter of the msdamo function must have 7 components. Normally you do not need to provide the values of the components
4 to 7, and the third one has to be the daily geomagnetic amplitude Ap instead of the 3-hourly Kp. If you modify the sw (switch) array (see MSIS-86 specific documentation), then you must provide 6 components of sf (2 to 7). The output values of MSIS also differs from all the other models. MSIS computes the atomic Nitrogen (N) number-density, while Jacchia models don't. The N number-density contributes to the density and molecular weight of the parameters wmol and dens, but the ad output array does not return its value. You may change it in the source code if you wish.

Download Models

All the libraries have the source code and all the necessary files to run a test sample. They are zip compressed. Please download and expand it at your project folder.

The FORTRAN versions come with the following files:
  • A sample main program (xxxx_test.for)
  • The thermospheric model library
  • A common library atsu.cpp or atsu.for for models based on Jacchia's 77, and
  • Three data files corresponding to the output of the Profile, Static and Dynamic functions.
Additionally, C versions bring
  • A common header atmodels.h for all the models, except MSIS.
  • A model specific header (like jac70.h, jac77.h or lafon.h, for instance).
  Model
FORTRAN C/C++
Download the Jacchia's 70 jac70f.zip jac70c.zip
Download the Jacchia-Roberts roberf.zip
roberc.zip
Download the Jacchia's 77 jac77f.zip jac77c.zip
Download the LaFontaine-Hughes lafonf.zip
lafonc.zip
Download the Jacchia-Lineberry (77 model) linebf.zip
linebc.zip
Download the MSIS-86
link msisc.zip

Solar Flux Data File

All the thermospheric models shown above depend on the mean and observed solar flux at 10.7 cm wave length at a given date. The flux varies with solar activity, which presents a 10.6 years cycle.
The Dominion Radio Astrophysical Observatory (DRAO, 2011) measures the flux daily. The data is available for download at the Observatory site.

During solar storms, the solar wind deflects the geomagnetic field, which causes a significant heating in the thermospheric molecules. This effect is named geomagnetic activity and is measured by the International Service of Geomagnetic Indices (ISGI, 2011), which also make it available for download.

Both the solar flux and geomagnetic activity are compiled and stored in a single file by the Space Physics Interactive Data Resources (SPDIR, 2011), branch of the National Geophysic Data Center. The data can be found at this link:


Nevertheless, the mean values of the solar flux can't be found anywhere.

So, in order to get a single computation of the thermospheric properties, one should find the daily solar flux, the geomagnetic indices and needs also to compute the mean flux, which depends on the thermospheric model. To make it easy to you, all this data and the program that generates it is available here for download - you may see it at next sections.

The data file stores the geomagnetic activity, the mean solar flux as well as the observed and adjusted daily fluxes. Jacchia (Jacchia 1964, 1977) didn't mention if the mean is performed on the observed or adjusted fluxes. To make things worse, there is some indirect evidence that the mean is "centered" on the day of computation, which means that the today's thermospheric temperature depends on the tomorrow's flux!!!! Even considering that this doesn't make sense, the file brings the "centered" mean values computed in two different ways: a gaussian mean, according to Jacchia's 77 model, and an arithmetic mean, adopted by the Jacchia's 70 and derived models. The gaussian mean uses a standard deviation of
70 days and is computed on +/- 3 sd interval, or 420 days plus the computed day. The arithmetic mean is performed in +/- 3 solar rotations (27 days each), or 163 days, including the centered day. Both were calculated with the observed flux.

Solar Flux Data File Download

The solar flux data file - flux_mean.dat - stores the following data in ASCII format:

1st record:
  • Modified Julian date (referred to 1950.0) of the first record (normally Jan, 1st, 1958). Integer, bytes 1-5
  • 1 blank space (column 6)
  • Modified Julian date of the last record. Integer, bytes 7-11
  • 36 blank spaces in columns 12 to 47 (this is necessary to allow direct record readings)
Remaining records
  • Modified Julian date (referred to 1950.0) of the record, bytes 1-5
  • Planetary 3 hour range index kp for 0 hour UTC, bytes 6-7
  • Planetary 3 hour range index kp for 3 hour UTC, bytes 8-9
  • Planetary 3 hour range index kp for 6 hour UTC, bytes 10-11
  • Planetary 3 hour range index kp for 9 hour UTC, bytes 12-13
  • Planetary 3 hour range index kp for 12 hour UTC, bytes 14-15
  • Planetary 3 hour range index kp for 15 hour UTC, bytes 16-17
  • Planetary 3 hour range index kp for 18 hour UTC, bytes 18-19
  • Planetary 3 hour range index kp for 21 hour UTC, bytes 20-21
  • Sum of the 8 values of kp, bytes 22-24
  • Planetary equivalent daily amplitude Ap, bytes 25-27
  • Observed daily solar flux at 10.7 cm, bytes 28-32
  • Daily solar flux adjusted to 1 Astronomical Unit, bytes 33-37
  • Gaussian mean of the solar flux, according to Jacchia (Jacchia, 1977), bytes 38-42
  • Arithmetic mean of the solar flux, according to Jacchia (Jacchia, 1964), bytes 43-47

Download the mean solar flux data file: 1958-1-1 to 2008-12-31 flux_mean.dat - flux_mean.zip

Solar Flux library

The flux data can be retrieved from the data archive by providing the date and the day fraction in seconds. The date shall be given in modified Julian date referred to 1950.0. Functions to convert the date (day, month and year) to and from modified Julian date and to compute the Greenwich sidereal time were provided in the package
. The library also have functions to retrieve the observed and mean data of a given day, and to provide the proper values to compute the thermospheric properties according to Jacchia's 70 and 77 models. You may want to see the header comments in the package in order to know how to use it. An example main program is also provided. Feel free to fit it to your needs.

The package contains:
  • A sample main program (soflud_test.for)
  • The function sfluxd, to retrieve the data of a given date from the mean solar flux data file (flux_mean.dat),
  • A function - sfdj70 - to retrieve the observed, mean and geomagnetic activity of a given date according to Jacchia's 70 model.
  • A function - sfdj77 - to retrieve the observed, mean and geomagnetic activity of a given date according to Jacchia's 77 model.
  • A function - gst - to calculate the Greenwich sidereal time (required by Jacchia's 77 model and its analytical versions).
  • A function - djm - to calculate the modified Julian date from a given date.
  • A function - djminv - to calculate the date given the modified Julian date.
  • The functions open_sfdf and close_sfdf to open and to close the solar flux data file (flux_mean.dat) - C/C++ only. The FORTRAN version of the package doesn't need these functions.

FORTRAN C/C++
Download the mean solar flux library
sflibf.zip sflibc.zip

Additionally a sample program (solflud_test) was included in the package, as well as its output data file.

Solar Flux Data File update program

Once you have successfully compiled and run the thermospheric model with solar flux data file, it is time to update it with fresh new values. You may update it once a year, since the observed flux and geomagnetic indexes are also available yearly. So, go to
ftp://ftp.ngdc.noaa.gov/STP/GEOMAGNETIC_DATA/INDICES/KP_AP/ and download the most recent data, probably something like 2011 or 2012. Since the average values are centered at the date, you need also the files for preceeding and succeding years. To compute the 2010 mean values, for instance, the 2009 and 2011 files are also required. These files shall be renamed so as to include a ".dat" extension. The program allows you to input the year of computation, as well as all the required information needed to run it properly. You should be advised that the solar flux update file program doesn't check the input data for consistency, and this means that you have to assure that the data is provided in right sequential order. After execution, a file yyyy_m.dat shall be created, where yyyy is the computation year. Now is up to you to append this file to the Solar Flux Data File (flux_mean.dat) using an ordinary text editor like Notepad or Notepad++. You shall assure that the modified julian date of the yyyy_m file first record follows the MJD of the flux_mean file last record. After that you have to change manually the first record of the file, in order to update the date of the last record (colunms 7 to 11) with the new value. Now you may save the new updated flux_mean.dat file.

The package provided bellow for download includes also the Windows console .exe file. The source code is available in C language only. Unless you want to modify it in order to include something new, I recommend you to simple run it when you need.

Download the Solar Flux Data File update program
solfluxmean.zip

So, that's all for now. Good luck and happy modeling.

Valdemir Carrara, Sep, 22, 2011.


References

Bowman, B. R.; Tobiska, W. K.; Marcos, F. A.; Valladares, C. The JB2006 empirical thermospheric density model. In: Journal of Atmospheric and Solar-Terrestrial Physics, 70, pp. 774-793, 2008a. Available in: http://sol.spacenvironment.net/~JB2006/pubs/JASTP_Bowman_JB2006_2008.pdf

Bowman, B. R.; Tobiska, W. K. Marcos, F. A.; Huang, C. Y.; Lin, C. S.; Burke, W. J. A New Empirical Thermospheric Density Model JB2008 Using New Solar and Geomagnetic Indices. AIAA/AAS Astrodynamics Specialist Conference. Proceedings. Honolulu, Aug. 2008b. Available in: http://sol.spacenvironment.net/~JB2008/pubs/AIAA_2008-6438_JB2008_Model.pdf

Carrara, V. Implementações de modelos atmosféricos para uso em propagadores de órbita e atitude. S. J. Campos, INPE, maio 1990 (INPE-5094-RPI/231).

CIRA 1986, Part I: Thermosphere Model, D. Rees (ed.), Adv. Space Res. Vol. 8, n. 5, n. 6, 1988.

DRAO – Dominion Radio Astrophysical Observatory – Herzberg Institute of Astrophysics. Solar Radio Monitoring Programme. http://www.drao-ofr.hia-iha.nrc-cnrc.gc.ca 2011.

Hedin, A. E. MSIS-86 Thermospheric Model. In: Jounal of Geophysical Research, Vol. 92, n. A5, p. 4649-4662, May, 1987.

ISGI – International Service of Geomagnetic Indices. Indices: Data http://isgi.cetp.ipsl.fr/. Access in 2011.

Jacchia, L. G. Static Diffusion Models of the Upper Atmosphere With Empirical Temperature Profiles. Cambridge, MA, SAO, 1964 (SAO Special Report 170).

Jacchia, L. G. Atmospheric models in the region from 110 to 2000 km. In: COMMITTEE ON SPACE RESEARCH (COSPAR). CIRA 1972. Berlim, AkademicVerlag, 1972. Part 3, p. 227-338.

Jacchia, L. G. Thermospheric Temperature, Density and Compostion: New Models. Cambridge,  Ma,  SAO, 1977. (SAO Special Report No 375).

Lafontaine, J.; Hughes, P. An Analytic Version of Jacchia's 1977 model atmosphere. Celestial Mechanics, 29(3-26) 1983.

Mattos, B. S. Modelamento da Atmosfera Superior da Terra. In: II COLÓQUIO DE DINÂMICA ORBITAL E RESSONÂNCIA GRAVITACIONAL, 1984. ATIBAIA-SP.

Mueller, A. C. Jacchia-Lineberry upper atmosphere density model. Huston NASA, 1982. (NASA-CR-167824).

Roberts Jr., C. E. An analytical model for upper atmospheric densities based upon Jacchia's 1970 models. Celestial Mechanics, 4(3/4):368-377, Dec. 1971.

SPDIR - Space Physics Interactive Data Resources - National Geophysical Data Center. Geomagnetic Indices Data Set. http://spidr.ngdc.noaa.gov/spidr/ 2011.