In-situ Sensing Hands-On Lesson (FORTRAN)





October 14, 2004



Overview




In this lesson you will develop a simple program illustrating how SPICE can be used to compute various kinds of geometry information applicable to the experiments carried out by an in-situ instrument flown on an interplanetary spacecraft.



References




This section provides a list of SPICE documents that are referred to in this lesson.

Of these documents, the ``Tutorials'' contains the highest level descriptions with the least number of details while the ``Required Reading'' documents contain much more detailed specifications. The most complete specifications are provided in the ``Headers'' -- the comments in the top section of the source file.

In some cases the lesson explanations also refer to the information provided in the meta-data area of the kernels used in the lesson examples. It is especially true in case of the FK and IK files, which often contain comprehensive descriptions of the frames, instrument FOVs, etc. Since both FK and IK are text kernels, the information provided in them can be viewed using any text editor, while the meta information provided in binary kernels -- SPKs and CKs -- can be viewed using ``commnt'' or ``spacit'' utility programs located in ``toolkit/exe'' of Toolkit installation tree.



Tutorials



The following SPICE tutorials are referred to by the explanation provided in this lesson:

   Name             Lesson steps/routines that it describes
   ---------------  -----------------------------------------
   Time             UTC to ET and SCLK to ET
   Loading Kernels  Loading SPICE kernels
   SCLK             SCLK to ET time conversion
   SPK              Computing positions and velocities
   Frames           Computing transformations between frames
These tutorials are available in printed form and as MS Office or PDF files from NAIF server at JPL:

   ftp://naif.jpl.nasa.gov/pub/naif/toolkit_docs/Tutorials


Required Reading Documents



The Required Reading documents are provided with the Toolkit and are located in ``toolkit/doc'' directory of the Toolkit installation trees.

   Name             Lesson steps/routines that it describes
   ---------------  -----------------------------------------
   time.req         UTC to ET time conversion
   kernel.req       Loading SPICE kernels
   sclk.req         SCLK to ET time conversion
   naif_ids.req     Body and reference frame names
   spk.req          Computing positions and velocities
Another very useful document, also distributed with the Toolkit, is ``Permuted Index'', called ``spicelib.idx'' for FORTRAN or ``cspice.idx'' for C and IDL located under ``doc'' directory in the Toolkit installation tree.

This text document provides an easy way to find what SPICE routine(s) performs a particular function of interest and the name of the source file that contains this function (this is especially useful for FORTRAN because some of the routines are entry points and, therefore, their name is different from the name of the source file in which they are located.)



Source Code Header Comments



The most detailed specification of a given SPICE FORTRAN or C routine is contained in the header section of its source code. The source code is distributed with the Toolkit and is located under ``toolkit/src/spicelib'' in FORTRAN and under ``cspice/src/cspice'' in C Toolkits.

For example the source code of the STR2ET/str2et_c routine is

   toolkit/src/spicelib/str2et.for
in the FORTRAN Toolkit and in

   cspice/src/cspice/str2et_c.c
in the C Toolkit.

Since some of the FORTRAN routines are entry points they are usually part of a source file that has different name. The ``Permuted Index'' document mentioned above can be used to locate the name of their source file.



Kernels Used




The following kernels are used in examples provided in this lesson:

   File Name                 Type Description
   ------------------------- ---- --------------------------
   naif0007.tls              LSK  Generic LSK
   cpck05Mar2004.tpc         PCK  Cassini project PCK
   cas00084.tsc              SCLK Cassini SCLK
   020514_SE_SAT105.bsp      SPK  Saturnian Satellite Ephemeris SPK
   030201AP_SK_SM546_T45.bsp SPK  Cassini Spacecraft SPK
   981005_PLTEPH-DE405S.bsp  SPK  Planetary Ephemeris SPK
   sat128.bsp                SPK  Saturnian Satellite Ephemeris SPK
   04135_04171pc_psiv2.bc    CK   Cassini Spacecraft CK
   cas_v37.tf                FK   Cassini FK
These kernels are available from the NAIF server at JPL:

   ftp://naif.jpl.nasa.gov/pub/naif/CASSINI/kernels


SPICE Routines Used




The example provided in this lesson uses the following SPICE routines:

   Name        Function that it performs
   ----------  ---------------------------------------------------
   FURNSH      Loads kernels, individually or listed in meta-kernel
   STR2ET      Converts UTC to ET
   SCS2E       Converts SCLK to ET
   SPKEZR      Computes states (position & velocity)
   SPKPOS      Computes positions
   VHAT        Find unit vector along a 3d vector
   SUBPT       Computes body-fixed coordinates of sub-observer point
   RECLAT      Converts rectangular coordinated to latitudinal
   VSUB        Subtracts 3d vectors
   PXFORM      Computes 3x3 matrix rotating vectors between frames
   MXV         Multiplies 3d vector by 3x3 matrix
   VPACK       Packs 3 number into a 3d vector
The most detailed documentation source for these routines are their headers.



Step-1: ``UTC to ET''







``UTC to ET'' Task Statement




Write a program that computes and prints Ephemeris Time (ET), expressed as the number of ephemeris seconds past J2000, that corresponds to ``2004-06-11T19:32:00'' UTC.



``UTC to ET'' Hints




Find out what SPICE kernel(s) is(are) needed to support this conversion. Look at the ``time.req'' and/or ``Time'' tutorial.

Find necessary kernel(s) on the NAIF's FTP site.

Find out what routine should be called to load necessary kernel(s). Look at the ``kernel.req'' and/or ``Loading Kernels'' tutorial.

