1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
|
SUBROUTINE sla_DH2E (AZ, EL, PHI, HA, DEC)
*+
* - - - - -
* D E 2 H
* - - - - -
*
* Horizon to equatorial coordinates: Az,El to HA,Dec
*
* (double precision)
*
* Given:
* AZ d azimuth
* EL d elevation
* PHI d observatory latitude
*
* Returned:
* HA d hour angle
* DEC d declination
*
* Notes:
*
* 1) All the arguments are angles in radians.
*
* 2) The sign convention for azimuth is north zero, east +pi/2.
*
* 3) HA is returned in the range +/-pi. Declination is returned
* in the range +/-pi/2.
*
* 4) The latitude is (in principle) geodetic. In critical
* applications, corrections for polar motion should be applied.
*
* 5) In some applications it will be important to specify the
* correct type of elevation in order to produce the required
* type of HA,Dec. In particular, it may be important to
* distinguish between the elevation as affected by refraction,
* which will yield the "observed" HA,Dec, and the elevation
* in vacuo, which will yield the "topocentric" HA,Dec. If the
* effects of diurnal aberration can be neglected, the
* topocentric HA,Dec may be used as an approximation to the
* "apparent" HA,Dec.
*
* 6) No range checking of arguments is done.
*
* 7) In applications which involve many such calculations, rather
* than calling the present routine it will be more efficient to
* use inline code, having previously computed fixed terms such
* as sine and cosine of latitude.
*
* P.T.Wallace Starlink 21 February 1996
*
* Copyright (C) 1996 Rutherford Appleton Laboratory
*-
IMPLICIT NONE
DOUBLE PRECISION AZ,EL,PHI,HA,DEC
DOUBLE PRECISION SA,CA,SE,CE,SP,CP,X,Y,Z,R
* Useful trig functions
SA=SIN(AZ)
CA=COS(AZ)
SE=SIN(EL)
CE=COS(EL)
SP=SIN(PHI)
CP=COS(PHI)
* HA,Dec as x,y,z
X=-CA*CE*SP+SE*CP
Y=-SA*CE
Z=CA*CE*CP+SE*SP
* To HA,Dec
R=SQRT(X*X+Y*Y)
IF (R.EQ.0D0) THEN
HA=0D0
ELSE
HA=ATAN2(Y,X)
END IF
DEC=ATAN2(Z,R)
END
|
