source: trunk/MagicSoft/slalib/refcoq.c@ 19550

Last change on this file since 19550 was 731, checked in by tbretz, 24 years ago
*** empty log message ***
  • Property svn:executable set to *
File size: 7.8 KB
Line 
1#include "slalib.h"
2#include "slamac.h"
3void slaRefcoq ( double tdk, double pmb, double rh, double wl,
4 double *refa, double *refb )
5/*
6** - - - - - - - - - -
7** s l a R e f c o q
8** - - - - - - - - - -
9**
10** Determine the constants A and B in the atmospheric refraction
11** model dZ = A tan Z + B tan^3 Z. This is a fast alternative
12** to the slaRefco routine - see notes.
13**
14** Z is the "observed" zenith distance (i.e. affected by refraction)
15** and dZ is what to add to Z to give the "topocentric" (i.e. in vacuo)
16** zenith distance.
17**
18** Given:
19** tdk double ambient temperature at the observer (deg K)
20** pmb double pressure at the observer (millibar)
21** rh double relative humidity at the observer (range 0-1)
22** wl double effective wavelength of the source (micrometre)
23**
24** Returned:
25** refa double* tan Z coefficient (radian)
26** refb double* tan^3 Z coefficient (radian)
27**
28** The radio refraction is chosen by specifying WL > 100 micrometres.
29**
30** Notes:
31**
32** 1 The model is an approximation, for moderate zenith distances,
33** to the predictions of the slaRefro routine. The approximation
34** is maintained across a range of conditions, and applies to
35** both optical/IR and radio.
36**
37** 2 The algorithm is a fast alternative to the slaRefco routine.
38** The latter calls the slaRefro routine itself: this involves
39** integrations through a model atmosphere, and is costly in
40** processor time. However, the model which is produced is precisely
41** correct for two zenith distance (45 degrees and about 76 degrees)
42** and at other zenith distances is limited in accuracy only by the
43** A tan Z + B tan^3 Z formulation itself. The present routine
44** is not as accurate, though it satisfies most practical
45** requirements.
46**
47** 3 The model omits the effects of (i) height above sea level (apart
48** from the reduced pressure itself), (ii) latitude (i.e. the
49** flattening of the Earth) and (iii) variations in tropospheric
50** lapse rate.
51**
52** The model was tested using the following range of conditions:
53**
54** lapse rates 0.0055, 0.0065, 0.0075 deg/metre
55** latitudes 0, 25, 50, 75 degrees
56** heights 0, 2500, 5000 metres ASL
57** pressures mean for height -10% to +5% in steps of 5%
58** temperatures -10 deg to +20 deg with respect to 280 deg at SL
59** relative humidity 0, 0.5, 1
60** wavelengths 0.4, 0.6, ... 2 micron, + radio
61** zenith distances 15, 45, 75 degrees
62**
63** The accuracy with respect to direct use of the slaRefro routine
64** was as follows:
65**
66** worst RMS
67**
68** optical/IR 62 mas 8 mas
69** radio 319 mas 49 mas
70**
71** For this particular set of conditions:
72**
73** lapse rate 0.0065 degK/metre
74** latitude 50 degrees
75** sea level
76** pressure 1005 mB
77** temperature 280.15 degK
78** humidity 80%
79** wavelength 5740 Angstroms
80**
81** the results were as follows:
82**
83** ZD slaRefro slaRefcoq Saastamoinen
84**
85** 10 10.27 10.27 10.27
86** 20 21.19 21.20 21.19
87** 30 33.61 33.61 33.60
88** 40 48.82 48.83 48.81
89** 45 58.16 58.18 58.16
90** 50 69.28 69.30 69.27
91** 55 82.97 82.99 82.95
92** 60 100.51 100.54 100.50
93** 65 124.23 124.26 124.20
94** 70 158.63 158.68 158.61
95** 72 177.32 177.37 177.31
96** 74 200.35 200.38 200.32
97** 76 229.45 229.43 229.42
98** 78 267.44 267.29 267.41
99** 80 319.13 318.55 319.10
100**
101** deg arcsec arcsec arcsec
102**
103** The values for Saastamoinen's formula (which includes terms
104** up to tan^5) are taken from Hohenkerk and Sinclair (1985).
105**
106** The results from the much slower but more accurate slaRefco
107** routine have not been included in the tabulation as they are
108** identical to those in the slaRefro column to the 0.01 arcsec
109** resolution used.
110**
111** 4 Outlandish input parameters are silently limited to mathematically
112** safe values. Zero pressure is permissible, and causes zeroes to
113** be returned.
114**
115** 5 The algorithm draws on several sources, as follows:
116**
117** a) The formula for the saturation vapour pressure of water as
118** a function of temperature and temperature is taken from
119** expressions A4.5-A4.7 of Gill (1982).
120**
121** b) The formula for the water vapour pressure, Given the
122** saturation pressure and the relative humidity, is from
123** Crane (1976), expression 2.5.5.
124**
125** c) The refractivity of air is a function of temperature,
126** total pressure, water-vapour pressure and, in the case
127** of optical/IR but not radio, wavelength. The formulae
128** for the two cases are developed from the Essen and Froome
129** expressions adopted in Resolution 1 of the 12th International
130** Geodesy Association General Assembly (1963).
131**
132** The above three items are as used in the slaRefro routine.
133**
134** d) The formula for beta, the ratio of the scale height of the
135** atmosphere to the geocentric distance of the observer, is
136** an adaption of expression 9 from Stone (1996). The
137** adaptations, arrived at empirically, consist of (i) a
138** small adjustment to the coefficient and (ii) a humidity
139** term for the radio case only.
140**
141** e) The formulae for the refraction constants as a function of
142** n-1 and beta are from Green (1987), expression 4.31.
143**
144** References:
145**
146** Crane, R.K., Meeks, M.L. (ed), "Refraction Effects in the Neutral
147** Atmosphere", Methods of Experimental Physics: Astrophysics 12B,
148** Academic Press, 1976.
149**
150** Gill, Adrian E., "Atmosphere-Ocean Dynamics", Academic Press, 1982.
151**
152** Hohenkerk, C.Y., & Sinclair, A.T., NAO Technical Note No. 63, 1985.
153**
154** International Geodesy Association General Assembly, Bulletin
155** Geodesique 70 p390, 1963.
156**
157** Stone, Ronald C., P.A.S.P. 108 1051-1058, 1996.
158**
159** Green, R.M., "Spherical Astronomy", Cambridge University Press, 1987.
160**
161** Last revision: 17 March 1999
162**
163** Copyright P.T.Wallace. All rights reserved.
164*/
165{
166 int optic;
167 double t, p, r,w, tdc, ps, pw, wlsq, gamma, beta;
168
169
170/* Decide whether optical/IR or radio case: switch at 100 microns. */
171 optic = ( wl <= 100.0 );
172
173/* Restrict parameters to safe values. */
174 t = gmax ( tdk, 100.0 );
175 t = gmin ( t, 500.0 );
176 p = gmax ( pmb, 0.0 );
177 p = gmin ( p, 10000.0 );
178 r = gmax ( rh, 0.0 );
179 r = gmin ( r, 1.0 );
180 w = gmax ( wl, 0.1 );
181 w = gmin ( w, 1e6 );
182
183/* Water vapour pressure at the observer. */
184 if ( p > 0.0 ) {
185 tdc = t - 273.15;
186 ps = pow ( 10.0, ( 0.7859 + 0.03477 * tdc ) /
187 ( 1.0 + 0.00412 * tdc ) ) *
188 ( 1.0 + p * ( 4.5e-6 + 6e-10 * tdc * tdc ) );
189 pw = r * ps / ( 1.0 - ( 1.0 - r ) * ps / p );
190 } else {
191 pw = 0.0;
192 }
193
194/* Refractive index minus 1 at the observer. */
195 if ( optic ) {
196 wlsq = wl * wl;
197 gamma = ( ( 77.532e-6 + ( 4.391e-7 + 3.57e-9 / wlsq ) / wlsq ) * p
198 - 11.2684e-6 * pw ) / t;
199 } else {
200 gamma = ( 77.624e-6 * p - ( 12.92e-6 - 0.371897 / t ) * pw ) / t;
201 }
202
203/* Formula for beta from Stone, with empirical adjustments. */
204 beta = 4.4474e-6 * t;
205 if ( !optic ) beta -= 0.0074 * pw * beta;
206
207/* Refraction constants from Green. */
208 *refa = gamma * ( 1.0 - beta );
209 *refb = - gamma * ( beta - gamma / 2.0 );
210}
Note: See TracBrowser for help on using the repository browser.