chgres_cube  1.9.0
All Data Structures Namespaces Files Functions Variables Pages
wam_climo_data.f90
Go to the documentation of this file.
1 
12 
13 !-----------------------------------------------------------------------
30 ! msise-00 01-feb-02
31 !
32  real :: pt1(50)
33  real :: pt2(50)
34  real :: pt3(50)
35  real :: pa1(50)
36  real :: pa2(50)
37  real :: pa3(50)
38  real :: pb1(50)
39  real :: pb2(50)
40  real :: pb3(50)
41  real :: pc1(50)
42  real :: pc2(50)
43  real :: pc3(50)
44  real :: pd1(50)
45  real :: pd2(50)
46  real :: pd3(50)
47  real :: pe1(50)
48  real :: pe2(50)
49  real :: pe3(50)
50  real :: pf1(50)
51  real :: pf2(50)
52  real :: pf3(50)
53  real :: pg1(50)
54  real :: pg2(50)
55  real :: pg3(50)
56  real :: ph1(50)
57  real :: ph2(50)
58  real :: ph3(50)
59  real :: pi1(50)
60  real :: pi2(50)
61  real :: pi3(50)
62  real :: pj1(50)
63  real :: pj2(50)
64  real :: pj3(50)
65  real :: pk1(50)
66  real :: pl1(50)
67  real :: pl2(50)
68  real :: pm1(50)
69  real :: pm2(50)
70  real :: pn1(50)
71  real :: pn2(50)
72  real :: po1(50)
73  real :: po2(50)
74  real :: pp1(50)
75  real :: pp2(50)
76  real :: pq1(50)
77  real :: pq2(50)
78  real :: pr1(50)
79  real :: pr2(50)
80  real :: ps1(50)
81  real :: ps2(50)
82  real :: pu1(50)
83  real :: pu2(50)
84  real :: pv1(50)
85  real :: pv2(50)
86  real :: pw1(50)
87  real :: pw2(50)
88  real :: px1(50)
89  real :: px2(50)
90  real :: py1(50)
91  real :: py2(50)
92  real :: pz1(50)
93  real :: pz2(50)
94  real :: paa1(50)
95  real :: paa2(50)
96 !
97  real :: ptm(10)
98  real :: pdm(10,8)
99 !
100  real :: pavgm(10)
101 !
102  character*4:: isdate(3)
103  character*4:: istime(2)
104  character*4:: name(2)
105 !
106  integer :: imr
107 !
108  real :: pr65(2,65)
109  real :: pr151(2,151)
110 
111  data imr/0/
112  data isdate/'01-f','eb-0','2 '/,istime/'15:4','9:27'/
113  data name/'msis','e-00'/
114 ! temperature
115  data pt1/ &
116  9.86573e-01, 1.62228e-02, 1.55270e-02,-1.04323e-01,-3.75801e-03,&
117  -1.18538e-03,-1.24043e-01, 4.56820e-03, 8.76018e-03,-1.36235e-01,&
118  -3.52427e-02, 8.84181e-03,-5.92127e-03,-8.61650e+00, 0.00000e+00,&
119  1.28492e-02, 0.00000e+00, 1.30096e+02, 1.04567e-02, 1.65686e-03,&
120  -5.53887e-06, 2.97810e-03, 0.00000e+00, 5.13122e-03, 8.66784e-02,&
121  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00,-7.27026e-06,&
122  0.00000e+00, 6.74494e+00, 4.93933e-03, 2.21656e-03, 2.50802e-03,&
123  0.00000e+00, 0.00000e+00,-2.08841e-02,-1.79873e+00, 1.45103e-03,&
124  2.81769e-04,-1.44703e-03,-5.16394e-05, 8.47001e-02, 1.70147e-01,&
125  5.72562e-03, 5.07493e-05, 4.36148e-03, 1.17863e-04, 4.74364e-03/
126  data pt2/ &
127  6.61278e-03, 4.34292e-05, 1.44373e-03, 2.41470e-05, 2.84426e-03,&
128  8.56560e-04, 2.04028e-03, 0.00000e+00,-3.15994e+03,-2.46423e-03,&
129  1.13843e-03, 4.20512e-04, 0.00000e+00,-9.77214e+01, 6.77794e-03,&
130  5.27499e-03, 1.14936e-03, 0.00000e+00,-6.61311e-03,-1.84255e-02,&
131  -1.96259e-02, 2.98618e+04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
132  6.44574e+02, 8.84668e-04, 5.05066e-04, 0.00000e+00, 4.02881e+03,&
133  -1.89503e-03, 0.00000e+00, 0.00000e+00, 8.21407e-04, 2.06780e-03,&
134  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
135  -1.20410e-02,-3.63963e-03, 9.92070e-05,-1.15284e-04,-6.33059e-05,&
136  -6.05545e-01, 8.34218e-03,-9.13036e+01, 3.71042e-04, 0.00000e+00/
137  data pt3/ &
138  4.19000e-04, 2.70928e-03, 3.31507e-03,-4.44508e-03,-4.96334e-03,&
139  -1.60449e-03, 3.95119e-03, 2.48924e-03, 5.09815e-04, 4.05302e-03,&
140  2.24076e-03, 0.00000e+00, 6.84256e-03, 4.66354e-04, 0.00000e+00,&
141  -3.68328e-04, 0.00000e+00, 0.00000e+00,-1.46870e+02, 0.00000e+00,&
142  0.00000e+00, 1.09501e-03, 4.65156e-04, 5.62583e-04, 3.21596e+00,&
143  6.43168e-04, 3.14860e-03, 3.40738e-03, 1.78481e-03, 9.62532e-04,&
144  5.58171e-04, 3.43731e+00,-2.33195e-01, 5.10289e-04, 0.00000e+00,&
145  0.00000e+00,-9.25347e+04, 0.00000e+00,-1.99639e-03, 0.00000e+00,&
146  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
147  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
148 ! he density
149  data pa1/ &
150  1.09979e+00,-4.88060e-02,-1.97501e-01,-9.10280e-02,-6.96558e-03,&
151  2.42136e-02, 3.91333e-01,-7.20068e-03,-3.22718e-02, 1.41508e+00,&
152  1.68194e-01, 1.85282e-02, 1.09384e-01,-7.24282e+00, 0.00000e+00,&
153  2.96377e-01,-4.97210e-02, 1.04114e+02,-8.61108e-02,-7.29177e-04,&
154  1.48998e-06, 1.08629e-03, 0.00000e+00, 0.00000e+00, 8.31090e-02,&
155  1.12818e-01,-5.75005e-02,-1.29919e-02,-1.78849e-02,-2.86343e-06,&
156  0.00000e+00,-1.51187e+02,-6.65902e-03, 0.00000e+00,-2.02069e-03,&
157  0.00000e+00, 0.00000e+00, 4.32264e-02,-2.80444e+01,-3.26789e-03,&
158  2.47461e-03, 0.00000e+00, 0.00000e+00, 9.82100e-02, 1.22714e-01,&
159  -3.96450e-02, 0.00000e+00,-2.76489e-03, 0.00000e+00, 1.87723e-03/
160  data pa2/ &
161  -8.09813e-03, 4.34428e-05,-7.70932e-03, 0.00000e+00,-2.28894e-03,&
162  -5.69070e-03,-5.22193e-03, 6.00692e-03,-7.80434e+03,-3.48336e-03,&
163  -6.38362e-03,-1.82190e-03, 0.00000e+00,-7.58976e+01,-2.17875e-02,&
164  -1.72524e-02,-9.06287e-03, 0.00000e+00, 2.44725e-02, 8.66040e-02,&
165  1.05712e-01, 3.02543e+04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
166  -6.01364e+03,-5.64668e-03,-2.54157e-03, 0.00000e+00, 3.15611e+02,&
167  -5.69158e-03, 0.00000e+00, 0.00000e+00,-4.47216e-03,-4.49523e-03,&
168  4.64428e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
169  4.51236e-02, 2.46520e-02, 6.17794e-03, 0.00000e+00, 0.00000e+00,&
170  -3.62944e-01,-4.80022e-02,-7.57230e+01,-1.99656e-03, 0.00000e+00/
171  data pa3/ &
172  -5.18780e-03,-1.73990e-02,-9.03485e-03, 7.48465e-03, 1.53267e-02,&
173  1.06296e-02, 1.18655e-02, 2.55569e-03, 1.69020e-03, 3.51936e-02,&
174  -1.81242e-02, 0.00000e+00,-1.00529e-01,-5.10574e-03, 0.00000e+00,&
175  2.10228e-03, 0.00000e+00, 0.00000e+00,-1.73255e+02, 5.07833e-01,&
176  -2.41408e-01, 8.75414e-03, 2.77527e-03,-8.90353e-05,-5.25148e+00,&
177  -5.83899e-03,-2.09122e-02,-9.63530e-03, 9.77164e-03, 4.07051e-03,&
178  2.53555e-04,-5.52875e+00,-3.55993e-01,-2.49231e-03, 0.00000e+00,&
179  0.00000e+00, 2.86026e+01, 0.00000e+00, 3.42722e-04, 0.00000e+00,&
180  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
181  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
182 ! o density
183  data pb1/ &
184  1.02315e+00,-1.59710e-01,-1.06630e-01,-1.77074e-02,-4.42726e-03,&
185  3.44803e-02, 4.45613e-02,-3.33751e-02,-5.73598e-02, 3.50360e-01,&
186  6.33053e-02, 2.16221e-02, 5.42577e-02,-5.74193e+00, 0.00000e+00,&
187  1.90891e-01,-1.39194e-02, 1.01102e+02, 8.16363e-02, 1.33717e-04,&
188  6.54403e-06, 3.10295e-03, 0.00000e+00, 0.00000e+00, 5.38205e-02,&
189  1.23910e-01,-1.39831e-02, 0.00000e+00, 0.00000e+00,-3.95915e-06,&
190  0.00000e+00,-7.14651e-01,-5.01027e-03, 0.00000e+00,-3.24756e-03,&
191  0.00000e+00, 0.00000e+00, 4.42173e-02,-1.31598e+01,-3.15626e-03,&
192  1.24574e-03,-1.47626e-03,-1.55461e-03, 6.40682e-02, 1.34898e-01,&
193  -2.42415e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00, 6.13666e-04/
194  data pb2/ &
195  -5.40373e-03, 2.61635e-05,-3.33012e-03, 0.00000e+00,-3.08101e-03,&
196  -2.42679e-03,-3.36086e-03, 0.00000e+00,-1.18979e+03,-5.04738e-02,&
197  -2.61547e-03,-1.03132e-03, 1.91583e-04,-8.38132e+01,-1.40517e-02,&
198  -1.14167e-02,-4.08012e-03, 1.73522e-04,-1.39644e-02,-6.64128e-02,&
199  -6.85152e-02,-1.34414e+04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
200  6.07916e+02,-4.12220e-03,-2.20996e-03, 0.00000e+00, 1.70277e+03,&
201  -4.63015e-03, 0.00000e+00, 0.00000e+00,-2.25360e-03,-2.96204e-03,&
202  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
203  3.92786e-02, 1.31186e-02,-1.78086e-03, 0.00000e+00, 0.00000e+00,&
204  -3.90083e-01,-2.84741e-02,-7.78400e+01,-1.02601e-03, 0.00000e+00/
205  data pb3/ &
206  -7.26485e-04,-5.42181e-03,-5.59305e-03, 1.22825e-02, 1.23868e-02,&
207  6.68835e-03,-1.03303e-02,-9.51903e-03, 2.70021e-04,-2.57084e-02,&
208  -1.32430e-02, 0.00000e+00,-3.81000e-02,-3.16810e-03, 0.00000e+00,&
209  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
210  0.00000e+00,-9.05762e-04,-2.14590e-03,-1.17824e-03, 3.66732e+00,&
211  -3.79729e-04,-6.13966e-03,-5.09082e-03,-1.96332e-03,-3.08280e-03,&
212  -9.75222e-04, 4.03315e+00,-2.52710e-01, 0.00000e+00, 0.00000e+00,&
213  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
214  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
215  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
216 ! n2 density
217  data pc1/ &
218  1.16112e+00, 0.00000e+00, 0.00000e+00, 3.33725e-02, 0.00000e+00,&
219  3.48637e-02,-5.44368e-03, 0.00000e+00,-6.73940e-02, 1.74754e-01,&
220  0.00000e+00, 0.00000e+00, 0.00000e+00, 1.74712e+02, 0.00000e+00,&
221  1.26733e-01, 0.00000e+00, 1.03154e+02, 5.52075e-02, 0.00000e+00,&
222  0.00000e+00, 8.13525e-04, 0.00000e+00, 0.00000e+00, 8.66784e-02,&
223  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
224  0.00000e+00,-2.50482e+01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
225  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-2.48894e-03,&
226  6.16053e-04,-5.79716e-04, 2.95482e-03, 8.47001e-02, 1.70147e-01,&
227  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
228  data pc2/ &
229  0.00000e+00, 2.47425e-05, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
230  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
231  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
232  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
233  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
234  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
235  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
236  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
237  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
238  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
239  data pc3/ &
240  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
241  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
242  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
243  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
244  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
245  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
246  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
247  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
248  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
249  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
250 ! tlb
251  data pd1/ &
252  9.44846e-01, 0.00000e+00, 0.00000e+00,-3.08617e-02, 0.00000e+00,&
253  -2.44019e-02, 6.48607e-03, 0.00000e+00, 3.08181e-02, 4.59392e-02,&
254  0.00000e+00, 0.00000e+00, 0.00000e+00, 1.74712e+02, 0.00000e+00,&
255  2.13260e-02, 0.00000e+00,-3.56958e+02, 0.00000e+00, 1.82278e-04,&
256  0.00000e+00, 3.07472e-04, 0.00000e+00, 0.00000e+00, 8.66784e-02,&
257  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
258  0.00000e+00, 0.00000e+00, 3.83054e-03, 0.00000e+00, 0.00000e+00,&
259  -1.93065e-03,-1.45090e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
260  0.00000e+00,-1.23493e-03, 1.36736e-03, 8.47001e-02, 1.70147e-01,&
261  3.71469e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
262  data pd2/ &
263  5.10250e-03, 2.47425e-05, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
264  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
265  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
266  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
267  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
268  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
269  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
270  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
271  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
272  0.00000e+00, 3.68756e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00/
273  data pd3/ &
274  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
275  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
276  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
277  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
278  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
279  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
280  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
281  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
282  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
283  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
284 ! o2 density
285  data pe1/ &
286  1.35580e+00, 1.44816e-01, 0.00000e+00, 6.07767e-02, 0.00000e+00,&
287  2.94777e-02, 7.46900e-02, 0.00000e+00,-9.23822e-02, 8.57342e-02,&
288  0.00000e+00, 0.00000e+00, 0.00000e+00, 2.38636e+01, 0.00000e+00,&
289  7.71653e-02, 0.00000e+00, 8.18751e+01, 1.87736e-02, 0.00000e+00,&
290  0.00000e+00, 1.49667e-02, 0.00000e+00, 0.00000e+00, 8.66784e-02,&
291  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
292  0.00000e+00,-3.67874e+02, 5.48158e-03, 0.00000e+00, 0.00000e+00,&
293  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
294  0.00000e+00, 0.00000e+00, 0.00000e+00, 8.47001e-02, 1.70147e-01,&
295  1.22631e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
296  data pe2/ &
297  8.17187e-03, 3.71617e-05, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
298  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
299  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-2.10826e-03,&
300  -3.13640e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
301  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
302  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
303  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
304  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
305  -7.35742e-02,-5.00266e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
306  0.00000e+00, 1.94965e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00/
307  data pe3/ &
308  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
309  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
310  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
311  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
312  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
313  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
314  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
315  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
316  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
317  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
318 ! ar density
319  data pf1/ &
320  1.04761e+00, 2.00165e-01, 2.37697e-01, 3.68552e-02, 0.00000e+00,&
321  3.57202e-02,-2.14075e-01, 0.00000e+00,-1.08018e-01,-3.73981e-01,&
322  0.00000e+00, 3.10022e-02,-1.16305e-03,-2.07596e+01, 0.00000e+00,&
323  8.64502e-02, 0.00000e+00, 9.74908e+01, 5.16707e-02, 0.00000e+00,&
324  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 8.66784e-02,&
325  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
326  0.00000e+00, 3.46193e+02, 1.34297e-02, 0.00000e+00, 0.00000e+00,&
327  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-3.48509e-03,&
328  -1.54689e-04, 0.00000e+00, 0.00000e+00, 8.47001e-02, 1.70147e-01,&
329  1.47753e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
330  data pf2/ &
331  1.89320e-02, 3.68181e-05, 1.32570e-02, 0.00000e+00, 0.00000e+00,&
332  3.59719e-03, 7.44328e-03,-1.00023e-03,-6.50528e+03, 0.00000e+00,&
333  1.03485e-02,-1.00983e-03,-4.06916e-03,-6.60864e+01,-1.71533e-02,&
334  1.10605e-02, 1.20300e-02,-5.20034e-03, 0.00000e+00, 0.00000e+00,&
335  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
336  -2.62769e+03, 7.13755e-03, 4.17999e-03, 0.00000e+00, 1.25910e+04,&
337  0.00000e+00, 0.00000e+00, 0.00000e+00,-2.23595e-03, 4.60217e-03,&
338  5.71794e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
339  -3.18353e-02,-2.35526e-02,-1.36189e-02, 0.00000e+00, 0.00000e+00,&
340  0.00000e+00, 2.03522e-02,-6.67837e+01,-1.09724e-03, 0.00000e+00/
341  data pf3/ &
342  -1.38821e-02, 1.60468e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
343  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 1.51574e-02,&
344  -5.44470e-04, 0.00000e+00, 7.28224e-02, 6.59413e-02, 0.00000e+00,&
345  -5.15692e-03, 0.00000e+00, 0.00000e+00,-3.70367e+03, 0.00000e+00,&
346  0.00000e+00, 1.36131e-02, 5.38153e-03, 0.00000e+00, 4.76285e+00,&
347  -1.75677e-02, 2.26301e-02, 0.00000e+00, 1.76631e-02, 4.77162e-03,&
348  0.00000e+00, 5.39354e+00, 0.00000e+00,-7.51710e-03, 0.00000e+00,&
349  0.00000e+00,-8.82736e+01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
350  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
351  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
352 ! h density
353  data pg1/ &
354  1.26376e+00,-2.14304e-01,-1.49984e-01, 2.30404e-01, 2.98237e-02,&
355  2.68673e-02, 2.96228e-01, 2.21900e-02,-2.07655e-02, 4.52506e-01,&
356  1.20105e-01, 3.24420e-02, 4.24816e-02,-9.14313e+00, 0.00000e+00,&
357  2.47178e-02,-2.88229e-02, 8.12805e+01, 5.10380e-02,-5.80611e-03,&
358  2.51236e-05,-1.24083e-02, 0.00000e+00, 0.00000e+00, 8.66784e-02,&
359  1.58727e-01,-3.48190e-02, 0.00000e+00, 0.00000e+00, 2.89885e-05,&
360  0.00000e+00, 1.53595e+02,-1.68604e-02, 0.00000e+00, 1.01015e-02,&
361  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.84552e-04,&
362  -1.22181e-03, 0.00000e+00, 0.00000e+00, 8.47001e-02, 1.70147e-01,&
363  -1.04927e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00,-5.91313e-03/
364  data pg2/ &
365  -2.30501e-02, 3.14758e-05, 0.00000e+00, 0.00000e+00, 1.26956e-02,&
366  8.35489e-03, 3.10513e-04, 0.00000e+00, 3.42119e+03,-2.45017e-03,&
367  -4.27154e-04, 5.45152e-04, 1.89896e-03, 2.89121e+01,-6.49973e-03,&
368  -1.93855e-02,-1.48492e-02, 0.00000e+00,-5.10576e-02, 7.87306e-02,&
369  9.51981e-02,-1.49422e+04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
370  2.65503e+02, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
371  0.00000e+00, 0.00000e+00, 0.00000e+00, 6.37110e-03, 3.24789e-04,&
372  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
373  6.14274e-02, 1.00376e-02,-8.41083e-04, 0.00000e+00, 0.00000e+00,&
374  0.00000e+00,-1.27099e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00/
375  data pg3/ &
376  -3.94077e-03,-1.28601e-02,-7.97616e-03, 0.00000e+00, 0.00000e+00,&
377  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
378  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
379  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
380  0.00000e+00,-6.71465e-03,-1.69799e-03, 1.93772e-03, 3.81140e+00,&
381  -7.79290e-03,-1.82589e-02,-1.25860e-02,-1.04311e-02,-3.02465e-03,&
382  2.43063e-03, 3.63237e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
383  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
384  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
385  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
386 ! n density
387  data ph1/ &
388  7.09557e+01,-3.26740e-01, 0.00000e+00,-5.16829e-01,-1.71664e-03,&
389  9.09310e-02,-6.71500e-01,-1.47771e-01,-9.27471e-02,-2.30862e-01,&
390  -1.56410e-01, 1.34455e-02,-1.19717e-01, 2.52151e+00, 0.00000e+00,&
391  -2.41582e-01, 5.92939e-02, 4.39756e+00, 9.15280e-02, 4.41292e-03,&
392  0.00000e+00, 8.66807e-03, 0.00000e+00, 0.00000e+00, 8.66784e-02,&
393  1.58727e-01, 9.74701e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
394  0.00000e+00, 6.70217e+01,-1.31660e-03, 0.00000e+00,-1.65317e-02,&
395  0.00000e+00, 0.00000e+00, 8.50247e-02, 2.77428e+01, 4.98658e-03,&
396  6.15115e-03, 9.50156e-03,-2.12723e-02, 8.47001e-02, 1.70147e-01,&
397  -2.38645e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00, 1.37380e-03/
398  data ph2/ &
399  -8.41918e-03, 2.80145e-05, 7.12383e-03, 0.00000e+00,-1.66209e-02,&
400  1.03533e-04,-1.68898e-02, 0.00000e+00, 3.64526e+03, 0.00000e+00,&
401  6.54077e-03, 3.69130e-04, 9.94419e-04, 8.42803e+01,-1.16124e-02,&
402  -7.74414e-03,-1.68844e-03, 1.42809e-03,-1.92955e-03, 1.17225e-01,&
403  -2.41512e-02, 1.50521e+04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
404  1.60261e+03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
405  0.00000e+00, 0.00000e+00, 0.00000e+00,-3.54403e-04,-1.87270e-02,&
406  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
407  2.76439e-02, 6.43207e-03,-3.54300e-02, 0.00000e+00, 0.00000e+00,&
408  0.00000e+00,-2.80221e-02, 8.11228e+01,-6.75255e-04, 0.00000e+00/
409  data ph3/ &
410  -1.05162e-02,-3.48292e-03,-6.97321e-03, 0.00000e+00, 0.00000e+00,&
411  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
412  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
413  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
414  0.00000e+00,-1.45546e-03,-1.31970e-02,-3.57751e-03,-1.09021e+00,&
415  -1.50181e-02,-7.12841e-03,-6.64590e-03,-3.52610e-03,-1.87773e-02,&
416  -2.22432e-03,-3.93895e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
417  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
418  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
419  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
420 ! hot o density
421  data pi1/ &
422  6.04050e-02, 1.57034e+00, 2.99387e-02, 0.00000e+00, 0.00000e+00,&
423  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-1.51018e+00,&
424  0.00000e+00, 0.00000e+00, 0.00000e+00,-8.61650e+00, 1.26454e-02,&
425  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
426  0.00000e+00, 5.50878e-03, 0.00000e+00, 0.00000e+00, 8.66784e-02,&
427  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
428  0.00000e+00, 0.00000e+00, 6.23881e-02, 0.00000e+00, 0.00000e+00,&
429  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
430  0.00000e+00, 0.00000e+00, 0.00000e+00, 8.47001e-02, 1.70147e-01,&
431  -9.45934e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
432  data pi2/ &
433  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
434  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
435  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
436  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
437  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
438  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
439  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
440  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
441  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
442  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
443  data pi3/ &
444  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
445  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
446  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
447  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
448  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
449  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
450  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
451  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
452  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
453  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
454 ! s param
455  data pj1/ &
456  9.56827e-01, 6.20637e-02, 3.18433e-02, 0.00000e+00, 0.00000e+00,&
457  3.94900e-02, 0.00000e+00, 0.00000e+00,-9.24882e-03,-7.94023e-03,&
458  0.00000e+00, 0.00000e+00, 0.00000e+00, 1.74712e+02, 0.00000e+00,&
459  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
460  0.00000e+00, 2.74677e-03, 0.00000e+00, 1.54951e-02, 8.66784e-02,&
461  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
462  0.00000e+00, 0.00000e+00, 0.00000e+00,-6.99007e-04, 0.00000e+00,&
463  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
464  0.00000e+00, 1.24362e-02,-5.28756e-03, 8.47001e-02, 1.70147e-01,&
465  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
466  data pj2/ &
467  0.00000e+00, 2.47425e-05, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
468  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
469  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
470  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
471  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
472  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
473  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
474  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
475  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
476  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
477  data pj3/ &
478  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
479  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
480  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
481  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
482  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
483  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
484  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
485  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
486  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
487  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
488 ! turbo
489  data pk1/ &
490  1.09930e+00, 3.90631e+00, 3.07165e+00, 9.86161e-01, 1.63536e+01,&
491  4.63830e+00, 1.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
492  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
493  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
494  0.00000e+00, 0.00000e+00, 1.28840e+00, 3.10302e-02, 1.18339e-01,&
495  1.00000e+00, 7.00000e-01, 1.15020e+00, 3.44689e+00, 1.28840e+00,&
496  1.00000e+00, 1.08738e+00, 1.22947e+00, 1.10016e+00, 7.34129e-01,&
497  1.15241e+00, 2.22784e+00, 7.95046e-01, 4.01612e+00, 4.47749e+00,&
498  1.23435e+02,-7.60535e-02, 1.68986e-06, 7.44294e-01, 1.03604e+00,&
499  1.72783e+02, 1.15020e+00, 3.44689e+00,-7.46230e-01, 9.49154e-01/
500 ! lower boundary
501  data ptm/ &
502  1.04130e+03, 3.86000e+02, 1.95000e+02, 1.66728e+01, 2.13000e+02,&
503  1.20000e+02, 2.40000e+02, 1.87000e+02,-2.00000e+00, 0.00000e+00/
504  data pdm/ &
505  2.45600e+07, 6.71072e-06, 1.00000e+02, 0.00000e+00, 1.10000e+02,&
506  1.00000e+01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
507  8.59400e+10, 1.00000e+00, 1.05000e+02,-8.00000e+00, 1.10000e+02,&
508  1.00000e+01, 9.00000e+01, 2.00000e+00, 0.00000e+00, 0.00000e+00,&
509  2.81000e+11, 0.00000e+00, 1.05000e+02, 2.80000e+01, 2.89500e+01,&
510  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
511  3.30000e+10, 2.68270e-01, 1.05000e+02, 1.00000e+00, 1.10000e+02,&
512  1.00000e+01, 1.10000e+02,-1.00000e+01, 0.00000e+00, 0.00000e+00,&
513  1.33000e+09, 1.19615e-02, 1.05000e+02, 0.00000e+00, 1.10000e+02,&
514  1.00000e+01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
515  1.76100e+05, 1.00000e+00, 9.50000e+01,-8.00000e+00, 1.10000e+02,&
516  1.00000e+01, 9.00000e+01, 2.00000e+00, 0.00000e+00, 0.00000e+00,&
517  1.00000e+07, 1.00000e+00, 1.05000e+02,-8.00000e+00, 1.10000e+02,&
518  1.00000e+01, 9.00000e+01, 2.00000e+00, 0.00000e+00, 0.00000e+00,&
519  1.00000e+06, 1.00000e+00, 1.05000e+02,-8.00000e+00, 5.50000e+02,&
520  7.60000e+01, 9.00000e+01, 2.00000e+00, 0.00000e+00, 4.00000e+03/
521 ! tn1(2)
522  data pl1/ &
523  1.00858e+00, 4.56011e-02,-2.22972e-02,-5.44388e-02, 5.23136e-04,&
524  -1.88849e-02, 5.23707e-02,-9.43646e-03, 6.31707e-03,-7.80460e-02,&
525  -4.88430e-02, 0.00000e+00, 0.00000e+00,-7.60250e+00, 0.00000e+00,&
526  -1.44635e-02,-1.76843e-02,-1.21517e+02, 2.85647e-02, 0.00000e+00,&
527  0.00000e+00, 6.31792e-04, 0.00000e+00, 5.77197e-03, 8.66784e-02,&
528  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
529  0.00000e+00,-8.90272e+03, 3.30611e-03, 3.02172e-03, 0.00000e+00,&
530  -2.13673e-03,-3.20910e-04, 0.00000e+00, 0.00000e+00, 2.76034e-03,&
531  2.82487e-03,-2.97592e-04,-4.21534e-03, 8.47001e-02, 1.70147e-01,&
532  8.96456e-03, 0.00000e+00,-1.08596e-02, 0.00000e+00, 0.00000e+00/
533  data pl2/ &
534  5.57917e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
535  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
536  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
537  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
538  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
539  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
540  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
541  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
542  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
543  0.00000e+00, 9.65405e-03, 0.00000e+00, 0.00000e+00, 2.00000e+00/
544 ! tn1(3)
545  data pm1/ &
546  9.39664e-01, 8.56514e-02,-6.79989e-03, 2.65929e-02,-4.74283e-03,&
547  1.21855e-02,-2.14905e-02, 6.49651e-03,-2.05477e-02,-4.24952e-02,&
548  0.00000e+00, 0.00000e+00, 0.00000e+00, 1.19148e+01, 0.00000e+00,&
549  1.18777e-02,-7.28230e-02,-8.15965e+01, 1.73887e-02, 0.00000e+00,&
550  0.00000e+00, 0.00000e+00,-1.44691e-02, 2.80259e-04, 8.66784e-02,&
551  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
552  0.00000e+00, 2.16584e+02, 3.18713e-03, 7.37479e-03, 0.00000e+00,&
553  -2.55018e-03,-3.92806e-03, 0.00000e+00, 0.00000e+00,-2.89757e-03,&
554  -1.33549e-03, 1.02661e-03, 3.53775e-04, 8.47001e-02, 1.70147e-01,&
555  -9.17497e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
556  data pm2/ &
557  3.56082e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
558  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
559  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
560  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
561  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
562  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
563  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
564  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
565  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
566  0.00000e+00,-1.00902e-02, 0.00000e+00, 0.00000e+00, 2.00000e+00/
567 ! tn1(4)
568  data pn1/ &
569  9.85982e-01,-4.55435e-02, 1.21106e-02, 2.04127e-02,-2.40836e-03,&
570  1.11383e-02,-4.51926e-02, 1.35074e-02,-6.54139e-03, 1.15275e-01,&
571  1.28247e-01, 0.00000e+00, 0.00000e+00,-5.30705e+00, 0.00000e+00,&
572  -3.79332e-02,-6.24741e-02, 7.71062e-01, 2.96315e-02, 0.00000e+00,&
573  0.00000e+00, 0.00000e+00, 6.81051e-03,-4.34767e-03, 8.66784e-02,&
574  1.58727e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
575  0.00000e+00, 1.07003e+01,-2.76907e-03, 4.32474e-04, 0.00000e+00,&
576  1.31497e-03,-6.47517e-04, 0.00000e+00,-2.20621e+01,-1.10804e-03,&
577  -8.09338e-04, 4.18184e-04, 4.29650e-03, 8.47001e-02, 1.70147e-01,&
578  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
579  data pn2/ &
580  -4.04337e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
581  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
582  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-9.52550e-04,&
583  8.56253e-04, 4.33114e-04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
584  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 1.21223e-03,&
585  2.38694e-04, 9.15245e-04, 1.28385e-03, 8.67668e-04,-5.61425e-06,&
586  1.04445e+00, 3.41112e+01, 0.00000e+00,-8.40704e-01,-2.39639e+02,&
587  7.06668e-01,-2.05873e+01,-3.63696e-01, 2.39245e+01, 0.00000e+00,&
588  -1.06657e-03,-7.67292e-04, 1.54534e-04, 0.00000e+00, 0.00000e+00,&
589  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
590 ! tn1(5) tn2(1)
591  data po1/ &
592  1.00320e+00, 3.83501e-02,-2.38983e-03, 2.83950e-03, 4.20956e-03,&
593  5.86619e-04, 2.19054e-02,-1.00946e-02,-3.50259e-03, 4.17392e-02,&
594  -8.44404e-03, 0.00000e+00, 0.00000e+00, 4.96949e+00, 0.00000e+00,&
595  -7.06478e-03,-1.46494e-02, 3.13258e+01,-1.86493e-03, 0.00000e+00,&
596  -1.67499e-02, 0.00000e+00, 0.00000e+00, 5.12686e-04, 8.66784e-02,&
597  1.58727e-01,-4.64167e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
598  4.37353e-03,-1.99069e+02, 0.00000e+00,-5.34884e-03, 0.00000e+00,&
599  1.62458e-03, 2.93016e-03, 2.67926e-03, 5.90449e+02, 0.00000e+00,&
600  0.00000e+00,-1.17266e-03,-3.58890e-04, 8.47001e-02, 1.70147e-01,&
601  0.00000e+00, 0.00000e+00, 1.38673e-02, 0.00000e+00, 0.00000e+00/
602  data po2/ &
603  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
604  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
605  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 1.60571e-03,&
606  6.28078e-04, 5.05469e-05, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
607  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-1.57829e-03,&
608  -4.00855e-04, 5.04077e-05,-1.39001e-03,-2.33406e-03,-4.81197e-04,&
609  1.46758e+00, 6.20332e+00, 0.00000e+00, 3.66476e-01,-6.19760e+01,&
610  3.09198e-01,-1.98999e+01, 0.00000e+00,-3.29933e+02, 0.00000e+00,&
611  -1.10080e-03,-9.39310e-05, 1.39638e-04, 0.00000e+00, 0.00000e+00,&
612  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
613 ! tn2(2)
614  data pp1/ &
615  9.81637e-01,-1.41317e-03, 3.87323e-02, 0.00000e+00, 0.00000e+00,&
616  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-3.58707e-02,&
617  -8.63658e-03, 0.00000e+00, 0.00000e+00,-2.02226e+00, 0.00000e+00,&
618  -8.69424e-03,-1.91397e-02, 8.76779e+01, 4.52188e-03, 0.00000e+00,&
619  2.23760e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
620  0.00000e+00,-7.07572e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
621  -4.11210e-03, 3.50060e+01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
622  0.00000e+00, 0.00000e+00,-8.36657e-03, 1.61347e+01, 0.00000e+00,&
623  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
624  0.00000e+00, 0.00000e+00,-1.45130e-02, 0.00000e+00, 0.00000e+00/
625  data pp2/ &
626  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
627  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
628  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 1.24152e-03,&
629  6.43365e-04, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
630  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 1.33255e-03,&
631  2.42657e-03, 1.60666e-03,-1.85728e-03,-1.46874e-03,-4.79163e-06,&
632  1.22464e+00, 3.53510e+01, 0.00000e+00, 4.49223e-01,-4.77466e+01,&
633  4.70681e-01, 8.41861e+00,-2.88198e-01, 1.67854e+02, 0.00000e+00,&
634  7.11493e-04, 6.05601e-04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
635  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
636 ! tn2(3)
637  data pq1/ &
638  1.00422e+00,-7.11212e-03, 5.24480e-03, 0.00000e+00, 0.00000e+00,&
639  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-5.28914e-02,&
640  -2.41301e-02, 0.00000e+00, 0.00000e+00,-2.12219e+01,-1.03830e-02,&
641  -3.28077e-03, 1.65727e-02, 1.68564e+00,-6.68154e-03, 0.00000e+00,&
642  1.45155e-02, 0.00000e+00, 8.42365e-03, 0.00000e+00, 0.00000e+00,&
643  0.00000e+00,-4.34645e-03, 0.00000e+00, 0.00000e+00, 2.16780e-02,&
644  0.00000e+00,-1.38459e+02, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
645  0.00000e+00, 0.00000e+00, 7.04573e-03,-4.73204e+01, 0.00000e+00,&
646  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
647  0.00000e+00, 0.00000e+00, 1.08767e-02, 0.00000e+00, 0.00000e+00/
648  data pq2/ &
649  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
650  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-8.08279e-03,&
651  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 5.21769e-04,&
652  -2.27387e-04, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
653  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 3.26769e-03,&
654  3.16901e-03, 4.60316e-04,-1.01431e-04, 1.02131e-03, 9.96601e-04,&
655  1.25707e+00, 2.50114e+01, 0.00000e+00, 4.24472e-01,-2.77655e+01,&
656  3.44625e-01, 2.75412e+01, 0.00000e+00, 7.94251e+02, 0.00000e+00,&
657  2.45835e-03, 1.38871e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
658  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
659 ! tn2(4) tn3(1)
660  data pr1/ &
661  1.01890e+00,-2.46603e-02, 1.00078e-02, 0.00000e+00, 0.00000e+00,&
662  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-6.70977e-02,&
663  -4.02286e-02, 0.00000e+00, 0.00000e+00,-2.29466e+01,-7.47019e-03,&
664  2.26580e-03, 2.63931e-02, 3.72625e+01,-6.39041e-03, 0.00000e+00,&
665  9.58383e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
666  0.00000e+00,-1.85291e-03, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
667  0.00000e+00, 1.39717e+02, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
668  0.00000e+00, 0.00000e+00, 9.19771e-03,-3.69121e+02, 0.00000e+00,&
669  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
670  0.00000e+00, 0.00000e+00,-1.57067e-02, 0.00000e+00, 0.00000e+00/
671  data pr2/ &
672  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
673  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-7.07265e-03,&
674  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-2.92953e-03,&
675  -2.77739e-03,-4.40092e-04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
676  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.47280e-03,&
677  2.95035e-04,-1.81246e-03, 2.81945e-03, 4.27296e-03, 9.78863e-04,&
678  1.40545e+00,-6.19173e+00, 0.00000e+00, 0.00000e+00,-7.93632e+01,&
679  4.44643e-01,-4.03085e+02, 0.00000e+00, 1.15603e+01, 0.00000e+00,&
680  2.25068e-03, 8.48557e-04,-2.98493e-04, 0.00000e+00, 0.00000e+00,&
681  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
682 ! tn3(2)
683  data ps1/ &
684  9.75801e-01, 3.80680e-02,-3.05198e-02, 0.00000e+00, 0.00000e+00,&
685  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 3.85575e-02,&
686  5.04057e-02, 0.00000e+00, 0.00000e+00,-1.76046e+02, 1.44594e-02,&
687  -1.48297e-03,-3.68560e-03, 3.02185e+01,-3.23338e-03, 0.00000e+00,&
688  1.53569e-02, 0.00000e+00,-1.15558e-02, 0.00000e+00, 0.00000e+00,&
689  0.00000e+00, 4.89620e-03, 0.00000e+00, 0.00000e+00,-1.00616e-02,&
690  -8.21324e-03,-1.57757e+02, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
691  0.00000e+00, 0.00000e+00, 6.63564e-03, 4.58410e+01, 0.00000e+00,&
692  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
693  0.00000e+00, 0.00000e+00,-2.51280e-02, 0.00000e+00, 0.00000e+00/
694  data ps2/ &
695  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
696  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 9.91215e-03,&
697  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-8.73148e-04,&
698  -1.29648e-03,-7.32026e-05, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
699  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-4.68110e-03,&
700  -4.66003e-03,-1.31567e-03,-7.39390e-04, 6.32499e-04,-4.65588e-04,&
701  -1.29785e+00,-1.57139e+02, 0.00000e+00, 2.58350e-01,-3.69453e+01,&
702  4.10672e-01, 9.78196e+00,-1.52064e-01,-3.85084e+03, 0.00000e+00,&
703  -8.52706e-04,-1.40945e-03,-7.26786e-04, 0.00000e+00, 0.00000e+00,&
704  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
705 ! tn3(3)
706  data pu1/ &
707  9.60722e-01, 7.03757e-02,-3.00266e-02, 0.00000e+00, 0.00000e+00,&
708  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.22671e-02,&
709  4.10423e-02, 0.00000e+00, 0.00000e+00,-1.63070e+02, 1.06073e-02,&
710  5.40747e-04, 7.79481e-03, 1.44908e+02, 1.51484e-04, 0.00000e+00,&
711  1.97547e-02, 0.00000e+00,-1.41844e-02, 0.00000e+00, 0.00000e+00,&
712  0.00000e+00, 5.77884e-03, 0.00000e+00, 0.00000e+00, 9.74319e-03,&
713  0.00000e+00,-2.88015e+03, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
714  0.00000e+00, 0.00000e+00,-4.44902e-03,-2.92760e+01, 0.00000e+00,&
715  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
716  0.00000e+00, 0.00000e+00, 2.34419e-02, 0.00000e+00, 0.00000e+00/
717  data pu2/ &
718  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
719  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 5.36685e-03,&
720  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-4.65325e-04,&
721  -5.50628e-04, 3.31465e-04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
722  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-2.06179e-03,&
723  -3.08575e-03,-7.93589e-04,-1.08629e-04, 5.95511e-04,-9.05050e-04,&
724  1.18997e+00, 4.15924e+01, 0.00000e+00,-4.72064e-01,-9.47150e+02,&
725  3.98723e-01, 1.98304e+01, 0.00000e+00, 3.73219e+03, 0.00000e+00,&
726  -1.50040e-03,-1.14933e-03,-1.56769e-04, 0.00000e+00, 0.00000e+00,&
727  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
728 ! tn3(4)
729  data pv1/ &
730  1.03123e+00,-7.05124e-02, 8.71615e-03, 0.00000e+00, 0.00000e+00,&
731  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-3.82621e-02,&
732  -9.80975e-03, 0.00000e+00, 0.00000e+00, 2.89286e+01, 9.57341e-03,&
733  0.00000e+00, 0.00000e+00, 8.66153e+01, 7.91938e-04, 0.00000e+00,&
734  0.00000e+00, 0.00000e+00, 4.68917e-03, 0.00000e+00, 0.00000e+00,&
735  0.00000e+00, 7.86638e-03, 0.00000e+00, 0.00000e+00, 9.90827e-03,&
736  0.00000e+00, 6.55573e+01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
737  0.00000e+00, 0.00000e+00, 0.00000e+00,-4.00200e+01, 0.00000e+00,&
738  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
739  0.00000e+00, 0.00000e+00, 7.07457e-03, 0.00000e+00, 0.00000e+00/
740  data pv2/ &
741  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
742  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 5.72268e-03,&
743  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-2.04970e-04,&
744  1.21560e-03,-8.05579e-06, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
745  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-2.49941e-03,&
746  -4.57256e-04,-1.59311e-04, 2.96481e-04,-1.77318e-03,-6.37918e-04,&
747  1.02395e+00, 1.28172e+01, 0.00000e+00, 1.49903e-01,-2.63818e+01,&
748  0.00000e+00, 4.70628e+01,-2.22139e-01, 4.82292e-02, 0.00000e+00,&
749  -8.67075e-04,-5.86479e-04, 5.32462e-04, 0.00000e+00, 0.00000e+00,&
750  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
751 ! tn3(5) surface temp tsl
752  data pw1/ &
753  1.00828e+00,-9.10404e-02,-2.26549e-02, 0.00000e+00, 0.00000e+00,&
754  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-2.32420e-02,&
755  -9.08925e-03, 0.00000e+00, 0.00000e+00, 3.36105e+01, 0.00000e+00,&
756  0.00000e+00, 0.00000e+00,-1.24957e+01,-5.87939e-03, 0.00000e+00,&
757  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
758  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
759  0.00000e+00, 2.79765e+01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
760  0.00000e+00, 0.00000e+00, 0.00000e+00, 2.01237e+03, 0.00000e+00,&
761  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
762  0.00000e+00, 0.00000e+00,-1.75553e-02, 0.00000e+00, 0.00000e+00/
763  data pw2/ &
764  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
765  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
766  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 3.29699e-03,&
767  1.26659e-03, 2.68402e-04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
768  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 1.17894e-03,&
769  1.48746e-03, 1.06478e-04, 1.34743e-04,-2.20939e-03,-6.23523e-04,&
770  6.36539e-01, 1.13621e+01, 0.00000e+00,-3.93777e-01, 2.38687e+03,&
771  0.00000e+00, 6.61865e+02,-1.21434e-01, 9.27608e+00, 0.00000e+00,&
772  1.68478e-04, 1.24892e-03, 1.71345e-03, 0.00000e+00, 0.00000e+00,&
773  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
774 ! tgn3(2) surface grad tslg
775  data px1/ &
776  1.57293e+00,-6.78400e-01, 6.47500e-01, 0.00000e+00, 0.00000e+00,&
777  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-7.62974e-02,&
778  -3.60423e-01, 0.00000e+00, 0.00000e+00, 1.28358e+02, 0.00000e+00,&
779  0.00000e+00, 0.00000e+00, 4.68038e+01, 0.00000e+00, 0.00000e+00,&
780  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
781  0.00000e+00,-1.67898e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
782  0.00000e+00, 2.90994e+04, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
783  0.00000e+00, 0.00000e+00, 0.00000e+00, 3.15706e+01, 0.00000e+00,&
784  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
785  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
786  data px2/ &
787  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
788  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
789  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
790  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
791  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
792  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
793  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
794  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
795  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
796  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
797 ! tgn2(1) tgn1(2)
798  data py1/ &
799  8.60028e-01, 3.77052e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
800  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-1.17570e+00,&
801  0.00000e+00, 0.00000e+00, 0.00000e+00, 7.77757e-03, 0.00000e+00,&
802  0.00000e+00, 0.00000e+00, 1.01024e+02, 0.00000e+00, 0.00000e+00,&
803  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
804  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
805  0.00000e+00, 6.54251e+02, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
806  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
807  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
808  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
809  data py2/ &
810  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
811  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
812  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
813  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
814  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,-1.56959e-02,&
815  1.91001e-02, 3.15971e-02, 1.00982e-02,-6.71565e-03, 2.57693e-03,&
816  1.38692e+00, 2.82132e-01, 0.00000e+00, 0.00000e+00, 3.81511e+02,&
817  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
818  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
819  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
820 ! tgn3(1) tgn2(2)
821  data pz1/ &
822  1.06029e+00,-5.25231e-02, 3.73034e-01, 0.00000e+00, 0.00000e+00,&
823  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 3.31072e-02,&
824  -3.88409e-01, 0.00000e+00, 0.00000e+00,-1.65295e+02,-2.13801e-01,&
825  -4.38916e-02,-3.22716e-01,-8.82393e+01, 1.18458e-01, 0.00000e+00,&
826  -4.35863e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
827  0.00000e+00,-1.19782e-01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
828  0.00000e+00, 2.62229e+01, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
829  0.00000e+00, 0.00000e+00, 0.00000e+00,-5.37443e+01, 0.00000e+00,&
830  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
831  0.00000e+00, 0.00000e+00,-4.55788e-01, 0.00000e+00, 0.00000e+00/
832  data pz2/ &
833  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
834  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
835  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 3.84009e-02,&
836  3.96733e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
837  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 5.05494e-02,&
838  7.39617e-02, 1.92200e-02,-8.46151e-03,-1.34244e-02, 1.96338e-02,&
839  1.50421e+00, 1.88368e+01, 0.00000e+00, 0.00000e+00,-5.13114e+01,&
840  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
841  5.11923e-02, 3.61225e-02, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
842  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 2.00000e+00/
843 ! semiannual mult sam
844  data paa1/ &
845  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00,&
846  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00,&
847  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00,&
848  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00,&
849  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00,&
850  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00,&
851  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00,&
852  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00,&
853  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00,&
854  1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00, 1.00000e+00/
855  data paa2/ &
856  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
857  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
858  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
859  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
860  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
861  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
862  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
863  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
864  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00,&
865  0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00, 0.00000e+00/
866 ! middle atmosphere averages
867  data pavgm/ &
868  2.61000e+02, 2.64000e+02, 2.29000e+02, 2.17000e+02, 2.17000e+02,&
869  2.23000e+02, 2.86760e+02,-2.93940e+00, 2.50000e+00, 0.00000e+00/
870  end module wam_gtd7bk_mod
871 
872 !-----------------------------------------------------------------------
879  module gettemp_mod
880 !
881  real :: tlb
882  real :: s
883  real :: db04
884  real :: db16
885  real :: db28
886  real :: db32
887  real :: db40
888  real :: db48
889  real :: db01
890  real :: za
891  real :: t0
892  real :: z0
893  real :: g0
894  real :: rl
895  real :: dd
896  real :: db14
897  real :: tr12
898 !
899  real :: tn1(5)
900  real :: tn2(4)
901  real :: tn3(5)
902  real :: tgn1(2)
903  real :: tgn2(2)
904  real :: tgn3(2)
905 !
906  real :: pt(150)
907  real :: pd(150,9)
908  real :: ps(150)
909  real :: pdl(25,2)
910  real :: ptl(100,4)
911  real :: pma(100,10)
912  real :: sam(100)
913 !
914  real :: sw(25)
915  real :: swc(25)
916 !
917  real :: dm04
918  real :: dm16
919  real :: dm28
920  real :: dm32
921  real :: dm40
922  real :: dm01
923  real :: dm14
924 !
925  real :: gsurf
926  real :: re
927 !
928  real :: tinfg
929  real :: tt(15)
930 !
931  real :: plg(9,4)
932  real :: ctloc
933  real :: stloc
934  real :: c2tloc
935  real :: s2tloc
936  real :: c3tloc
937  real :: s3tloc
938  real :: day
939  real :: df
940  real :: dfa
941  real :: apd
942  real :: apdf
943  real :: apt(4)
944  real :: xlong
945 !
946  integer :: isw
947  integer :: iyr
948 !
949  end module gettemp_mod
950 
951 ! ----------------------------------------------------------------------
967  subroutine gettemp(iday,nday,xlat,nlat,pr,np,temp,n_o,n_o2,n_n2)
968  implicit none
969  integer, intent(in) :: nday ! number of days
970  integer, intent(in) :: nlat ! number of latitudes
971  integer, intent(in) :: np ! number of pressure layer
972  real, intent(in) :: pr(np) ! pressure in mb
973  real, intent(in) :: xlat(nlat) ! latitude in degree
974  integer, intent(in) :: iday(nday) ! calender day
975  real, intent(out) :: temp(np,nlat,nday) ! temperature
976  real, intent(out) :: n_o(np,nlat,nday) ! number density of o
977  real, intent(out) :: n_o2(np,nlat,nday) ! number density of o2
978  real, intent(out) :: n_n2(np,nlat,nday) ! number density of n2
979  real :: alt(np,nlat,nday) ! altitude in km
980  real :: d(9),t(2),sw(25),ap(7),ut,xlong,xlst,f107, &
981  f107a
982  integer :: k,il,ip
983 ! set magnetic index average value
984  data ap/7*9./
985 ! set swich 7,8,10,14 zero to avoid diurnal changes in output temperatu
986 ! swich #7 is for diurnal,#8 is for semidiurnal,# 10 is for all ut/longi
987 ! effect,#14 is for terdiurnal
988  data sw/1.,1.,1.,1.,1.,1.,0.,0.,1.,0.,1.,1.,1.,0.,1.,1.,1.,1.,1., &
989  1.,1.,1.,1.,1.,1./
990 ! set 10.7cm flux be average value
991  f107=150.
992  f107a=150.
993 ! turn on swich
994  call tselec(sw)
995 ! set longitude, ut, local time , it should not make difference to outpu
996  ut=0.
997  xlong=0.
998  xlst=ut/3600.+xlong/15.
999 ! calculate temperature for each lat,pres level,day
1000  do k=1,nday
1001  do il=1,nlat
1002  do ip=1,np
1003  call ghp7(iday(k),ut,alt(ip,il,k),xlat(il),xlong,xlst,f107a,f107, &
1004  ap,d,t,pr(ip))
1005  temp(ip,il,k)=t(2)
1006  n_o(ip,il,k)=d(2)
1007  n_o2(ip,il,k)=d(4)
1008  n_n2(ip,il,k)=d(3)
1009  enddo
1010  enddo
1011  enddo
1012  end subroutine gettemp
1013 
1014 !-----------------------------------------------------------------------
1143  subroutine gtd7(iyd,sec,alt,glat,glong,stl,f107a,f107,ap,mass,d,t)
1145  use gettemp_mod, only: dd, tn1,tn2,tn3,tgn1,tgn2,tgn3, &
1146  pt,pd,ps,pdl,ptl,pma,sam,sw,isw,&
1147  dm28,gsurf,re
1148 
1149  dimension d(9),t(2),ap(7),ds(9),ts(2)
1150  dimension zn3(5),zn2(4),sv(25)
1151 !
1152  save
1153  data mn3/5/,zn3/32.5,20.,15.,10.,0./
1154  data mn2/4/,zn2/72.5,55.,45.,32.5/
1155  data zmix/62.5/,alast/99999./,mssl/-999/
1156  data sv/25*1./
1157 ! ==== assign common/parm7/
1158  pt(1:50) =pt1(1:50); pt(51:100) =pt2(1:50); pt(101:150) =pt3(1:50)
1159  pd(1:50,1)=pa1(1:50); pd(51:100,1)=pa2(1:50); pd(101:150,1)=pa3(1:50)
1160  pd(1:50,2)=pb1(1:50); pd(51:100,2)=pb2(1:50); pd(101:150,2)=pb3(1:50)
1161  pd(1:50,3)=pc1(1:50); pd(51:100,3)=pc2(1:50); pd(101:150,3)=pc3(1:50)
1162  pd(1:50,4)=pd1(1:50); pd(51:100,4)=pd2(1:50); pd(101:150,4)=pd3(1:50)
1163  pd(1:50,5)=pe1(1:50); pd(51:100,5)=pe2(1:50); pd(101:150,5)=pe3(1:50)
1164  pd(1:50,6)=pf1(1:50); pd(51:100,6)=pf2(1:50); pd(101:150,6)=pf3(1:50)
1165  pd(1:50,7)=pg1(1:50); pd(51:100,7)=pg2(1:50); pd(101:150,7)=pg3(1:50)
1166  pd(1:50,8)=ph1(1:50); pd(51:100,8)=ph2(1:50); pd(101:150,8)=ph3(1:50)
1167  pd(1:50,9)=pi1(1:50); pd(51:100,9)=pi2(1:50); pd(101:150,9)=pi3(1:50)
1168  ps(1:50) =pj1(1:50); ps(51:100) =pj2(1:50); ps(101:150) =pj3(1:50)
1169  pdl(1:25,1)=pk1(1:25); pdl(1:25,2)=pk1(26:50)
1170  ptl(1:50,1)=pl1(1:50); ptl(51:100,1)=pl2(1:50)
1171  ptl(1:50,2)=pm1(1:50); ptl(51:100,2)=pm2(1:50)
1172  ptl(1:50,3)=pn1(1:50); ptl(51:100,3)=pn2(1:50)
1173  ptl(1:50,4)=po1(1:50); ptl(51:100,4)=po2(1:50)
1174  pma(1:50,1)=pp1(1:50); pma(51:100,1)=pp2(1:50)
1175  pma(1:50,2)=pq1(1:50); pma(51:100,2)=pq2(1:50)
1176  pma(1:50,3)=pr1(1:50); pma(51:100,3)=pr2(1:50)
1177  pma(1:50,4)=ps1(1:50); pma(51:100,4)=ps2(1:50)
1178  pma(1:50,5)=pu1(1:50); pma(51:100,5)=pu2(1:50)
1179  pma(1:50,6)=pv1(1:50); pma(51:100,6)=pv2(1:50)
1180  pma(1:50,7)=pw1(1:50); pma(51:100,7)=pw2(1:50)
1181  pma(1:50,8)=px1(1:50); pma(51:100,8)=px2(1:50)
1182  pma(1:50,9)=py1(1:50); pma(51:100,9)=py2(1:50)
1183  pma(1:50,10)=pz1(1:50); pma(51:100,10)=pz2(1:50)
1184  sam(1:50)=paa1(1:50); sam(51:100)=paa2(1:50)
1185 !
1186  if(isw.ne.64999) call tselec(sv)
1187 !
1188 ! test for changed input
1189  v1=vtst7(iyd,sec,glat,glong,stl,f107a,f107,ap,1)
1190 ! latitude variation of gravity (none for sw(2)=0)
1191  xlat=glat
1192  if(sw(2).eq.0) xlat=45.
1193  call glatf(xlat,gsurf,re)
1194 !
1195  xmm=pdm(5,3)
1196 !
1197 ! thermosphere/mesosphere (above zn2(1))
1198  altt=amax1(alt,zn2(1))
1199  mss=mass
1200 ! only calculate n2 in thermosphere if alt in mixed region
1201  if(alt.lt.zmix.and.mass.gt.0) mss=28
1202 ! only calculate thermosphere if input parameters changed
1203 ! or altitude above zn2(1) in mesosphere
1204  if(v1.eq.1..or.alt.gt.zn2(1).or.alast.gt.zn2(1).or.mss.ne.mssl) then
1205  call gts7(iyd,sec,altt,glat,glong,stl,f107a,f107,ap,mss,ds,ts)
1206  dm28m=dm28
1207 ! metric adjustment
1208  if(imr.eq.1) dm28m=dm28*1.e6
1209  mssl=mss
1210  endif
1211  t(1)=ts(1)
1212  t(2)=ts(2)
1213  if(alt.ge.zn2(1)) then
1214  do 5 j=1,9
1215  d(j)=ds(j)
1216  5 continue
1217  goto 10
1218  endif
1219 !
1220 ! lower mesosphere/upper stratosphere [between zn3(1) and zn2(1)]
1221 ! temperature at nodes and gradients at end nodes
1222 ! inverse temperature a linear function of spherical harmonics
1223 ! only calculate nodes if input changed
1224  if(v1.eq.1..or.alast.ge.zn2(1)) then
1225  tgn2(1)=tgn1(2)
1226  tn2(1)=tn1(5)
1227  tn2(2)=pma(1,1)*pavgm(1)/(1.-sw(20)*glob7s(pma(1,1)))
1228  tn2(3)=pma(1,2)*pavgm(2)/(1.-sw(20)*glob7s(pma(1,2)))
1229  tn2(4)=pma(1,3)*pavgm(3)/(1.-sw(20)*sw(22)*glob7s(pma(1,3)))
1230  tgn2(2)=pavgm(9)*pma(1,10)*(1.+sw(20)*sw(22)*glob7s(pma(1,10))) &
1231  *tn2(4)*tn2(4)/(pma(1,3)*pavgm(3))**2
1232  tn3(1)=tn2(4)
1233  endif
1234  if(alt.ge.zn3(1)) goto 6
1235 !
1236 ! lower stratosphere and troposphere [below zn3(1)]
1237 ! temperature at nodes and gradients at end nodes
1238 ! inverse temperature a linear function of spherical harmonics
1239 ! only calculate nodes if input changed
1240  if(v1.eq.1..or.alast.ge.zn3(1)) then
1241  tgn3(1)=tgn2(2)
1242  tn3(2)=pma(1,4)*pavgm(4)/(1.-sw(22)*glob7s(pma(1,4)))
1243  tn3(3)=pma(1,5)*pavgm(5)/(1.-sw(22)*glob7s(pma(1,5)))
1244  tn3(4)=pma(1,6)*pavgm(6)/(1.-sw(22)*glob7s(pma(1,6)))
1245  tn3(5)=pma(1,7)*pavgm(7)/(1.-sw(22)*glob7s(pma(1,7)))
1246  tgn3(2)=pma(1,8)*pavgm(8)*(1.+sw(22)*glob7s(pma(1,8))) &
1247  *tn3(5)*tn3(5)/(pma(1,7)*pavgm(7))**2
1248  endif
1249  6 continue
1250  if(mass.eq.0) goto 50
1251 ! linear transition to full mixing below zn2(1)
1252  dmc=0
1253  if(alt.gt.zmix) dmc=1.-(zn2(1)-alt)/(zn2(1)-zmix)
1254  dz28=ds(3)
1255 ! ***** n2 density ****
1256  dmr=ds(3)/dm28m-1.
1257  d(3)=densm(alt,dm28m,xmm,tz,mn3,zn3,tn3,tgn3,mn2,zn2,tn2,tgn2)
1258  d(3)=d(3)*(1.+dmr*dmc)
1259 ! ***** he density ****
1260  d(1)=0
1261  if(mass.ne.4.and.mass.ne.48) goto 204
1262  dmr=ds(1)/(dz28*pdm(2,1))-1.
1263  d(1)=d(3)*pdm(2,1)*(1.+dmr*dmc)
1264  204 continue
1265 ! **** o density ****
1266  d(2)=0
1267  d(9)=0
1268 ! ***** o2 density ****
1269  d(4)=0
1270  if(mass.ne.32.and.mass.ne.48) goto 232
1271  dmr=ds(4)/(dz28*pdm(2,4))-1.
1272  d(4)=d(3)*pdm(2,4)*(1.+dmr*dmc)
1273  232 continue
1274 ! ***** ar density ****
1275  d(5)=0
1276  if(mass.ne.40.and.mass.ne.48) goto 240
1277  dmr=ds(5)/(dz28*pdm(2,5))-1.
1278  d(5)=d(3)*pdm(2,5)*(1.+dmr*dmc)
1279  240 continue
1280 ! ***** hydrogen density ****
1281  d(7)=0
1282 ! ***** atomic nitrogen density ****
1283  d(8)=0
1284 !
1285 ! total mass density
1286 !
1287  if(mass.eq.48) then
1288  d(6) = 1.66e-24*(4.*d(1)+16.*d(2)+28.*d(3)+32.*d(4)+40.*d(5)+ &
1289  d(7)+14.*d(8))
1290  if(imr.eq.1) d(6)=d(6)/1000.
1291  endif
1292  t(2)=tz
1293  10 continue
1294  goto 90
1295  50 continue
1296  dd=densm(alt,1.,0.,tz,mn3,zn3,tn3,tgn3,mn2,zn2,tn2,tgn2)
1297  t(2)=tz
1298  90 continue
1299  alast=alt
1300  return
1301  end subroutine gtd7
1302 
1303 !-----------------------------------------------------------------------
1356  subroutine gtd7d(iyd,sec,alt,glat,glong,stl,f107a,f107,ap,mass,d,t)
1357  use wam_gtd7bk_mod, only: imr
1358 !
1359  dimension d(9),t(2),ap(7)
1360  call gtd7(iyd,sec,alt,glat,glong,stl,f107a,f107,ap,mass,d,t)
1361 ! total mass density
1362 !
1363  if(mass.eq.48) then
1364  d(6) = 1.66e-24*(4.*d(1)+16.*d(2)+28.*d(3)+32.*d(4)+40.*d(5)+ &
1365  d(7)+14.*d(8)+16.*d(9))
1366  if(imr.eq.1) d(6)=d(6)/1000.
1367  endif
1368  return
1369  end subroutine gtd7d
1370 
1371 !-----------------------------------------------------------------------
1411  subroutine ghp7(iyd,sec,alt,glat,glong,stl,f107a,f107,ap,d,t,press)
1412  use gettemp_mod, only: gsurf,re
1413  use wam_gtd7bk_mod, only: imr
1414 !
1415  dimension d(9),t(2),ap(7)
1416  save
1417  data bm/1.3806e-19/,rgas/831.4/
1418  data test/.00043/,ltest/12/
1419 !! bm=1.3806e-19; rgas=831.4
1420 !! test=.00043; ltest=12
1421  pl=alog10(press)
1422 ! initial altitude estimate
1423  if(pl.ge.-5.) then
1424  if(pl.gt.2.5) zi=18.06*(3.00-pl)
1425  if(pl.gt..75.and.pl.le.2.5) zi=14.98*(3.08-pl)
1426  if(pl.gt.-1..and.pl.le..75) zi=17.8*(2.72-pl)
1427  if(pl.gt.-2..and.pl.le.-1.) zi=14.28*(3.64-pl)
1428  if(pl.gt.-4..and.pl.le.-2.) zi=12.72*(4.32-pl)
1429  if(pl.le.-4.) zi=25.3*(.11-pl)
1430  iday=mod(iyd,1000)
1431  cl=glat/90.
1432  cl2=cl*cl
1433  if(iday.lt.182) cd=1.-iday/91.25
1434  if(iday.ge.182) cd=iday/91.25-3.
1435  ca=0
1436  if(pl.gt.-1.11.and.pl.le.-.23) ca=1.0
1437  if(pl.gt.-.23) ca=(2.79-pl)/(2.79+.23)
1438  if(pl.le.-1.11.and.pl.gt.-3.) ca=(-2.93-pl)/(-2.93+1.11)
1439  z=zi-4.87*cl*cd*ca-1.64*cl2*ca+.31*ca*cl
1440  endif
1441  if(pl.lt.-5.) z=22.*(pl+4.)**2+110
1442 ! iteration loop
1443  l=0
1444  10 continue
1445  l=l+1
1446  call gtd7(iyd,sec,z,glat,glong,stl,f107a,f107,ap,48,d,t)
1447  xn=d(1)+d(2)+d(3)+d(4)+d(5)+d(7)+d(8)
1448  p=bm*xn*t(2)
1449  if(imr.eq.1) p=p*1.e-6
1450  diff=pl-alog10(p)
1451  if(abs(diff).lt.test .or. l.eq.ltest) goto 20
1452  xm=d(6)/xn/1.66e-24
1453  if(imr.eq.1) xm = xm*1.e3
1454  g=gsurf/(1.+z/re)**2
1455  sh=rgas*t(2)/(xm*g)
1456 ! new altitude estimate using scale height
1457  if(l.lt.6) then
1458  z=z-sh*diff*2.302
1459  else
1460  z=z-sh*diff
1461  endif
1462  goto 10
1463  20 continue
1464  if(l.eq.ltest) write(6,100) press,diff
1465  100 format(1x,29hghp7 not converging for press, 1pe12.2,e12.2)
1466  alt=z
1467  return
1468  end subroutine ghp7
1469 
1470 !-----------------------------------------------------------------------
1478  subroutine glatf(lat,gv,reff)
1479  real lat
1480  save
1481  data dgtr/1.74533e-2/
1482 !! dgtr=1.74533e-2
1483  c2 = cos(2.*dgtr*lat)
1484  gv = 980.616*(1.-.0026373*c2)
1485  reff = 2.*gv/(3.085462e-6 + 2.27e-9*c2)*1.e-5
1486  return
1487  end subroutine glatf
1488 
1489 !-----------------------------------------------------------------------
1517  function vtst7(iyd,sec,glat,glong,stl,f107a,f107,ap,ic)
1518  use gettemp_mod, only: sw,swc
1519 !
1520  dimension ap(7),iydl(2),secl(2),glatl(2),gll(2),stll(2)
1521  dimension fal(2),fl(2),apl(7,2),swl(25,2),swcl(25,2)
1522  save
1523  data iydl/2*-999/,secl/2*-999./,glatl/2*-999./,gll/2*-999./
1524  data stll/2*-999./,fal/2*-999./,fl/2*-999./,apl/14*-999./
1525  data swl/50*-999./,swcl/50*-999./
1526  vtst7=0
1527  if(iyd.ne.iydl(ic)) goto 10
1528  if(sec.ne.secl(ic)) goto 10
1529  if(glat.ne.glatl(ic)) goto 10
1530  if(glong.ne.gll(ic)) goto 10
1531  if(stl.ne.stll(ic)) goto 10
1532  if(f107a.ne.fal(ic)) goto 10
1533  if(f107.ne.fl(ic)) goto 10
1534  do 5 i=1,7
1535  if(ap(i).ne.apl(i,ic)) goto 10
1536  5 end do
1537  do 7 i=1,25
1538  if(sw(i).ne.swl(i,ic)) goto 10
1539  if(swc(i).ne.swcl(i,ic)) goto 10
1540  7 end do
1541  goto 20
1542  10 continue
1543  vtst7=1
1544  iydl(ic)=iyd
1545  secl(ic)=sec
1546  glatl(ic)=glat
1547  gll(ic)=glong
1548  stll(ic)=stl
1549  fal(ic)=f107a
1550  fl(ic)=f107
1551  do 15 i=1,7
1552  apl(i,ic)=ap(i)
1553  15 end do
1554  do 16 i=1,25
1555  swl(i,ic)=sw(i)
1556  swcl(i,ic)=swc(i)
1557  16 end do
1558  20 continue
1559  return
1560  end function vtst7
1561 
1562 !-----------------------------------------------------------------------
1605  subroutine gts7(iyd,sec,alt,glat,glong,stl,f107a,f107,ap,mass,d,t)
1606  use wam_gtd7bk_mod, only: &
1607  ptm,pdm, &
1608  imr
1609 
1610  use gettemp_mod, only: tlb,s,db04,db16,db28,db32,db40,db48,db01,za,t0, &
1611  z0,g0,rl,dd,db14,tr12,tn1,tgn1, &
1612  pt,pd,ps,pdl,ptl,pma,sw, &
1614 
1615  dimension zn1(5),alpha(9)
1616  dimension d(9),t(2),mt(11),ap(*),altl(8)
1617  save
1618  data mt/48,0,4,16,28,32,40,1,49,14,17/
1619  data altl/200.,300.,160.,250.,240.,450.,320.,450./
1620  data mn1/5/,zn1/120.,110.,100.,90.,72.5/
1621  data dgtr/1.74533e-2/,dr/1.72142e-2/,alast/-999./
1622  data alpha/-0.38,0.,0.,0.,0.17,0.,-0.38,0.,0./
1623 ! test for changed input
1624  v2=vtst7(iyd,sec,glat,glong,stl,f107a,f107,ap,2)
1625 !
1626  yrd=iyd
1627  za=pdl(16,2)
1628  zn1(1)=za
1629  do 2 j=1,9
1630  d(j)=0.
1631  2 end do
1632 ! tinf variations not important below za or zn1(1)
1633  if(alt.gt.zn1(1)) then
1634  if(v2.eq.1..or.alast.le.zn1(1)) tinf=ptm(1)*pt(1) &
1635  *(1.+sw(16)*globe7(yrd,sec,glat,glong,stl,f107a,f107,ap,pt))
1636  else
1637  tinf=ptm(1)*pt(1)
1638  endif
1639  t(1)=tinf
1640 ! gradient variations not important below zn1(5)
1641  if(alt.gt.zn1(5)) then
1642  if(v2.eq.1.or.alast.le.zn1(5)) g0=ptm(4)*ps(1) &
1643  *(1.+sw(19)*globe7(yrd,sec,glat,glong,stl,f107a,f107,ap,ps))
1644  else
1645  g0=ptm(4)*ps(1)
1646  endif
1647 ! calculate these temperatures only if input changed
1648  if(v2.eq.1. .or. alt.lt.300.) &
1649  tlb=ptm(2)*(1.+sw(17)*globe7(yrd,sec,glat,glong,stl, &
1650  f107a,f107,ap,pd(1,4)))*pd(1,4)
1651  s=g0/(tinf-tlb)
1652 ! lower thermosphere temp variations not significant for
1653 ! density above 300 km
1654  if(alt.lt.300.) then
1655  if(v2.eq.1..or.alast.ge.300.) then
1656  tn1(2)=ptm(7)*ptl(1,1)/(1.-sw(18)*glob7s(ptl(1,1)))
1657  tn1(3)=ptm(3)*ptl(1,2)/(1.-sw(18)*glob7s(ptl(1,2)))
1658  tn1(4)=ptm(8)*ptl(1,3)/(1.-sw(18)*glob7s(ptl(1,3)))
1659  tn1(5)=ptm(5)*ptl(1,4)/(1.-sw(18)*sw(20)*glob7s(ptl(1,4)))
1660  tgn1(2)=ptm(9)*pma(1,9)*(1.+sw(18)*sw(20)*glob7s(pma(1,9))) &
1661  *tn1(5)*tn1(5)/(ptm(5)*ptl(1,4))**2
1662  endif
1663  else
1664  tn1(2)=ptm(7)*ptl(1,1)
1665  tn1(3)=ptm(3)*ptl(1,2)
1666  tn1(4)=ptm(8)*ptl(1,3)
1667  tn1(5)=ptm(5)*ptl(1,4)
1668  tgn1(2)=ptm(9)*pma(1,9) &
1669  *tn1(5)*tn1(5)/(ptm(5)*ptl(1,4))**2
1670  endif
1671 !
1672  z0=zn1(4)
1673  t0=tn1(4)
1674  tr12=1.
1675 !
1676  if(mass.eq.0) go to 50
1677 ! n2 variation factor at zlb
1678  g28=sw(21)*globe7(yrd,sec,glat,glong,stl,f107a,f107, &
1679  ap,pd(1,3))
1680  day=amod(yrd,1000.)
1681 ! variation of turbopause height
1682  zhf=pdl(25,2) &
1683  *(1.+sw(5)*pdl(25,1)*sin(dgtr*glat)*cos(dr*(day-pt(14))))
1684  yrd=iyd
1685  t(1)=tinf
1686  xmm=pdm(5,3)
1687  z=alt
1688 !
1689  do 10 j = 1,11
1690  if(mass.eq.mt(j)) go to 15
1691  10 end do
1692  write(6,100) mass
1693  go to 90
1694  15 if(z.gt.altl(6).and.mass.ne.28.and.mass.ne.48) go to 17
1695 !
1696 ! **** n2 density ****
1697 !
1698 ! diffusive density at zlb
1699  db28 = pdm(1,3)*exp(g28)*pd(1,3)
1700 ! diffusive density at alt
1701  d(3)=densu(z,db28,tinf,tlb, 28.,alpha(3),t(2),ptm(6),s,mn1,zn1, &
1702  tn1,tgn1)
1703  dd=d(3)
1704 ! turbopause
1705  zh28=pdm(3,3)*zhf
1706  zhm28=pdm(4,3)*pdl(6,2)
1707  xmd=28.-xmm
1708 ! mixed density at zlb
1709  b28=densu(zh28,db28,tinf,tlb,xmd,alpha(3)-1.,tz,ptm(6),s,mn1, &
1710  zn1,tn1,tgn1)
1711  if(z.gt.altl(3).or.sw(15).eq.0.) go to 17
1712 ! mixed density at alt
1713  dm28=densu(z,b28,tinf,tlb,xmm,alpha(3),tz,ptm(6),s,mn1, &
1714  zn1,tn1,tgn1)
1715 ! net density at alt
1716  d(3)=dnet(d(3),dm28,zhm28,xmm,28.)
1717  17 continue
1718  go to (20,50,20,25,90,35,40,45,25,48,46), j
1719  20 continue
1720 !
1721 ! **** he density ****
1722 !
1723 ! density variation factor at zlb
1724  g4 = sw(21)*globe7(yrd,sec,glat,glong,stl,f107a,f107,ap,pd(1,1))
1725 ! diffusive density at zlb
1726  db04 = pdm(1,1)*exp(g4)*pd(1,1)
1727 ! diffusive density at alt
1728  d(1)=densu(z,db04,tinf,tlb, 4.,alpha(1),t(2),ptm(6),s,mn1,zn1, &
1729  tn1,tgn1)
1730  dd=d(1)
1731  if(z.gt.altl(1).or.sw(15).eq.0.) go to 24
1732 ! turbopause
1733  zh04=pdm(3,1)
1734 ! mixed density at zlb
1735  b04=densu(zh04,db04,tinf,tlb,4.-xmm,alpha(1)-1., &
1736  t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1737 ! mixed density at alt
1738  dm04=densu(z,b04,tinf,tlb,xmm,0.,t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1739  zhm04=zhm28
1740 ! net density at alt
1741  d(1)=dnet(d(1),dm04,zhm04,xmm,4.)
1742 ! correction to specified mixing ratio at ground
1743  rl=alog(b28*pdm(2,1)/b04)
1744  zc04=pdm(5,1)*pdl(1,2)
1745  hc04=pdm(6,1)*pdl(2,2)
1746 ! net density corrected at alt
1747  d(1)=d(1)*ccor(z,rl,hc04,zc04)
1748  24 continue
1749  if(mass.ne.48) go to 90
1750  25 continue
1751 !
1752 ! **** o density ****
1753 !
1754 ! density variation factor at zlb
1755  g16= sw(21)*globe7(yrd,sec,glat,glong,stl,f107a,f107,ap,pd(1,2))
1756 ! diffusive density at zlb
1757  db16 = pdm(1,2)*exp(g16)*pd(1,2)
1758 ! diffusive density at alt
1759  d(2)=densu(z,db16,tinf,tlb, 16.,alpha(2),t(2),ptm(6),s,mn1, &
1760  zn1,tn1,tgn1)
1761  dd=d(2)
1762  if(z.gt.altl(2).or.sw(15).eq.0.) go to 34
1763 ! corrected from pdm(3,1) to pdm(3,2) 12/2/85
1764 ! turbopause
1765  zh16=pdm(3,2)
1766 ! mixed density at zlb
1767  b16=densu(zh16,db16,tinf,tlb,16-xmm,alpha(2)-1., &
1768  t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1769 ! mixed density at alt
1770  dm16=densu(z,b16,tinf,tlb,xmm,0.,t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1771  zhm16=zhm28
1772 ! net density at alt
1773  d(2)=dnet(d(2),dm16,zhm16,xmm,16.)
1774 ! 3/16/99 change form to match o2 departure from diff equil near 150
1775 ! km and add dependence on f10.7
1776 ! rl=alog(b28*pdm(2,2)*abs(pdl(17,2))/b16)
1777  rl=pdm(2,2)*pdl(17,2)*(1.+sw(1)*pdl(24,1)*(f107a-150.))
1778  hc16=pdm(6,2)*pdl(4,2)
1779  zc16=pdm(5,2)*pdl(3,2)
1780  hc216=pdm(6,2)*pdl(5,2)
1781  d(2)=d(2)*ccor2(z,rl,hc16,zc16,hc216)
1782 ! chemistry correction
1783  hcc16=pdm(8,2)*pdl(14,2)
1784  zcc16=pdm(7,2)*pdl(13,2)
1785  rc16=pdm(4,2)*pdl(15,2)
1786 ! net density corrected at alt
1787  d(2)=d(2)*ccor(z,rc16,hcc16,zcc16)
1788  34 continue
1789  if(mass.ne.48.and.mass.ne.49) go to 90
1790  35 continue
1791 !
1792 ! **** o2 density ****
1793 !
1794 ! density variation factor at zlb
1795  g32= sw(21)*globe7(yrd,sec,glat,glong,stl,f107a,f107,ap,pd(1,5))
1796 ! diffusive density at zlb
1797  db32 = pdm(1,4)*exp(g32)*pd(1,5)
1798 ! diffusive density at alt
1799  d(4)=densu(z,db32,tinf,tlb, 32.,alpha(4),t(2),ptm(6),s,mn1, &
1800  zn1,tn1,tgn1)
1801  if(mass.eq.49) then
1802  dd=dd+2.*d(4)
1803  else
1804  dd=d(4)
1805  endif
1806  if(sw(15).eq.0.) go to 39
1807  if(z.gt.altl(4)) go to 38
1808 ! turbopause
1809  zh32=pdm(3,4)
1810 ! mixed density at zlb
1811  b32=densu(zh32,db32,tinf,tlb,32.-xmm,alpha(4)-1., &
1812  t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1813 ! mixed density at alt
1814  dm32=densu(z,b32,tinf,tlb,xmm,0.,t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1815  zhm32=zhm28
1816 ! net density at alt
1817  d(4)=dnet(d(4),dm32,zhm32,xmm,32.)
1818 ! correction to specified mixing ratio at ground
1819  rl=alog(b28*pdm(2,4)/b32)
1820  hc32=pdm(6,4)*pdl(8,2)
1821  zc32=pdm(5,4)*pdl(7,2)
1822  d(4)=d(4)*ccor(z,rl,hc32,zc32)
1823  38 continue
1824 ! correction for general departure from diffusive equilibrium above
1825  hcc32=pdm(8,4)*pdl(23,2)
1826  hcc232=pdm(8,4)*pdl(23,1)
1827  zcc32=pdm(7,4)*pdl(22,2)
1828  rc32=pdm(4,4)*pdl(24,2)*(1.+sw(1)*pdl(24,1)*(f107a-150.))
1829 ! net density corrected at alt
1830  d(4)=d(4)*ccor2(z,rc32,hcc32,zcc32,hcc232)
1831  39 continue
1832  if(mass.ne.48) go to 90
1833  40 continue
1834 !
1835 ! **** ar density ****
1836 !
1837 ! density variation factor at zlb
1838  g40= sw(21)*globe7(yrd,sec,glat,glong,stl,f107a,f107,ap,pd(1,6))
1839 ! diffusive density at zlb
1840  db40 = pdm(1,5)*exp(g40)*pd(1,6)
1841 ! diffusive density at alt
1842  d(5)=densu(z,db40,tinf,tlb, 40.,alpha(5),t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1843  dd=d(5)
1844  if(z.gt.altl(5).or.sw(15).eq.0.) go to 44
1845 ! turbopause
1846  zh40=pdm(3,5)
1847 ! mixed density at zlb
1848  b40=densu(zh40,db40,tinf,tlb,40.-xmm,alpha(5)-1., &
1849  t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1850 ! mixed density at alt
1851  dm40=densu(z,b40,tinf,tlb,xmm,0.,t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1852  zhm40=zhm28
1853 ! net density at alt
1854  d(5)=dnet(d(5),dm40,zhm40,xmm,40.)
1855 ! correction to specified mixing ratio at ground
1856  rl=alog(b28*pdm(2,5)/b40)
1857  hc40=pdm(6,5)*pdl(10,2)
1858  zc40=pdm(5,5)*pdl(9,2)
1859 ! net density corrected at alt
1860  d(5)=d(5)*ccor(z,rl,hc40,zc40)
1861  44 continue
1862  if(mass.ne.48) go to 90
1863  45 continue
1864 !
1865 ! **** hydrogen density ****
1866 !
1867 ! density variation factor at zlb
1868  g1 = sw(21)*globe7(yrd,sec,glat,glong,stl,f107a,f107,ap,pd(1,7))
1869 ! diffusive density at zlb
1870  db01 = pdm(1,6)*exp(g1)*pd(1,7)
1871 ! diffusive density at alt
1872  d(7)=densu(z,db01,tinf,tlb,1.,alpha(7),t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1873  dd=d(7)
1874  if(z.gt.altl(7).or.sw(15).eq.0.) go to 47
1875 ! turbopause
1876  zh01=pdm(3,6)
1877 ! mixed density at zlb
1878  b01=densu(zh01,db01,tinf,tlb,1.-xmm,alpha(7)-1., &
1879  t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1880 ! mixed density at alt
1881  dm01=densu(z,b01,tinf,tlb,xmm,0.,t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1882  zhm01=zhm28
1883 ! net density at alt
1884  d(7)=dnet(d(7),dm01,zhm01,xmm,1.)
1885 ! correction to specified mixing ratio at ground
1886  rl=alog(b28*pdm(2,6)*abs(pdl(18,2))/b01)
1887  hc01=pdm(6,6)*pdl(12,2)
1888  zc01=pdm(5,6)*pdl(11,2)
1889  d(7)=d(7)*ccor(z,rl,hc01,zc01)
1890 ! chemistry correction
1891  hcc01=pdm(8,6)*pdl(20,2)
1892  zcc01=pdm(7,6)*pdl(19,2)
1893  rc01=pdm(4,6)*pdl(21,2)
1894 ! net density corrected at alt
1895  d(7)=d(7)*ccor(z,rc01,hcc01,zcc01)
1896  47 continue
1897  if(mass.ne.48) go to 90
1898  48 continue
1899 !
1900 ! **** atomic nitrogen density ****
1901 !
1902 ! density variation factor at zlb
1903  g14 = sw(21)*globe7(yrd,sec,glat,glong,stl,f107a,f107,ap,pd(1,8))
1904 ! diffusive density at zlb
1905  db14 = pdm(1,7)*exp(g14)*pd(1,8)
1906 ! diffusive density at alt
1907  d(8)=densu(z,db14,tinf,tlb,14.,alpha(8),t(2),ptm(6),s,mn1, &
1908  zn1,tn1,tgn1)
1909  dd=d(8)
1910  if(z.gt.altl(8).or.sw(15).eq.0.) go to 49
1911 ! turbopause
1912  zh14=pdm(3,7)
1913 ! mixed density at zlb
1914  b14=densu(zh14,db14,tinf,tlb,14.-xmm,alpha(8)-1., &
1915  t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1916 ! mixed density at alt
1917  dm14=densu(z,b14,tinf,tlb,xmm,0.,t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1918  zhm14=zhm28
1919 ! net density at alt
1920  d(8)=dnet(d(8),dm14,zhm14,xmm,14.)
1921 ! correction to specified mixing ratio at ground
1922  rl=alog(b28*pdm(2,7)*abs(pdl(3,1))/b14)
1923  hc14=pdm(6,7)*pdl(2,1)
1924  zc14=pdm(5,7)*pdl(1,1)
1925  d(8)=d(8)*ccor(z,rl,hc14,zc14)
1926 ! chemistry correction
1927  hcc14=pdm(8,7)*pdl(5,1)
1928  zcc14=pdm(7,7)*pdl(4,1)
1929  rc14=pdm(4,7)*pdl(6,1)
1930 ! net density corrected at alt
1931  d(8)=d(8)*ccor(z,rc14,hcc14,zcc14)
1932  49 continue
1933  if(mass.ne.48) go to 90
1934  46 continue
1935 !
1936 ! **** anomalous oxygen density ****
1937 !
1938  g16h = sw(21)*globe7(yrd,sec,glat,glong,stl,f107a,f107,ap,pd(1,9))
1939  db16h = pdm(1,8)*exp(g16h)*pd(1,9)
1940  tho=pdm(10,8)*pdl(7,1)
1941  dd=densu(z,db16h,tho,tho,16.,alpha(9),t2,ptm(6),s,mn1, &
1942  zn1,tn1,tgn1)
1943  zsht=pdm(6,8)
1944  zmho=pdm(5,8)
1945  zsho=scalh(zmho,16.,tho)
1946  d(9)=dd*exp(-zsht/zsho*(exp(-(z-zmho)/zsht)-1.))
1947  if(mass.ne.48) go to 90
1948 !
1949 ! total mass density
1950 !
1951  d(6) = 1.66e-24*(4.*d(1)+16.*d(2)+28.*d(3)+32.*d(4)+40.*d(5)+ &
1952  d(7)+14.*d(8))
1953  db48=1.66e-24*(4.*db04+16.*db16+28.*db28+32.*db32+40.*db40+db01+ &
1954  14.*db14)
1955  go to 90
1956 ! temperature at altitude
1957  50 continue
1958  z=abs(alt)
1959  ddum = densu(z,1., tinf,tlb,0.,0.,t(2),ptm(6),s,mn1,zn1,tn1,tgn1)
1960  90 continue
1961 ! adjust densities from cgs to kgm
1962  if(imr.eq.1) then
1963  do 95 i=1,9
1964  d(i)=d(i)*1.e6
1965  95 continue
1966  d(6)=d(6)/1000.
1967  endif
1968  alast=alt
1969  return
1970  100 format(1x,'mass', i5, ' not valid')
1971  end subroutine gts7
1972 
1973 !-----------------------------------------------------------------------
1977  subroutine meters(meter)
1978  use wam_gtd7bk_mod, only: imr
1979 !
1980  logical meter
1981  save
1982  imr=0
1983  if(meter) imr=1
1984  end subroutine meters
1985 
1986 !-----------------------------------------------------------------------
1994  function scalh(alt,xm,temp)
1995  use gettemp_mod, only: gsurf,re
1996 !
1997  save
1998  data rgas/831.4/
1999  g=gsurf/(1.+alt/re)**2
2000  scalh=rgas*temp/(g*xm)
2001  return
2002  end function scalh
2003 
2004 !-----------------------------------------------------------------------
2028  function globe7(yrd,sec,lat,long,tloc,f107a,f107,ap,p)
2029  use gettemp_mod, only: tinf=>tinfg,t=>tt, &
2030  sw,swc,isw, &
2033  real lat, long
2034  dimension p(*),sv(25),ap(*)
2035 !---- functions ------
2036 ! 3hr magnetic activity functions
2037 ! eq. a24d
2038  g0(a)=(a-4.+(p(26)-1.)*(a-4.+(exp(-abs(p(25))*(a-4.))-1.)/ &
2039  abs(p(25))))
2040 ! eq. a24c
2041  sumex(ex)=1.+(1.-ex**19)/(1.-ex)*ex**(.5)
2042 ! eq. a24a
2043  sg0(ex)=(g0(ap(2))+(g0(ap(3))*ex+g0(ap(4))*ex*ex+g0(ap(5))*ex**3 &
2044  +(g0(ap(6))*ex**4+g0(ap(7))*ex**12)*(1.-ex**8)/(1.-ex)) &
2045  )/sumex(ex)
2046 !---------------------
2047  save
2048  data dgtr/1.74533e-2/,dr/1.72142e-2/, xl/1000./,tll/1000./
2049  data sw9/1./,dayl/-1./,p14/-1000./,p18/-1000./,p32/-1000./
2050  data hr/.2618/,sr/7.2722e-5/,sv/25*1./,nsw/14/,p39/-1000./
2051  if(isw.ne.64999) call tselec(sv)
2052  do 10 j=1,14
2053  t(j)=0
2054  10 end do
2055  if(sw(9).gt.0) sw9=1.
2056  if(sw(9).lt.0) sw9=-1.
2057  iyr = nint(yrd/1000.)
2058  day = yrd - iyr*1000.
2059  xlong=long
2060 ! eq. a22 (remainder of code)
2061  if(xl.eq.lat) go to 15
2062 ! calculate legendre polynomials
2063  c = sin(lat*dgtr)
2064  s = cos(lat*dgtr)
2065  c2 = c*c
2066  c4 = c2*c2
2067  s2 = s*s
2068  plg(2,1) = c
2069  plg(3,1) = 0.5*(3.*c2 -1.)
2070  plg(4,1) = 0.5*(5.*c*c2-3.*c)
2071  plg(5,1) = (35.*c4 - 30.*c2 + 3.)/8.
2072  plg(6,1) = (63.*c2*c2*c - 70.*c2*c + 15.*c)/8.
2073  plg(7,1) = (11.*c*plg(6,1) - 5.*plg(5,1))/6.
2074 ! plg(8,1) = (13.*c*plg(7,1) - 6.*plg(6,1))/7.
2075  plg(2,2) = s
2076  plg(3,2) = 3.*c*s
2077  plg(4,2) = 1.5*(5.*c2-1.)*s
2078  plg(5,2) = 2.5*(7.*c2*c-3.*c)*s
2079  plg(6,2) = 1.875*(21.*c4 - 14.*c2 +1.)*s
2080  plg(7,2) = (11.*c*plg(6,2)-6.*plg(5,2))/5.
2081 ! plg(8,2) = (13.*c*plg(7,2)-7.*plg(6,2))/6.
2082 ! plg(9,2) = (15.*c*plg(8,2)-8.*plg(7,2))/7.
2083  plg(3,3) = 3.*s2
2084  plg(4,3) = 15.*s2*c
2085  plg(5,3) = 7.5*(7.*c2 -1.)*s2
2086  plg(6,3) = 3.*c*plg(5,3)-2.*plg(4,3)
2087  plg(7,3)=(11.*c*plg(6,3)-7.*plg(5,3))/4.
2088  plg(8,3)=(13.*c*plg(7,3)-8.*plg(6,3))/5.
2089  plg(4,4) = 15.*s2*s
2090  plg(5,4) = 105.*s2*s*c
2091  plg(6,4)=(9.*c*plg(5,4)-7.*plg(4,4))/2.
2092  plg(7,4)=(11.*c*plg(6,4)-8.*plg(5,4))/3.
2093  xl=lat
2094  15 continue
2095  if(tll.eq.tloc) go to 16
2096  if(sw(7).eq.0.and.sw(8).eq.0.and.sw(14).eq.0) goto 16
2097  stloc = sin(hr*tloc)
2098  ctloc = cos(hr*tloc)
2099  s2tloc = sin(2.*hr*tloc)
2100  c2tloc = cos(2.*hr*tloc)
2101  s3tloc = sin(3.*hr*tloc)
2102  c3tloc = cos(3.*hr*tloc)
2103  tll = tloc
2104  16 continue
2105  if(day.ne.dayl.or.p(14).ne.p14) cd14=cos(dr*(day-p(14)))
2106  if(day.ne.dayl.or.p(18).ne.p18) cd18=cos(2.*dr*(day-p(18)))
2107  if(day.ne.dayl.or.p(32).ne.p32) cd32=cos(dr*(day-p(32)))
2108  if(day.ne.dayl.or.p(39).ne.p39) cd39=cos(2.*dr*(day-p(39)))
2109  dayl = day
2110  p14 = p(14)
2111  p18 = p(18)
2112  p32 = p(32)
2113  p39 = p(39)
2114 ! f10.7 effect
2115  df = f107 - f107a
2116  dfa=f107a-150.
2117  t(1) = p(20)*df*(1.+p(60)*dfa) + p(21)*df*df + p(22)*dfa &
2118  + p(30)*dfa**2
2119  f1 = 1. + (p(48)*dfa +p(20)*df+p(21)*df*df)*swc(1)
2120  f2 = 1. + (p(50)*dfa+p(20)*df+p(21)*df*df)*swc(1)
2121 ! time independent
2122  t(2) = &
2123  (p(2)*plg(3,1) + p(3)*plg(5,1)+p(23)*plg(7,1)) &
2124  +(p(15)*plg(3,1))*dfa*swc(1) &
2125  +p(27)*plg(2,1)
2126 ! symmetrical annual
2127  t(3) = &
2128  (p(19) )*cd32
2129 ! symmetrical semiannual
2130  t(4) = &
2131  (p(16)+p(17)*plg(3,1))*cd18
2132 ! asymmetrical annual
2133  t(5) = f1* &
2134  (p(10)*plg(2,1)+p(11)*plg(4,1))*cd14
2135 ! asymmetrical semiannual
2136  t(6) = p(38)*plg(2,1)*cd39
2137 ! diurnal
2138  if(sw(7).eq.0) goto 200
2139  t71 = (p(12)*plg(3,2))*cd14*swc(5)
2140  t72 = (p(13)*plg(3,2))*cd14*swc(5)
2141  t(7) = f2* &
2142  ((p(4)*plg(2,2) + p(5)*plg(4,2) + p(28)*plg(6,2) &
2143  + t71)*ctloc &
2144  + (p(7)*plg(2,2) + p(8)*plg(4,2) +p(29)*plg(6,2) &
2145  + t72)*stloc)
2146  200 continue
2147 ! semidiurnal
2148  if(sw(8).eq.0) goto 210
2149  t81 = (p(24)*plg(4,3)+p(36)*plg(6,3))*cd14*swc(5)
2150  t82 = (p(34)*plg(4,3)+p(37)*plg(6,3))*cd14*swc(5)
2151  t(8) = f2* &
2152  ((p(6)*plg(3,3) + p(42)*plg(5,3) + t81)*c2tloc &
2153  +(p(9)*plg(3,3) + p(43)*plg(5,3) + t82)*s2tloc)
2154  210 continue
2155 ! terdiurnal
2156  if(sw(14).eq.0) goto 220
2157  t(14) = f2* &
2158  ((p(40)*plg(4,4)+(p(94)*plg(5,4)+p(47)*plg(7,4))*cd14*swc(5))* &
2159  s3tloc &
2160  +(p(41)*plg(4,4)+(p(95)*plg(5,4)+p(49)*plg(7,4))*cd14*swc(5))* &
2161  c3tloc)
2162  220 continue
2163 ! magnetic activity based on daily ap
2164 
2165  if(sw9.eq.-1.) go to 30
2166  apd=(ap(1)-4.)
2167  p44=p(44)
2168  p45=p(45)
2169  if(p44.lt.0) p44=1.e-5
2170  apdf = apd+(p45-1.)*(apd+(exp(-p44 *apd)-1.)/p44)
2171  if(sw(9).eq.0) goto 40
2172  t(9)=apdf*(p(33)+p(46)*plg(3,1)+p(35)*plg(5,1)+ &
2173  (p(101)*plg(2,1)+p(102)*plg(4,1)+p(103)*plg(6,1))*cd14*swc(5)+ &
2174  (p(122)*plg(2,2)+p(123)*plg(4,2)+p(124)*plg(6,2))*swc(7)* &
2175  cos(hr*(tloc-p(125))))
2176  go to 40
2177  30 continue
2178  if(p(52).eq.0) go to 40
2179  exp1 = exp(-10800.*abs(p(52))/(1.+p(139)*(45.-abs(lat))))
2180  if(exp1.gt..99999) exp1=.99999
2181  if(p(25).lt.1.e-4) p(25)=1.e-4
2182  apt(1)=sg0(exp1)
2183 ! apt(2)=sg2(exp1)
2184 ! apt(3)=sg0(exp2)
2185 ! apt(4)=sg2(exp2)
2186  if(sw(9).eq.0) goto 40
2187  t(9) = apt(1)*(p(51)+p(97)*plg(3,1)+p(55)*plg(5,1)+ &
2188  (p(126)*plg(2,1)+p(127)*plg(4,1)+p(128)*plg(6,1))*cd14*swc(5)+ &
2189  (p(129)*plg(2,2)+p(130)*plg(4,2)+p(131)*plg(6,2))*swc(7)* &
2190  cos(hr*(tloc-p(132))))
2191  40 continue
2192  if(sw(10).eq.0.or.long.le.-1000.) go to 49
2193 ! longitudinal
2194  if(sw(11).eq.0) goto 230
2195  t(11)= (1.+p(81)*dfa*swc(1))* &
2196  ((p(65)*plg(3,2)+p(66)*plg(5,2)+p(67)*plg(7,2) &
2197  +p(104)*plg(2,2)+p(105)*plg(4,2)+p(106)*plg(6,2) &
2198  +swc(5)*(p(110)*plg(2,2)+p(111)*plg(4,2)+p(112)*plg(6,2))*cd14)* &
2199  cos(dgtr*long) &
2200  +(p(91)*plg(3,2)+p(92)*plg(5,2)+p(93)*plg(7,2) &
2201  +p(107)*plg(2,2)+p(108)*plg(4,2)+p(109)*plg(6,2) &
2202  +swc(5)*(p(113)*plg(2,2)+p(114)*plg(4,2)+p(115)*plg(6,2))*cd14)* &
2203  sin(dgtr*long))
2204  230 continue
2205 ! ut and mixed ut,longitude
2206  if(sw(12).eq.0) goto 240
2207  t(12)=(1.+p(96)*plg(2,1))*(1.+p(82)*dfa*swc(1))* &
2208  (1.+p(120)*plg(2,1)*swc(5)*cd14)* &
2209  ((p(69)*plg(2,1)+p(70)*plg(4,1)+p(71)*plg(6,1))* &
2210  cos(sr*(sec-p(72))))
2211  t(12)=t(12)+swc(11)* &
2212  (p(77)*plg(4,3)+p(78)*plg(6,3)+p(79)*plg(8,3))* &
2213  cos(sr*(sec-p(80))+2.*dgtr*long)*(1.+p(138)*dfa*swc(1))
2214  240 continue
2215 ! ut,longitude magnetic activity
2216  if(sw(13).eq.0) goto 48
2217  if(sw9.eq.-1.) go to 45
2218  t(13)= apdf*swc(11)*(1.+p(121)*plg(2,1))* &
2219  ((p( 61)*plg(3,2)+p( 62)*plg(5,2)+p( 63)*plg(7,2))* &
2220  cos(dgtr*(long-p( 64)))) &
2221  +apdf*swc(11)*swc(5)* &
2222  (p(116)*plg(2,2)+p(117)*plg(4,2)+p(118)*plg(6,2))* &
2223  cd14*cos(dgtr*(long-p(119))) &
2224  + apdf*swc(12)* &
2225  (p( 84)*plg(2,1)+p( 85)*plg(4,1)+p( 86)*plg(6,1))* &
2226  cos(sr*(sec-p( 76)))
2227  goto 48
2228  45 continue
2229  if(p(52).eq.0) goto 48
2230  t(13)=apt(1)*swc(11)*(1.+p(133)*plg(2,1))* &
2231  ((p(53)*plg(3,2)+p(99)*plg(5,2)+p(68)*plg(7,2))* &
2232  cos(dgtr*(long-p(98)))) &
2233  +apt(1)*swc(11)*swc(5)* &
2234  (p(134)*plg(2,2)+p(135)*plg(4,2)+p(136)*plg(6,2))* &
2235  cd14*cos(dgtr*(long-p(137))) &
2236  +apt(1)*swc(12)* &
2237  (p(56)*plg(2,1)+p(57)*plg(4,1)+p(58)*plg(6,1))* &
2238  cos(sr*(sec-p(59)))
2239  48 continue
2240 ! parms not used: 83, 90,100,140-150
2241  49 continue
2242  tinf=p(31)
2243  do i = 1,nsw
2244  tinf = tinf + abs(sw(i))*t(i)
2245  enddo
2246  globe7 = tinf
2247  return
2248  end function globe7
2249 
2250 !-----------------------------------------------------------------------
2262  subroutine tselec(sv)
2263  use gettemp_mod, only: sw,swc,isw
2264 !
2265  dimension sv(*),sav(25),svv(*)
2266  save
2267  do 100 i = 1,25
2268  sav(i)=sv(i)
2269  sw(i)=amod(sv(i),2.)
2270  if(abs(sv(i)).eq.1.or.abs(sv(i)).eq.2.) then
2271  swc(i)=1.
2272  else
2273  swc(i)=0.
2274  endif
2275  100 end do
2276  isw=64999
2277  return
2278  entry tretrv(svv)
2279  do 200 i=1,25
2280  svv(i)=sav(i)
2281  200 end do
2282  end subroutine tselec
2283 
2284 !-----------------------------------------------------------------------
2291  function glob7s(p)
2292  use mpi
2294  day,dfa,apdf,apt,long=>xlong,sw,swc
2295  dimension p(*),t(14)
2296 !
2297  integer :: rc
2298  save
2299  data dr/1.72142e-2/,dgtr/1.74533e-2/,pset/2./
2300  data dayl/-1./,p32,p18,p14,p39/4*-1000./
2301  if(p(100).eq.0) p(100)=pset
2302  if(p(100).ne.pset) then
2303  write(6,900) pset,p(100)
2304  900 format(1x,'FATAL ERROR: Wrong parameter set for glob7s',3f10.1)
2305  call mpi_abort(mpi_comm_world, 999, rc)
2306  endif
2307  do 10 j=1,14
2308  t(j)=0.
2309  10 end do
2310  if(day.ne.dayl.or.p32.ne.p(32)) cd32=cos(dr*(day-p(32)))
2311  if(day.ne.dayl.or.p18.ne.p(18)) cd18=cos(2.*dr*(day-p(18)))
2312  if(day.ne.dayl.or.p14.ne.p(14)) cd14=cos(dr*(day-p(14)))
2313  if(day.ne.dayl.or.p39.ne.p(39)) cd39=cos(2.*dr*(day-p(39)))
2314  dayl=day
2315  p32=p(32)
2316  p18=p(18)
2317  p14=p(14)
2318  p39=p(39)
2319 !
2320 ! f10.7
2321  t(1)=p(22)*dfa
2322 ! time independent
2323  t(2)=p(2)*plg(3,1)+p(3)*plg(5,1)+p(23)*plg(7,1) &
2324  +p(27)*plg(2,1)+p(15)*plg(4,1)+p(60)*plg(6,1)
2325 ! symmetrical annual
2326  t(3)=(p(19)+p(48)*plg(3,1)+p(30)*plg(5,1))*cd32
2327 ! symmetrical semiannual
2328  t(4)=(p(16)+p(17)*plg(3,1)+p(31)*plg(5,1))*cd18
2329 ! asymmetrical annual
2330  t(5)=(p(10)*plg(2,1)+p(11)*plg(4,1)+p(21)*plg(6,1))*cd14
2331 ! asymmetrical semiannual
2332  t(6)=(p(38)*plg(2,1))*cd39
2333 ! diurnal
2334  if(sw(7).eq.0) goto 200
2335  t71 = p(12)*plg(3,2)*cd14*swc(5)
2336  t72 = p(13)*plg(3,2)*cd14*swc(5)
2337  t(7) = &
2338  ((p(4)*plg(2,2) + p(5)*plg(4,2) &
2339  + t71)*ctloc &
2340  + (p(7)*plg(2,2) + p(8)*plg(4,2) &
2341  + t72)*stloc)
2342  200 continue
2343 ! semidiurnal
2344  if(sw(8).eq.0) goto 210
2345  t81 = (p(24)*plg(4,3)+p(36)*plg(6,3))*cd14*swc(5)
2346  t82 = (p(34)*plg(4,3)+p(37)*plg(6,3))*cd14*swc(5)
2347  t(8) = &
2348  ((p(6)*plg(3,3) + p(42)*plg(5,3) + t81)*c2tloc &
2349  +(p(9)*plg(3,3) + p(43)*plg(5,3) + t82)*s2tloc)
2350  210 continue
2351 ! terdiurnal
2352  if(sw(14).eq.0) goto 220
2353  t(14) = p(40)*plg(4,4)*s3tloc +p(41)*plg(4,4)*c3tloc
2354  220 continue
2355 ! magnetic activity
2356  if(sw(9).eq.0) goto 40
2357  if(sw(9).eq.1) &
2358  t(9)=apdf*(p(33)+p(46)*plg(3,1)*swc(2))
2359  if(sw(9).eq.-1) &
2360  t(9)=(p(51)*apt(1)+p(97)*plg(3,1)*apt(1)*swc(2))
2361  40 continue
2362  if(sw(10).eq.0.or.sw(11).eq.0.or.long.le.-1000.) go to 49
2363 ! longitudinal
2364  t(11)= (1.+plg(2,1)*(p(81)*swc(5)*cos(dr*(day-p(82))) &
2365  +p(86)*swc(6)*cos(2.*dr*(day-p(87)))) &
2366  +p(84)*swc(3)*cos(dr*(day-p(85))) &
2367  +p(88)*swc(4)*cos(2.*dr*(day-p(89)))) &
2368  *((p(65)*plg(3,2)+p(66)*plg(5,2)+p(67)*plg(7,2) &
2369  +p(75)*plg(2,2)+p(76)*plg(4,2)+p(77)*plg(6,2) &
2370  )*cos(dgtr*long) &
2371  +(p(91)*plg(3,2)+p(92)*plg(5,2)+p(93)*plg(7,2) &
2372  +p(78)*plg(2,2)+p(79)*plg(4,2)+p(80)*plg(6,2) &
2373  )*sin(dgtr*long))
2374  49 continue
2375  tt=0.
2376  do i=1,14
2377  tt=tt+abs(sw(i))*t(i)
2378  enddo
2379  glob7s=tt
2380  return
2381  end function glob7s
2382 
2383 !--------------------------------------------------------------------
2403  function densu(alt,dlb,tinf,tlb,xm,alpha,tz,zlb,s2, &
2404  mn1,zn1,tn1,tgn1)
2405  use gettemp_mod, only: gsurf,re
2406 !
2407  dimension zn1(mn1),tn1(mn1),tgn1(2),xs(5),ys(5),y2out(5)
2408 !function
2409  zeta(zz,zl)=(zz-zl)*(re+zl)/(re+zz)
2410  save
2411  data rgas/831.4/
2412 !! rgas=831.4
2413 !cccccwrite(6,*) 'db',alt,dlb,tinf,tlb,xm,alpha,zlb,s2,mn1,zn1,tn1
2414  densu=1.
2415 ! joining altitude of bates and spline
2416  za=zn1(1)
2417  z=amax1(alt,za)
2418 ! geopotential altitude difference from zlb
2419  zg2=zeta(z,zlb)
2420 ! bates temperature
2421  tt=tinf-(tinf-tlb)*exp(-s2*zg2)
2422  ta=tt
2423  tz=tt
2424  densu=tz
2425  if(alt.ge.za) go to 10
2426 !
2427 ! calculate temperature below za
2428 ! temperature gradient at za from bates profile
2429  dta=(tinf-ta)*s2*((re+zlb)/(re+za))**2
2430  tgn1(1)=dta
2431  tn1(1)=ta
2432  z=amax1(alt,zn1(mn1))
2433  mn=mn1
2434  z1=zn1(1)
2435  z2=zn1(mn)
2436  t1=tn1(1)
2437  t2=tn1(mn)
2438 ! geopotental difference from z1
2439  zg=zeta(z,z1)
2440  zgdif=zeta(z2,z1)
2441 ! set up spline nodes
2442  do 20 k=1,mn
2443  xs(k)=zeta(zn1(k),z1)/zgdif
2444  ys(k)=1./tn1(k)
2445  20 end do
2446 ! end node derivatives
2447  yd1=-tgn1(1)/(t1*t1)*zgdif
2448  yd2=-tgn1(2)/(t2*t2)*zgdif*((re+z2)/(re+z1))**2
2449 ! calculate spline coefficients
2450  call spline(xs,ys,mn,yd1,yd2,y2out)
2451  x=zg/zgdif
2452  call splint(xs,ys,y2out,mn,x,y)
2453 ! temperature at altitude
2454  tz=1./y
2455  densu=tz
2456  10 if(xm.eq.0.) go to 50
2457 !
2458 ! calculate density above za
2459  glb=gsurf/(1.+zlb/re)**2
2460  gamma=xm*glb/(s2*rgas*tinf)
2461  expl=exp(-s2*gamma*zg2)
2462  if(expl.gt.50.or.tt.le.0.) then
2463  expl=50.
2464  endif
2465 ! density at altitude
2466  densa=dlb*(tlb/tt)**(1.+alpha+gamma)*expl
2467  densu=densa
2468  if(alt.ge.za) go to 50
2469 !
2470 ! calculate density below za
2471  glb=gsurf/(1.+z1/re)**2
2472  gamm=xm*glb*zgdif/rgas
2473 ! integrate spline temperatures
2474  call splini(xs,ys,y2out,mn,x,yi)
2475  expl=gamm*yi
2476  if(expl.gt.50..or.tz.le.0.) then
2477  expl=50.
2478  endif
2479 ! density at altitude
2480  densu=densu*(t1/tz)**(1.+alpha)*exp(-expl)
2481  50 continue
2482  return
2483  end function densu
2484 
2485 !--------------------------------------------------------------------
2503  function densm(alt,d0,xm,tz,mn3,zn3,tn3,tgn3,mn2,zn2,tn2,tgn2)
2504  use gettemp_mod, only: gsurf,re
2505 !
2506  dimension zn3(mn3),tn3(mn3),tgn3(2),xs(10),ys(10),y2out(10)
2507  dimension zn2(mn2),tn2(mn2),tgn2(2)
2508 ! function
2509  zeta(zz,zl)=(zz-zl)*(re+zl)/(re+zz)
2510  save
2511  data rgas/831.4/
2512 !! rgas=831.4
2513  densm=d0
2514  if(alt.gt.zn2(1)) goto 50
2515 ! stratosphere/mesosphere temperature
2516  z=amax1(alt,zn2(mn2))
2517  mn=mn2
2518  z1=zn2(1)
2519  z2=zn2(mn)
2520  t1=tn2(1)
2521  t2=tn2(mn)
2522  zg=zeta(z,z1)
2523  zgdif=zeta(z2,z1)
2524 ! set up spline nodes
2525  do 210 k=1,mn
2526  xs(k)=zeta(zn2(k),z1)/zgdif
2527  ys(k)=1./tn2(k)
2528  210 end do
2529  yd1=-tgn2(1)/(t1*t1)*zgdif
2530  yd2=-tgn2(2)/(t2*t2)*zgdif*((re+z2)/(re+z1))**2
2531 ! calculate spline coefficients
2532  call spline(xs,ys,mn,yd1,yd2,y2out)
2533  x=zg/zgdif
2534  call splint(xs,ys,y2out,mn,x,y)
2535 ! temperature at altitude
2536  tz=1./y
2537  if(xm.eq.0.) go to 20
2538 !
2539 ! calculate stratosphere/mesosphere density
2540  glb=gsurf/(1.+z1/re)**2
2541  gamm=xm*glb*zgdif/rgas
2542 ! integrate temperature profile
2543  call splini(xs,ys,y2out,mn,x,yi)
2544  expl=gamm*yi
2545  if(expl.gt.50.) expl=50.
2546 ! density at altitude
2547  densm=densm*(t1/tz)*exp(-expl)
2548  20 continue
2549  if(alt.gt.zn3(1)) goto 50
2550 !
2551 ! troposphere/stratosphere temperature
2552  z=alt
2553  mn=mn3
2554  z1=zn3(1)
2555  z2=zn3(mn)
2556  t1=tn3(1)
2557  t2=tn3(mn)
2558  zg=zeta(z,z1)
2559  zgdif=zeta(z2,z1)
2560 ! set up spline nodes
2561  do 220 k=1,mn
2562  xs(k)=zeta(zn3(k),z1)/zgdif
2563  ys(k)=1./tn3(k)
2564  220 end do
2565  yd1=-tgn3(1)/(t1*t1)*zgdif
2566  yd2=-tgn3(2)/(t2*t2)*zgdif*((re+z2)/(re+z1))**2
2567 ! calculate spline coefficients
2568  call spline(xs,ys,mn,yd1,yd2,y2out)
2569  x=zg/zgdif
2570  call splint(xs,ys,y2out,mn,x,y)
2571 ! temperature at altitude
2572  tz=1./y
2573  if(xm.eq.0.) go to 30
2574 !
2575 ! calculate tropospheric/stratosphere density
2576 !
2577  glb=gsurf/(1.+z1/re)**2
2578  gamm=xm*glb*zgdif/rgas
2579 ! integrate temperature profile
2580  call splini(xs,ys,y2out,mn,x,yi)
2581  expl=gamm*yi
2582  if(expl.gt.50.) expl=50.
2583 ! density at altitude
2584  densm=densm*(t1/tz)*exp(-expl)
2585  30 continue
2586  50 continue
2587  if(xm.eq.0) densm=tz
2588  return
2589  end function densm
2590 
2591 !-----------------------------------------------------------------------
2604  subroutine spline(x,y,n,yp1,ypn,y2)
2605  parameter(nmax=100)
2606  dimension x(n),y(n),y2(n),u(nmax)
2607  save
2608  if(yp1.gt..99e30) then
2609  y2(1)=0
2610  u(1)=0
2611  else
2612  y2(1)=-.5
2613  u(1)=(3./(x(2)-x(1)))*((y(2)-y(1))/(x(2)-x(1))-yp1)
2614  endif
2615  do 11 i=2,n-1
2616  sig=(x(i)-x(i-1))/(x(i+1)-x(i-1))
2617  p=sig*y2(i-1)+2.
2618  y2(i)=(sig-1.)/p
2619  u(i)=(6.*((y(i+1)-y(i))/(x(i+1)-x(i))-(y(i)-y(i-1)) &
2620  /(x(i)-x(i-1)))/(x(i+1)-x(i-1))-sig*u(i-1))/p
2621  11 end do
2622  if(ypn.gt..99e30) then
2623  qn=0
2624  un=0
2625  else
2626  qn=.5
2627  un=(3./(x(n)-x(n-1)))*(ypn-(y(n)-y(n-1))/(x(n)-x(n-1)))
2628  endif
2629  y2(n)=(un-qn*u(n-1))/(qn*y2(n-1)+1.)
2630  do 12 k=n-1,1,-1
2631  y2(k)=y2(k)*y2(k+1)+u(k)
2632  12 end do
2633  return
2634  end subroutine spline
2635 
2636 !-----------------------------------------------------------------------
2648  subroutine splint(xa,ya,y2a,n,x,y)
2649  dimension xa(n),ya(n),y2a(n)
2650  save
2651  klo=1
2652  khi=n
2653  1 continue
2654  if(khi-klo.gt.1) then
2655  k=(khi+klo)/2
2656  if(xa(k).gt.x) then
2657  khi=k
2658  else
2659  klo=k
2660  endif
2661  goto 1
2662  endif
2663  h=xa(khi)-xa(klo)
2664  if(h.eq.0) write(6,*) 'bad xa input to splint'
2665  a=(xa(khi)-x)/h
2666  b=(x-xa(klo))/h
2667  y=a*ya(klo)+b*ya(khi)+ &
2668  ((a*a*a-a)*y2a(klo)+(b*b*b-b)*y2a(khi))*h*h/6.
2669  return
2670  end subroutine splint
2671 
2672 !-----------------------------------------------------------------------
2683  subroutine splini(xa,ya,y2a,n,x,yi)
2684  dimension xa(n),ya(n),y2a(n)
2685  save
2686  yi=0
2687  klo=1
2688  khi=2
2689  1 continue
2690  if(x.gt.xa(klo).and.khi.le.n) then
2691  xx=x
2692  if(khi.lt.n) xx=amin1(x,xa(khi))
2693  h=xa(khi)-xa(klo)
2694  a=(xa(khi)-xx)/h
2695  b=(xx-xa(klo))/h
2696  a2=a*a
2697  b2=b*b
2698  yi=yi+((1.-a2)*ya(klo)/2.+b2*ya(khi)/2.+ &
2699  ((-(1.+a2*a2)/4.+a2/2.)*y2a(klo)+ &
2700  (b2*b2/4.-b2/2.)*y2a(khi))*h*h/6.)*h
2701  klo=klo+1
2702  khi=khi+1
2703  goto 1
2704  endif
2705  return
2706  end subroutine splini
2707 
2708 !-----------------------------------------------------------------------
2719  function dnet(dd,dm,zhm,xmm,xm)
2720  save
2721  a=zhm/(xmm-xm)
2722  if(dm.gt.0.and.dd.gt.0) goto 5
2723  write(6,*) 'dnet log error',dm,dd,xm
2724  if(dd.eq.0.and.dm.eq.0) dd=1.
2725  if(dm.eq.0) goto 10
2726  if(dd.eq.0) goto 20
2727  5 continue
2728  ylog=a*alog(dm/dd)
2729  if(ylog.lt.-10.) go to 10
2730  if(ylog.gt.10.) go to 20
2731  dnet=dd*(1.+exp(ylog))**(1/a)
2732  go to 50
2733  10 continue
2734  dnet=dd
2735  go to 50
2736  20 continue
2737  dnet=dm
2738  go to 50
2739  50 continue
2740  return
2741  end function dnet
2742 
2743 !-----------------------------------------------------------------------
2753  function ccor(alt, r,h1,zh)
2754  save
2755  e=(alt-zh)/h1
2756  if(e.gt.70.) go to 20
2757  if(e.lt.-70.) go to 10
2758  ex=exp(e)
2759  ccor=r/(1.+ex)
2760  go to 50
2761  10 ccor=r
2762  go to 50
2763  20 ccor=0.
2764  go to 50
2765  50 continue
2766  ccor=exp(ccor)
2767  return
2768  end function ccor
2769 
2770 !-----------------------------------------------------------------------
2781  function ccor2(alt, r,h1,zh,h2)
2782  e1=(alt-zh)/h1
2783  e2=(alt-zh)/h2
2784  if(e1.gt.70. .or. e2.gt.70.) go to 20
2785  if(e1.lt.-70. .and. e2.lt.-70) go to 10
2786  ex1=exp(e1)
2787  ex2=exp(e2)
2788  ccor2=r/(1.+.5*(ex1+ex2))
2789  go to 50
2790  10 ccor2=r
2791  go to 50
2792  20 ccor2=0.
2793  go to 50
2794  50 continue
2795  ccor2=exp(ccor2)
2796  return
2797  end function ccor2
real, dimension(50) pa3
block space data for he denisity
real, dimension(25) sw
weighting
real tlb
labeled temperature
real day
day in a year
function densm(alt, d0, xm, tz, mn3, zn3, tn3, tgn3, mn2, zn2, tn2, tgn2)
Calculate temperature and density profiles for lower atmos.
real, dimension(50) px2
block space data for tgn3(2) surface grad tslg
real, dimension(50) ph3
block space data for n density
real, dimension(50) pe3
block space data for o2 density
subroutine gts7(iyd, sec, alt, glat, glong, stl, f107a, f107, ap, mass, d, t)
Thermospheric portion of nrlmsise-00.
real, dimension(50) pl1
block space data for tn1(2)
real, dimension(50) py2
block space data for tgn2(1) tgn1(2)
function globe7(yrd, sec, lat, long, tloc, f107a, f107, ap, p)
Calculate g(l) function for upper thermosphere parameters
function ccor(alt, r, h1, zh)
Chemistry/dissociation correction.
real z0
initial height
real dm01
mixed density at alt01
real, dimension(10) ptm
block space data for lower boundary
subroutine meters(meter)
Convert outputs to kg & meters if meter true.
real, dimension(50) pg1
block space data for h density
real, dimension(50) pg2
block space data for h density
real, dimension(50) pw2
block space data for tn3(5) surface temperature tsl
real c3tloc
cosine of 3 time location
real, dimension(2) tgn2
temperature gradient at node 2 (~stratosphere)
real, dimension(50) pf3
block space data for ar density
real, dimension(50) pc3
block space data for n2 density
real, dimension(10) pavgm
block space data for middle atmosphere averages
real, dimension(50) pj2
block space data for s param
real, dimension(100, 4) ptl
upper temperature
function dnet(dd, dm, zhm, xmm, xm)
Turbopause correction.
real dm04
mixed density at alt04
real, dimension(50) paa1
block space data for semiannual mult sam
real, dimension(2, 151) pr151
define pressures
real, dimension(50) pb2
block space data for o density
subroutine tselec(sv)
Set switches.
Use moduke for common blocks.
character *4, dimension(2) istime
define time
real, dimension(50) pz1
block space data for tgn3(1) tgn2(2)
real, dimension(50) po1
block space data for tn1(5) tn2(1)
real, dimension(25) swc
weighting
real, dimension(50) py1
block space data for tgn2(1) tgn1(2)
real gsurf
surface gravitation force at given latitude
subroutine ghp7(iyd, sec, alt, glat, glong, stl, f107a, f107, ap, d, t, press)
Find altitude of pressure surface (press) from gtd7.
real dm16
mixed density at alt16
real, dimension(50) pm1
block space data for tn1(3)
real db01
diffusive density at zlb for g01
real xlong
a given longitude
integer isw
indix for sw
real dm40
mixed density at alt40
function glob7s(p)
Version of globe for lower atmosphere.
character *4, dimension(3) isdate
define date
real, dimension(50) pt3
block space data for temperature
real, dimension(50) pg3
block space data for h density
real, dimension(5) tn3
temperature at node 3 (~troposphere)
real, dimension(50) pl2
block space data for tn1(2)
real, dimension(50) po2
block space data for tn1(5) tn2(1)
real, dimension(50) pk1
block space data for turbo
real, dimension(50) pb1
block space data for o density
real apdf
the same as apd
real, dimension(150) pt
temperature
real, dimension(50) ps2
block space data for tn3(2)
real dd
diffusive density at alt
real, dimension(100, 10) pma
middle and low temperature
subroutine splint(xa, ya, y2a, n, x, y)
Calculate cubic spline interp value.
real ctloc
cosine of the location
real, dimension(50) pe1
block space data for o2 density
real, dimension(50) pc1
block space data for n2 density
function vtst7(iyd, sec, glat, glong, stl, f107a, f107, ap, ic)
Test variable condition.
real db16
diffusive density at zlb for g18
subroutine glatf(lat, gv, reff)
Calculate latitude variable.
real, dimension(50) pj3
block space data for s param
real, dimension(50) pu1
block space data for tn3(3)
real, dimension(50) pn1
block space data for tn1(4)
real dm28
mixed density at alt28
real, dimension(50) paa2
block space data for semiannual mult sam
real g0
initial gradient variations
real stloc
sine of the location
real dm14
mixed density at alt14
subroutine spline(x, y, n, yp1, ypn, y2)
Calculate 2nd derivatives of cubic spline interp function.
real, dimension(50) pp1
block space data for tn2(2)
character *4, dimension(2) name
define data name
real, dimension(50) pi2
block space data for hot o density
real, dimension(25, 2) pdl
turbo
real dfa
the difference to reference value
Use moduke for blockdata gtd7bk.
real, dimension(50) ps1
block space data for tn3(2)
real apd
parameter calcumate for magnetic activity
real, dimension(100) sam
semiannual mult sam
integer iyr
integer for a given year
real dm32
mixed density at alt32
real, dimension(50) pt2
block space data for temperature
real, dimension(4) apt
daily magnetic activity
real, dimension(50) ph1
block space data for n density
integer imr
define version
real, dimension(50) pf1
block space data for ar density
real, dimension(50) pz2
block space data for tgn3(1) tgn2(2)
real, dimension(50) pj1
block space data for s param
real s
scale inverse to temperature difference
real, dimension(50) pn2
block space data for tn1(4)
real, dimension(15) tt
referenced temperature
real, dimension(5) tn1
temperature at node 1 (~mesosphere)
real re
referenced height related to gsurf
real tinfg
startinf referenced point for tt
real tr12
try factor 1 or 2
real db14
diffusive density at zlb for g14
real, dimension(4) tn2
temperature at node 2 (~stratosphere)
function densu(alt, dlb, tinf, tlb, xm, alpha, tz, zlb, s2, mn1, zn1, tn1, tgn1)
Calculate temperature and density profiles.
real, dimension(10, 8) pdm
block space data for lower boundary
subroutine gettemp(iday, nday, xlat, nlat, pr, np, temp, n_o, n_o2, n_n2)
Entry routine to get WAM needed temperature and composition profiles.
real s3tloc
sine of 3 time location
real, dimension(50) pp2
block space data for tn2(2)
real za
joining altitude of bates and spline
real, dimension(150, 9) pd
he density
real s2tloc
sine of 2 time location
real, dimension(2) tgn3
temperature gradient at node 3 (~troposphere)
real, dimension(50) pc2
block space data for n2 density
real t0
initial temperature
real, dimension(50) pd3
block space data for tlb
real, dimension(50) pb3
block space data for o density
real, dimension(2) tgn1
temperature gradient at node 1 (~mesosphere)
real db40
diffusive density at zlb for g40
real db48
diffusive density at zlb for g48
real, dimension(150) ps
s parameter
real, dimension(50) pv2
block space data for tn3(4)
subroutine gtd7d(iyd, sec, alt, glat, glong, stl, f107a, f107, ap, mass, d, t)
The nrlmsise-00 subroutine gtd7d.
real, dimension(50) pi1
block space data for hot o density
real, dimension(50) pq1
block space data for tn2(3)
real db28
diffusive density at zlb for g28
real, dimension(50) pm2
block space data for tn1(3)
real, dimension(50) px1
block space data for tgn3(2) surface grad tslg
real, dimension(50) pd2
block space data for tlb
function ccor2(alt, r, h1, zh, h2)
O and O2 chemistry/dissociation correction.
real, dimension(50) pr1
block space data for tn2(4) tn3(1)
real, dimension(50) pt1
block space data for temperature
function scalh(alt, xm, temp)
Calculate scale height (km)
real, dimension(50) pq2
block space data for tn2(3)
subroutine gtd7(iyd, sec, alt, glat, glong, stl, f107a, f107, ap, mass, d, t)
The nrlmsise-00 subroutine gtd7.
real, dimension(50) pi3
block space data for hot o density
real rl
correction to specified mixing ratio at ground
real, dimension(50) pa2
block space data for he denisity
real, dimension(2, 65) pr65
define pressures
real db04
diffusive density at zlb for g4
real, dimension(50) pa1
block space data for he denisity
real, dimension(50) pd1
block space data for tlb
subroutine splini(xa, ya, y2a, n, x, yi)
Integrate cubic spline function.
real c2tloc
cosine of 2 time location
real db32
diffusive density at zlb for g32
real, dimension(9, 4) plg
Legendre polynomial points.
real, dimension(50) pe2
block space data for o2 density
real, dimension(50) pr2
block space data for tn2(4) tn3(1)
real, dimension(50) pu2
block space data for tn3(3)
real, dimension(50) pv1
block space data for tn3(4)
real, dimension(50) pf2
block space data for ar density
real, dimension(50) ph2
block space data for n density
real, dimension(50) pw1
block space data for tn3(5) surface temperature tsl
real df
the difference of f10.7 effect