Program Listing for File spice.cc¶
↰ Return to documentation for file (interfaces/spice.cc)
#include "lupnt/interfaces/spice.h"
#include <cspice/SpiceUsr.h>
#include <cspice/SpiceZfc.h>
#include <string.h>
#include <cstdlib>
#include <ctime>
#include <filesystem>
#include <fstream>
#include <iostream>
#include <mutex>
#include <optional>
#include <regex>
#include <sstream>
#include <stdexcept>
#include "lupnt/conversions/coordinate_conversions.h"
#include "lupnt/conversions/time_conversions.h"
#include "lupnt/core/error.h"
#include "lupnt/core/file.h"
#include "lupnt/core/logger.h"
#include "lupnt/environment/body.h"
#include "lupnt/interfaces/spice_cheby.h"
#include "lupnt/numerics/math_utils.h"
namespace lupnt {
namespace spice {
static segment_t* cheby_s;
static long cheby_n;
bool spice_loaded = false;
void CheckSpiceFailure(const std::string& context) {
if (!failed_c()) return;
SpiceChar short_msg[1841];
SpiceChar long_msg[1841];
getmsg_c("SHORT", sizeof(short_msg), short_msg);
getmsg_c("LONG", sizeof(long_msg), long_msg);
reset_c();
throw std::runtime_error(fmt::format("{} failed: {} {}", context, short_msg, long_msg));
}
void SetSpiceReturnMode() {
SpiceChar action[] = "RETURN";
erract_c("SET", 0, action);
}
namespace {
constexpr const char* kNaifGenericKernelsUrl
= "https://naif.jpl.nasa.gov/pub/naif/generic_kernels";
bool IsKernelDownloadDisabled() {
const char* val = std::getenv("LUPNT_SKIP_SPICE_KERNEL_DOWNLOAD");
if (val == nullptr) return false;
std::string s(val);
return !(s.empty() || s == "0" || s == "false" || s == "FALSE" || s == "OFF" || s == "off");
}
bool RunShellCommand(const std::string& cmd) { return std::system(cmd.c_str()) == 0; }
bool DownloadToFile(const std::string& url, const std::filesystem::path& dest_path) {
if (IsKernelDownloadDisabled()) return false;
std::filesystem::path tmp_path = dest_path;
tmp_path += ".part";
std::error_code ec;
std::filesystem::remove(tmp_path, ec);
std::string cmd
= fmt::format("curl -fsSL --connect-timeout 10 --max-time 600 -o \"{}\" \"{}\"",
tmp_path.string(), url);
bool ok = RunShellCommand(cmd) && std::filesystem::exists(tmp_path)
&& std::filesystem::file_size(tmp_path) > 0;
if (!ok) {
std::filesystem::remove(tmp_path, ec);
return false;
}
std::filesystem::rename(tmp_path, dest_path, ec);
if (ec) {
// `rename` can fail across filesystems; fall back to copy + remove.
ec.clear();
std::filesystem::copy_file(tmp_path, dest_path,
std::filesystem::copy_options::overwrite_existing, ec);
std::filesystem::remove(tmp_path, ec);
}
return !ec;
}
std::optional<std::string> FetchTextFromUrl(const std::string& url) {
if (IsKernelDownloadDisabled()) return std::nullopt;
auto tmp_path
= std::filesystem::temp_directory_path()
/ fmt::format("lupnt_naif_listing_{}.html", static_cast<long>(std::time(nullptr)));
std::string cmd
= fmt::format("curl -fsSL --connect-timeout 10 --max-time 60 -o \"{}\" \"{}\"",
tmp_path.string(), url);
std::error_code ec;
if (!RunShellCommand(cmd) || !std::filesystem::exists(tmp_path)) {
std::filesystem::remove(tmp_path, ec);
return std::nullopt;
}
std::ifstream ifs(tmp_path);
std::stringstream ss;
ss << ifs.rdbuf();
ifs.close();
std::filesystem::remove(tmp_path, ec);
return ss.str();
}
std::optional<std::string> FindLatestFilenameInListing(const std::string& listing_url,
const std::regex& pattern) {
auto contents = FetchTextFromUrl(listing_url);
if (!contents.has_value()) return std::nullopt;
std::string latest;
for (auto it = std::sregex_iterator(contents->begin(), contents->end(), pattern);
it != std::sregex_iterator(); ++it) {
const std::string match = it->str();
if (match > latest) latest = match;
}
if (latest.empty()) return std::nullopt;
return latest;
}
std::optional<std::string> FindLatestFilenameLocally(const std::filesystem::path& dir,
const std::regex& pattern) {
std::string latest;
std::error_code ec;
for (const auto& entry : std::filesystem::directory_iterator(dir, ec)) {
if (!entry.is_regular_file()) continue;
std::string name = entry.path().filename().string();
if (std::regex_match(name, pattern) && name > latest) latest = name;
}
if (latest.empty()) return std::nullopt;
return latest;
}
bool RefreshAndLoadKernel(const std::string& filename, const std::string& url,
const std::filesystem::path& kernel_dir) {
std::filesystem::path local_path = kernel_dir / filename;
bool had_local_copy = std::filesystem::exists(local_path);
if (DownloadToFile(url, local_path)) {
Logger::Info(fmt::format("SPICE: downloaded latest NAIF kernel '{}'", filename), "Spice");
} else if (had_local_copy) {
Logger::Debug(
fmt::format("SPICE: could not refresh '{}' from NAIF; using cached local copy",
filename),
"Spice");
} else {
Logger::Warn(fmt::format("SPICE: could not download '{}' from NAIF ({}) and no local "
"copy is available; this kernel will not be loaded",
filename, url),
"Spice");
return false;
}
furnsh_c(filename.c_str());
CheckSpiceFailure(fmt::format("Loading {}", filename));
return true;
}
} // namespace
void LoadSpiceKernel(void) {
if (spice_loaded) return;
SpiceInt kcount;
ktotal_c("ALL", &kcount);
if (kcount > 0) return;
SetSpiceReturnMode();
std::string orig_dir = std::filesystem::current_path().string();
std::filesystem::path kernel_dir = GetCspiceKernelDir();
std::filesystem::current_path(kernel_dir);
furnsh_c("naif0012.tls"); // leap seconds
CheckSpiceFailure("Loading naif0012.tls");
furnsh_c("de440.bsp"); // planetary ephemeris
CheckSpiceFailure("Loading de440.bsp");
furnsh_c("pck00011.tpc"); // planetary constants
CheckSpiceFailure("Loading pck00011.tpc");
// ----------------------------------------------------------------
// High-accuracy Earth and lunar orientation kernels
//
// The default text PCK ("pck00011.tpc") only provides the low-accuracy
// IAU_EARTH / IAU_MOON body-fixed frames (errors of roughly 0.06 deg /
// 6 km for Earth, and up to ~0.005 deg / ~155 m for the Moon).
// High-accuracy orientation requires binary PCKs that define the
// ITRF93 (Earth) and MOON_PA_DExxx / MOON_ME_DExxx (Moon) frames; see
// the NAIF tutorial "'High Accuracy' Orientation and Body-fixed Frames
// for the Moon and Earth" (April 2023).
//
// NAIF continually republishes these kernels as new orientation data
// becomes available, so during initialization LuPNT looks up and
// downloads the *current latest* files directly from the NAIF generic
// kernels server, caching them under the local kernel directory and
// falling back to the cached copy (or older bundled kernels) if the
// server cannot be reached. Set LUPNT_SKIP_SPICE_KERNEL_DOWNLOAD=1 to
// force fully-offline operation using only locally-cached kernels.
// ----------------------------------------------------------------
const std::string kPckUrl = std::string(kNaifGenericKernelsUrl) + "/pck/";
const std::string kLunarFkUrl = std::string(kNaifGenericKernelsUrl) + "/fk/satellites/";
// --- Earth --------------------------------------------------------
// NAIF publishes the latest reconstructed + short-term-predicted
// high-accuracy Earth orientation under a stable file name that never
// changes ("earth_latest_high_prec.bpc"), making it trivial to always
// fetch the most recent data.
const std::string kEarthLatestPck = "earth_latest_high_prec.bpc";
bool earth_loaded
= RefreshAndLoadKernel(kEarthLatestPck, kPckUrl + kEarthLatestPck, kernel_dir);
if (!earth_loaded) {
// Fall back to whichever (older, hard-coded) Earth orientation
// kernels may already be cached locally from a previous LuPNT
// release, so that LuPNT keeps working fully offline.
for (const char* fallback :
{"earth_200101_990825_predict.bpc", "earth_000101_241014_240722.bpc"}) {
if (std::filesystem::exists(fallback)) {
furnsh_c(fallback);
CheckSpiceFailure(fmt::format("Loading {}", fallback));
}
}
}
// --- Moon ---------------------------------------------------------
// Unlike Earth, NAIF does not publish the latest lunar orientation
// kernels under a stable file name: a new binary PCK + frame kernel
// (FK) pair is released -- named after the underlying JPL DE/LE
// ephemeris model and release date -- each time a new ephemeris is
// produced (e.g. "moon_pa_de440_200625.bpc" + "moon_de440_250416.tf").
// Determine the current latest pair by inspecting the NAIF directory
// listings: file names embed zero-padded DE version numbers and
// dates, so the lexicographically-largest match is also the most
// recent. If the listings cannot be retrieved, fall back to the
// newest matching kernel already cached locally.
static const std::regex kMoonPckPattern(R"(moon_pa_de\d+_[0-9]{6}\.bpc)");
static const std::regex kMoonFkPattern(R"(moon_de\d+_[0-9]{6}\.tf)");
auto moon_pck_name = FindLatestFilenameInListing(kPckUrl, kMoonPckPattern);
if (!moon_pck_name.has_value())
moon_pck_name = FindLatestFilenameLocally(kernel_dir, kMoonPckPattern);
auto moon_fk_name = FindLatestFilenameInListing(kLunarFkUrl, kMoonFkPattern);
if (!moon_fk_name.has_value())
moon_fk_name = FindLatestFilenameLocally(kernel_dir, kMoonFkPattern);
// Load the binary PCK before the frame kernel: the FK connects the
// MOON_PA_DExxx frame name to the orientation data provided by the PCK.
if (moon_pck_name.has_value()) {
RefreshAndLoadKernel(*moon_pck_name, kPckUrl + *moon_pck_name, kernel_dir);
} else {
Logger::Warn(
"SPICE: could not determine the latest lunar orientation binary PCK from "
"NAIF; high-accuracy lunar orientation (MOON_PA / MOON_ME) will be "
"unavailable",
"Spice");
}
if (moon_fk_name.has_value()) {
RefreshAndLoadKernel(*moon_fk_name, kLunarFkUrl + *moon_fk_name, kernel_dir);
} else {
Logger::Warn(
"SPICE: could not determine the latest lunar frame kernel (FK) from NAIF; "
"high-accuracy lunar orientation (MOON_PA / MOON_ME) will be unavailable",
"Spice");
}
// Make the high-accuracy MOON_PA frame the default lunar body-fixed
// frame for the legacy (pre-N0062) SPICE APIs (LSPCN, ET2LST, ILLUM,
// SRFXPT, SUBPT, SUBSOL, ...).
const std::string kMoonAssocPa = "moon_assoc_pa.tf";
RefreshAndLoadKernel(kMoonAssocPa, kLunarFkUrl + kMoonAssocPa, kernel_dir);
// Mars
if (std::filesystem::exists("mars097.bsp")) {
furnsh_c("mars097.bsp");
CheckSpiceFailure("Loading mars097.bsp");
}
// Load Chebyshev coefficients
cheby_s = spk_extract("de440.bsp", &cheby_n);
if (cheby_s == nullptr) {
throw std::runtime_error(
"Could not load SPK file - Please Download the SPK file. See "
"data/ephemeris/readme.md for instructions");
}
std::filesystem::current_path(orig_dir);
spice_loaded = true;
return;
}
void LoadSpiceKernel(const std::filesystem::path& filepath) {
LoadSpiceKernel();
SetSpiceReturnMode();
furnsh_c(filepath.string().c_str());
CheckSpiceFailure(fmt::format("Loading {}", filepath.string()));
}
int GetNaifId(const std::string& body_name) {
LoadSpiceKernel();
SpiceInt code;
SpiceBoolean found;
#pragma omp critical
{
bodn2c_c(body_name.c_str(), &code, &found);
}
CheckSpiceFailure(fmt::format("Resolving NAIF body {}", body_name));
if (!found) throw std::runtime_error(fmt::format("NAIF body not found: {}", body_name));
return static_cast<int>(code);
}
std::string GetNaifName(int naif_id) {
LoadSpiceKernel();
SpiceChar name[256];
SpiceBoolean found;
#pragma omp critical
{
bodc2n_c(static_cast<SpiceInt>(naif_id), sizeof(name), name, &found);
}
CheckSpiceFailure(fmt::format("Resolving NAIF body id {}", naif_id));
if (!found) return std::to_string(naif_id);
return std::string(name);
}
bool HasNaifBody(const std::string& body_name) {
LoadSpiceKernel();
SpiceInt code;
SpiceBoolean found;
#pragma omp critical
{
bodn2c_c(body_name.c_str(), &code, &found);
}
CheckSpiceFailure(fmt::format("Resolving NAIF body {}", body_name));
return static_cast<bool>(found);
}
bool HasFrame(const std::string& frame_name) {
LoadSpiceKernel();
SpiceInt frame_code = 0;
#pragma omp critical
{
namfrm_c(frame_name.c_str(), &frame_code);
}
CheckSpiceFailure(fmt::format("Resolving NAIF frame {}", frame_name));
return frame_code != 0;
}
void ExtractPckCoeffs() {
int handle;
SpiceInt pck_handle;
ConstSpiceChar* pck_file = "../data/ephemeris/moon_pa_de440_200625.bpc";
double t_tdb = 8000.0;
int body = 301;
double descr[5];
char ident[40];
int found;
SpiceDouble record[120];
int rsize;
int pdeg;
// double ra;
// double dec;
// double w;
// double lambda;
static char bref[32];
double eulang[6];
// The fundamental quantities defined by PCK orientation models are actually
// Euler angles, not matrices. These Euler angles, which we call ``RA, DEC,
// and W,'' are related to the transformation operator returned from
// pxform_c by the equation rotate = [ W ] [ Pi/2 - DEC ] [ Pi/2 + RA ]
// 3 1 3
// To directly retrieve these angles, use the call:
// bodeul_( &body, &t_tdb, &ra, &dec, &w, &lambda );
// std::cout << "body: " << body << std::endl;
// std::cout << "t_tdb: " << t_tdb << std::endl;
// std::cout << "ra: " << ra << std::endl;
// std::cout << "dec: " << dec << std::endl;
// std::cout << "w: " << w << std::endl;
// std::cout << " " << std::endl;
pckeul_(&body, &t_tdb, &found, bref, eulang, (ftnlen)32);
std::cout << "found:" << found << std::endl;
std::cout << "phi: " << eulang[0] << std::endl;
std::cout << "delta: " << eulang[1] << std::endl;
std::cout << " " << std::endl;
pcklof_c(pck_file, &pck_handle); // load the PCK file
pcksfs_(&body, &t_tdb, &handle, descr, ident, &found, (ftnlen)40);
std::cout << "pck handle: :" << pck_handle << std::endl;
std::cout << "handle: :" << handle << std::endl;
std::cout << "descr: " << &descr << std::endl;
std::cout << "ident: " << &ident << std::endl;
std::cout << "found:" << found << std::endl;
if (found) {
pckr02_(&handle, descr, &t_tdb, record);
rsize = record[1];
pdeg = (rsize - 2) / 3 - 1;
std::cout << "Polynomial Size:" << rsize << std::endl;
std::cout << "Polynomial Degree:" << pdeg << std::endl;
}
// extract coefficients from the CSPICE PCK file
// pcksfs_(body, t_tdb, &handle, descr, ident, found, (ftnlen)40);
// pckr02_c(handle, target)
}
Vec3d GetBodyPosSpice(Real t_tdb, BodyId obs, BodyId target, const std::string& refFrame,
const std::string& abCorrection) {
if (!spice_loaded) {
LoadSpiceKernel();
}
std::string targ_str = std::to_string((int)target);
std::string obs_str = std::to_string((int)obs);
SpiceDouble ptarg[3];
SpiceDouble et = t_tdb.val();
const char* targ = strcpy(new char[targ_str.length() + 1], targ_str.c_str());
const char* ref = strcpy(new char[refFrame.length() + 1], refFrame.c_str());
const char* abcorr = strcpy(new char[abCorrection.length() + 1], abCorrection.c_str());
const char* obs_spice = strcpy(new char[obs_str.length() + 1], obs_str.c_str());
SpiceDouble lt;
// void spkpos_c(ConstSpiceChar * targ, SpiceDouble t_tdb, ConstSpiceChar *
// ref, ConstSpiceChar * abcorr,
// ConstSpiceChar * obs, SpiceDouble ptarg[3], SpiceDouble *
// lt)
spkpos_c(targ, et, ref, abcorr, obs_spice, ptarg, <);
Vec3d r;
for (int i = 0; i < 3; i++) r(i) = ptarg[i];
return r;
}
Vec6d GetBodyPosVelSpice(Real t_tdb, BodyId obs, BodyId target, const std::string& refFrame,
const std::string& abCorrection) {
if (!spice_loaded) LoadSpiceKernel();
SpiceDouble starg[6];
Vec6d rv;
std::string targ_str = std::to_string((int)target);
std::string obs_str = std::to_string((int)obs);
SpiceDouble et = t_tdb.val();
const char* targ = strcpy(new char[targ_str.length() + 1], targ_str.c_str());
const char* ref = strcpy(new char[refFrame.length() + 1], refFrame.c_str());
const char* abcorr = strcpy(new char[abCorrection.length() + 1], abCorrection.c_str());
const char* obs_spice = strcpy(new char[obs_str.length() + 1], obs_str.c_str());
SpiceDouble lt;
// void spkez_c ( SpiceInt targ,
// SpiceDouble et,
// ConstSpiceChar *ref,
// ConstSpiceChar *abcorr,
// SpiceInt obs,
// SpiceDouble starg[6],
// SpiceDouble *lt )
#pragma omp critical
{
spkezr_c(targ, et, ref, abcorr, obs_spice, starg, <);
}
for (int i = 0; i < 6; i++) rv(i) = starg[i];
return rv;
}
GroundStationSpiceData GetGroundStationDataSpice(Real t_tdb, const std::string& station_name,
BodyId center, const std::string& refFrame,
const std::string& abCorrection) {
int station_id = GetNaifId(station_name);
GroundStationSpiceData data
= GetGroundStationDataSpice(t_tdb, station_id, center, refFrame, abCorrection);
data.name = station_name;
return data;
}
GroundStationSpiceData GetGroundStationDataSpice(Real t_tdb, int station_id, BodyId center,
const std::string& refFrame,
const std::string& abCorrection) {
LoadSpiceKernel();
SetSpiceReturnMode();
SpiceDouble et = t_tdb.val();
SpiceDouble pos_km[3];
SpiceDouble lt;
std::string station = std::to_string(station_id);
std::string center_id = std::to_string(static_cast<int>(center));
#pragma omp critical
{
spkpos_c(station.c_str(), et, refFrame.c_str(), abCorrection.c_str(), center_id.c_str(),
pos_km, <);
}
CheckSpiceFailure(fmt::format("Fetching ground station {} from SPICE", station_id));
GroundStationSpiceData data;
data.naif_id = station_id;
data.name = GetNaifName(station_id);
data.body_id = center;
data.frame = refFrame;
for (int i = 0; i < 3; ++i) data.position_m(i) = pos_km[i] * M_KM;
BodyData body_data = GetBodyData(center);
Cart3 pos(data.position_m.cast<Real>(), body_data.fixed_frame);
LatLonAlt lla = CartToLatLonAlt(pos, body_data.R, body_data.flattening);
data.latitude_deg = lla(0).val();
data.longitude_deg = lla(1).val();
data.altitude_m = lla(2).val();
return data;
}
Mat6d GetFrameConversionMat(Real t_tdb, const std::string& from_frame,
const std::string& to_frame) {
if (!spice_loaded) LoadSpiceKernel();
SpiceDouble et_spice = (SpiceDouble)t_tdb.val();
double xform[6][6];
Mat6d M_rot;
const char* from_frame_char = strcpy(new char[from_frame.length() + 1], from_frame.c_str());
const char* to_frame_char = strcpy(new char[to_frame.length() + 1], to_frame.c_str());
#pragma omp critical
{
sxform_c(from_frame_char, to_frame_char, et_spice, xform);
}
for (int i = 0; i < 6; i++) {
for (int j = 0; j < 6; j++) {
M_rot(i, j) = xform[i][j];
}
}
return M_rot;
}
Vec3d GetPlanetOrientation(BodyId id, Real t_tdb) {
if (!spice_loaded) LoadSpiceKernel();
std::string to_frame = "J2000";
char bodyname[36];
SpiceBoolean found;
#pragma omp critical
{
bodc2n_c((SpiceInt)id, 36, bodyname, &found);
}
if (!found) {
throw std::runtime_error("Invalid planet ID");
}
std::string from_frame = "IAU_" + std::string(bodyname);
SpiceDouble rotmat[3][3];
SpiceDouble et_spice = (SpiceDouble)t_tdb.val();
const char* from_frame_char = strcpy(new char[from_frame.length() + 1], from_frame.c_str());
const char* to_frame_char = strcpy(new char[to_frame.length() + 1], to_frame.c_str());
#pragma omp critical
{
pxform_c(from_frame_char, to_frame_char, et_spice, rotmat);
}
SpiceDouble psi, theta, phi;
#pragma omp critical
m2eul_c(rotmat, 3, 1, 3, &psi, &theta, &phi);
Vec3d angles;
double W = -(phi + PI); // R_z(-W)
double delta0 = PI / 2 - theta; // R_x(PI/2 - delta0)
double alpha0 = PI / 2 - psi; // R_z(-alpha0)
angles << alpha0, delta0, W;
return angles;
}
Real StringToTdb(const std::string& str) {
if (!spice_loaded) LoadSpiceKernel();
SpiceDouble t_tdb;
str2et_c(str.c_str(), &t_tdb);
return t_tdb;
}
Real StringToTai(const std::string& str) {
if (!spice_loaded) LoadSpiceKernel();
Real t_tdb = StringToTdb(str);
Real t_tai = spice::ConvertTime(t_tdb, Time::TDB, Time::TAI);
return t_tai;
}
std::string TDBtoStringUTC(Real t_tdb, int prec = 3) {
if (!spice_loaded) LoadSpiceKernel();
SpiceDouble et = t_tdb.val();
SpiceChar str[100];
et2utc_c(et, "C", prec, 100, str);
return std::string(str);
}
std::string TAItoStringUTC(Real t_tai, int prec = 3) {
if (!spice_loaded) LoadSpiceKernel();
Real et_tdb = spice::ConvertTime(t_tai, Time::TAI, Time::TDB);
std::string str = TDBtoStringUTC(et_tdb, prec);
return str;
}
Real ConvertTime(Real t, Time from, Time to) {
if (!spice_loaded) LoadSpiceKernel();
if (from == to) return t;
SpiceDouble t_in = t.val();
SpiceDouble t_out_spice;
#pragma omp critical
{
t_out_spice
= unitim_c(t_in, time_to_string.at(from).c_str(), time_to_string.at(to).c_str());
}
double offset = t_out_spice - t_in; // offset in seconds
Real t_out = t + offset; // this is to convert to real
return t_out;
}
MatX6 GetBodyPosVel(const VecX& t_tdb, BodyId center, BodyId target) {
if (!spice_loaded) LoadSpiceKernel();
MatX6 retState(t_tdb.size(), 6);
for (int i = 0; i < t_tdb.size(); i++)
retState.row(i) = GetBodyPosVel(t_tdb(i), center, target).transpose();
return retState;
}
Vec6 GetBodyPosVelBase(const Real t_tdb, BodyId center, BodyId target) {
if (!spice_loaded) LoadSpiceKernel();
if (center == target) return Vec6::Zero();
for (int i = 0; i < cheby_n; i++) {
if (cheby_s[i].center == (int)center && cheby_s[i].target == (int)target)
return cheby_posvel_ad(t_tdb, cheby_s[i].seg, cheby_s[i].len);
if (cheby_s[i].center == (int)target && cheby_s[i].target == (int)center)
return -cheby_posvel_ad(t_tdb, cheby_s[i].seg, cheby_s[i].len);
}
LUPNT_CHECK(false, "Chebyshev coefficients not found", "SpiceInterface");
return Vec6::Zero();
}
Vec6 GetBodyPosVel(const Real t_tdb, BodyId center, BodyId target) {
if (!spice_loaded) LoadSpiceKernel();
if (center == target) return Vec6::Zero();
Vec6 rv = Vec6::Zero();
if (target == BodyId::EARTH) {
rv += GetBodyPosVelBase(t_tdb, BodyId::EMB, BodyId::EARTH);
target = BodyId::EMB;
}
if (center == BodyId::EARTH) {
rv -= GetBodyPosVelBase(t_tdb, BodyId::EMB, BodyId::EARTH);
center = BodyId::EMB;
}
if (target == BodyId::MOON) {
rv += GetBodyPosVelBase(t_tdb, BodyId::EMB, BodyId::MOON);
target = BodyId::EMB;
}
if (center == BodyId::MOON) {
rv -= GetBodyPosVelBase(t_tdb, BodyId::EMB, BodyId::MOON);
center = BodyId::EMB;
}
if (target == BodyId::MERCURY) {
rv += GetBodyPosVelBase(t_tdb, BodyId::MERCURY_BARYCENTER, BodyId::MERCURY);
target = BodyId::MERCURY_BARYCENTER;
}
if (center == BodyId::MERCURY) {
rv -= GetBodyPosVelBase(t_tdb, BodyId::MERCURY_BARYCENTER, BodyId::MERCURY);
center = BodyId::MERCURY_BARYCENTER;
}
if (target == BodyId::VENUS) {
rv += GetBodyPosVelBase(t_tdb, BodyId::VENUS_BARYCENTER, BodyId::VENUS);
target = BodyId::VENUS_BARYCENTER;
}
if (center == BodyId::VENUS) {
rv -= GetBodyPosVelBase(t_tdb, BodyId::VENUS_BARYCENTER, BodyId::VENUS);
center = BodyId::VENUS_BARYCENTER;
}
rv += GetBodyPosVelBase(t_tdb, center, target);
return rv;
}
} // namespace spice
} // namespace lupnt