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
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
|
SUBROUTINE sla_PLANTE (DATE, ELONG, PHI, JFORM, EPOCH,
: ORBINC, ANODE, PERIH, AORQ, E,
: AORL, DM, RA, DEC, R, JSTAT)
*+
* - - - - - - -
* P L A N T E
* - - - - - - -
*
* Topocentric apparent RA,Dec of a Solar-System object whose
* heliocentric orbital elements are known.
*
* Given:
* DATE d MJD of observation (JD - 2400000.5)
* ELONG d observer's east longitude (radians)
* PHI d observer's geodetic latitude (radians)
* JFORM i choice of element set (1-3; Note 4)
* EPOCH d epoch of elements (TT MJD)
* ORBINC d inclination (radians)
* ANODE d longitude of the ascending node (radians)
* PERIH d longitude or argument of perihelion (radians)
* AORQ d mean distance or perihelion distance (AU)
* E d eccentricity
* AORL d mean anomaly or longitude (radians, JFORM=1,2 only)
* DM d daily motion (radians, JFORM=1 only )
*
* Returned:
* RA,DEC d RA, Dec (topocentric apparent, radians)
* R d distance from observer (AU)
* JSTAT i status: 0 = OK
* -1 = illegal JFORM
* -2 = illegal E
* -3 = illegal AORQ
* -4 = illegal DM
* -5 = numerical error
*
* Notes:
*
* 1 DATE is the instant for which the prediction is required. It is
* in the TT timescale (formerly Ephemeris Time, ET) and is a
* Modified Julian Date (JD-2400000.5).
*
* 2 The longitude and latitude allow correction for geocentric
* parallax. This is usually a small effect, but can become
* important for Earth-crossing asteroids. Geocentric positions
* can be generated by appropriate use of routines sla_EVP and
* sla_PLANEL.
*
* 3 The elements are with respect to the J2000 ecliptic and equinox.
*
* 4 Three different element-format options are available:
*
* Option JFORM=1, suitable for the major planets:
*
* EPOCH = epoch of elements (TT MJD)
* ORBINC = inclination i (radians)
* ANODE = longitude of the ascending node, big omega (radians)
* PERIH = longitude of perihelion, curly pi (radians)
* AORQ = mean distance, a (AU)
* E = eccentricity, e
* AORL = mean longitude L (radians)
* DM = daily motion (radians)
*
* Option JFORM=2, suitable for minor planets:
*
* EPOCH = epoch of elements (TT MJD)
* ORBINC = inclination i (radians)
* ANODE = longitude of the ascending node, big omega (radians)
* PERIH = argument of perihelion, little omega (radians)
* AORQ = mean distance, a (AU)
* E = eccentricity, e
* AORL = mean anomaly M (radians)
*
* Option JFORM=3, suitable for comets:
*
* EPOCH = epoch of perihelion (TT MJD)
* ORBINC = inclination i (radians)
* ANODE = longitude of the ascending node, big omega (radians)
* PERIH = argument of perihelion, little omega (radians)
* AORQ = perihelion distance, q (AU)
* E = eccentricity, e
*
* 5 Unused elements (DM for JFORM=2, AORL and DM for JFORM=3) are
* not accessed.
*
* Called: sla_GMST, sla_DT, sla_EPJ, sla_PVOBS, sla_PRENUT,
* sla_PLANEL, sla_DMXV, sla_DCC2S, sla_DRANRM
*
* P.T.Wallace Starlink 17 March 1999
*
* Copyright (C) 1999 Rutherford Appleton Laboratory
*-
IMPLICIT NONE
DOUBLE PRECISION DATE,ELONG,PHI
INTEGER JFORM
DOUBLE PRECISION EPOCH,ORBINC,ANODE,PERIH,AORQ,E,
: AORL,DM,RA,DEC,R
INTEGER JSTAT
* Light time for unit distance (sec)
DOUBLE PRECISION TAU
PARAMETER (TAU=499.004782D0)
INTEGER I
DOUBLE PRECISION DVB(3),DPB(3),VSG(6),VSP(6),V(6),RMAT(3,3),
: VGP(6),STL,VGO(6),DX,DY,DZ,D,TL
DOUBLE PRECISION sla_GMST,sla_DT,sla_EPJ,sla_DRANRM
* Sun to geocentre (J2000).
CALL sla_EVP(DATE,2000D0,DVB,DPB,VSG(4),VSG)
* Sun to planet (J2000).
CALL sla_PLANEL(DATE,JFORM,EPOCH,ORBINC,ANODE,PERIH,AORQ,
: E,AORL,DM,VSP,JSTAT)
* Geocentre to planet (J2000).
DO I=1,6
V(I)=VSP(I)-VSG(I)
END DO
* Precession and nutation to date.
CALL sla_PRENUT(2000D0,DATE,RMAT)
CALL sla_DMXV(RMAT,V,VGP)
CALL sla_DMXV(RMAT,V(4),VGP(4))
* Geocentre to observer (date).
STL=sla_GMST(DATE-sla_DT(sla_EPJ(DATE))/86400D0)+ELONG
CALL sla_PVOBS(PHI,0D0,STL,VGO)
* Observer to planet (date).
DO I=1,6
V(I)=VGP(I)-VGO(I)
END DO
* Geometric distance (AU).
DX=V(1)
DY=V(2)
DZ=V(3)
D=SQRT(DX*DX+DY*DY+DZ*DZ)
* Light time (sec).
TL=TAU*D
* Correct position for planetary aberration
DO I=1,3
V(I)=V(I)-TL*V(I+3)
END DO
* To RA,Dec.
CALL sla_DCC2S(V,RA,DEC)
RA=sla_DRANRM(RA)
R=D
END
|
