| 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 | }
|
|---|