Find the ``loader'' routine calling sequence specification. Look at the ``time.req'' and that routine's source code header. This routine may be an entry point, in which case there will be no source file with the same name. To find out in which source file this entry point is, search for its name in the ``Permuted Index''.

Find the routine(s) used to convert time between UTC and ET. Look at the ``time.req'' and/or ``Time'' tutorial.

Find the ``converter'' routine(s) calling sequence specification. Look in the ``time.req'' and the routine's source code header.

Put all calls together in a program, add variable declarations (the routine header's ``Declarations'' and ``Examples'' sections are a good place to look for declaration specification and examples) and output print statements. Compile it and link it against SPICELIB.



``UTC to ET'' Solution Steps




Only one kernel file is needed to support this conversion -- an LSK file ``naif0007.tls''.

As any other SPICE kernel this file can be loaded by the FURNSH routine (entry in ``toolkit/src/spicelib/keeper.f''). For that, the name of the file can put be provided as a sole argument of this routine:

   CHARACTER*(256)      LSKFLE
   ...
   LSKFLE = 'naif0007.tls'
   ...
   CALL FURNSH( LSKFLE )
or it can be listed in a meta-kernel:

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        )
   \begintext
the name of which, let's call it ``spice_example.mk'', can be then provided as a sole argument of the FURNSH routine:

   CHARACTER*(256)      MKFILE
   ...
   MKFILE = 'spice_example.mk'
   ...
   CALL FURNSH( MKFILE )
While the second option seems to involve a bit more work -- it requires making an extra file -- it is a much better way to go if you plan to load more kernels as you extend the program. With the meta-kernel approach simply adding more kernels to the list in KERNEL_TO_LOAD without changing the program code will accomplish that.

The highest level SPICE time routine converting UTC to ET is STR2ET (``toolkit/src/spicelib/str2et.f'').

It has two arguments -- input time string representing UTC in a variety of formats (see STR2ET header's section ``Particulars'' for the complete description of input time formats) and output DP number of ET seconds past J2000. A call to STR2ET converting a given UTC to ET could look like this:

   CHARACTER*(32)       UTC
   DOUBLE PRECISION     ET
   ...
   UTC = '2004-06-11T19:32:00'
   CALL STR2ET( UTC, ET )
By combining FURNSH and STR2ET calls and required declarations and by adding a simple print statement, one would get a complete program that prints ET for the given UTC epoch.

The program's source code then needs to be compiled and linked against SPICELIB. Assuming that the program was saved in a file called ``spice_example.f'', this can be done with the following command on a Sun workstation (the command below assumes that the library file(s) ``spicelib.a'' is located in current directory, which may not be the case):

   > f77 -C -u -o spice_example spice_example.f spicelib.a
When you run the program's executable, ``spice_example'', it will produce the following output (the output below was generated by this program compiled with g77 on a PC running Linux; your output may differ slightly in its format and numeric precision):

   > ./spice_example
    UTC  = 2004-06-11T19:32:00
    ET   =   140254384.
   >


``UTC to ET'' Code




Program ``spice_example.f'':

   IMPLICIT NONE
 
   CHARACTER*(256)      MKFILE
   CHARACTER*(32)       UTC
   DOUBLE PRECISION     ET
 
   MKFILE = 'spice_example.mk'
   CALL FURNSH( MKFILE )
 
   UTC = '2004-06-11T19:32:00'
   CALL STR2ET( UTC, ET )
 
   WRITE(*,*) 'UTC  = ', UTC
   WRITE(*,*) 'ET   = ', ET
 
   END
Meta-kernel file ``spice_example.mk'':

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        )
   \begintext


Step-2: ``SCLK to ET''







``SCLK to ET'' Task Statement




Extend the program from Step-1 to compute and print ET for the following CASSINI on-board clock epoch ``1465674964.105''.



``SCLK to ET'' Hints




Find out what additional (to those already loaded in Step-1) SPICE kernel(s) is(are) needed to support SCLK to ET conversion. Look at the ``sclk.req'' and/or ``SCLK'' tutorial.

Find necessary kernel(s) on the NAIF's FTP site.

Modify the program or meta-kernel to load this(these) kernels.

Find the routine(s) that is(are) used to convert time between SCLK and ET. Look at the ``sclk.req'' and/or ``Time'' and ``SCLK'' tutorials.

Find the ``converter'' routine's calling sequence specification. Look in the ``sclk.req'' and the routine's source code header.

Look at ``naif_ids.req'' and the comments in the additional kernel(s) that you have loaded for information on proper values of input arguments of this routine.

Add calls to the ``converter'' routine(s), necessary variable declarations (the routine header's ``Declarations'' and ``Examples'' sections are a good place to look for declaration specification and examples), and output print statements to the program. Re-compile and re-link it against SPICELIB.



``SCLK to ET'' Solution Steps




A CASSINI SCLK file is needed additionally to the LSK file loaded in the Step-1 to support this conversion.

No code change is needed in the loading portion of the program if a meta-kernel approach was used in the Step-1. The program will load the file if it will be added to the list of kernels in the KERNELS_TO_LOAD variable:

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        'kernels/sclk/cas00084.tsc'
                        )
   \begintext
The highest level SPICE routine converting SCLK to ET is SCS2E (``toolkit/src/spicelib/scs2e.f'').

It has three arguments -- NAIF ID for CASSINI s/c (-82 as described by ``naif_ids.req'' document), input time string representing CASSINI SCLK, and output DP number of ET seconds past J2000. A call to STR2ET converting given SCLK to ET could look like this:

   CHARACTER*(32)       SCLK
   INTEGER              SCID
   ...
   SCID = -82
   SCLK = '1465674964.105'
   CALL SCS2E( SCID, SCLK, ET )
