-
Notifications
You must be signed in to change notification settings - Fork 3
/
Copy pathmsis_gtd8d.F90
147 lines (136 loc) · 6.04 KB
/
msis_gtd8d.F90
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
!#######################################################################
! MSIS® (NRL-SOF-014-1) SOFTWARE
! NRLMSIS® empirical atmospheric model software. Use is governed by the
! Open Source Academic Research License Agreement contained in the file
! nrlmsis2.1_license.txt, which is part of this software package. BY
! USING OR MODIFYING THIS SOFTWARE, YOU ARE AGREEING TO THE TERMS AND
! CONDITIONS OF THE LICENSE.
!#######################################################################
!!! ===========================================================================
!!! NRLMSIS 2.1:
!!! Neutral atmosphere empirical model from the surface to lower exosphere
!!! ===========================================================================
!!!
!!! GTD8D: Legacy wrapper with input and output arguments used in NRLMSISE-00
!
! PREREQUISITES:
! Must first run MSISINIT to load parameters and set switches. The
! MSISCALC subroutine (called by this wrapper) checks for initialization
! and does a default initialization if necessary. This self-initialization
! will be removed in future versions.
!
! CALLING SEQUENCE:
! CALL GTD8D(IYD, SEC, ALT, GLAT, GLONG, STL, F107A, F107, AP, MASS, D, T)
!
! INPUT VARIABLES:
! IYD Year and day as YYDDD (day of year from 1 to 365 (or 366))
! (Year is ignored in current model)
! SEC Universal time (seconds)
! ALT Geodetic altitude (km)
! GLAT Geodetic latitude (deg)
! GLONG Geodetic longitude (deg)
! STL Local solar time (Ignored; calculated from SEC and GLONG)
! F107A 81 day average, centered on input time, of F10.7 solar activity
! index
! F107 Daily F10.7 for previous day
! AP Geomagnetic activity index array:
! (1) Daily Ap
! (2) 3 hr ap index for current time
! (3) 3 hr ap index for 3 hrs before current time
! (4) 3 hr ap index for 6 hrs before current time
! (5) 3 hr ap index for 9 hrs before current time
! (6) Average of eight 3 hr ap indices from 12 to 33 hrs
! prior to current time
! (7) Average of eight 3 hr ap indices from 36 to 57 hrs
! prior to current time
! AP(2:7) are only used when switch_legacy(9) = -1.0 in MSISINIT
! MASS Mass number (Ignored in 2.0)
!
! NOTES ON INPUT VARIABLES:
! - If lzalt_type = .false. in the MSISINIT call, then the ALT input
! argument is treated as geopotential height.
! - The STL input argument is ignored in NRLMSIS 2.0. Instead, local time
! is computed from universal time and longitude.
! - F107 and F107A values are the 10.7 cm radio flux at the Sun-Earth
! distance, not the radio flux at 1 AU.
! - The MASS input argument is ignored in NRLMSIS 2.0; species to be
! calculated are set in MSISINIT.
!
! OUTPUT VARIABLES:
! D(1) He number density (cm-3)
! D(2) O number density (cm-3)
! D(3) N2 number density (cm-3)
! D(4) O2 number density (cm-3)
! D(5) Ar number density (cm-3)
! D(6) Total mass density (g/cm3)
! D(7) H number density (cm-3)
! D(8) N number density (cm-3)
! D(9) Anomalous oxygen number density (cm-3)
! D(10) NO number density (cm-3)
! T(1) Exospheric temperature (K)
! T(2) Temperature at altitude (K)
!
! NOTES ON OUTPUT VARIABLES:
! - Missing density values are returned as 9.999e-38
! - Species included in mass density calculation are set in MSISINIT
!
!!! =========================================================================
!==================================================================================================
! GTD8D: Legacy wrapper
!==================================================================================================
subroutine gtd8d(iyd,sec,alt,glat,glong,stl,f107a,f107,ap,mass,d,t)
use msis_constants, only : rp, dmissing
use msis_init, only : msisinit
use msis_calc, only : msiscalc
implicit none
! MSIS Legacy subroutine arguments
integer, intent(in) :: iyd
real(4), intent(in) :: sec
real(4), intent(in) :: alt
real(4), intent(in) :: glat
real(4), intent(in) :: glong
real(4), intent(in) :: stl
real(4), intent(in) :: f107a
real(4), intent(in) :: f107
real(4), intent(in) :: ap(7)
integer, intent(in) :: mass
real(4), intent(out) :: d(10), t(2)
! MSIS 1.97 subroutine arguments
real(kind=rp) :: xday
real(kind=rp) :: xutsec
real(kind=rp) :: xalt
real(kind=rp) :: xlat
real(kind=rp) :: xlon
real(kind=rp) :: xsfluxavg, xsflux
real(kind=rp) :: xap(1:7)
real(kind=rp) :: xtn
real(kind=rp) :: xdn(1:10)
real(kind=rp) :: xtex
! Convert the legacy input arguments to the new interface values and precision
xday = mod(iyd,1000)
xutsec = sec
xalt = alt
xlat = glat
xlon = glong
xsfluxavg = f107a
xsflux = f107
xap = ap
! Call the new subroutine
call msiscalc(xday,xutsec,xalt,xlat,xlon,xsfluxavg,xsflux,xap,xtn,xdn,tex=xtex)
! Convert the output arguments to the legacy format (mks to cgs, re-order species)
t(1) = sngl(xtex) ! Expospheric temperature
t(2) = sngl(xtn) ! Temperature at altitude
where (xdn .ne. dmissing) xdn = xdn*1d-6
if (xdn(1) .ne. dmissing) xdn(1) = xdn(1)*1e3_rp
d(1) = sngl(xdn(5)) ! [He]
d(2) = sngl(xdn(4)) ! [O]
d(3) = sngl(xdn(2)) ! [N2]
d(4) = sngl(xdn(3)) ! [O2]
d(5) = sngl(xdn(7)) ! [Ar]
d(6) = sngl(xdn(1)) ! Mass density
d(7) = sngl(xdn(6)) ! [H]
d(8) = sngl(xdn(8)) ! [N]
d(9) = sngl(xdn(9)) ! [Anomalous O]
d(10) = sngl(xdn(10)) ! [NO]
return
end subroutine gtd8d