coord.f

Go to the documentation of this file.

      SUBROUTINE coord(X,GM,Y)
C  VERSION OF 4/1/85
C  PURPOSE
C    TRANSFORM ORBITAL ELEMENT A,E,I,CAPW,W,EA TO CARTESIAN COORD
C    X,Y,Z,XD,ZD,YD.
C  INPUT
C    X(1)   = A, SEMI-MAJOR AXIS
C     (2)   = E, ECCENTRICITY
C     (3)   = I, INCLINATION (RAD)
C     (4)   = CAPW, LONGITUDE OF ASCENDING NODE (RAD)
C     (5)   = W, ARGUMENT OF PERIAPSIS (RAD)
C     (6)   = EA, ECCENTRIC ANOMALY (RAD)
C    GM     = GRAVITATIONAL CONSTANT * MASS OF PLANET (KM**3/SEC**2)
C  OUTPUT
C    Y(1)   = X (KM)
C     (2)   = Y (KM)
C     (3)   = Z (KM)
C     (4)   = XD (KM/SEC)
C     (5)   = YD (KM/SEC)
C     (6)   = ZD (KM/SEC)
C  CALL SUBROUTINES
C    NONE
C  REFERENCES
C    JPL EM 312/87-153, 20 APRIL 1987
C    ORBITAL MOTION, A. E. ROY, 1978, PAGE 102
C  ANALYSIS
C    J. H. KWOK - JPL
C  PROGRAMMER
C    J. H. KWOK - JPL
C  PROGRAM MODIFICATIONS
C    NONE
C
      IMPLICIT DOUBLE PRECISION (a-h,o-z)
      dimension x(6),y(6),td(8)
      DATA one/1.d0/
      cc=dcos(x(4))
      sc=dsin(x(4))
      cw=dcos(x(5))
      sw=dsin(x(5))
      ci=dcos(x(3))
      si=dsin(x(3))
      td(1)=cc*cw-sc*sw*ci
      td(2)=sc*cw+cc*sw*ci
      td(3)=sw*si
      td(4)=-cc*sw-sc*cw*ci
      td(5)=-sc*sw+cc*cw*ci
      td(6)=cw*si
      td(7)=dsqrt(one-x(2)**2)*x(1)
      td(8)=dsqrt(gm/x(1)**3)
      ce=dcos(x(6))
      se=dsin(x(6))
      r=x(1)*(one-x(2)*ce)
      DO 10 i=1,3
      y(i)=x(1)*td(i)*(ce-x(2))+td(7)*td(i+3)*se
   10 y(i+3)=x(1)*td(8)*(td(7)*ce*td(i+3)-x(1)*se*td(i))/r
      RETURN
      END