By adding the SCS2E call, required declarations and a simple print statement, one would get a complete program that prints ET for the given SCLK epoch.

The program's source code then needs to be re-compiled and re-linked against SPICELIB. It can be done using the same compile command as in Step-1:

   > f77 -C -u -o spice_example spice_example.f spicelib.a
When you run the program's executable, ``spice_example'', it will produce the following output (the output below was generated by this program compiled with g77 on a PC running Linux; your output may differ slightly in its format and numeric precision):

   > ./spice_example
    UTC  = 2004-06-11T19:32:00
    ET   =   140254384.
    SCLK = 1465674964.105
    ET   =   140254384.
   >


``SCLK to ET'' Code




Program ``spice_example.f'':

   IMPLICIT NONE
 
   CHARACTER*(256)      MKFILE
   CHARACTER*(32)       UTC
   CHARACTER*(32)       SCLK
   DOUBLE PRECISION     ET
   INTEGER              SCID
 
   MKFILE = 'spice_example.mk'
   CALL FURNSH( MKFILE )
 
   UTC = '2004-06-11T19:32:00'
   CALL STR2ET( UTC, ET )
 
   WRITE(*,*) 'UTC  = ', UTC
   WRITE(*,*) 'ET   = ', ET
 
   SCID = -82
   SCLK = '1465674964.105'
   CALL SCS2E( SCID, SCLK, ET )
 
   WRITE(*,*) 'SCLK = ', SCLK
   WRITE(*,*) 'ET   = ', ET
 
   END
Meta-kernel file ``spice_example.mk'':

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        'kernels/sclk/cas00084.tsc'
                        )
   \begintext


Step-3: ``Spacecraft State''







``Spacecraft State'' Task Statement




Extend the program from Step-2 to compute geometric state -- position and velocity -- of the CASSINI spacecraft with respect to the Sun in Ecliptic frame at the epoch specified by SCLK time from Step-2.



``Spacecraft State'' Hints




Find out what additional (to those already loaded in Steps-1&2) SPICE kernel(s) is(are) needed to support state computation. Look at the ``spk.req'' and/or ``SPK'' tutorial.

Find necessary kernel(s) on the NAIF's FTP site.

Verify that the kernels contain enough data to compute the state of interest. Use ``brief'' utility program located under ``toolkit/exe'' directory for that.

Modify the meta-kernel to load this(these) kernels.

Find the routine(s) that is(are) used to compute states. Look at the ``spk.req'' and/or ``SPK'' tutorial presentation.

Find the the routine(s) calling sequence specification. Look in the ``spk.req'' and the routine's source code header.

Look at the ``naif_ids.req'' and ``frames.req'' and the routine(s) header ``Inputs'' and ``Particulars'' sections to determine proper values of the input arguments of this routine.

Add calls to the routine(s), necessary variable declarations and output print statements to the program. Re-compile and re-link it against SPICELIB.



``Spacecraft State'' Solution Steps




A CASSINI spacecraft trajectory SPK and generic planetary ephemeris SPK files are needed to support computation of the state of interest.

The file names can be added to the meta-kernel to get them loaded into the program:

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        'kernels/sclk/cas00084.tsc'
                        'kernels/spk/020514_SE_SAT105.bsp'
                        'kernels/spk/030201AP_SK_SM546_T45.bsp'
                        'kernels/spk/981005_PLTEPH-DE405S.bsp'
                        'kernels/spk/sat128.bsp'
                        )
   \begintext
The highest level SPICE routine computing states is SPKEZR (``toolkit/src/spicelib/spkezr.f'').

We are interested in computing CASSINI position and velocity with respect to the Sun, therefore the target and observer names should be set to 'CASSINI' and 'Sun' (both names can be found in ``naif_ids.req'').

The state should be in ecliptic frame, therefore the name of the frame in which the state should be computed is 'ECLIPJ2000' (see ``frames.req'' document.)

Since the geometric position is sought for, the ABCORR argument of the routine should be set to 'NONE' (see aberration correction discussion in the (``toolkit/src/spicelib/spkezr.f''). header).

Putting it all together, we get:

   CHARACTER*(32)       TARGET
   CHARACTER*(32)       FRAME
   CHARACTER*(32)       CORRTN
   CHARACTER*(32)       OBSERV
   DOUBLE PRECISION     STATE  ( 6 )
   DOUBLE PRECISION     LT
   ...
   TARGET = 'CASSINI'
   FRAME  = 'ECLIPJ2000'
   CORRTN = 'NONE'
   OBSERV = 'SUN'
   CALL SPKEZR( TARGET, ET, FRAME, CORRTN, OBSERV, STATE, LT )
The updated program with added calls, required declarations and a simple print statements produces the following output (the output below was generated by this program compiled with g77 on a PC running Linux; your output may differ slightly in its format and numeric precision):

   > ./spice_example
    UTC  = 2004-06-11T19:32:00
    ET   =   140254384.
    SCLK = 1465674964.105
    ET   =   140254384.
    STATE =  -376599062.  1.29448778E+09 -7064853.05 -5.16422619
    0.801718909   0.0406030574
   >


``Spacecraft State'' Code




