FACT++  1.0
palRefz.c
Go to the documentation of this file.
1 /*
2 *+
3 * Name:
4 * palRefz
5 
6 * Purpose:
7 * Adjust unrefracted zenith distance
8 
9 * Language:
10 * Starlink ANSI C
11 
12 * Type of Module:
13 * Library routine
14 
15 * Invocation:
16 * void palRefz ( double zu, double refa, double refb, double *zr );
17 
18 * Arguments:
19 * zu = double (Given)
20 * Unrefracted zenith distance of the source (radians)
21 * refa = double (Given)
22 * tan Z coefficient (radians)
23 * refb = double (Given)
24 * tan**3 Z coefficient (radian)
25 * zr = double * (Returned)
26 * Refracted zenith distance (radians)
27 
28 * Description:
29 * Adjust an unrefracted zenith distance to include the effect of
30 * atmospheric refraction, using the simple A tan Z + B tan**3 Z
31 * model (plus special handling for large ZDs).
32 
33 * Authors:
34 * PTW: Patrick T. Wallace
35 * TIMJ: Tim Jenness (JAC, Hawaii)
36 * {enter_new_authors_here}
37 
38 * Notes:
39 * - This routine applies the adjustment for refraction in the
40 * opposite sense to the usual one - it takes an unrefracted
41 * (in vacuo) position and produces an observed (refracted)
42 * position, whereas the A tan Z + B tan**3 Z model strictly
43 * applies to the case where an observed position is to have the
44 * refraction removed. The unrefracted to refracted case is
45 * harder, and requires an inverted form of the text-book
46 * refraction models; the formula used here is based on the
47 * Newton-Raphson method. For the utmost numerical consistency
48 * with the refracted to unrefracted model, two iterations are
49 * carried out, achieving agreement at the 1D-11 arcseconds level
50 * for a ZD of 80 degrees. The inherent accuracy of the model
51 * is, of course, far worse than this - see the documentation for
52 * palRefco for more information.
53 *
54 * - At ZD 83 degrees, the rapidly-worsening A tan Z + B tan^3 Z
55 * model is abandoned and an empirical formula takes over. For
56 * optical/IR wavelengths, over a wide range of observer heights and
57 * corresponding temperatures and pressures, the following levels of
58 * accuracy (arcsec, worst case) are achieved, relative to numerical
59 * integration through a model atmosphere:
60 *
61 * ZR error
62 *
63 * 80 0.7
64 * 81 1.3
65 * 82 2.4
66 * 83 4.7
67 * 84 6.2
68 * 85 6.4
69 * 86 8
70 * 87 10
71 * 88 15
72 * 89 30
73 * 90 60
74 * 91 150 } relevant only to
75 * 92 400 } high-elevation sites
76 *
77 * For radio wavelengths the errors are typically 50% larger than
78 * the optical figures and by ZD 85 deg are twice as bad, worsening
79 * rapidly below that. To maintain 1 arcsec accuracy down to ZD=85
80 * at the Green Bank site, Condon (2004) has suggested amplifying
81 * the amount of refraction predicted by palRefz below 10.8 deg
82 * elevation by the factor (1+0.00195*(10.8-E_t)), where E_t is the
83 * unrefracted elevation in degrees.
84 *
85 * The high-ZD model is scaled to match the normal model at the
86 * transition point; there is no glitch.
87 *
88 * - Beyond 93 deg zenith distance, the refraction is held at its
89 * 93 deg value.
90 *
91 * - See also the routine palRefv, which performs the adjustment in
92 * Cartesian Az/El coordinates, and with the emphasis on speed
93 * rather than numerical accuracy.
94 
95 * References:
96 * Condon,J.J., Refraction Corrections for the GBT, PTCS/PN/35.2,
97 * NRAO Green Bank, 2004.
98 
99 * History:
100 * 2012-08-24 (TIMJ):
101 * Initial version, ported directly from Fortran SLA
102 * Adapted with permission from the Fortran SLALIB library.
103 * {enter_further_changes_here}
104 
105 * Copyright:
106 * Copyright (C) 2004 Rutherford Appleton Laboratory
107 * Copyright (C) 2012 Science and Technology Facilities Council.
108 * All Rights Reserved.
109 
110 * Licence:
111 * This program is free software; you can redistribute it and/or
112 * modify it under the terms of the GNU General Public License as
113 * published by the Free Software Foundation; either version 3 of
114 * the License, or (at your option) any later version.
115 *
116 * This program is distributed in the hope that it will be
117 * useful, but WITHOUT ANY WARRANTY; without even the implied
118 * warranty of MERCHANTABILITY or FITNESS FOR A PARTICULAR
119 * PURPOSE. See the GNU General Public License for more details.
120 *
121 * You should have received a copy of the GNU General Public License
122 * along with this program; if not, write to the Free Software
123 * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston,
124 * MA 02110-1301, USA.
125 
126 * Bugs:
127 * {note_any_bugs_here}
128 *-
129 */
130 
131 #include <math.h>
132 
133 #include "pal.h"
134 #include "palmac.h"
135 
136 void palRefz ( double zu, double refa, double refb, double *zr ) {
137 
138  /* Constants */
139 
140  /* Largest usable ZD (deg) */
141  const double D93 = 93.0;
142 
143  /* ZD at which one model hands over to the other (radians) */
144  const double Z83 = 83.0 * PAL__DD2R;
145 
146  /* coefficients for high ZD model (used beyond ZD 83 deg) */
147  const double C1 = +0.55445;
148  const double C2 = -0.01133;
149  const double C3 = +0.00202;
150  const double C4 = +0.28385;
151  const double C5 = +0.02390;
152 
153  /* High-ZD-model prefiction (deg) for that point */
154  const double REF83 = (C1+C2*7.0+C3*49.0)/(1.0+C4*7.0+C5*49.0);
155 
156  double zu1,zl,s,c,t,tsq,tcu,ref,e,e2;
157 
158  /* perform calculations for zu or 83 deg, whichever is smaller */
159  zu1 = DMIN(zu,Z83);
160 
161  /* functions of ZD */
162  zl = zu1;
163  s = sin(zl);
164  c = cos(zl);
165  t = s/c;
166  tsq = t*t;
167  tcu = t*tsq;
168 
169  /* refracted zd (mathematically to better than 1 mas at 70 deg) */
170  zl = zl-(refa*t+refb*tcu)/(1.0+(refa+3.0*refb*tsq)/(c*c));
171 
172  /* further iteration */
173  s = sin(zl);
174  c = cos(zl);
175  t = s/c;
176  tsq = t*t;
177  tcu = t*tsq;
178  ref = zu1-zl+
179  (zl-zu1+refa*t+refb*tcu)/(1.0+(refa+3.0*refb*tsq)/(c*c));
180 
181  /* special handling for large zu */
182  if (zu > zu1) {
183  e = 90.0-DMIN(D93,zu*PAL__DR2D);
184  e2 = e*e;
185  ref = (ref/REF83)*(C1+C2*e+C3*e2)/(1.0+C4*e+C5*e2);
186  }
187 
188  /* return refracted zd */
189  *zr = zu-ref;
190 
191 }
#define DMIN(A, B)
Definition: palmac.h:129
static const double PAL__DD2R
Definition: palmac.h:72
static const double PAL__DR2D
Definition: palmac.h:81
TT t
Definition: test_client.c:26
void palRefz(double zu, double refa, double refb, double *zr)
Definition: palRefz.c:136