1 | #include "slalib.h"
|
---|
2 | #include "slamac.h"
|
---|
3 | void slaRefz ( double zu, double refa, double refb, double *zr )
|
---|
4 | /*
|
---|
5 | ** - - - - - - - -
|
---|
6 | ** s l a R e f z
|
---|
7 | ** - - - - - - - -
|
---|
8 | **
|
---|
9 | ** Adjust an unrefracted zenith distance to include the effect of
|
---|
10 | ** atmospheric refraction, using the simple A tan z + B tan^3 z
|
---|
11 | ** model.
|
---|
12 | **
|
---|
13 | ** Given:
|
---|
14 | ** zu double unrefracted zenith distance of the source (radian)
|
---|
15 | ** refa double A: tan z coefficient (radian)
|
---|
16 | ** refb double B: tan^3 z coefficient (radian)
|
---|
17 | **
|
---|
18 | ** Returned:
|
---|
19 | ** *zr double refracted zenith distance (radian)
|
---|
20 | **
|
---|
21 | ** Notes:
|
---|
22 | **
|
---|
23 | ** 1 This routine applies the adjustment for refraction in the
|
---|
24 | ** opposite sense to the usual one - it takes an unrefracted
|
---|
25 | ** (in vacuo) position and produces an observed (refracted)
|
---|
26 | ** position, whereas the A tan Z + B tan^3 Z model strictly
|
---|
27 | ** applies to the case where an observed position is to have the
|
---|
28 | ** refraction removed. The unrefracted to refracted case is
|
---|
29 | ** harder, and requires an inverted form of the text-book
|
---|
30 | ** refraction models; the formula used here is based on the
|
---|
31 | ** Newton-Raphson method. For the utmost numerical consistency
|
---|
32 | ** with the refracted to unrefracted model, two iterations are
|
---|
33 | ** carried out, achieving agreement at the 1D-11 arcseconds level
|
---|
34 | ** for a ZD of 80 degrees. The inherent accuracy of the model
|
---|
35 | ** is, of course, far worse than this - see the documentation for
|
---|
36 | ** slaRefco for more information.
|
---|
37 | **
|
---|
38 | ** 2 At ZD 83 degrees, the rapidly-worsening A tan Z + B tan^3 Z
|
---|
39 | ** model is abandoned and an empirical formula takes over. Over a
|
---|
40 | ** wide range of observer heights and corresponding temperatures and
|
---|
41 | ** pressures, the following levels of accuracy (arcsec) are achieved,
|
---|
42 | ** relative to numerical integration through a model atmosphere:
|
---|
43 | **
|
---|
44 | ** ZR error
|
---|
45 | **
|
---|
46 | ** 80 0.4
|
---|
47 | ** 81 0.8
|
---|
48 | ** 82 1.5
|
---|
49 | ** 83 3.2
|
---|
50 | ** 84 4.9
|
---|
51 | ** 85 5.8
|
---|
52 | ** 86 6.1
|
---|
53 | ** 87 7.1
|
---|
54 | ** 88 10
|
---|
55 | ** 89 20
|
---|
56 | ** 90 40
|
---|
57 | ** 91 100 } relevant only to
|
---|
58 | ** 92 200 } high-elevation sites
|
---|
59 | **
|
---|
60 | ** The high-ZD model is scaled to match the normal model at the
|
---|
61 | ** transition point; there is no glitch.
|
---|
62 | **
|
---|
63 | **
|
---|
64 | ** 3 Beyond 93 deg zenith distance, the refraction is held at its
|
---|
65 | ** 93 deg value.
|
---|
66 | **
|
---|
67 | ** 4 See also the routine slaRefv, which performs the adjustment in
|
---|
68 | ** Cartesian Az/El coordinates, and with the emphasis on speed
|
---|
69 | ** rather than numerical accuracy.
|
---|
70 | **
|
---|
71 | ** Defined in slamac.h: DR2D
|
---|
72 | **
|
---|
73 | ** Last revision: 4 June 1997
|
---|
74 | **
|
---|
75 | ** Copyright P.T.Wallace. All rights reserved.
|
---|
76 | */
|
---|
77 | {
|
---|
78 | double zu1, zl, s, c, t, tsq, tcu, ref, e, e2;
|
---|
79 |
|
---|
80 | /* Coefficients for high ZD model (used beyond ZD 83 deg */
|
---|
81 | const double c1 = 0.55445,
|
---|
82 | c2 = -0.01133,
|
---|
83 | c3 = 0.00202,
|
---|
84 | c4 = 0.28385,
|
---|
85 | c5 = 0.02390;
|
---|
86 |
|
---|
87 | /* Largest usable ZD (deg) */
|
---|
88 | const double z93 = 93.0;
|
---|
89 |
|
---|
90 | /* ZD at which one model hands over to the other (radians) */
|
---|
91 | const double z83 = 83.0 / DR2D;
|
---|
92 |
|
---|
93 | /* High-ZD-model prediction (deg) for that point */
|
---|
94 | const double ref83 = ( c1 + c2 * 7.0 + c3 * 49.0 ) /
|
---|
95 | ( 1.0 + c4 * 7.0 + c5 * 49.0 );
|
---|
96 |
|
---|
97 |
|
---|
98 | /* Perform calculations for zu or 83 deg, whichever is smaller */
|
---|
99 | zu1 = gmin ( zu, z83 );
|
---|
100 |
|
---|
101 | /* Functions of ZD */
|
---|
102 | zl = zu1;
|
---|
103 | s = sin ( zl );
|
---|
104 | c = cos ( zl );
|
---|
105 | t = s / c;
|
---|
106 | tsq = t * t;
|
---|
107 | tcu = t * tsq;
|
---|
108 |
|
---|
109 | /* Refracted ZD (mathematically to better than 1mas at 70 deg) */
|
---|
110 | zl -= ( refa * t + refb * tcu )
|
---|
111 | / ( 1.0 + ( refa + 3.0 * refb * tsq ) / ( c * c ) );
|
---|
112 |
|
---|
113 | /* Further iteration */
|
---|
114 | s = sin ( zl );
|
---|
115 | c = cos ( zl );
|
---|
116 | t = s / c;
|
---|
117 | tsq = t * t;
|
---|
118 | tcu = t * tsq;
|
---|
119 | ref = zu1 - zl + ( zl - zu1 + refa * t + refb * tcu )
|
---|
120 | / ( 1.0 + ( refa + 3.0 * refb * tsq ) / ( c * c ) );
|
---|
121 |
|
---|
122 | /* Special handling for large zu */
|
---|
123 | if ( zu > zu1 ) {
|
---|
124 | e = 90.0 - gmin ( z93, zu * DR2D );
|
---|
125 | e2 = e * e;
|
---|
126 | ref = ( ref / ref83 ) * ( c1 + c2 * e + c3 * e2 ) /
|
---|
127 | ( 1.0 + c4 * e + c5 * e2 );
|
---|
128 | }
|
---|
129 |
|
---|
130 | /* Refracted ZD */
|
---|
131 | *zr = zu - ref;
|
---|
132 | }
|
---|