Program ``spice_example.f'':

   IMPLICIT NONE
 
   CHARACTER*(256)      MKFILE
   CHARACTER*(32)       UTC
   CHARACTER*(32)       SCLK
   CHARACTER*(32)       TARGET
   CHARACTER*(32)       FRAME
   CHARACTER*(32)       CORRTN
   CHARACTER*(32)       OBSERV
   DOUBLE PRECISION     ET
   DOUBLE PRECISION     STATE  ( 6 )
   DOUBLE PRECISION     LT
   INTEGER              SCID
 
   MKFILE = 'spice_example.mk'
   CALL FURNSH( MKFILE )
 
   UTC = '2004-06-11T19:32:00'
   CALL STR2ET( UTC, ET )
 
   WRITE(*,*) 'UTC  = ', UTC
   WRITE(*,*) 'ET   = ', ET
 
   SCID = -82
   SCLK = '1465674964.105'
   CALL SCS2E( SCID, SCLK, ET )
 
   WRITE(*,*) 'SCLK = ', SCLK
   WRITE(*,*) 'ET   = ', ET
 
   TARGET = 'CASSINI'
   FRAME  = 'ECLIPJ2000'
   CORRTN = 'NONE'
   OBSERV = 'SUN'
   CALL SPKEZR( TARGET, ET, FRAME, CORRTN, OBSERV, STATE, LT )
   WRITE(*,*) 'STATE = ', STATE
 
   END
Meta-kernel file ``spice_example.mk'':

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        'kernels/sclk/cas00084.tsc'
                        'kernels/spk/020514_SE_SAT105.bsp'
                        'kernels/spk/030201AP_SK_SM546_T45.bsp'
                        'kernels/spk/981005_PLTEPH-DE405S.bsp'
                        'kernels/spk/sat128.bsp'
                        )
   \begintext


Step-4: ``Sun Direction''







``Sun Direction'' Task Statement




Extend the program from Step-3 to compute apparent direction of the Sun in the INMS frame at the epoch specified by SCLK time from Step-2.



``Sun Direction'' Hints




Find out what additional (to those already loaded in previous steps) SPICE kernel(s) is(are) needed to support the direction computation, knowing that they should provide the s/c and instrument frame orientation. Get these kernels from the NAIF's FTP site.

Verify that the orientation data in the kernels have adequate coverage to support computation of the direction of interest. Use ``ckbrief'' utility program located under ``toolkit/exe'' directory for that.

Modify the meta-kernel to load this(these) kernels.

Determine which input arguments should used in the call to SPKPOS in order to compute this direction (which is the position portion of the output state). Look through the Frames Kernel find the name of the frame to used.

Add calls to the routine(s), necessary variable declarations and output print statements to the program. Re-compile and re-link it against SPICELIB.



``Sun Direction'' Solution Steps




A CASSINI spacecraft orientation CK file, providing s/c orientation with respect to an inertial frame, and CASSINI FK file, providing orientation of the INMS frame with respect to the s/ frame, are needed additionally to already loaded kernels to support computation of this direction.

The file names can be added to the meta-kernel to get them loaded into the program:

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        'kernels/sclk/cas00084.tsc'
                        'kernels/spk/020514_SE_SAT105.bsp'
                        'kernels/spk/030201AP_SK_SM546_T45.bsp'
                        'kernels/spk/981005_PLTEPH-DE405S.bsp'
                        'kernels/spk/sat128.bsp'
                        'kernels/ck/04135_04171pc_psiv2.bc'
                        'kernels/fk/cas_v37.tf'
                        )
   \begintext
The same highest level SPICE routine computing positions, SPKPOS, can be used to compute this direction.

Since this is the direction of the Sun as seen from the s/c, the target argument should be set to 'Sun' and the observer argument should be set to 'CASSINI'. The name of the INMS frame is 'CASSINI_INMS', the definition and description of this frame are provided in the CASSINI FK file, ``cassini_v02.tf''.

Since the apparent, or ``as seen'', position is sought for, the ABCORR argument of the routine should be set to 'LT+S' (see aberration correction discussion in the (``toolkit/src/spicelib/spkpos.f'') header).

If desired, the position can then be turned into a unit vector using VHAT routine (``toolkit/src/spicelib/vhat.f''). Putting it all together, we get:

   DOUBLE PRECISION     SUNDIR ( 3 )
   ...
   TARGET = 'SUN'
   FRAME  = 'CASSINI_INMS'
   CORRTN = 'LT+S'
   OBSERV = 'CASSINI'
   CALL SPKPOS( TARGET, ET, FRAME, CORRTN, OBSERV, SUNDIR, LT )
   CALL VHAT  ( SUNDIR, SUNDIR )
The updated program with added calls, required declarations and a simple print statements produces the following output (the output below was generated by this program compiled with g77 on a PC running Linux; your output may differ slightly in its format and numeric precision):

   > ./spice_example
    UTC  = 2004-06-11T19:32:00
    ET   =   140254384.
    SCLK = 1465674964.105
    ET   =   140254384.
    STATE =  -376599062.  1.29448778E+09 -7064853.05 -5.16422619
    0.801718909   0.0406030574
    SUNDIR =  -0.290204017  0.881631188  0.372166733
   >


``Sun Direction'' Code




