1 | #include "slalib.h"
|
---|
2 | #include "slamac.h"
|
---|
3 | void slaAopqk ( double rap, double dap, double aoprms[14],
|
---|
4 | double *aob, double *zob, double *hob,
|
---|
5 | double *dob, double *rob )
|
---|
6 | /*
|
---|
7 | ** - - - - - - - - -
|
---|
8 | ** s l a A o p q k
|
---|
9 | ** - - - - - - - - -
|
---|
10 | **
|
---|
11 | ** Quick apparent to observed place (but see note 8, below, for
|
---|
12 | ** remarks about speed).
|
---|
13 | **
|
---|
14 | ** Given:
|
---|
15 | ** rap double geocentric apparent right ascension
|
---|
16 | ** dap double geocentric apparent declination
|
---|
17 | ** aoprms double[14] star-independent apparent-to-observed parameters:
|
---|
18 | **
|
---|
19 | ** (0) geodetic latitude (radians)
|
---|
20 | ** (1,2) sine and cosine of geodetic latitude
|
---|
21 | ** (3) magnitude of diurnal aberration vector
|
---|
22 | ** (4) height (hm)
|
---|
23 | ** (5) ambient temperature (t)
|
---|
24 | ** (6) pressure (p)
|
---|
25 | ** (7) relative humidity (rh)
|
---|
26 | ** (8) wavelength (wl)
|
---|
27 | ** (9) lapse rate (tlr)
|
---|
28 | ** (10,11) refraction constants A and B (radians)
|
---|
29 | ** (12) longitude + eqn of equinoxes + sidereal DUT (radians)
|
---|
30 | ** (13) local apparent sidereal time (radians)
|
---|
31 | **
|
---|
32 | ** Returned:
|
---|
33 | ** *aob double observed azimuth (radians: N=0,E=90)
|
---|
34 | ** *zob double observed zenith distance (radians)
|
---|
35 | ** *hob double observed hour angle (radians)
|
---|
36 | ** *dob double observed declination (radians)
|
---|
37 | ** *rob double observed right ascension (radians)
|
---|
38 | **
|
---|
39 | ** Notes:
|
---|
40 | **
|
---|
41 | ** 1) This routine returns zenith distance rather than elevation
|
---|
42 | ** in order to reflect the fact that no allowance is made for
|
---|
43 | ** depression of the horizon.
|
---|
44 | **
|
---|
45 | ** 2) The accuracy of the result is limited by the corrections for
|
---|
46 | ** refraction. Providing the meteorological parameters are
|
---|
47 | ** known accurately and there are no gross local effects, the
|
---|
48 | ** observed RA,Dec predicted by this routine should be within
|
---|
49 | ** about 0.1 arcsec for a zenith distance of less than 70 degrees.
|
---|
50 | ** Even at a topocentric zenith distance of 90 degrees, the
|
---|
51 | ** accuracy in elevation should be better than 1 arcmin; useful
|
---|
52 | ** results are available for a further 3 degrees, beyond which
|
---|
53 | ** the slaRefro routine returns a fixed value of the refraction.
|
---|
54 | ** The complementary routines slaAop (or slaAopqk) and slaOap
|
---|
55 | ** (or slaOapqk) are self-consistent to better than 1 micro-
|
---|
56 | ** arcsecond all over the celestial sphere.
|
---|
57 | **
|
---|
58 | ** 3) It is advisable to take great care with units, as even
|
---|
59 | ** unlikely values of the input parameters are accepted and
|
---|
60 | ** processed in accordance with the models used.
|
---|
61 | **
|
---|
62 | ** 4) "Apparent" place means the geocentric apparent right ascension
|
---|
63 | ** and declination, which is obtained from a catalogue mean place
|
---|
64 | ** by allowing for space motion, parallax, precession, nutation,
|
---|
65 | ** annual aberration, and the Sun's gravitational lens effect. For
|
---|
66 | ** star positions in the FK5 system (i.e. J2000), these effects can
|
---|
67 | ** be applied by means of the slaMap etc routines. Starting from
|
---|
68 | ** other mean place systems, additional transformations will be
|
---|
69 | ** needed; for example, FK4 (i.e. B1950) mean places would first
|
---|
70 | ** have to be converted to FK5, which can be done with the
|
---|
71 | ** slaFk425 etc routines.
|
---|
72 | **
|
---|
73 | ** 5) "Observed" Az,El means the position that would be seen by a
|
---|
74 | ** perfect theodolite located at the observer. This is obtained
|
---|
75 | ** from the geocentric apparent RA,Dec by allowing for Earth
|
---|
76 | ** orientation and diurnal aberration, rotating from equator
|
---|
77 | ** to horizon coordinates, and then adjusting for refraction.
|
---|
78 | ** The HA,Dec is obtained by rotating back into equatorial
|
---|
79 | ** coordinates, using the geodetic latitude corrected for polar
|
---|
80 | ** motion, and is the position that would be seen by a perfect
|
---|
81 | ** equatorial located at the observer and with its polar axis
|
---|
82 | ** aligned to the Earth's axis of rotation (n.b. not to the
|
---|
83 | ** refracted pole). Finally, the RA is obtained by subtracting
|
---|
84 | ** the HA from the local apparent ST.
|
---|
85 | **
|
---|
86 | ** 6) To predict the required setting of a real telescope, the
|
---|
87 | ** observed place produced by this routine would have to be
|
---|
88 | ** adjusted for the tilt of the azimuth or polar axis of the
|
---|
89 | ** mounting (with appropriate corrections for mount flexures),
|
---|
90 | ** for non-perpendicularity between the mounting axes, for the
|
---|
91 | ** position of the rotator axis and the pointing axis relative
|
---|
92 | ** to it, for tube flexure, for gear and encoder errors, and
|
---|
93 | ** finally for encoder zero points. Some telescopes would, of
|
---|
94 | ** course, exhibit other properties which would need to be
|
---|
95 | ** accounted for at the appropriate point in the sequence.
|
---|
96 | **
|
---|
97 | ** 7) The star-independent apparent-to-observed-place parameters
|
---|
98 | ** in aoprms may be computed by means of the slaAoppa routine.
|
---|
99 | ** If nothing has changed significantly except the time, the
|
---|
100 | ** slaAoppat routine may be used to perform the requisite
|
---|
101 | ** partial recomputation of aoprms.
|
---|
102 | **
|
---|
103 | ** 8) At zenith distances beyond about 76 degrees, the need for
|
---|
104 | ** special care with the corrections for refraction causes a
|
---|
105 | ** marked increase in execution time. Moreover, the effect
|
---|
106 | ** gets worse with increasing zenith distance. Adroit
|
---|
107 | ** programming in the calling application may allow the
|
---|
108 | ** problem to be reduced. Prepare an alternative aoprms array,
|
---|
109 | ** computed for zero air-pressure; this will disable the
|
---|
110 | ** refraction corrections and cause rapid execution. Using
|
---|
111 | ** this aoprms array, a preliminary call to the present routine
|
---|
112 | ** will, depending on the application, produce a rough position
|
---|
113 | ** which may be enough to establish whether the full, slow
|
---|
114 | ** calculation (using the real aoprms array) is worthwhile.
|
---|
115 | ** For example, there would be no need for the full calculation
|
---|
116 | ** if the preliminary call had already established that the
|
---|
117 | ** source was well below the elevation limits for a particular
|
---|
118 | ** telescope.
|
---|
119 | **
|
---|
120 | ** 9) The azimuths etc produced by the present routine are with
|
---|
121 | ** respect to the celestial pole. Corrections to the terrestrial
|
---|
122 | ** pole can be computed using slaPolmo.
|
---|
123 | **
|
---|
124 | ** Called: slaDcs2c, slaRefz, slaRefro, slaDcc2s, slaDranrm
|
---|
125 | **
|
---|
126 | ** Last revision: 22 February 1996
|
---|
127 | **
|
---|
128 | ** Copyright P.T.Wallace. All rights reserved.
|
---|
129 | */
|
---|
130 | {
|
---|
131 | /*
|
---|
132 | ** Breakpoint for fast/slow refraction algorithm:
|
---|
133 | ** ZD greater than arctan(4), (see slaRefco routine)
|
---|
134 | ** or vector z less than cosine(arctan(z)) = 1/sqrt(17)
|
---|
135 | */
|
---|
136 | static double zbreak = 0.242535625;
|
---|
137 |
|
---|
138 | int i;
|
---|
139 |
|
---|
140 | double sphi, cphi, st, v[3], xhd, yhd, zhd, diurab, f,
|
---|
141 | xhdt, yhdt, zhdt, xaet, yaet, zaet, azobs,
|
---|
142 | zdt, refa, refb, zdobs, dzd, dref, ce,
|
---|
143 | xaeo, yaeo, zaeo, hmobs, dcobs, raobs;
|
---|
144 |
|
---|
145 | /* Sin, cos of latitude */
|
---|
146 | sphi = aoprms[1];
|
---|
147 | cphi = aoprms[2];
|
---|
148 |
|
---|
149 | /* Local apparent sidereal time */
|
---|
150 | st = aoprms[13];
|
---|
151 |
|
---|
152 | /* Apparent RA,Dec to Cartesian -HA,Dec */
|
---|
153 | slaDcs2c(rap - st, dap, v);
|
---|
154 | xhd = v[0];
|
---|
155 | yhd = v[1];
|
---|
156 | zhd = v[2];
|
---|
157 |
|
---|
158 | /* Diurnal aberration */
|
---|
159 | diurab = aoprms[3];
|
---|
160 | f = 1.0 - diurab * yhd;
|
---|
161 | xhdt = f * xhd;
|
---|
162 | yhdt = f * ( yhd + diurab );
|
---|
163 | zhdt = f * zhd;
|
---|
164 |
|
---|
165 | /* Cartesian -HA,Dec to Cartesian az,el (S=0,E=90) */
|
---|
166 | xaet = sphi * xhdt - cphi * zhdt;
|
---|
167 | yaet = yhdt;
|
---|
168 | zaet = cphi * xhdt + sphi * zhdt;
|
---|
169 |
|
---|
170 | /* Azimuth (N=0,E=90) */
|
---|
171 | azobs = ( (xaet == 0.0) && (yaet == 0.0) ) ?
|
---|
172 | 0.0 : atan2 ( yaet, -xaet );
|
---|
173 |
|
---|
174 | /* Topocentric zenith distance */
|
---|
175 | zdt = atan2 ( sqrt ( xaet * xaet + yaet * yaet ), zaet );
|
---|
176 |
|
---|
177 | /*
|
---|
178 | ** Refraction
|
---|
179 | ** ----------
|
---|
180 | */
|
---|
181 |
|
---|
182 | /* Fast algorithm using two constant model */
|
---|
183 | refa = aoprms[10];
|
---|
184 | refb = aoprms[11];
|
---|
185 | slaRefz ( zdt, refa, refb, &zdobs );
|
---|
186 |
|
---|
187 | /* Large zenith distance? */
|
---|
188 | if ( cos ( zdobs ) < zbreak ) {
|
---|
189 |
|
---|
190 | /* Yes: use rigorous algorithm */
|
---|
191 |
|
---|
192 | /* Initialize loop (maximum of 10 iterations) */
|
---|
193 | i = 1;
|
---|
194 | do {
|
---|
195 |
|
---|
196 | /* Compute refraction using current estimate of observed ZD */
|
---|
197 | slaRefro ( zdobs, aoprms[4], aoprms[5], aoprms[6],
|
---|
198 | aoprms[7], aoprms[8], aoprms[0],
|
---|
199 | aoprms[9], 1e-8, &dref );
|
---|
200 |
|
---|
201 | /* Remaining discrepancy */
|
---|
202 | dzd = zdobs + dref - zdt;
|
---|
203 |
|
---|
204 | /* Update the estimate */
|
---|
205 | zdobs -= dzd;
|
---|
206 |
|
---|
207 | /* Increment the iteration counter */
|
---|
208 | i++;
|
---|
209 |
|
---|
210 | } while ( fabs ( dzd ) > 1e-10 && i <= 10 );
|
---|
211 | }
|
---|
212 |
|
---|
213 | /* To Cartesian az/ZD */
|
---|
214 | ce = sin ( zdobs );
|
---|
215 | xaeo = -cos ( azobs ) * ce;
|
---|
216 | yaeo = sin ( azobs ) * ce;
|
---|
217 | zaeo = cos ( zdobs );
|
---|
218 |
|
---|
219 | /* Cartesian az/ZD to Cartesian -HA,Dec */
|
---|
220 | v[0] = sphi * xaeo + cphi * zaeo;
|
---|
221 | v[1] = yaeo;
|
---|
222 | v[2] = -cphi * xaeo + sphi * zaeo;
|
---|
223 |
|
---|
224 | /* To spherical -HA,dec */
|
---|
225 | slaDcc2s ( v, &hmobs, &dcobs );
|
---|
226 |
|
---|
227 | /* Right ascension */
|
---|
228 | raobs = slaDranrm ( st + hmobs );
|
---|
229 |
|
---|
230 | /* Return the results */
|
---|
231 | *aob = azobs;
|
---|
232 | *zob = zdobs;
|
---|
233 | *hob = -hmobs;
|
---|
234 | *dob = dcobs;
|
---|
235 | *rob = raobs;
|
---|
236 | }
|
---|