PCL
pcl::JPLEphemeris Class Reference

JPL planetary ephemeris. More...

#include <JPLEphemeris.h>

Public Types

using constant = EphemerisConstant
 
using constant_list = EphemerisConstantList
 
using item_index = JPLEphemerisItem::value_type
 

Public Member Functions

 JPLEphemeris (const String &filePath)
 
void AddData (const String &filePath)
 
void ComputeState (Vector &r, double jd1, double jd2, item_index i) const
 
void ComputeState (Vector &r, TimePoint t, item_index i) const
 
void ComputeState (Vector &r, Vector &v, double jd1, double jd2, item_index i) const
 
void ComputeState (Vector &r, Vector &v, TimePoint t, item_index i) const
 
const constant_listConstants () const
 
double ConstantValue (const IsoString &name) const
 
int DENumber () const
 
double EndJD () const
 
double EphemerisEndJD () const
 
double EphemerisStartJD () const
 
bool IsItemAvailable (int i) const
 
int LENumber () const
 
double StartJD () const
 
IsoString Summary () const
 
bool Test (const String &filePath, bool verbose=true) const
 

Static Public Member Functions

static int ComponentsForItem (int i)
 

Detailed Description

This class implements JPL DE/LE planetary ephemerides computed from original files in ASCII format. On the PixInsight/PCL platform, this is a utility class used to generate truncated binary ephemeris files in XEPH format by reinterpolation with Chebyshev polynomials from the original ephemeris data.

Source JPL ephemeris ASCII files, documentation, and original test and utility programs and examples, are available at:

ftp://ssd.jpl.nasa.gov/pub/eph/planets/

See also
EphemerisFile

Definition at line 160 of file JPLEphemeris.h.

Member Typedef Documentation

◆ constant

Represents a numerical integration ephemeris constant defined by its name and value.

Definition at line 174 of file JPLEphemeris.h.

◆ constant_list

Represents a list of integration ephemeris constants.

Definition at line 179 of file JPLEphemeris.h.

◆ item_index

using pcl::JPLEphemeris::item_index = JPLEphemerisItem::value_type

Represents an ephemeris item. See the JPLEphemerisItem namespace for a list of the available items.

Definition at line 168 of file JPLEphemeris.h.

Constructor & Destructor Documentation

◆ JPLEphemeris()

pcl::JPLEphemeris::JPLEphemeris ( const String filePath)

Constructs a &JPLEphemeris instance initialized from the specified DE/LE ephemeris header file in ASCII format.

In the event of errors or invalid header data, this constructor throws the appropriate Error exceptions.

Member Function Documentation

◆ AddData()

void pcl::JPLEphemeris::AddData ( const String filePath)

Adds new ephemeris Chebyshev coefficients from the specified JPL DE/LE ephemeris data file in ASCII format.

If this is not the first call to this member function for this object, the specified data file must cover a time span contiguous to the existing coefficients, that is, you must call this function successively with a list of consecutive ephemeris data files. For example, for DE/LE 432, you can call this function in sequence for the following file names:

ascp01550.432, ascp01650.432, ascp01750.432, ..., ascp02550.432

to cover the entire ephemeris time span from JDE 2287184.5 to 2688976.5. You can also call this function to load a smaller time span; for example, calling it with the two file names:

ascp02050.432, ascp02150.432

will cover a limited range from JDE 2469776.5 to 2542864.5.

◆ ComponentsForItem()

static int pcl::JPLEphemeris::ComponentsForItem ( int  i)
inlinestatic

Returns the number of components calculated for the specified ephemeris item:

3 for the position and velocity of all major solar system bodies: Sun, Mercury to Pluto, and Moon.

2 for nutation angles.

3 for lunar libration angles and lunar mantle velocity.

1 for the difference TT-TDB.

Definition at line 366 of file JPLEphemeris.h.

◆ ComputeState() [1/4]

void pcl::JPLEphemeris::ComputeState ( Vector r,
double  jd1,
double  jd2,
item_index  i 
) const
inline

Computes a state vector.

Parameters
[out]rReference to a vector where the components of the computed state will be stored.
jd1,jd2The requested time point in the TDB time scale, equal to jd1+jd2.
iIndex of the requested ephemeris item.

Rectangular position coordinates are provided in au, except for the geocentric Moon, whose position is provided in km.

Angles are provided in radians.

TT-TDB differences are provided in seconds.

Definition at line 389 of file JPLEphemeris.h.

References pcl::GenericVector< T >::Begin().

Referenced by ComputeState().

◆ ComputeState() [2/4]

void pcl::JPLEphemeris::ComputeState ( Vector r,
TimePoint  t,
item_index  i 
) const
inline

Computes a state vector for the specified time coordinate t.

Calling this member function is equivalent to:

ComputeState( r, t.JDI(), t.JDF(), i );
void ComputeState(Vector &r, double jd1, double jd2, item_index i) const
Definition: JPLEphemeris.h:389

Definition at line 477 of file JPLEphemeris.h.

References ComputeState(), pcl::TimePoint::JDF(), and pcl::TimePoint::JDI().

◆ ComputeState() [3/4]

void pcl::JPLEphemeris::ComputeState ( Vector r,
Vector v,
double  jd1,
double  jd2,
item_index  i 
) const
inline

Computes a state vector and its first derivative.

Parameters
[out]rReference to a vector where the components of the computed state will be stored.
[out]vReference to a vector where the components of the computed first derivative will be stored.
jd1,jd2The requested time point in the TDB time scale, equal to jd1+jd2.
iIndex of the requested ephemeris item.

