C $Id: ic_gci_to_geo.f,v 1.2 1998/07/24 21:54:29 asc Exp $
CCCC
C
C  IC_GCI_TO_GEO - return a transformation matrix
C 
C  PURPOSE:  Calculate the transformation matrix from GCI 
C            coordinates to GEO coordinates at a given date and time.
C
C  UNIT TYPE:  SUBROUTINE
C
C  INVOCATION METHOD:  CALL IC_GCI_TO_GEO (orb_pos_time, 
C                                          transform_matrix,
C                                          deriv_matrix,
C                                          velocity_req)
C
C  ARGUMENT LIST:
C
C  NAME	                  TYPE   USE  DESCRIPTION
C  ----                   ----   ---  -----------
C  ORB_POS_TIME(2)        I*4    I    TIME OF ORB. VECTOR, YEAR-DAY-MILLI OF DAY
C  TRANSFORM_MATRIX(3,3)  R*8    O    TRANSFORMATION MATRIX 
C  DERIV_MATRIX(3,3)      R*8    O    TIME DERIVATIVE OF TRANSFORMATION 
C                                      MATRIX.
C  VELOCITY_REQ           I*4    I    FLAG TO SHOW WHETHER THE DERIV_MATRIX
C                                      IS NEEDED.  IF =0, THEN DERIV_MATRIX
C                                      WILL BE ZERO-FILLED.  IF >0, THEN 
C                                      DERIV_MATRIX WILL HAVE NON-ZERO 
C                                      VALUES IN IT.
C
C  FILE/RECORD REFERENCES:  NONE
C
C  EXTERNAL VARIABLES:  NONE
C
C  EXTERNAL REFERENCES:
C	IC_GRNWCH_SIDEREAL - Returns the Greenwich sidereal time in radians
C	MULMAT - Multiplies two matrices
C	IC_CONV_MATRIX - Returns the conversion matrix to rotate
C                        from mean of date to true of date
C
C  ABNORMAL TERMINATION CONDITIONS, ERROR MESSAGES:  NONE
C
C  ASSUMPTIONS, CONSTRAINTS, RESTRICTIONS:  NONE
C
C  DEVELOPMENT HISTORY
C
C  AUTHOR	CHANGE ID	RELEASE	  DATE	    DESCRIPTION OF CHANGE
C  ------	---------	-------   ----	    ---------------------
C  J. LUBELCZYK                 B1R1      11/21/90  INITIAL PDL
C  J. Lubelczyk			B1R1      12/10/90  CODING
C  J. Lubelczyk                           09/18/91  Updated to return true
C                                                   of date trans matrix
C  J. LUBELCZYK ICCR #83, CCR #'S 130, 137 11/91    B3 update
C  C. RAYMOND   CCR 462                   07/12/93  Added code needed for
C                                                    velocity 
C                                                    transformations.
C  B. SAMUELSON SPOF PORT       NONE      04/11/94  Change: required to port
C     (CSC)     (see notes)                          icss routines to SPOF
C
C  NOTES:
C
C     1.  The changes recorded under ID SPOF-PORT are required to make the
C         ICSS coordinate conversion routines, originally developed under
C         VAX-VMS 5.4 run on the UNIX-based workstations of the SPOF.  (Sun
C         SPARCstations and DEC DECstations).  The changes are as follows:
C            a.  Delete references to ICSS_INC
C            b.  Define Message texts or files to correspond to the messages
C                ICSS_SUCCESSFUL, etc. which are embedded in the error handling.
C            c.  Remove references to the NAG routines F01CRF and F01CKF (matrix
C                transposition and matrix multiplication routines).
C     2.  In addition, to successfully run the software packages, copies of the
C         Solar/Lunar/Planetary (SLP) file and timing coefficients file (TCC)
C         must be ported onto the SPOF.
C
CCCC
C
C  PDL:
C
C  call IC_GRNWCH_SIDEREAL to get greenwich sidereal time in radians
C
C  calculate the sin and cos of the greenwich sidereal time
C
C  Fill the mean_matrix using the sin and cos
C   of the greenwich mean sidereal time
C
C  call IC_CONV_MATRIX to get the conversion matrix, CMATRIX, 
C                      to rotate to true of date
C
C  call MULMAT to multiple MEAN_MATRIX by CMATRIX to get the transform matrix
C
C  IF (velocity_req .gt. 0 )
C
C  THEN
C     
C     LOAD D_MATRIX with values for derivative
C
C     CALL MULMAT to post-multiply D_MATRIX by CMATRIX to get final DERIV_
C                                                                   MATRIX
C  
C  ELSE
C
C     ZERO out DERIV_MATRIX
C
C  ENDIF 
C
C  RETURN
C
CCCC
	subroutine IC_GCI_TO_GEO (ORB_POS_TIME, TRANSFORM_MATRIX,
	1                         DERIV_MATRIX, VELOCITY_REQ )
