Table of contents
CSPICE_SXFORM returns the state transformation matrix from one
frame to another at a specified epoch.
Given:
from the scalar string name of a reference frame in which a state is
known.
help, from
STRING = Scalar
to the scalar string name of a reference frame in which it is
desired to represent the state.
help, to
STRING = Scalar
et the double precision scalar or N-vector of epochs in ephemeris
seconds past the epoch of J2000 (TDB) at which the state
transformation matrix should be evaluated.
help, et
DOUBLE = Scalar or DOUBLE = Array[N]
the call:
cspice_sxform, from, to, et, xform
returns:
xform a double precision, 6x6 or 6x6xN array state transformation
matrix that transforms states from the reference frame `from' to
the frame `to' at epoch `et'.
help, xform
DOUBLE = Array[6,6] or DOUBLE = Array[6,6,N]
Recall the difference in CSPICE and IDL matrix representations.
None.
Any numerical results shown for this example may differ between
platforms as the results depend on the SPICE kernels used as input
and the machine specific arithmetic implementation.
1) Suppose you have geodetic coordinates of a station on the
surface of Earth and that you need the inertial (J2000)
state of this station. The following code fragment
illustrates how to transform the geodetic state of the
station to a J2000 state.
Use the meta-kernel shown below to load the required SPICE
kernels.
KPL/MK
File name: sxform_ex1.tm
This meta-kernel is intended to support operation of SPICE
example programs. The kernels shown here should not be
assumed to contain adequate or correct versions of data
required by SPICE-based user applications.
In order for an application to use this meta-kernel, the
kernels referenced here must be present in the user's
current working directory.
The names and contents of the kernels referenced
by this meta-kernel are as follows:
File name Contents
--------- --------
pck00008.tpc Planet orientation and
radii
naif0009.tls Leapseconds
\begindata
KERNELS_TO_LOAD = ( 'pck00008.tpc',
'naif0009.tls' )
\begintext
End of meta-kernel
Example code begins here.
PRO sxform_ex1
;;
;; Load the PCK and LSK kernels.
;;
cspice_furnsh, 'sxform_ex1.tm'
;;
;; Define a geodetic longitude, latitude, altitude
;; coordinate set. These coordinates are defined in the
;; non-inertial, earth fixed frame "IAU_EARTH".
;;
lon = 118.25d * cspice_rpd()
lat = 34.05d * cspice_rpd()
alt = 0.d
;;
;; Define a UTC time of interest. Convert the 'utc' string
;; to ephemeris time J2000.
;;
utc = 'January 1, 1990'
cspice_str2et, utc, et
;;
;; Retrieve the equatorial and polar axis of the earth (body 399).
;;
cspice_bodvrd, 'EARTH', 'RADII', 3, abc
equatr = abc[0]
polar = abc[2]
;;
;; Calculate the flattening factor for earth.
;;
f = ( equatr - polar ) / equatr
;;
;; Calculate the Cartesian coordinates on earth for the
;; location at 'lon', 'lat', 'alt'.
;;
cspice_georec, lon, lat, alt, equatr, f, estate
;;
;; cspice_georec returned the position vector of the geodetic
;; coordinates, but we want the state vector. Since it is a fixed
;; location referenced in the "IAU_EARTH" frame, the location has
;; no velocity. We need to extend estate to a 6-vector, the final
;; three elements with value 0.d.
;;
estate = [ temporary(estate) , dblarr(3) ]
;;
;; Retrieve the transformation matrix from "IAU_EARTH"
;; to "J2000" at epoch 'et'.
;;
cspice_sxform, 'IAU_EARTH', 'J2000', et, xform
;;
;; Transform the Cartesian state vector from "IAU_EARTH"
;; to "J2000."
;;
jstate = transpose(xform) # estate
cspice_et2utc, et, 'C', 3, utcstr
print, 'Epoch : ', utcstr
print, FORMAT='("Position in J2000 frame (km): ",3F11.4)', $
jstate[0:2]
print, FORMAT='("Velocity in J2000 frame (km/s): ",3F11.4)', $
jstate[3:5]
print
;;
;; Return the state transformation matrices from "IAU_EARTH"
;; to "J2000" approximately every three months for the time
;; interval February 1, 1990 to February 1, 1991 (UTC).
;;
;;
;; Define the time bounds for the time interval,
;; 1 year, convert to ephemeris time J2000.
;;
utc_bounds = [ '1 Feb 1990', '1 Feb 1991' ]
cspice_str2et, utc_bounds, et_bounds
;;
;; Step in units of three months. 1 years -> 4 steps.
;;
step = (et_bounds[1] - et_bounds[0]) / 4.d
;;
;; Create an array of 4 ephemeris times starting at
;; et_bounds[0] in intervals of 'step'.
;;
et = dindgen(4)*step + et_bounds[0]
;;
;; Convert the 4-vector of 'et' to an array of corresponding
;; transformation matrices (dimensions (6,6,4) ).
;;
cspice_sxform, 'IAU_EARTH', 'J2000', et, xform
;;
;; Transform the Cartesian state vector from "IAU_EARTH"
;; to "J2000".
;;
cspice_et2utc, et, 'C', 3, utcstr
for i = 0, 3 do begin
jstate = transpose( xform(*,*,i) ) # estate
print, 'Epoch : ', utcstr(i)
print, FORMAT='("Position in J2000 frame (km): ",3F11.4)', $
jstate[0:2]
print, FORMAT='("Velocity in J2000 frame (km/s): ",3F11.4)', $
jstate[3:5]
print
endfor
;;
;; It's always good form to unload kernels after use,
;; particularly in IDL due to data persistence.
;;
cspice_kclear
END
When this program was executed on a Mac/Intel/IDL8.x/64-bit
platform, the output was:
Epoch : 1990 JAN 01 00:00:00.000
Position in J2000 frame (km): -4131.4630 -3308.3707 3547.0215
Velocity in J2000 frame (km/s): 0.2412 -0.3010 0.0002
Epoch : 1990 FEB 01 00:00:00.000
Position in J2000 frame (km): -1876.4713 -4947.4745 3549.2275
Velocity in J2000 frame (km/s): 0.3608 -0.1366 0.0003
Epoch : 1990 MAY 03 06:00:00.249
Position in J2000 frame (km): 1875.1001 4945.4239 3552.8083
Velocity in J2000 frame (km/s): -0.3606 0.1370 -0.0003
Epoch : 1990 AUG 02 12:00:00.502
Position in J2000 frame (km): -1887.0731 -4943.3818 3549.3093
Velocity in J2000 frame (km/s): 0.3605 -0.1374 0.0003
Epoch : 1990 NOV 01 18:00:00.752
Position in J2000 frame (km): 1886.0424 4941.3201 3552.7263
Velocity in J2000 frame (km/s): -0.3603 0.1378 -0.0003
This routine provides the user level interface for computing
state transformations from one reference frame to another.
Note that the reference frames may be inertial or non-inertial.
However, the user must take care that sufficient SPICE kernel
information is loaded to provide a complete state transformation
path from the `from' frame to the `to' frame.
To perform a transformation of a state vector from one reference
to another:
either returning DOUBLE ARRAY [1,6]
to_state = xform ## from_state
or the classic IDL format returning DOUBLE ARRAY[6]
to_state = transpose(xform) # from_state
1) If sufficient information has not been supplied via loaded
SPICE kernels to compute the transformation between the two
frames, an error is signaled by a routine in the call tree of
this routine.
2) If either frame `from' or `to' is not recognized, the error
SPICE(UNKNOWNFRAME) is signaled by a routine in the call tree
of this routine.
3) If any of the input arguments, `from', `to' or `et', is
undefined, an error is signaled by the IDL error handling
system.
4) If any of the input arguments, `from', `to' or `et', is not of
the expected type, or it does not have the expected dimensions
and size, an error is signaled by the Icy interface.
5) If the output argument `xform' is not a named variable, an
error is signaled by the Icy interface.
None.
None.
ICY.REQ
ROTATION.REQ
FRAMES.REQ
None.
J. Diaz del Rio (ODC Space)
E.D. Wright (JPL)
-Icy Version 1.1.2, 10-AUG-2021 (JDR)
Updated -I/O section to comply with NAIF standard.
Edited the -Examples section to comply with NAIF standard. Added
example's problem statement and meta-kernel. Reformatted
example's output and added call to cspice_kclear.
Added -Parameters, -Exceptions, -Files, -Restrictions,
-Literature_References and -Author_and_Institution sections, and
completed -Particulars section.
Removed reference to the routine's corresponding CSPICE header from
-Abstract section.
Added arguments' type and size information in the -I/O section.
-Icy Version 1.1.1, 18-MAY-2005 (EDW)
Corrected the expression for the flattening factor,
from:
f = ( polar - equatr ) / equatr
to
f = ( equatr - polar ) / equatr
-Icy Version 1.1.0, 26-OCT-2004 (EDW)
Added capability to process vector `et' as
input returning a 6x6xN `xform' array.
-Icy Version 1.0.0, 16-JUN-2003 (EDW)
Find a state transformation matrix
|