aboutsummaryrefslogtreecommitdiff
path: root/math/slalib/precss.f.sav
blob: c18707dc8222ab88d40341164c566009f8e86d98 (plain) (blame)
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
SUBROUTINE slPRCS (SYSTEM, EP0, EP1, RA, DC)
*+
*     - - - - - - -
*      P R C E
*     - - - - - - -
*
*  Precession - either FK4 (Bessel-Newcomb, pre IAU 1976) or
*  FK5 (Fricke, post IAU 1976) as required.
*
*  Given:
*     SYSTEM     int    precession to be applied: 1 = FK4 or 2 = FK5
*     EP0,EP1    dp     starting and ending epoch
*     RA,DC      dp     RA,Dec, mean equator & equinox of epoch EP0
*
*  Returned:
*     RA,DC      dp     RA,Dec, mean equator & equinox of epoch EP1
*
*  Called:    slDA2P, slPRBN, slPREC, slDS2C,
*             slDMXV, slDC2S
*
*  Notes:
*
*     1)  Lowercase characters in SYSTEM are acceptable.
*
*     2)  The epochs are Besselian if SYSTEM=FK4 and Julian if FK5.
*         For example, to precess coordinates in the old system from
*         equinox 1900.0 to 1950.0 the call would be:
*             CALL slPRCS (1, 1900D0, 1950D0, RA, DC)
*
*     3)  This routine will NOT correctly convert between the old and
*         the new systems - for example conversion from B1950 to J2000.
*         For these purposes see slFK45, slFK54, slF45Z and
*         slF54Z.
*
*     4)  If an invalid SYSTEM is supplied, values of -99D0,-99D0 will
*         be returned for both RA and DC.
*
*  P.T.Wallace   Starlink   20 April 1990
*-

      IMPLICIT NONE

      INTEGER SYSTEM
      DOUBLE PRECISION EP0,EP1,RA,DC

      DOUBLE PRECISION PM(3,3),V1(3),V2(3)

      DOUBLE PRECISION slDA2P


*  Convert to uppercase and validate SYSTEM
      IF (SYSTEM.NE.1.AND.SYSTEM.NE.2) THEN
         RA=-99D0
         DC=-99D0
      ELSE

*     Generate appropriate precession matrix
         IF (SYSTEM.EQ.1) THEN
            CALL slPRBN(EP0,EP1,PM)
         ELSE
            CALL slPREC(EP0,EP1,PM)
         END IF

*     Convert RA,Dec to x,y,z
         CALL slDS2C(RA,DC,V1)

*     Precess
         CALL slDMXV(PM,V1,V2)

*     Back to RA,Dec
         CALL slDC2S(V2,RA,DC)
         RA=slDA2P(RA)

      END IF

      END