Program ``spice_example.f'':

   IMPLICIT NONE
 
   CHARACTER*(256)      MKFILE
   CHARACTER*(32)       UTC
   CHARACTER*(32)       SCLK
   CHARACTER*(32)       TARGET
   CHARACTER*(32)       FRAME
   CHARACTER*(32)       CORRTN
   CHARACTER*(32)       OBSERV
   DOUBLE PRECISION     ET
   DOUBLE PRECISION     STATE  ( 6 )
   DOUBLE PRECISION     SUNDIR ( 3 )
   DOUBLE PRECISION     LT
   INTEGER              SCID
 
   MKFILE = 'spice_example.mk'
   CALL FURNSH( MKFILE )
 
   UTC = '2004-06-11T19:32:00'
   CALL STR2ET( UTC, ET )
 
   WRITE(*,*) 'UTC  = ', UTC
   WRITE(*,*) 'ET   = ', ET
 
   SCID = -82
   SCLK = '1465674964.105'
   CALL SCS2E( SCID, SCLK, ET )
 
   WRITE(*,*) 'SCLK = ', SCLK
   WRITE(*,*) 'ET   = ', ET
 
   TARGET = 'CASSINI'
   FRAME  = 'ECLIPJ2000'
   CORRTN = 'NONE'
   OBSERV = 'SUN'
   CALL SPKEZR( TARGET, ET, FRAME, CORRTN, OBSERV, STATE, LT )
   WRITE(*,*) 'STATE = ', STATE
 
   TARGET = 'SUN'
   FRAME  = 'CASSINI_INMS'
   CORRTN = 'LT+S'
   OBSERV = 'CASSINI'
   CALL SPKPOS( TARGET, ET, FRAME, CORRTN, OBSERV, SUNDIR, LT )
   CALL VHAT  ( SUNDIR, SUNDIR )
   WRITE(*,*) 'SUNDIR = ', SUNDIR
   END
Meta-kernel file ``spice_example.mk'':

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        'kernels/sclk/cas00084.tsc'
                        'kernels/spk/020514_SE_SAT105.bsp'
                        'kernels/spk/030201AP_SK_SM546_T45.bsp'
                        'kernels/spk/981005_PLTEPH-DE405S.bsp'
                        'kernels/spk/sat128.bsp'
                        'kernels/ck/04135_04171pc_psiv2.bc'
                        'kernels/fk/cas_v37.tf'
                        )
   \begintext


Step-5: ``Sub-Spacecraft Point''







``Sub-Spacecraft Point'' Task Statement




Extend the program from Step-4 to compute planetocentric longitude and and latitude of the sub-spacecraft point and the direction from the spacecraft to that point in the INMS frame.



``Sub-Spacecraft Point'' Hints




Find the SPICE routine that computes sub-observer point coordinates. Use ``Permuted Index'' or ``subpt'' cookbook program for that.

Refer to the routine's header to determine which additional kernels should be loaded to support this direction computation. Get these kernels from the NAIF's FTP site. Modify the meta-kernel to load this(these) kernels.

Determine which input arguments should used in the call to this routine. Use the routine's header for that.

Convert Cartesian vector returned by this routine to latitudinal coordinates. Use ``Permuted Index'' to find the routine that does this conversion. Refer to the routine's header for input/output argument specifications.