Rectangular position and velocity coordinates are provided in au and au/day, respectively, except for the geocentric Moon, whose position and velocity are provided in km and km/day, respectively.

Angles and their variations are provided in radians and radians/day.

TT-TDB differences and their variations are provided in seconds and seconds/day.

Definition at line 438 of file JPLEphemeris.h.

References pcl::GenericVector< T >::Begin().

◆ ComputeState() [4/4]

void pcl::JPLEphemeris::ComputeState ( Vector r,
Vector v,
TimePoint  t,
item_index  i 
) const
inline

Computes a state vector and its first derivative for the specified time coordinate t.

Calling this member function is equivalent to:

ComputeState( r, v, t.JDI(), t.JDF(), i );

Definition at line 490 of file JPLEphemeris.h.

References ComputeState(), pcl::TimePoint::JDF(), and pcl::TimePoint::JDI().

◆ Constants()

const constant_list& pcl::JPLEphemeris::Constants ( ) const
inline

Returns a reference to the list of numerical integration constants used by this DE/LE ephemeris. The integration constants are loaded from an ASCII header file by the class constructor.

The returned list is sorted by constant name in ascending order.

Definition at line 321 of file JPLEphemeris.h.

◆ ConstantValue()

double pcl::JPLEphemeris::ConstantValue ( const IsoString name) const
inline

Returns the value of an integration constant given by its name.

If no integration constant is available with the specified name (case-insensitive), this function throws an Error exception.

Definition at line 332 of file JPLEphemeris.h.

References pcl::Array< T, A >::Begin(), pcl::BinarySearch(), and pcl::Array< T, A >::End().

Referenced by DENumber(), and LENumber().

◆ DENumber()

int pcl::JPLEphemeris::DENumber ( ) const
inline

Returns the DE (Development Ephemeris) number for the data loaded by this object. The returned number is the value of the "DENUM" integration constant.

Definition at line 245 of file JPLEphemeris.h.

References ConstantValue(), and pcl::TruncInt().

◆ EndJD()

double pcl::JPLEphemeris::EndJD ( ) const
inline

Returns a Julian date corresponding to the last date covered by the ephemeris data loaded in this object.

The value returned by this function is the largest time point for which this object can compute ephemerides. The effective ephemeris time span is defined by the data loaded in successive calls to AddData().

Definition at line 309 of file JPLEphemeris.h.

◆ EphemerisEndJD()

double pcl::JPLEphemeris::EphemerisEndJD ( ) const
inline

Returns a Julian date corresponding to the last date covered by the JPL ephemeris for which this object has been initialized.

The value returned by this function is for informative purposes only. Note that the ephemeris time span does not have to be the same as the effective time span covered by the data loaded in this object. See the AddData() member function for more information on ephemeris data.

Definition at line 283 of file JPLEphemeris.h.

◆ EphemerisStartJD()

double pcl::JPLEphemeris::EphemerisStartJD ( ) const
inline

Returns a Julian date corresponding to the first date covered by the JPL ephemeris for which this object has been initialized.

The value returned by this function is for informative purposes only. Note that the ephemeris time span does not have to be the same as the effective time span covered by the data loaded in this object. See the AddData() member function for more information on ephemeris data.

Definition at line 269 of file JPLEphemeris.h.

◆ IsItemAvailable()

bool pcl::JPLEphemeris::IsItemAvailable ( int  i) const
inline

Returns true if the specified ephemeris item is available in the data loaded by this object. See the JPLEphemerisItem namespace for information on ephemeris items.

Definition at line 346 of file JPLEphemeris.h.

◆ LENumber()

int pcl::JPLEphemeris::LENumber ( ) const
inline

Returns the LE (Lunar Ephemeris) number for the data loaded by this object. The returned number is the value of the "LENUM" integration constant.

Definition at line 255 of file JPLEphemeris.h.

References ConstantValue(), and pcl::TruncInt().

◆ StartJD()

double pcl::JPLEphemeris::StartJD ( ) const
inline

Returns a Julian date corresponding to the first date covered by the ephemeris data loaded in this object.

The value returned by this function is the smallest time point for which this object can compute ephemerides. The effective ephemeris time span is defined by the data loaded in successive calls to AddData().

Definition at line 296 of file JPLEphemeris.h.

◆ Summary()

IsoString pcl::JPLEphemeris::Summary ( ) const

Generates a summary of the ephemeris loaded by this object, giving information on general ephemeris properties and a list of integration constants.

The generated information is organized as a sequence of text lines terminated with UNIX end-of-line characters.

◆ Test()

bool pcl::JPLEphemeris::Test ( const String filePath,
bool  verbose = true 
) const

Verifies validity of the calculated ephemerides for a standard JPL DE/LE test file in ASCII format. Returns true iff the test was performed and it was successful.

A test is successful if the computed coordinates and components do not differ from the values provided by the test file by more than 1e-13 in absolute value. If this member function returns true, you can be confident that the ephemerides being calculated with this object are correct.

The specified file must be a standard "testpo" JPL ephemeris file corresponding to the DE/LE ephemeris loaded by this object. For example, for DE/LE 432, the required test file is "testpo.432".

The test is only performed within the time span covered by this object, and unavailable target and center items are ignored. See the AddData() member function for more information.

If the verbose argument is true, the function will output diagnostics information to stdout.


The documentation for this class was generated from the following file: