Ephemerides.DAF_RECORD_LENGTH
— ConstantEphemerides.FTPSTR
— ConstantEphemerides.SPK_SEGMENTLIST_MAPPING
— ConstantSPK_SEGMENT_MAPPING
A dictionary mapping SPK segment types to the field index of the SPKSegmentList
.
Ephemerides.TCB_SEGMENTS
— ConstantTCB_SEGMENTS
List of the SPK/PCK segment types for which the time argument is expressed in the TCB scale.
Ephemerides.TDB_SEGMENTS
— ConstantTDB_SEGMENTS
List of the SPK/PCK segment types for which the time argument is expressed in the TDB scale.
Ephemerides.AbstractEphemRecord
— TypeAbstractEphemRecord
Abstract type for ephemeris segment records.
Ephemerides.AbstractSPKCache
— TypeAbstractSPKCache
Abstract type for all SPK segment type caches.
Ephemerides.AbstractSPKHeader
— TypeAbstractSPKHeader
Abstract type for all SPK segment type headers.
Ephemerides.AbstractSPKSegment
— TypeAbstractSPKSegment
Abstract type for all SPK segment types.
Ephemerides.DAF
— TypeDAF
Container to hold the information of NAIF's Double precision Array File (DAF).
Fields
filepath
–String
system filepath of the DAFarray
–Vector{UInt8}
binary content of the DAFheader
–DAFHeader
file record of the DAFcomment
–String
text within the DAF comment areaftype
–Int
file type, equals 1 for SPK and 2 for PCKdesc
– DAF PCK/SPK segment descriptorsseglist
–SPKSegmentList
list of the SPK/PCK segments within the DAF
References
See Also
See also DAFHeader
, Ephemerides.SPKSegmentList
and EphemerisProvider
Ephemerides.DAF
— MethodDAF(filename::String)
Parse a DAF file and retrieve its content.
This function does not initialise the segment types. That operation can be done in a second moment through the initialise_segments!
function.
Only DAF files with identification word equal to "DAF/SPK" or "DAF/PCK" are accepted, otherwise an error is thrown.
See Also
See also initialise_segments!
.
Ephemerides.DAFHeader
— TypeDAFHeader
The DAF header, or file record, is the first physical record in a DAF and stores general information about the content of the file.
Fields
nd
–Int32
number of double components in each array summaryni
–Int32
number of integer components in each array summaryfwd
–Int32
record number of initial summary recordbwd
–Int32
record number of final summary recordffa
–Int32
first free address of the filename
–String
internal name of the filelend
–Bool
true if the file was generated in little endian
References
See Also
See also DAF
and EphemerisProvider
Ephemerides.DAFHeader
— MethodDAFHeader(record::Vector{UInt8})
Parse the header (file record) of the DAF file, i.e., the first physical record in a DAF which contains global information about the file.
References
See Also
See also DAF
, parse_daf_comment
and parse_daf_summaries
Ephemerides.DAFSegmentDescriptor
— TypeDAFSegmentDescriptor
A container object to store both SPK and PCK descriptors information.
Fields
segtype
–Int32
SPK/PCK segment typetstart
–Float64
initial segment type, in TDB seconds since J2000.0tend
–Float64
final segment type, in TDB seconds since J2000.0tid
–Int32
target object NAIF IDcid
–Int32
center object NAIF IDaxesid
–Int32
reference axes ID. Defaults to -1 for PCKsiaa
–Int32
initial array addressfaa
–Int32
final array address
References
Ephemerides.DAFSegmentDescriptor
— MethodDAFSegmentDescriptor(summary::Vector{UInt8}, head::DAFHeader, isspk::Bool)
Generate an SPK or PCK descriptor by parsing a DAF summary.
Ephemerides.EphemRecordPCK
— TypeEphemRecordPCK <: AbstractEphemRecord
Store the PCK metadata relative to a given (target, center) axes pair.
Fields
target
–Int
target axes NAIF IDcenter
–Int
center axes NAIF IDt_start
– start times of each sub-window, in TDB seconds since J2000t_end
– final times of each sub-window, in TDB seconds since J2000
Ephemerides.EphemRecordSPK
— TypeEphemRecordSPK <: AbstractEphemRecord
Store the SPK metadata relative to a given (target, center) objects pair.
Fields
target
–Int
target object NAIF IDcenter
–Int
center object NAIF IDaxes
–Int
reference axes IDt_start
– start times of each sub-window, in TDB seconds since J2000t_end
– final times of each sub-window, in TDB seconds since J2000
Ephemerides.EphemerisProvider
— TypeEphemerisProvider(file::String)
EphemerisProvider(files::Vector{String})
Create an EphemerisProvider
instance by loading a single or multiple binary ephemeris kernel files specified by files
. Currently, only NAIF Double precision Array File (DAF) kernels (i.e., SPK and PCK) are accepted.
Example
julia> eph = EphemerisProvider("PATH_TO_KERNEL")
EphemerisProvider([...])
julia> eph = EphemerisProvider(["PATH_TO_KERNEL_1", "PATH_TO_KERNEL_2"])
EphemerisProvider([])
Ephemerides.InterpCache
— TypeInterpCache{T}
Ephemerides.InterpCache
— MethodInterpCache{T}(len::Int, buffsize::Int)
Create an InterpCache
instance of type T
with len
buffers, each of size buffsize
.
Ephemerides.SPKLink
— TypeSPKLink
A link object to create a mapping between DAFSegmentDescriptor
and its actual location within an EphemerisProvider
object.
Fields
desc
–DAFSegmentDescriptor
for the segment associated to this linkfid
–Int
index of the DAF containg the link data.lid
–Int
field number in theSPKSegmentList
for this segment type.eid
–Int
index of the inner segment list that stores this SPK segment.fct
–Int
1 or -1 depending on whether the (from, to) directions must be reversed.
See Also
See also SPKLinkTable
, SPKSegmentList
and add_spklinks!
.
Ephemerides.SPKLinkTable
— TypeSPKLinkTable
Dictionary object providing all the SPKLink
available between a set of (from, to) objects
Ephemerides.SPKSegmentCache1
— TypeSPKSegmentCache1 <: AbstractSPKCache
Cache instance for SPK segments of type 1 and 21. The fields contained within this cache are taken from the FORTRAN NAIF's SPICE implementation for type 1 SPK segments.
Fields
tl
– Reference epoch of the difference line.g
– Stepsize function vector.refpos
– Reference position vector.refvel
– Reference velocity vector.dt
– Modified Divided Difference arrays, with size (maxdim, 3)kqmax
– Maximum integration order plus 1.kq
– Integration order array.id
– Index of the currently loaded logical record.fc
– Buffer for the MDA computations.wc
– Buffer for the MDA computations.w
– Buffer for the MDA computations.vct
– Buffer for the MDA computations.
Ephemerides.SPKSegmentCache1
— MethodSPKSegmentCache1(head::SPKSegmentHeader1)
Initialise the cache for an SPK segment of type 1.
Ephemerides.SPKSegmentCache18
— TypeSPKSegmentCache18 <: AbstractSPKCache
Cache instance for SPK segments of type 18.
Fields
p
– Vector storing indexes of the first and last points as well as the window size.epochs
– Epochs of the interpolating points.states
– Matrix storing the states of the interpolating points.buff
– Buffers to compute the interpolating polynomials.
Ephemerides.SPKSegmentCache19
— TypeSPKSegmentCache19 <: AbstractSPKCache
Cache instance for SPK segments of type 19.
Fields
minihead
– Header with the mini-segment properties.minidata
– Cache for the mini-segment.id
– Index of the currently loaded mini-segment.
Ephemerides.SPKSegmentCache19
— MethodSPKSegmentCache19(daf::DAF, head::SPKSegmentHeader2)
Initialise the cache for an SPK segment of type 19.
Ephemerides.SPKSegmentCache2
— TypeSPKSegmentCache2 <: AbstractSPKCache
Cache instance for SPK segments of type 2 and 3.
Fields
A
– Chebyshev's polynomial coefficients, with size (ncomp, order)p
– Stores the record mid point and radius and scale factorbuff
– Stores the buffers for the Chebyshev polynomialsid
– Index of the currently loaded logical record
Ephemerides.SPKSegmentCache2
— MethodSPKSegmentCache14(head::SPKSegmentHeader14)
Initialise the cache for an SPK segment of type 14.
Ephemerides.SPKSegmentCache2
— MethodSPKSegmentCache2(head::SPKSegmentHeader2)
Initialise the cache for an SPK segment of type 2 and 3.
Ephemerides.SPKSegmentCache20
— TypeSPKSegmentCache20 <: AbstractSPKCache
Cache instance for SPK segments of type 20.
Fields
id
– Index of the currently loaded logical recordp
– Stores the record position constantsA
– Chebyshev's polynomial coefficients, with size (ncomp, order)buff
– Stores the buffers for the Chebyshev polynomials
Ephemerides.SPKSegmentCache20
— MethodSPKSegmentCache20(head::SPKSegmentHeader2)
Initialise the cache for an SPK segment of type 20.
Ephemerides.SPKSegmentCache5
— TypeSPKSegmentCache5 <: AbstractSPKCache
Cache instance for SPK segments of type 5.
Fields
c1
– Twobody propagation cache for the left state.c2
– Twobody propagation cache for the right state.epochs
– Epochs associated to the two states.id
– Index of the currently loaded logical record.
Ephemerides.SPKSegmentCache5
— MethodSPKSegmentCache5(head::SPKSegmentHeader5)
Initialise the cache for an SPK segment of type 5.
Ephemerides.SPKSegmentCache8
— TypeSPKSegmentCache8 <: AbstractSPKCache
Cache instance for SPK segments of type 8 and 12.
Fields
states
– Matrix storing the states of the interpolating points.buff
– Buffers to compute the interpolating polynomials.id
– Index of the currently loaded logical record.
Ephemerides.SPKSegmentCache8
— MethodSPKSegmentCache8(head::SPKSegmentHeader2)
Initialise the cache for an SPK segment of type 8 and 12.
Ephemerides.SPKSegmentCache9
— TypeSPKSegmentCache9 <: AbstractSPKCache
Cache instance for SPK segments of type 9 and 13.
Fields
epochs
– Epochs of the interpolating points.states
– Matrix storing the states of the interpolating points.buff
– Buffers to compute the interpolating polynomials.id
– Index of the currently loaded logical record.
Ephemerides.SPKSegmentCache9
— MethodSPKSegmentCache9(head::SPKSegmentHeader2)
Initialise the cache for an SPK segment of type 9 and 13.
Ephemerides.SPKSegmentHeader1
— TypeSPKSegmentHeader1 <: AbstractSPKHeader
Header instance for SPK segments of type 1 and 21.
Fields
n
–Int
number of records in the segmentndirs
–Int
number of directory epochsepochs
– Storage for directory epochs or epochs (when ndirs = 0)iaa
-Int
initial segment file addressetid
–Int
initial address for the epoch table (after all the MDA records)recsize
-Int
Number of double numbers stored in each MDA recordmaxdim
-Int
MDA dimension (fixed to 15 for type 1)
Ephemerides.SPKSegmentHeader1
— MethodSPKSegmentHeader1(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 1.
Ephemerides.SPKSegmentHeader14
— TypeSPKSegmentHeader14 <: AbstractSPKHeader
Header instance for SPK segments of type 14.
Fields
order
–Int
interpolating polynomial degreen
–Int
number of packets in the segmentndirs
–Int
number of epoch directoriesepochs
– Storage for directory epochs or epochs (when ndirs = 0)etid
–Int
initial address for the epoch table (after all the state data)ptid
–Int
initial address for the packet table (after the constants)pktsize
–Int
size of each data packet excluding the packet information area.pktoff
–Int
offset of the packet data from the packet startncomp
–Int
number of states coefficients (= 6 for SPK 14)N
–Int
number of polynomial coefficients
Ephemerides.SPKSegmentHeader14
— MethodSPKSegmentHeader14(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 14.
Ephemerides.SPKSegmentHeader15
— TypeSPKSegmentHeader15 <: AbstractSPKHeader
Header instance for SPK segments of type 15.
Fields
epoch
– Epoch of periapsistp
– Trajectory pole, i.e., vector parallel to the angular momentum of the orbitpv
– Central body north pole unit vectorpa
– Periapsis unit vector at epochp
– Semi-latus rectumecc
– Eccentricityj2f
– J2 processing flagvj2
– J2 validation flag, true if the orbit shape is compliant with J2 pertubations.GM
– Central body gravitational constant (km³/s²)J2
– Central body J2R
– Central body radius (km)dmdt
– Mean anomaly rate of change (rad/s)kn
– Gain factor for the regression of the nodeskp
– Gain factor for the precession of the pericenter
Ephemerides.SPKSegmentHeader15
— MethodSPKSegmentHeader15(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 15.
Ephemerides.SPKSegmentHeader17
— TypeSPKSegmentHeader17 <: AbstractSPKHeader
Header instance for SPK segments of type 17.
Fields
epoch
: epoch of periapsis (s)sma
: semi-major axis (km)h
: H term of the equinoctial elementsk
: K term of the equinoctial elementslon
: mean longitude at epoch (rad)p
: P term of the equinoctial elementsq
: Q term of the equinoctial elementsdlpdt
: rate of longitude of the periapse (rad/s)dmldt
: mean longitude rate (mean motion rate), (rad/s)dnodedt
: longitude of the ascending node rate (rad/s)ra
: equatorial pole right ascension (rad)de
: equatorial pole declination (rad)R
: Rotation matrix from planetary equator to inertial reference frame
Ephemerides.SPKSegmentHeader17
— MethodSPKSegmentHeader17(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 17.
Ephemerides.SPKSegmentHeader18
— TypeSPKSegmentHeader18 <: AbstractSPKHeader
Header instance for SPK segments of type 18.
Fields
n
–Int
number of states in the segmentndirs
–Int
number of epoch directoriesepochs
– Storage for directory epochs or epochs (when ndirs = 0)iaa
-Int
initial segment file addressetid
–Int
initial address for the epoch table (after all the state data)order
–Int
interpolating polynomial degreeN
–Int
group sizesubtype
–Int
type 18 subtype, either 0 (Hermite) or 1 (Lagrange)packetsize
–Int
packet size for each point, either 12 (Hermite) or 6 (Lagrange)
Ephemerides.SPKSegmentHeader18
— MethodSPKSegmentHeader18(n::Int)
Create a default header for type 18 segments, whose epoch array has length n
.
Ephemerides.SPKSegmentHeader19
— TypeSPKSegmentHeader19 <: AbstractSPKHeader
Header instance for SPK segments of type 19.
Fields
n
–Int
number of states in the segment.ndirs
–Int
number of epoch directories.times
– Storage for interval directories or start times (when ndirs = 0).iaa
-Int
initial segment file address.etid
–Int
byte address for the interval table (after all the minisegment data).ptid
–Int
byte for the pointer table.usefirst
–Bool
boundary flag, true if the preceding segment should be used.type
–Int
either type 18 or 19.
Ephemerides.SPKSegmentHeader19
— MethodSPKSegmentHeader19(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 18 and 19.
Ephemerides.SPKSegmentHeader2
— TypeSPKSegmentHeader2 <: AbstractSPKHeader
Header instance for SPK segments of type 2 and 3.
Fields
tstart
–Float64
initial epoch of the first record, in seconds since J2000tlen
–Float64
interval length covered by each record, in secondsorder
–Int
polynomial orderN
–Int
number of coefficients in each windown
–Int
number of records in the segmentrecsize
–Int
byte size of each logical recordncomp
–Int
number of vector componentsiaa
–Int
initial segment file addresstype
–Int
SPK segment type, either 2 or 2
Ephemerides.SPKSegmentHeader2
— MethodSPKSegmentHeader2(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 2 and 3.
Ephemerides.SPKSegmentHeader20
— TypeSPKSegmentHeader20 <: AbstractSPKHeader
Header instance for SPK segments of type 20.
Fields
dscale
–Float64
length conversion factortscale
–Float64
time conversion factortstart
–Float64
initial epoch of the first record, in seconds since J2000tlen
–Float64
interval length covered by each record, in secondsrecsize
–Int
byte size of each logical recordorder
–Int
polynomial orderN
–Int
number of coefficients in each windown
–Int
number of records in the segmentiaa
–Int
initial segment file address
Ephemerides.SPKSegmentHeader20
— MethodSPKSegmentHeader20(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 20
Ephemerides.SPKSegmentHeader5
— TypeSPKSegmentHeader5 <: AbstractSPKHeader
Header instance for SPK segments of type 5.
Fields
GM
–Float64
Gravitational constantn
–Int
number of statesndirs
–Int
number of epoch directoriesetid
–Int
initial address for the epoch table (after all the state data)epochs
– Storage for directory epochs or epochs (when ndirs = 0)iaa
-Int
initial segment file address
Ephemerides.SPKSegmentHeader5
— MethodSPKSegmentHeader5(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 5.
Ephemerides.SPKSegmentHeader8
— TypeSPKSegmentHeader8 <: AbstractSPKHeader
Header instance for SPK segments of type 8 and 12.
Fields
tstart
–Float64
segment starting epoch, in TDB seconds since J2000tlen
–Float64
interval length, in secondsorder
–Int
interpolating polynomial degreeN
–Int
group sizen
–Int
number of states in the segmentiaa
-Int
initial segment file addressiseven
–Bool
true for even group sizetype
–Int
SPK type (either 8 or 12)
Ephemerides.SPKSegmentHeader8
— MethodSPKSegmentHeader8(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 8 and 12.
Ephemerides.SPKSegmentHeader9
— TypeSPKSegmentHeader9 <: AbstractSPKHeader
Header instance for SPK segments of type 9 and 13.
Fields
n
–Int
number of states in the segmentndirs
–Int
number of epoch directoriesepochs
– Storage for directory epochs or epochs (when ndirs = 0)iaa
-Int
initial segment file addressetid
–Int
initial address for the epoch table (after all the state data)order
–Int
interpolating polynomial degreeN
–Int
group sizeiseven
–Bool
true for even group sizetype
–Int
SPK type (either 9 or 13)
Ephemerides.SPKSegmentHeader9
— MethodSPKSegmentHeader9(daf::DAF, desc::DAFSegmentDescriptor)
Create the segment header for an SPK segment of type 9 and 13.
Ephemerides.SPKSegmentList
— TypeSPKSegmentList
A container object to efficiently store all the different SPK segments that are contained within a single DAF file.
SPKSegmentList()
Initialises an empty SPKSegmentList
object.
See also
See also Ephemerides.add_segment!
Ephemerides.SPKSegmentType1
— TypeSPKSegmentType1 <: AbstractSPKSegment
Segment instance for SPK segments of type 1 and 21, which contain Modified Difference Arrays (MDA). This data type is normally used for spacecraft whose ephemerides are produced by JPL's principal trajectory integrator DPTRAJ.
Fields
head
– Segment headercache
– Segment cache
References
Ephemerides.SPKSegmentType1
— MethodSPKSegmentType1(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 1.
Ephemerides.SPKSegmentType14
— TypeSPKSegmentType14 <: AbstractSPKSegment
Segment instance for SPK segments of type 14.
Fields
head
– Segment headercache
– Segment cache
References
Ephemerides.SPKSegmentType14
— MethodSPKSegmentType14(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 14.
Ephemerides.SPKSegmentType15
— TypeSPKSegmentType15 <: AbstractSPKSegment
Segment instance for SPK segments of type 15.
Fields
head
– Segment headercache
– Segment cache
References
Ephemerides.SPKSegmentType15
— MethodSPKSegmentType15(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 15.
Ephemerides.SPKSegmentType17
— TypeSPKSegmentType17 <: AbstractSPKSegment
Segment instance for SPK segments of type 17.
Fields
head
– Segment header
SPK segments of type 17 do not require a cache because they do not extract any additional coefficients at runtime.
References
Ephemerides.SPKSegmentType17
— MethodSPKSegmentType17(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 17.
Ephemerides.SPKSegmentType19
— TypeSPKSegmentType19 <: AbstractSPKSegment
Segment instance for SPK segments of type 18 and 19. Type 18 segments are treated as special cases of a type 19 with a single mini-segment.
Fields
head
– Segment headercache
– Segment cache
References
Ephemerides.SPKSegmentType19
— MethodSPKSegmentType19(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 19.
Ephemerides.SPKSegmentType2
— TypeSPKSegmentType2 <: AbstractSPKSegment
Segment instance for SPK segments of type 2 and 3, which contain Chebyshev polynomial coefficients for the position and/or state of the body as function of time. This data type is normally used for planet barycenters, and for satellites whose ephemerides are integrated.
Fields
head
– Segment headercache
– Segment cache
References
Ephemerides.SPKSegmentType2
— MethodSPKSegmentType2(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 2.
Ephemerides.SPKSegmentType20
— TypeSPKSegmentType2 <: AbstractSPKSegment
Segment instance for SPK segments of type 20, which contain Chebyshev polynomial coefficients for the position and/or state of the body as function of time. This data type is normally used for planet barycenters, and for satellites whose ephemerides are integrated.
Fields
head
– Segment headercache
– Segment cache
References
Ephemerides.SPKSegmentType20
— MethodSPKSegmentType20(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 2.
Ephemerides.SPKSegmentType5
— TypeSPKSegmentType5 <: AbstractSPKSegment
Segment instance for SPK segments of type 5.
Fields
head
– Segment headercache
– Segment cache
References
Ephemerides.SPKSegmentType5
— MethodSPKSegmentType2(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 5.
Ephemerides.SPKSegmentType8
— TypeSPKSegmentType8 <: AbstractSPKSegment
Segment instance for SPK segments of type 8 and 12.
Fields
head
– Segment headercache
– Segment cache
References
Ephemerides.SPKSegmentType8
— MethodSPKSegmentType8(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 8 and 12.
Ephemerides.SPKSegmentType9
— TypeSPKSegmentType9 <: AbstractSPKSegment
Segment instance for SPK segments of type 9 and 13.
Fields
head
– Segment headercache
– Segment cache
References
Ephemerides.SPKSegmentType9
— MethodSPKSegmentType9(daf::DAF, desc::DAFSegmentDescriptor)
Create the object representing an SPK segment of type 9 and 13.
Ephemerides.TwoBodyUniversalCache
— TypeTwoBodyUniversalCache
A container to store precomputed quantities required for the two-body propagation based on universal variables. Only the quantities that depend on the initial state are computed.
Ephemerides.add_segment!
— Methodadd_segment!(list::SPKSegmentList, spk::AbstractSPKSegment)
Add the SPK segment to the proper vector within the given Ephemerides.SPKSegmentList
list
Ephemerides.add_spklinks!
— Methodadd_spklinks!(table::SPKLinkTable, daf::DAF, fid::Int)
Insert in the input SPKLinkTable
all the SPK or PCK links associated to the segment descriptors of the input DAF.
Ephemerides.analyse_timespan
— Methodanalyse_timespan(records)
Analyse a set of AbstractEphemRecord
, returning the minimum and maximum covered times, in TDB seconds since J2000, together with a continuity parameter.
References
- CALCEPH C++ library
See Also
See also ephem_spk_timespan
and ephem_pck_timespan
.
Ephemerides.array
— Methodget_array(daf::DAF)
Return the byte content of the DAF file.
Ephemerides.axes
— Methodaxes(desc::DAFSegmentDescriptor)
Return the NAIF integer code for the reference axes. It is valid only for SPK files and defaults to -1 for PCKs.
Ephemerides.cache
— Methodcache(spk::AbstractSPKSegment)
Return the segment cache data.
Ephemerides.center
— Methodcenter(desc::DAFSegmentDescriptor)
Return the NAIF integer code for the reference object or axes for SPK and PCK, respectively.
Ephemerides.chebyshev
— Functionchebyshev(cache::InterpCache, cₖ, t::Number, idx::Int, N::Int, ibuff::Int)
Evaluate a sum of Cheybyshev polynomials of the first kind at t
using a recursive algorithm. It simultenously evalutes the 3 state components. idx
is the index of the starting row (in 0-based notation) in the matrix of coefficients cₖ
. ibuff
is the index of the first free buffer.
x
is a re-work of the actual ascissa value that lies between [-1, 1]
See Also
See also ∂chebyshev
, ∂²chebyshev
and ∂³chebyshev
Ephemerides.comment
— Methodget_comment(daf::DAF)
Return the comment written in the DAF comment section.
Ephemerides.compute_mda_position
— Methodcompute_mda_position(cache::SPKSegmentCache1, Δ::Number)
Ephemerides.compute_mda_velocity
— Methodcompute_mda_velocity(cache::SPKSegmentCache1, Δ::Number)
Ephemerides.create_linktables
— Methodcreate_linktables(dafs::Vector{DAF})
Create the SPK and PCK SPKLinkTable
for all the segments stores in the input DAFs.
Ephemerides.create_spk_segment
— Methodcreate_spk_segment(daf::DAF, desc::DAFSegmentDescriptor)
Initialise an SPK segment according to the segment type defined in the DAFSegmentDescriptor
desc
.
Ephemerides.descriptor
— Methoddescriptor(link::SPKLink)
Return the SPK/PCK segment descriptor associated to this link.
Ephemerides.descriptors
— Methodget_descriptors(daf::DAF)
Return the SPK/PCK segment descriptors contained in the DAF.
Ephemerides.element_id
— Methodelement_id(link::SPKLink)
Return the segment index in the inner SPK/PCK segment list.
Ephemerides.endian
— Methodendian(head::DAFHeader)
Return true
if the DAF file is in little-endian.
Ephemerides.endian
— Methodendian(daf::DAF)
Return true
if the DAF is in little-endian.
Ephemerides.ephem_get_axes
— Methodephem_get_axes(eph::EphemerisProvider)
Return a list of Frame IDs representing axes with available orientation data.
Ephemerides.ephem_get_points
— Methodephem_get_points(eph::EphemerisProvider)
Return a list of NAIFIds representing bodies with available ephemeris data.
Ephemerides.ephem_pck_records
— Methodephem_pck_records(eph::EphemerisProvider)
Return a list of EphemRecordPCK
storing metadata relative to each (target, center) axes pairs in the loaded PCK kernels.
Ephemerides.ephem_pck_timespan
— Methodephem_pck_timespan(eph::EphemerisProvider)
Return the minimum and maximum time available in the PCK kernels loaded within eph
, in TDB seconds since J2000, together with a continuity parameter defined as follows:
- 0 no PCK data is available.
- 1 the quantities of all axes are available for any time between the first and last time.
- 2 the quantities of some axes are available on discontinuous time intervals between the first and last time.
- 3 the quantities of each axis are available on a continuous time interval between the first and the last time, but not available for any time between the first and last time.
References
- CALCEPH C++ library
See Also
See also ephem_spk_timespan
.
Ephemerides.ephem_rotation12
— Methodephem_rotation12(eph::EphemerisProvider, from::Int, to::Int, time::Number)
Compute the 12-elements orientation angles of one set of axes (to) relative to another (from) at time
, expressed in TDB/TCB seconds since J2000, in accordance with the kernel timescale.
For the orientation angles, it is not possible to automatically compute the reverse transformation , i.e., if the orientation of PA440 is defined with respect to the ICRF, it is not possible to compute the rotation from the PA440 to the ICRF with this routine.
Ephemerides.ephem_rotation3
— Methodephem_rotation3(eph::EphemerisProvider, from::Int, to::Int, time::Number)
Compute the 3-elements orientation angles of one set of axes (to) relative to another (from) at time
, expressed in TDB/TCB seconds since J2000, in accordance with the kernel timescale.
For the orientation angles, it is not possible to automatically compute the reverse transformation , i.e., if the orientation of PA440 is defined with respect to the ICRF, it is not possible to compute the rotation from the PA440 to the ICRF with this routine.
Ephemerides.ephem_rotation6
— Methodephem_rotation6(eph::EphemerisProvider, from::Int, to::Int, time::Number)
Compute the 6-elements orientation angles of one set of axes (to) relative to another (from) at time
, expressed in TDB/TCB seconds since J2000, in accordance with the kernel timescale.
For the orientation angles, it is not possible to automatically compute the reverse transformation , i.e., if the orientation of PA440 is defined with respect to the ICRF, it is not possible to compute the rotation from the PA440 to the ICRF with this routine.
Ephemerides.ephem_rotation9
— Methodephem_rotation9(eph::EphemerisProvider, from::Int, to::Int, time::Number)
Compute the 9-elements orientation angles of one set of axes (to) relative to another (from) at time
, expressed in TDB/TCB seconds since J2000, in accordance with the kernel timescale.
For the orientation angles, it is not possible to automatically compute the reverse transformation , i.e., if the orientation of PA440 is defined with respect to the ICRF, it is not possible to compute the rotation from the PA440 to the ICRF with this routine.
Ephemerides.ephem_spk_records
— Methodephem_spk_records(eph::EphemerisProvider)
Return a list of EphemRecordSPK
storing metadata relative to each (target, center) object pairs in the loaded SPK kernels.
Ephemerides.ephem_spk_timespan
— Methodephem_spk_timespan(eph::EphemerisProvider)
Return the minimum and maximum time available in the SPK kernels loaded within eph
, in TDB seconds since J2000, together with a continuity parameter defined as follows:
- 0 no SPK data is available.
- 1 the quantities of all bodies are available for any time between the first and last time.
- 2 the quantities of some bodies are available on discontinuous time intervals between the first and last time.
- 3 the quantities of each body are available on a continuous time interval between the first and the last time, but not available for any time between the first and last time.
References
- CALCEPH C++ library
See Also
See also ephem_pck_timespan
.
Ephemerides.ephem_timescale_id
— Methodephem_timescale_id(eph::EphemerisProvider)
Retrieve a timescale ID associated with the ephemeris handler eph
. It returns 1 for Barycentric Dynamical Time (TDB) and 2 for Barycentric Coordinate Time (TCB).
Ephemeris providers with mixed timescales are not supported. An error is thrown if in the ephemeris handler some segments are defined in TDB and some other segments in TCB.
Ephemerides.ephem_vector12
— Methodephem_vector12(eph::EphemerisProvider, from::Int, to::Int, time::Number)
Compute the 12-elements state of one body (to) relative to another (from) at time
, expressed in TDB/TCB seconds since J2000, in accordance with the kernel timescale.
Ephemerides.ephem_vector3
— Methodephem_vector3(eph::EphemerisProvider, from::Int, to::Int, time::Number)
Compute the 3-elements state of one body (to) relative to another (from) at time
, expressed in TDB/TCB seconds since J2000, in accordance with the kernel timescale.
Ephemerides.ephem_vector6
— Methodephem_vector6(eph::EphemerisProvider, from::Int, to::Int, time::Number)
Compute the 6-elements state of one body (to) relative to another (from) at time
, expressed in TDB/TCB seconds since J2000, in accordance with the kernel timescale.
Ephemerides.ephem_vector9
— Methodephem_vector9(eph::EphemerisProvider, from::Int, to::Int, time::Number)
Compute the 9-elements state of one body (to) relative to another (from) at time
, expressed in TDB/TCB seconds since J2000, in accordance with the kernel timescale.
Ephemerides.factor
— Methodfactor(link::SPKLink)
Return the direction multiplicative factor.
Ephemerides.file_id
— Methodfile_id(link::SPKLink)
Return the DAF file index.
Ephemerides.filename
— Methodfilename(head::DAFHeader)
Return the internal description of the DAF.
Ephemerides.filepath
— Methodfilepath(daf::DAF)
Return the system path of the DAF.
Ephemerides.final_address
— Methodfinal_address(desc::DAFSegmentDescriptor)
Return the final address of the segment array in teh DAF.
Ephemerides.final_record
— Methodfinal_record(head::DAFHeader)
Return the record number of the final summary record in the DAF
Ephemerides.final_record
— Methodfinal_record(daf::DAF)
Return the record number of the final summary record in the DAF.
Ephemerides.final_time
— Methodfinal_time(desc::DAFSegmentDescriptor)
Return the final epoch of the interval for which ephemeris data are contained in the segment, in seconds since J2000.0
Ephemerides.final_time
— Methodfinal_time(link::SPKLink)
Return the final epoch of the interval for which ephemeris data are contained in the segment associated to this link, in seconds since J2000.0
Ephemerides.final_times
— Methodfinal_times(record::AbstractEphemRecord)
Recover the final times of each sub-window in which the ephemeris data of record
is defined, expressed in TDB seconds since J2000
Ephemerides.find_logical_record
— Methodfind_logical_record(daf::DAF, head::SPKSegmentHeader1, time::Number)
Ephemerides.find_logical_record
— Methodfind_logical_record(daf::DAF, head::SPKSegmentHeader14, time::Number)
Ephemerides.find_logical_record
— Methodfind_logical_record(daf::DAF, head::SPKSegmentHeader18, time::Number)
Ephemerides.find_logical_record
— Methodfind_logical_record(daf::DAF, head::SPKSegmentHeader5, time::Number)
Ephemerides.find_logical_record
— Methodfind_logical_record(daf::DAF, head::SPKSegmentHeader1, time::Number)
Ephemerides.find_logical_record
— Methodfind_logical_record(head::SPKSegmentHeader2, time::Number)
Ephemerides.find_logical_record
— Methodfind_logical_record(head::SPKSegmentHeader20, time::Number)
Ephemerides.find_logical_record
— Methodfind_logical_record(head::SPKSegmentHeader8, time::Number)
Ephemerides.find_minisegment
— Methodfind_minirecord(daf::DAF, head::SPKSegmentHeader19, time::Number)
Ephemerides.free_address
— Methodfree_address(head::DAFHeader)
Return the first free address in the file, i.e., the address at which the first element of the next array is to be added.
Ephemerides.free_address
— Methodfree_address(daf::DAF)
Return the first free address in the file, i.e., the address at which the first element of the next array is to be added.
Ephemerides.get_buffer
— Methodget_buffer(c::InterpCache, idx::int, x::Number)
Return the idx
-th buffer from the corresponding DiffCache
depending on the type of x
.
Ephemerides.get_coefficients!
— Methodget_coefficients!(daf::DAF, head::SPKSegmentHeader1, cache::SPKSegmentCache1, index::Int)
Ephemerides.get_coefficients!
— Methodget_coefficients!(daf::DAF, head::SPKSegmentHeader14, cache::SPKSegmentCache14, index::Int)
Ephemerides.get_coefficients!
— Methodget_coefficients!(daf::DAF, head::SPKSegmentHeader18, cache::SPKSegmentCache18, first::Int, last::Int)
Ephemerides.get_coefficients!
— Methodget_coefficients!(daf::DAF, head::SPKSegmentHeader2, cache::SPKSegmentCache2, index::Int)
Ephemerides.get_coefficients!
— Methodget_coefficients!(daf::DAF, head::SPKSegmentHeader20, cache::SPKSegmentCache20, index::Int)
Ephemerides.get_coefficients!
— Methodget_coefficients!(daf::DAF, head::SPKSegmentHeader5, cache::SPKSegmentCache5, index::Int)
Ephemerides.get_coefficients!
— Methodget_coefficients!(daf::DAF, head::SPKSegmentHeader8, cache::SPKSegmentCache8, index::Int)
Ephemerides.get_coefficients!
— Methodget_coefficients!(daf::DAF, head::SPKSegmentHeader9, cache::SPKSegmentCache9, index::Int)
Ephemerides.get_daf
— Methoddaf(eph::EphemerisProvider, id::Int)
Return the DAF
file in the ephemeris provider at index id
.
Ephemerides.get_daf
— Methoddaf(eph::EphemerisProvider)
Return the DAF
files stored in the ephemeris provider.
Ephemerides.get_record
— Methodget_record(array::Vector{UInt8}, index::Integer)
Retrieve a whole DAF record at position index
.
Ephemerides.get_segment
— Methodget_segment(list::SPKSegmentList, lid::Int, eid::Int)
Return the segment contained in the lid
list at index eid
.
Ephemerides.get_segment_boundaries
— Methodget_segment_boundaries(desclist::Vector{DAFSegmentDescriptor})
Parse all the segment descriptors of a given (center, target) pair and return a set of initial and final times, in TDB seconds since J2000, representing all the time sub-windows in which the ephemeris data for this pair is defined.
Ephemerides.header
— Methodget_header(daf::DAF)
Return the DAFHeader
header of the DAF.
Ephemerides.header
— Methodheader(spk::AbstractSPKSegment)
Return the segment header.
Ephemerides.hermite
— Methodhermite(cache::InterpCache, states, epochs, x, idx::Int, N::Int)
Evalute a Hermite polynomial at x
using a recursive algorithm. This function handles unequally-spaced polynomials, where the coefficients in states
are interpolated at epochs
.
Ephemerides.hermite
— Methodhermite(cache::InterpCache, states, x, idx::Int, N::Int, Δt::Number)
Evalute a Hermite polynomial at x
using a recursive algorithm. This function is valid only for equally-spaced polynomials. idx
is the index of the desired state, N
is the number of coefficients of the polynomial and Δt
is the length of the polynomial interval.
x
is a re-work of the actual ascissa value that starts at 1
See Also
Ephemerides.initial_address
— Methodinitial_address(desc::DAFSegmentDescriptor)
Return the initial address of the segment array in the DAF.
Ephemerides.initial_record
— Methodinitial_record(head::DAFHeader)
Return the record number of the initial summary record in the DAF
Ephemerides.initial_record
— Methodinitial_record(daf::DAF)
Return the record number of the initial summary record in the DAF.
Ephemerides.initial_time
— Methodinitial_time(desc::DAFSegmentDescriptor)
Return the initial epoch of the interval for which ephemeris data are contained in the segment, in seconds since J2000.0
Ephemerides.initial_time
— Methodinitial_time(link::SPKLink)
Return the initial epoch of the interval for which ephemeris data are contained in the segment associated to this link, in seconds since J2000.0
Ephemerides.initial_times
— Methodinitial_times(record::AbstractEphemRecord)
Recover the initial times of each sub-window in which the ephemeris data of record
is defined, expressed in TDB seconds since J2000
Ephemerides.initialise_segments!
— Methodinitialise_segments!(daf::DAF)
Fill the Ephemerides.SPKSegmentList
by initialising the SPK/PCK segments associated to all the descriptors stores within the DAF.
See Also
See also DAF
and create_spk_segment
Ephemerides.is_little_endian
— Methodis_little_endian(array::Vector{UInt8})
Return true if the array corresponds to the string indicating a little-endian format.
Ephemerides.is_pck
— Methodis_pck(daf::DAF)
Return true
if the DAF stores PCK data.
Ephemerides.is_spk
— Methodis_spk(daf::DAF)
Return true
if the DAF stores SPK data.
Ephemerides.lagrange
— Methodlagrange(cache::InterpCache, states, epochs, x, idx::Int, N::Int)
Recursively evaluate a Lagrange polynomial at x
by using Neville's algorithm. This function handles unequally-spaced polynomials, where the coefficients in states
are interpolated at epochs
.
Ephemerides.lagrange
— Methodlagrange(cache::InterpCache, states, x, idx::Int, N::Int)
Recursively evaluate a Lagrange polynomial at x
by using Neville's algorithm. This function is valid only for equally-spaced polynomials. idx
is the index of the desired state and N
is the number of coefficients of the polynomial.
x
is a re-work of the actual ascissa value that starts at 1
See Also
See also ∂lagrange
and ∂²lagrange
.
Ephemerides.list_id
— Methodlist_id(link::SPKLink)
Return the index of the list containing the segments of the given SPK/PCK type.
Ephemerides.load_minisegment!
— Methodload_minisegment!(daf::DAF, head::SPKSegmentHeader19, cache::SPKSegmentCache19, index::Int)
Ephemerides.normalise_time
— Methodnormalise_time(cache::SPKSegmentCache2, time::Number)
Transform time
in an interval between [-1, 1] for compliance with Chebyshev polynomials.
Ephemerides.normalise_time
— Methodnormalise_time(head::SPKSegmentHeader20, time::Number, index::Int)
Ephemerides.normalise_time
— Methodnormalise_time(head::SPKSegmentHeader8, time::Number, index::Int)
Returned a normalised time that starts at 1 at the beginning of the interval.
Ephemerides.parse_daf_comment
— Methodparse_daf_comment(array::Vector{UInt8}, header::DAFHeader)
Retrieve the comment section of a binary DAF.
References
See Also
See also DAF
, DAFHeader
and parse_daf_summaries
Ephemerides.parse_daf_summaries
— Methodparse_daf_summaries(array::Vector{UInt8}, head::DAFHeader)
Parse the DAF binary content and retrieve all the summary records.
References
See Also
See also DAF
, DAFHeader
and parse_daf_comment
Ephemerides.parse_pck_segment_descriptor
— Methodparse_pck_segment_descriptor(summary::Vector{UInt8}, lend::Bool)
Create a DAFSegmentDescriptor
object by parsing a binary PCK segment descriptor. A default value of -1 is used to fill the reference frame field. The target and center fields are used for the actual target and center axes.
References
Ephemerides.parse_spk_segment_descriptor
— Methodparse_spk_segment_descriptor(summary::Vector{UInt8}, lend::Bool)
Create a DAFSegmentDescriptor
object by parsing a binary SPK segment descriptor.
References
Ephemerides.pck_links
— Methodpck_links(eph::EphemerisProvider)
Return the SPKLinkTable
for the PCK segments.
Ephemerides.propagate_twobody
— Methodpropagate_twobody(cache::TwoBodyUniversalCache, Δt::Number)
Propagate the state vector in cache
of Δt
using the universal variables formulation for Kepler's Equation and the Lagrange coefficients f and g.
This routine is valid for any type of orbit and uses a bisection method to find the root of the universal variables Kepler's equation. The algorithm has been directly taken from the SPICE toolkit prob2b.f
.
References
Ephemerides.reset_indexes!
— Methodreset_indexes!(cache::SPKSegmentCache18)
Reset the cache indexes to force the coefficients reload.
Ephemerides.reverse_link
— Methodreverse_link(link::SPKLink)
Reverse the sign, i.e. change the sign of the multiplicative factor, of the link.
Ephemerides.segment_list
— Methodget_segment_list(daf::DAF)
Return the Ephemerides.SPKSegmentList
list of segments stored in the DAF.
Ephemerides.segment_type
— Methodsegment_type(desc::DAFSegmentDescriptor)
Return the SPK/PCK segment type.
Ephemerides.spk_field
— Methodspk_field(spk::AbstractSPKSegment)
Return the field number in the Ephemerides.SPKSegmentList
associated to the given SPK segment type.
Ephemerides.spk_links
— Methodspk_links(eph::EphemerisProvider)
Return the [SPKLinkTable
] for the SPK segments.
Ephemerides.stumpff
— Methodstumpff(x::Number, p::AbstractVector)
Compute Stumpff's functions from C₀ up to C₃ at x
.
This routine uses the trigonometrical expressions of the functions when the absolute value of x
is greater or equal to 1. If that is not the case, the C₂ and C₃ functions are computed from a truncated expression of the Maclaurin series at order 11, which guarantees a higher precision and avoid overflow errors when x
is null.
References
Ephemerides.summary_size
— Methodsummary_size(head::DAFHeader)
Compute the size of a single summary record of a DAF file, in bytes.
Ephemerides.target
— Methodtarget(desc::DAFSegmentDescriptor)
Return the NAIF integer code for the target object or axes for SPK and PCK, respectively.
Ephemerides.update_cache!
— Methodupdate_cache!(c::TwoBodyUniversalCache)
Update the precomputed values in the cache using the position and velocity in c
.
Ephemerides.update_header!
— Methodupdate_header!(head::SPKSegmentHeader18, daf::DAF, iaa, faa, type)
Update the header of type 18 segments.
Ephemerides.∂chebyshev
— Function∂chebyshev(cache::InterpCache, cₖ, t::Number, idx::Int, N::Int, Δt, ibuff=1)
Evaluate a sum of Cheybyshev polynomials of the first kind and its derivative at t
using a recursive algorithm. It simultenously evalutes the 3 state components. idx
is the index of the starting row (in 0-based notation) in the matrix of coefficients cₖ
. ibuff
is the index of the first free buffer.
x
is a re-work of the actual ascissa value that lies between [-1, 1]
See Also
See also chebyshev
, ∂²chebyshev
and ∂³chebyshev
Ephemerides.∂hermite
— Method∂hermite(cache::InterpCache, states, epochs, x, idx::Int, N::Int)
Evalute a Hermite polynomial and its derivative at x
using a recursive algorithm. This function handles unequally-spaced polynomials, where the coefficients in states
are interpolated at epochs
.
Ephemerides.∂hermite
— Method∂hermite(cache::InterpCache, states, x, idx::Int, N::Int, Δt::Number)
Evalute a Hermite polynomial and its derivative at x
using a recursive algorithm. This function is valid only for equally-spaced polynomials. idx
is the index of the desired state, N
is the number of coefficients of the polynomial and Δt
is the length of the polynomial interval.
x
is a re-work of the actual ascissa value that starts at 1
See Also
Ephemerides.∂lagrange
— Method∂lagrange(cache::InterpCache, states, epochs, x, idx::Int, N::Int)
Recursively evaluate a Lagrange polynomial and its derivative at x
by using Neville's algorithm. This function handles unequally-spaced polynomials, where the coefficients in states
are interpolated at epochs
.
Ephemerides.∂lagrange
— Method∂lagrange(cache::InterpCache, states, x, idx::Int, N::Int, Δt::Number)
Recursively evaluate a Lagrange polynomial and its derivative at x
by using Neville's algorithm. This function is valid only for equally-spaced polynomials. idx
is the index of the desired state, N
is the number of coefficients of the polynomial and Δt
is the length of the polynomial interval
x
is a re-work of the actual ascissa value that starts at 1
See Also
See also lagrange
and ∂²lagrange
Ephemerides.∂²chebyshev
— Function∂²chebyshev(cache::InterpCache, cₖ, t::Number, idx::Int, N::Int, Δt, ibuff=1)
Evaluate a sum of Cheybyshev polynomials of the first kind and its two derivatives at t
using a recursive algorithm. It simultenously evalutes the 3 state components. idx
is the index of the starting row (in 0-based notation) in the matrix of coefficients cₖ
. ibuff
is the index of the first free buffer.
x
is a re-work of the actual ascissa value that lies between [-1, 1]
See Also
See also chebyshev
, ∂chebyshev
and ∂³chebyshev
Ephemerides.∂²hermite
— Method∂²hermite(cache::InterpCache, states, epochs, x, idx::Int, N::Int)
Evalute a Hermite polynomial and its two derivatives at x
using a recursive algorithm. This function handles unequally-spaced polynomials, where the coefficients in states
are interpolated at epochs
.
Ephemerides.∂²hermite
— Method∂²hermite(cache::InterpCache, states, x, idx::Int, N::Int, Δt::Number)
Evalute a Hermite polynomial and its two derivatives at x
using a recursive algorithm. This function is valid only for equally-spaced polynomials. idx
is the index of the desired state, N
is the number of coefficients of the polynomial and Δt
is the length of the polynomial interval.
x
is a re-work of the actual ascissa value that starts at 1
See Also
Ephemerides.∂²lagrange
— Method∂²lagrange(cache::InterpCache, states, epochs, x, idx::Int, N::Int)
Recursively evaluate a Lagrange polynomial and its two derivatives at x
by using Neville's algorithm. This function handles unequally-spaced polynomials, where the coefficients in states
are interpolated at epochs
.
See Also
Ephemerides.∂²lagrange
— Method∂²lagrange(cache::InterpCache, states, x, idx::Int, N::Int, Δt::Number)
Recursively evaluate a Lagrange polynomial and its two derivatives at x
by using Neville's algorithm. This function is valid only for equally-spaced polynomials. idx
is the index of the desired state, N
is the number of coefficients of the polynomial and Δt
is the length of the polynomial interval
x
is a re-work of the actual ascissa value that starts at 1
See Also
Ephemerides.∂³chebyshev
— Function∂³chebyshev(cache::InterpCache, cₖ, t::Number, idx::Int, N::Int, Δt, ibuff=1)
Evaluate a sum of Cheybyshev polynomials of the first kind and its three derivatives at t
using a recursive algorithm. It simultenously evalutes the 3 state components. idx
is the index of the starting row (in 0-based notation) in the matrix of coefficients cₖ
. ibuff
is the index of the first free buffer.
x
is a re-work of the actual ascissa value that lies between [-1, 1]
See Also
See also chebyshev
, ∂chebyshev
and ∂²chebyshev
Ephemerides.∂³hermite
— Method∂³hermite(cache::InterpCache, states, epochs, x, idx::Int, N::Int)
Evalute a Hermite polynomial and its three derivatives at x
using a recursive algorithm. This function handles unequally-spaced polynomials, where the coefficients in states
are interpolated at epochs
.
Ephemerides.∂³hermite
— Method∂³hermite(cache::InterpCache, states, x, idx::Int, N::Int, Δt::Number)
Evalute a Hermite polynomial and its three derivatives at x
using a recursive algorithm. This function is valid only for equally-spaced polynomials. idx
is the index of the desired state, N
is the number of coefficients of the polynomial and Δt
is the length of the polynomial interval.
x
is a re-work of the actual ascissa value that starts at 1
See Also
Ephemerides.∫chebyshev
— Method∫chebyshev(cache::InterpCache, cₖ, t::Number, N::Int, Δt, tlen, p₀)
Evaluate the integral of a sum of Cheybyshev polynomials of the first kind using a recursive algorithm. It simultenously evalutes the 3 state components. This function simultaneously computes both the Chebyshev polynomials as well as their integrals.
Ephemerides.∫chebyshev
— Method∫chebyshev(cache::InterpCache, cₖ, t::Number, N::Int, Δt, tlen, p₀)
Evaluate the integral of a sum of Cheybyshev polynomials of the first kind using a recursive algorithm. It simultenously evalutes the 3 state components. It assumes the Chebyshev polynomials up to degree N have already been computed and are stored in the buffer with index ibuff
. tlen
is the size of the record interval, Δt
is the timescale factor, and p₀
is a vector containing the position coefficients at the midpoint (i.e., when the integral is evaluated at t = 0).
x
is a re-work of the actual ascissa value that lies between [-1, 1]