Given that the direction from the spacecraft to the sub-spacecraft point is the difference between the sub-point position vector (which you already have) and the spacecraft position vector (which you don't have), compute the spacecraft position vector in the same frame (Phoebe body-fixed frame) using SPKPOS and subtract the two vectors. Use ``frames.req'' and/or ``Frames'' tutorial to find the name of the frame to be used in SPKPOS call.

Since the difference vector is computed in the Phoebe body-fixed frame, it should be transformed into the instrument frame get the direction we are looking for. Refer to ``frames.req'' and/or ``Frames'' tutorial to determine the name of the routine computing transformations and use it to compute transformation from Phoebe body-fixed to the INMS frame.

Using ``Permuted Index'' find the routine that multiplies 3x3 matrix by 3d vector and use it to rotate the difference vector to the instrument frame.

Add calls to the routine(s), necessary variable declarations and output print statements to the program. Re-compile and re-link it against SPICELIB.



``Sub-Spacecraft Point'' Solution Steps




The SUBPT routine (``toolkit/src/spicelib/subpt.f'') can be used to compute sub-observer point with a single call. To determine this point as the closest point on the Phoebe ellipsoid, the METHOD argument has to be set to 'Near point'. For our case the TARGET is 'PHOEBE' and the observer is 'CASSINI'.

Since the s/c is close to Phoebe, light time does not need to be taken into account and, therefore, the ABCORR argument can be set to 'NONE'.

In order for SUBPT to compute the nearest point location, a PCK file containing Phoebe radii has to be loaded into the program (see ``Files'' section of the routine's header.) All other files required for this computation are already being loaded by the program. With PCK file name added to it, the updated meta-kernel will look like this:

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        'kernels/sclk/cas00084.tsc'
                        'kernels/spk/020514_SE_SAT105.bsp'
                        'kernels/spk/030201AP_SK_SM546_T45.bsp'
                        'kernels/spk/981005_PLTEPH-DE405S.bsp'
                        'kernels/spk/sat128.bsp'
                        'kernels/ck/04135_04171pc_psiv2.bc'
                        'kernels/fk/cas_v37.tf'
                        'kernels/pck/cpck05Mar2004.tpc'
                        )
   \begintext
The sub-spacecraft point Cartesian vector can be converted to planetocentric radius, longitude and latitude using the RECLAT routine (``toolkit/src/spicelib/reclat.f'').

The other vector needed to compute direction from the spacecraft to the sub-spacecraft point -- the spacecraft position with respect to Phoebe -- can be computed using a call to SPKPOS. Note that both positions -- one of the sub-spacecraft point and the other of the spacecraft -- must be in the same frame in order to find the difference between them. Therefore, the frame in the SPKPOS call should be 'IAU_PHOEBE' which is SPICE's built-in name for Phoebe body-fixed frame. The other arguments will be 'CASSINI' (TARG), 'PHOEBE' (OBS), and 'NONE' (ABCORR), as the correction is also not essential in this case.

The VSUB routine (``toolkit/src/spicelib/vsub.f'') can be used to subtract the spacecraft position from the sub-spacecraft point position.

The computed difference vector has to be rotated from the body-fixed frame to the instrument frame. The name of the routine that computes 3x3 matrices rotating vectors from one frame to another is PXFORM (``toolkit/src/spicelib/pxform.f''). In our case the FROM argument should be set to 'IAU_PHOEBE' and the TO argument should be set to 'CASSINI_INMS' The difference vector should be then multiplied by this matrix using MXV routine (``toolkit/src/spicelib/mxv.f'') to rotate it to the instrument frame. Then, if desired, it can be unitized using the VHAT routine.

For output the longitude and latitude angles returned by RECLAT in radians can be converted to degrees by multiplying by DPR function (``toolkit/src/spicelib/dpr.f''). Putting it all together, we get:

   CHARACTER*(32)       METHOD
   CHARACTER*(32)       FROMFR
   CHARACTER*(32)       TOFR
   DOUBLE PRECISION     SPOINT ( 3 )
   DOUBLE PRECISION     ALT
   DOUBLE PRECISION     SRAD
   DOUBLE PRECISION     SLON
   DOUBLE PRECISION     SLAT
   DOUBLE PRECISION     SCPOS  ( 3 )
   DOUBLE PRECISION     SBPDIR ( 3 )
   DOUBLE PRECISION     M2IMAT ( 3, 3 )
   DOUBLE PRECISION     DPR
 
   ...
 
   METHOD = 'NEAR POINT'
   TARGET = 'PHOEBE'
   CORRTN = 'NONE'
   OBSERV = 'CASSINI'
   CALL SUBPT ( METHOD, TARGET, ET, CORRTN, OBSERV, SPOINT, ALT )
 
   CALL RECLAT( SPOINT, SRAD, SLON, SLAT )
 
   TARGET = 'CASSINI'
   FRAME  = 'IAU_PHOEBE'
   CORRTN = 'NONE'
   OBSERV = 'PHOEBE'
   CALL SPKPOS( TARGET, ET, FRAME, CORRTN, OBSERV, SCPOS, LT )
 
   CALL VSUB  ( SPOINT, SCPOS, SBPDIR )
 
   FROMFR = 'IAU_PHOEBE'
   TOFR   = 'CASSINI_INMS'
   CALL PXFORM( FROMFR, TOFR, ET, M2IMAT )
 
   CALL MXV   ( M2IMAT, SBPDIR, SBPDIR )
   CALL VHAT  ( SBPDIR, SBPDIR )
 
   WRITE(*,*) 'ALT = ', ALT
   WRITE(*,*) 'LON = ', SLON*DPR()
   WRITE(*,*) 'LAT = ', SLAT*DPR()
   WRITE(*,*) 'SBPDIR = ', SBPDIR
The updated program with added calls, required declarations and a simple print statements produces the following output (the output below was generated by this program compiled with g77 on a PC running Linux; your output may differ slightly in its format and numeric precision):

   > ./spice_example
    UTC  = 2004-06-11T19:32:00
    ET   =   140254384.
    SCLK = 1465674964.105
    ET   =   140254384.
    STATE =  -376599062.  1.29448778E+09 -7064853.05 -5.16422619
    0.801718909  0.0406030574
    SUNDIR =  -0.290204017  0.881631188  0.372166733
    ALT =   2084.11288
    LON =   23.4231584
    LAT =   3.70979716
    SBPDIR =  -0.000776207056 -0.999873199 -0.015905456
   >


``Sub-Spacecraft Point'' Code




Program ``spice_example.f'':

   IMPLICIT NONE
 
   CHARACTER*(256)      MKFILE
   CHARACTER*(32)       UTC
   CHARACTER*(32)       SCLK
   CHARACTER*(32)       TARGET
   CHARACTER*(32)       FRAME
   CHARACTER*(32)       CORRTN
   CHARACTER*(32)       OBSERV
   CHARACTER*(32)       METHOD
   CHARACTER*(32)       FROMFR
   CHARACTER*(32)       TOFR
   DOUBLE PRECISION     ET
   DOUBLE PRECISION     STATE  ( 6 )
   DOUBLE PRECISION     SUNDIR ( 3 )
   DOUBLE PRECISION     LT
   DOUBLE PRECISION     SPOINT ( 3 )
   DOUBLE PRECISION     ALT
   DOUBLE PRECISION     SRAD
   DOUBLE PRECISION     SLON
   DOUBLE PRECISION     SLAT
   DOUBLE PRECISION     SCPOS  ( 3 )
   DOUBLE PRECISION     SBPDIR ( 3 )
   DOUBLE PRECISION     M2IMAT ( 3, 3 )
   INTEGER              SCID
 
   DOUBLE PRECISION     DPR
 
 
   MKFILE = 'spice_example.mk'
   CALL FURNSH( MKFILE )
 
   UTC = '2004-06-11T19:32:00'
   CALL STR2ET( UTC, ET )
 
   WRITE(*,*) 'UTC  = ', UTC
   WRITE(*,*) 'ET   = ', ET
 
   SCID = -82
   SCLK = '1465674964.105'
   CALL SCS2E( SCID, SCLK, ET )
 
   WRITE(*,*) 'SCLK = ', SCLK
   WRITE(*,*) 'ET   = ', ET
 
   TARGET = 'CASSINI'
   FRAME  = 'ECLIPJ2000'
   CORRTN = 'NONE'
   OBSERV = 'SUN'
   CALL SPKEZR( TARGET, ET, FRAME, CORRTN, OBSERV, STATE, LT )
   WRITE(*,*) 'STATE = ', STATE
 
   TARGET = 'SUN'
   FRAME  = 'CASSINI_INMS'
   CORRTN = 'LT+S'
   OBSERV = 'CASSINI'
   CALL SPKPOS( TARGET, ET, FRAME, CORRTN, OBSERV, SUNDIR, LT )
   CALL VHAT  ( SUNDIR, SUNDIR )
   WRITE(*,*) 'SUNDIR = ', SUNDIR
 
   METHOD = 'NEAR POINT'
   TARGET = 'PHOEBE'
   CORRTN = 'NONE'
   OBSERV = 'CASSINI'
   CALL SUBPT ( METHOD, TARGET, ET, CORRTN, OBSERV, SPOINT, ALT )
 
   CALL RECLAT( SPOINT, SRAD, SLON, SLAT )
 
   TARGET = 'CASSINI'
   FRAME  = 'IAU_PHOEBE'
   CORRTN = 'NONE'
   OBSERV = 'PHOEBE'
   CALL SPKPOS( TARGET, ET, FRAME, CORRTN, OBSERV, SCPOS, LT )
 
   CALL VSUB  ( SPOINT, SCPOS, SBPDIR )
 
   FROMFR = 'IAU_PHOEBE'
   TOFR   = 'CASSINI_INMS'
   CALL PXFORM( FROMFR, TOFR, ET, M2IMAT )
 
   CALL MXV   ( M2IMAT, SBPDIR, SBPDIR )
   CALL VHAT  ( SBPDIR, SBPDIR )
 
   WRITE(*,*) 'ALT = ', ALT
   WRITE(*,*) 'LON = ', SLON*DPR()
   WRITE(*,*) 'LAT = ', SLAT*DPR()
   WRITE(*,*) 'SBPDIR = ', SBPDIR
   END
Meta-kernel file ``spice_example.mk'':

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        'kernels/sclk/cas00084.tsc'
                        'kernels/spk/020514_SE_SAT105.bsp'
                        'kernels/spk/030201AP_SK_SM546_T45.bsp'
                        'kernels/spk/981005_PLTEPH-DE405S.bsp'
                        'kernels/spk/sat128.bsp'
                        'kernels/ck/04135_04171pc_psiv2.bc'
                        'kernels/fk/cas_v37.tf'
                        'kernels/pck/cpck05Mar2004.tpc'
                        )
   \begintext