C
	implicit none
C
C*  Calling Parameters
C
	integer*4   ORB_POS_TIME(2)  !time of orb vector, YYYYDDD, milli of day
	real*8	    TRANSFORM_MATRIX(3,3) !transformation matrix
	real*8      DERIV_MATRIX(3,3) !time derivative of TRANSFORM_MATRIX
        integer*4   VELOCITY_REQ      !flag to show value of DERIV_MATRIX
C
C*  Other Variables
C
	real*8	GRNWCH_SIDEREAL_TIME  !greenwich mean sidereal time in degrees
	real*8	SIN_GST		      !double precision sin of GST in degrees
	real*8	COS_GST		      !double precision cos of GST in degrees
	real*8	MEAN_MATRIX(3,3)      !mean of date transform matrix
	real*8	CMATRIX(3,3)	      !The conversion matrix to rotate from 
C                                     !mean of date to true of date
	real*8  D_MATRIX(3,3)         !Time derivative of matrix used to go
C                                     !from true of date inertial to 
C                                     !geographic coordinates (MEAN_MATRIX)
        real*8  EARTH_ROT_RATE        !Rotation rate of the Earth, hard-
C                                     !coded to J2000 value.
	integer*4   i,j               !Do-loop indices
C
	data EARTH_ROT_RATE/7.2921151467D-05/ !Units=radians/second.
C
C*  Begin executable code
C
C*  call IC_GRNWCH_SIDEREAL to get greenwich mean sidereal time in radians
C
	call IC_GRNWCH_SIDEREAL (ORB_POS_TIME, GRNWCH_SIDEREAL_TIME)
C
C*  calculate the sin and cos of the greenwich mean sidereal time
C
	SIN_GST = dsin(GRNWCH_SIDEREAL_TIME)
	COS_GST = dcos(GRNWCH_SIDEREAL_TIME)
C
C*  Fill the mean of date transformation matrix using the sin and cos
C*   of the greenwich mean sidereal time
C
	MEAN_MATRIX(1,1) = COS_GST
	MEAN_MATRIX(2,1) = -SIN_GST
	MEAN_MATRIX(3,1) = 0.0d0
	MEAN_MATRIX(1,2) = SIN_GST
	MEAN_MATRIX(2,2) = COS_GST
	MEAN_MATRIX(3,2) = 0.0d0
	MEAN_MATRIX(1,3) = 0.0d0
	MEAN_MATRIX(2,3) = 0.0d0
	MEAN_MATRIX(3,3) = 1.0d0
C
	CALL IC_CONV_MATRIX (ORB_POS_TIME, CMATRIX)
C
	CALL MULMAT (MEAN_MATRIX, CMATRIX, 3, 3, 3, TRANSFORM_MATRIX)
C
	IF (VELOCITY_REQ .GT. 0 ) THEN
           D_MATRIX(1,1) = -SIN_GST*EARTH_ROT_RATE
	   D_MATRIX(1,2) =  COS_GST*EARTH_ROT_RATE
	   D_MATRIX(1,3) =  0.0d0
	   D_MATRIX(2,1) = -D_MATRIX(1,2)
	   D_MATRIX(2,2) =  D_MATRIX(1,1)
	   D_MATRIX(2,3) =  0.0d0
	   D_MATRIX(3,1) =  0.0d0
	   D_MATRIX(3,2) =  0.0d0
	   D_MATRIX(3,3) =  0.0d0
	   CALL MULMAT (D_MATRIX, CMATRIX, 3, 3, 3, DERIV_MATRIX)
	ELSE
	   DO 100 I=1,3
	   DO 100 J=1,3
	      DERIV_MATRIX(I,J)=0.0d0
100	   CONTINUE
	ENDIF
C
	RETURN
C
	END