Step-6: ``Spacecraft Velocity''







``Spacecraft Velocity'' Task Statement




Extend the program from Step-5 to compute the spacecraft velocity with respect to Phoebe in the INMS frame.



``Spacecraft Velocity'' Hints




Compute velocity of the spacecraft with respect to Phoebe in some inertial frame, for example J2000. Recall that velocity is the last three components of the state vector returned by SPKEZR.

Since the velocity vector is computed in the inertial frame, it should be rotated to the instrument frame. Look at the previous step the routine that compute necessary rotation and rotate vectors.

Add calls to the routine(s), necessary variable declarations and output print statements to the program. Re-compile and re-link it against SPICELIB.



``Spacecraft Velocity'' Solution Steps




All kernels required for computations in this step are already being loaded by the program, therefore, the meta-kernel does not need to be changed.

The spacecraft velocity vector is the last three components of the state returned by SPKEZR. To compute velocity of CASSINI with respect to Phoebe in the J2000 inertial frame the SPKEZR arguments should be set to 'CASSINI' (TARG), 'PHOEBE' (OBS), 'J2000' (REF) and 'NONE' (ABCORR).

For convenience the velocity can be copied from the output state in to a 3d vector using the VPACK routine (``toolkit/src/spicelib/vpack.f'').

The computed velocity vector has to be rotated from the J2000 frame to the instrument frame. The PXFORM routine used in the previous step can be used to compute the rotation matrix needed for that. In this case the frame name arguments should be set to 'J2000' (FROM) and 'CASSINI_INMS' (TO).

As in the previous step the difference vector should be then multiplied by this matrix using the MXV routine. Then, if desired, it can be unitized using the VHAT routine.

Putting it all together, we get:

   DOUBLE PRECISION     SCVDIR ( 3 )
   DOUBLE PRECISION     J2IMAT( 3, 3 )
 
   ...
 
   TARGET = 'CASSINI'
   FRAME  = 'J2000'
   CORRTN = 'NONE'
   OBSERV = 'PHOEBE'
   CALL SPKEZR( TARGET, ET, FRAME, CORRTN, OBSERV, STATE, LT )
   CALL VPACK ( STATE(4), STATE(5), STATE(6), SCVDIR )
 
   FROMFR = 'J2000'
   TOFR   = 'CASSINI_INMS'
   CALL PXFORM( FROMFR, TOFR, ET, J2IMAT )
 
   CALL MXV   ( J2IMAT, SCVDIR, SCVDIR )
   CALL VHAT  ( SCVDIR, SCVDIR )
 
   WRITE(*,*) 'SCVDIR = ', SCVDIR
The updated program with added calls, required declarations and a simple print statements produces the following output (the output below was generated by this program compiled with g77 on a PC running Linux; your output may differ slightly in its format and numeric precision):

   > ./spice_example
    UTC  = 2004-06-11T19:32:00
    ET   =   140254384.
    SCLK = 1465674964.105
    ET   =   140254384.
    STATE =  -376599062.  1.29448778E+09 -7064853.05 -5.16422619
    0.801718909  0.0406030574
    SUNDIR =  -0.290204017  0.881631188  0.372166733
    ALT =   2084.11288
    LON =   23.4231584
    LAT =   3.70979716
    SBPDIR =  -0.000776207056 -0.999873199 -0.015905456
    SCVDIR =   0.395784873 -0.292807673  0.870412546
   >
Note that computing the spacecraft velocity in the instrument frame by a single call to SPKEZR by specifying 'CASSINI_INMS' in the REF argument will produce an incorrect result. Such computation will take into account the spacecraft angular velocity from the CK files, which should not be considered in this case.



``Spacecraft Velocity'' Code Program ``spice_example.f'':




   IMPLICIT NONE
 
   CHARACTER*(256)      MKFILE
   CHARACTER*(32)       UTC
   CHARACTER*(32)       SCLK
   CHARACTER*(32)       TARGET
   CHARACTER*(32)       FRAME
   CHARACTER*(32)       CORRTN
   CHARACTER*(32)       OBSERV
   CHARACTER*(32)       METHOD
   CHARACTER*(32)       FROMFR
   CHARACTER*(32)       TOFR
   DOUBLE PRECISION     ET
   DOUBLE PRECISION     STATE  ( 6 )
   DOUBLE PRECISION     SUNDIR ( 3 )
   DOUBLE PRECISION     LT
   DOUBLE PRECISION     SPOINT ( 3 )
   DOUBLE PRECISION     ALT
   DOUBLE PRECISION     SRAD
   DOUBLE PRECISION     SLON
   DOUBLE PRECISION     SLAT
   DOUBLE PRECISION     SCPOS  ( 3 )
   DOUBLE PRECISION     SBPDIR ( 3 )
   DOUBLE PRECISION     M2IMAT ( 3, 3 )
   DOUBLE PRECISION     SCVDIR ( 3 )
   DOUBLE PRECISION     J2IMAT ( 3, 3 )
   INTEGER              SCID
 
   DOUBLE PRECISION     DPR
 
 
   MKFILE = 'spice_example.mk'
   CALL FURNSH( MKFILE )
 
   UTC = '2004-06-11T19:32:00'
   CALL STR2ET( UTC, ET )
 
   WRITE(*,*) 'UTC  = ', UTC
   WRITE(*,*) 'ET   = ', ET
 
   SCID = -82
   SCLK = '1465674964.105'
   CALL SCS2E( SCID, SCLK, ET )
 
   WRITE(*,*) 'SCLK = ', SCLK
   WRITE(*,*) 'ET   = ', ET
 
   TARGET = 'CASSINI'
   FRAME  = 'ECLIPJ2000'
   CORRTN = 'NONE'
   OBSERV = 'SUN'
   CALL SPKEZR( TARGET, ET, FRAME, CORRTN, OBSERV, STATE, LT )
   WRITE(*,*) 'STATE = ', STATE
 
   TARGET = 'SUN'
   FRAME  = 'CASSINI_INMS'
   CORRTN = 'LT+S'
   OBSERV = 'CASSINI'
   CALL SPKPOS( TARGET, ET, FRAME, CORRTN, OBSERV, SUNDIR, LT )
   CALL VHAT  ( SUNDIR, SUNDIR )
   WRITE(*,*) 'SUNDIR = ', SUNDIR
 
   METHOD = 'NEAR POINT'
   TARGET = 'PHOEBE'
   CORRTN = 'NONE'
   OBSERV = 'CASSINI'
   CALL SUBPT ( METHOD, TARGET, ET, CORRTN, OBSERV, SPOINT, ALT )
 
   CALL RECLAT( SPOINT, SRAD, SLON, SLAT )
 
   TARGET = 'CASSINI'
   FRAME  = 'IAU_PHOEBE'
   CORRTN = 'NONE'
   OBSERV = 'PHOEBE'
   CALL SPKPOS( TARGET, ET, FRAME, CORRTN, OBSERV, SCPOS, LT )
 
   CALL VSUB  ( SPOINT, SCPOS, SBPDIR )
 
   FROMFR = 'IAU_PHOEBE'
   TOFR   = 'CASSINI_INMS'
   CALL PXFORM( FROMFR, TOFR, ET, M2IMAT )
 
   CALL MXV   ( M2IMAT, SBPDIR, SBPDIR )
   CALL VHAT  ( SBPDIR, SBPDIR )
 
   WRITE(*,*) 'ALT = ', ALT
   WRITE(*,*) 'LON = ', SLON*DPR()
   WRITE(*,*) 'LAT = ', SLAT*DPR()
   WRITE(*,*) 'SBPDIR = ', SBPDIR
 
   TARGET = 'CASSINI'
   FRAME  = 'J2000'
   CORRTN = 'NONE'
   OBSERV = 'PHOEBE'
   CALL SPKEZR( TARGET, ET, FRAME, CORRTN, OBSERV, STATE, LT )
   CALL VPACK ( STATE(4), STATE(5), STATE(6), SCVDIR )
 
   FROMFR = 'J2000'
   TOFR   = 'CASSINI_INMS'
   CALL PXFORM( FROMFR, TOFR, ET, J2IMAT )
 
   CALL MXV   ( J2IMAT, SCVDIR, SCVDIR )
   CALL VHAT  ( SCVDIR, SCVDIR )
 
   WRITE(*,*) 'SCVDIR = ', SCVDIR
 
   END
Meta-kernel file ``spice_example.mk'':

   \begindata
      KERNELS_TO_LOAD = (
                        'kernels/lsk/naif0007.tls'
                        'kernels/sclk/cas00084.tsc'
                        'kernels/spk/020514_SE_SAT105.bsp'
                        'kernels/spk/030201AP_SK_SM546_T45.bsp'
                        'kernels/spk/981005_PLTEPH-DE405S.bsp'
                        'kernels/spk/sat128.bsp'
                        'kernels/ck/04135_04171pc_psiv2.bc'
                        'kernels/fk/cas_v37.tf'
                        'kernels/pck/cpck05Mar2004.tpc'
                        )
   \begintext