/***************************************************************************************** * * * OpenSpace * * * * Copyright (c) 2014-2022 * * * * Permission is hereby granted, free of charge, to any person obtaining a copy of this * * software and associated documentation files (the "Software"), to deal in the Software * * without restriction, including without limitation the rights to use, copy, modify, * * merge, publish, distribute, sublicense, and/or sell copies of the Software, and to * * permit persons to whom the Software is furnished to do so, subject to the following * * conditions: * * * * The above copyright notice and this permission notice shall be included in all copies * * or substantial portions of the Software. * * * * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR IMPLIED, * * INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, FITNESS FOR A * * PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT * * HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF * * CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE * * OR THE USE OR OTHER DEALINGS IN THE SOFTWARE. * ****************************************************************************************/ #include #include #include #include #include #include #include #include #include #include #include #include "SpiceUsr.h" #include "SpiceZpr.h" namespace { constexpr const char* _loggerCat = "SpiceManager"; // The value comes from // http://naif.jpl.nasa.gov/pub/naif/toolkit_docs/C/cspice/getmsg_c.html // as the maximum message length constexpr const unsigned SpiceErrorBufferSize = 1841; const char* toString(openspace::SpiceManager::FieldOfViewMethod m) { using SM = openspace::SpiceManager; switch (m) { case SM::FieldOfViewMethod::Ellipsoid: return "ELLIPSOID"; case SM::FieldOfViewMethod::Point: return "POINT"; default: throw ghoul::MissingCaseException(); } } const char* toString(openspace::SpiceManager::TerminatorType t) { using SM = openspace::SpiceManager; switch (t) { case SM::TerminatorType::Umbral: return "UMBRAL"; case SM::TerminatorType::Penumbral: return "PENUMBRAL"; default: throw ghoul::MissingCaseException(); } } } // namespace #include "spicemanager_lua.inl" namespace openspace { SpiceManager* SpiceManager::_instance = nullptr; SpiceManager::SpiceException::SpiceException(std::string msg) : ghoul::RuntimeError(std::move(msg), "Spice") { ghoul_assert( SpiceManager::ref().exceptionHandling() == SpiceManager::UseException::Yes, "No exceptions should be thrown when UseException is No" ); } SpiceManager::AberrationCorrection::AberrationCorrection(Type t, Direction d) : type(t) , direction(d) {} SpiceManager::AberrationCorrection::AberrationCorrection(const std::string& identifier) { const static std::map> Mapping = { { "NONE" , { Type::None, Direction::Reception } }, { "LT" , { Type::LightTime, Direction::Reception } }, { "LT+S" , { Type::LightTimeStellar, Direction::Reception } }, { "CN" , { Type::ConvergedNewtonian, Direction::Reception } }, { "CN+S" , { Type::ConvergedNewtonianStellar, Direction::Reception } }, { "XLT" , { Type::LightTime, Direction::Transmission } }, { "XLT+S" , { Type::LightTimeStellar, Direction::Transmission } }, { "XCN" , { Type::ConvergedNewtonian, Direction::Transmission } }, { "XCN+S" , { Type::ConvergedNewtonianStellar, Direction::Transmission } } }; auto it = Mapping.find(identifier); ghoul_assert(!identifier.empty(), "Identifier may not be empty"); ghoul_assert(it != Mapping.end(), fmt::format("Invalid identifer '{}'", identifier)); type = it->second.first; direction = it->second.second; } SpiceManager::AberrationCorrection::operator const char*() const { switch (type) { case Type::None: return "NONE"; case Type::LightTime: return (direction == Direction::Reception) ? "LT" : "XLT"; case Type::LightTimeStellar: return (direction == Direction::Reception) ? "LT+S" : "XLT+S"; case Type::ConvergedNewtonian: return (direction == Direction::Reception) ? "CN" : "XCN"; case Type::ConvergedNewtonianStellar: return (direction == Direction::Reception) ? "CN+S" : "XCN+S"; default: throw ghoul::MissingCaseException(); } } SpiceManager::FieldOfViewMethod SpiceManager::fieldOfViewMethodFromString( const std::string& method) { const static std::map Mapping = { { "ELLIPSOID", FieldOfViewMethod::Ellipsoid }, { "POINT", FieldOfViewMethod::Point } }; ghoul_assert(!method.empty(), "Method must not be empty"); return Mapping.at(method); } SpiceManager::TerminatorType SpiceManager::terminatorTypeFromString( const std::string& type) { const static std::map Mapping = { { "UMBRAL", TerminatorType::Umbral }, { "PENUMBRAL", TerminatorType::Penumbral } }; ghoul_assert(!type.empty(), "Type must not be empty"); return Mapping.at(type); } SpiceManager::SpiceManager() { // Set the SPICE library to not exit the program if an error occurs erract_c("SET", 0, const_cast("REPORT")); // NOLINT // But we do not want SPICE to print the errors, we will fetch them ourselves errprt_c("SET", 0, const_cast("NONE")); // NOLINT } SpiceManager::~SpiceManager() { for (const KernelInformation& i : _loadedKernels) { unload_c(i.path.c_str()); } // Set values back to default erract_c("SET", 0, const_cast("DEFAULT")); // NOLINT errprt_c("SET", 0, const_cast("DEFAULT")); // NOLINT } void SpiceManager::initialize() { ghoul_assert(!isInitialized(), "SpiceManager is already initialized"); _instance = new SpiceManager; } void SpiceManager::deinitialize() { ghoul_assert(isInitialized(), "SpiceManager is not initialized"); delete _instance; _instance = nullptr; } bool SpiceManager::isInitialized() { return _instance != nullptr; } SpiceManager& SpiceManager::ref() { ghoul_assert(isInitialized(), "SpiceManager is not initialized"); return *_instance; } // This method checks if one of the previous SPICE methods has failed. If it has, an // exception with the SPICE error message is thrown // If an error occurred, true is returned, otherwise, false void throwSpiceError(const std::string& errorMessage) { if (openspace::SpiceManager::ref().exceptionHandling()) { char buffer[SpiceErrorBufferSize]; getmsg_c("LONG", SpiceErrorBufferSize, buffer); reset_c(); throw openspace::SpiceManager::SpiceException(errorMessage + ": " + buffer); } else { reset_c(); } } SpiceManager::KernelHandle SpiceManager::loadKernel(std::string filePath) { ghoul_assert(!filePath.empty(), "Empty file path"); ghoul_assert( std::filesystem::is_regular_file(filePath), fmt::format("File '{}' ({}) does not exist", filePath, absPath(filePath)) ); ghoul_assert( std::filesystem::is_directory(std::filesystem::path(filePath).parent_path()), fmt::format( "File {} exists, but directory {} does not", absPath(filePath), std::filesystem::path(filePath).parent_path() ) ); std::filesystem::path path = absPath(std::move(filePath)); const auto it = std::find_if( _loadedKernels.begin(), _loadedKernels.end(), [path](const KernelInformation& info) { return info.path == path; } ); if (it != _loadedKernels.end()) { it->refCount++; return it->id; } // We need to set the current directory as meta-kernels are usually defined relative // to the directory they reside in. The directory change is not necessary for regular // kernels std::filesystem::path currentDirectory = std::filesystem::current_path(); std::filesystem::path p = path.parent_path(); std::filesystem::current_path(p); LINFO(fmt::format("Loading SPICE kernel {}", path)); // Load the kernel furnsh_c(path.string().c_str()); // Reset the current directory to the previous one std::filesystem::current_path(currentDirectory); if (failed_c()) { throwSpiceError("Kernel loading"); } std::filesystem::path fileExtension = path.extension(); if (fileExtension == ".bc" || fileExtension == ".BC") { findCkCoverage(path.string()); // binary ck kernel } else if (fileExtension == ".bsp" || fileExtension == ".BSP") { findSpkCoverage(path.string()); // binary spk kernel } KernelHandle kernelId = ++_lastAssignedKernel; ghoul_assert(kernelId != 0, fmt::format("Kernel Handle wrapped around to 0")); _loadedKernels.push_back({ path.string(), kernelId, 1 }); return kernelId; } void SpiceManager::unloadKernel(KernelHandle kernelId) { ghoul_assert(kernelId <= _lastAssignedKernel, "Invalid unassigned kernel"); ghoul_assert(kernelId != KernelHandle(0), "Invalid zero handle"); const auto it = std::find_if( _loadedKernels.begin(), _loadedKernels.end(), [&kernelId](const KernelInformation& info) { return info.id == kernelId; } ); if (it != _loadedKernels.end()) { // If there was only one part interested in the kernel, we can unload it if (it->refCount == 1) { // No need to check for errors as we do not allow empty path names LINFO(fmt::format("Unloading SPICE kernel {}", it->path)); unload_c(it->path.c_str()); _loadedKernels.erase(it); } // Otherwise, we hold on to it, but reduce the reference counter by 1 else { it->refCount--; LDEBUG(fmt::format("Reducing reference counter to: {}", it->refCount)); } } } void SpiceManager::unloadKernel(std::string filePath) { ghoul_assert(!filePath.empty(), "Empty filename"); std::filesystem::path path = absPath(std::move(filePath)); const auto it = std::find_if( _loadedKernels.begin(), _loadedKernels.end(), [&path](const KernelInformation& info) { return info.path == path; } ); if (it == _loadedKernels.end()) { if (_useExceptions) { throw SpiceException( fmt::format("{} did not correspond to a loaded kernel", path) ); } else { return; } } else { // If there was only one part interested in the kernel, we can unload it if (it->refCount == 1) { LINFO(fmt::format("Unloading SPICE kernel {}", path)); unload_c(path.string().c_str()); _loadedKernels.erase(it); } else { // Otherwise, we hold on to it, but reduce the reference counter by 1 it->refCount--; LDEBUG(fmt::format("Reducing reference counter to: {}", it->refCount)); } } } bool SpiceManager::hasSpkCoverage(const std::string& target, double et) const { ghoul_assert(!target.empty(), "Empty target"); const int id = naifId(target); // SOLAR SYSTEM BARYCENTER special case, implicitly included by Spice if (id == 0) { return true; } const auto it = _spkIntervals.find(id); if (it != _spkIntervals.end()) { const std::vector>& intervalVector = it->second; for (const std::pair& vecElement : intervalVector) { if ((vecElement.first < et) && (vecElement.second > et)) { return true; } } } return false; } std::vector> SpiceManager::spkCoverage( const std::string& target) const { ghoul_assert(!target.empty(), "Empty target"); const int id = naifId(target); const auto it = _spkIntervals.find(id); if (it != _spkIntervals.end()) { return it->second; } else { std::vector> emptyList; return emptyList; } } bool SpiceManager::hasCkCoverage(const std::string& frame, double et) const { ghoul_assert(!frame.empty(), "Empty target"); const int id = frameId(frame); const auto it = _ckIntervals.find(id); if (it != _ckIntervals.end()) { const std::vector>& intervalVector = it->second; for (const std::pair& i : intervalVector) { if ((i.first < et) && (i.second > et)) { return true; } } } return false; } std::vector> SpiceManager::ckCoverage( const std::string& target) const { ghoul_assert(!target.empty(), "Empty target"); int id = naifId(target); const auto it = _ckIntervals.find(id); if (it != _ckIntervals.end()) { return it->second; } else { id *= 1000; const auto it2 = _ckIntervals.find(id); if (it2 != _ckIntervals.end()) { return it2->second; } else { std::vector> emptyList; return emptyList; } } } std::vector> SpiceManager::spiceBodies( bool builtInFrames) const { std::vector> bodies; constexpr const int Frnmln = 33; static SpiceInt idsetBuffer[SPICE_CELL_CTRLSZ + 8192]; static SpiceCell idset = { SPICE_INT, 0, 8192, 0, SPICETRUE, SPICEFALSE, SPICEFALSE, &idsetBuffer, &(idsetBuffer[SPICE_CELL_CTRLSZ]) }; SpiceChar frname[Frnmln]; for (SpiceInt i = 1; i <= 6; i++) { if (i < 6) { if (builtInFrames) { bltfrm_c(i, &idset); } else { kplfrm_c(i, &idset); } } else { if (builtInFrames) { bltfrm_c(SPICE_FRMTYP_ALL, &idset); } else { kplfrm_c(SPICE_FRMTYP_ALL, &idset); } } for (SpiceInt j = 0; j < card_c(&idset); j++) { frmnam_c( (reinterpret_cast(idset.data))[j], Frnmln, frname ); bodies.push_back( std::make_pair( static_cast(reinterpret_cast(idset.data)[j]), frname ) ); } } return bodies; } bool SpiceManager::hasValue(int naifId, const std::string& item) const { return bodfnd_c(naifId, item.c_str()); } bool SpiceManager::hasValue(const std::string& body, const std::string& item) const { ghoul_assert(!body.empty(), "Empty body"); ghoul_assert(!item.empty(), "Empty item"); int id = naifId(body); return hasValue(id, item); } int SpiceManager::naifId(const std::string& body) const { ghoul_assert(!body.empty(), "Empty body"); SpiceBoolean success; SpiceInt id; bods2c_c(body.c_str(), &id, &success); if (!success && _useExceptions) { throw SpiceException(fmt::format("Could not find NAIF ID of body '{}'", body)); } return id; } bool SpiceManager::hasNaifId(const std::string& body) const { ghoul_assert(!body.empty(), "Empty body"); SpiceBoolean success; SpiceInt id; bods2c_c(body.c_str(), &id, &success); reset_c(); return success; } int SpiceManager::frameId(const std::string& frame) const { ghoul_assert(!frame.empty(), "Empty frame"); SpiceInt id; namfrm_c(frame.c_str(), &id); if (id == 0 && _useExceptions) { throw SpiceException(fmt::format("Could not find NAIF ID of frame '{}'", frame)); } return id; } bool SpiceManager::hasFrameId(const std::string& frame) const { ghoul_assert(!frame.empty(), "Empty frame"); SpiceInt id; namfrm_c(frame.c_str(), &id); return id != 0; } void getValueInternal(const std::string& body, const std::string& value, int size, double* v) { ghoul_assert(!body.empty(), "Empty body"); ghoul_assert(!value.empty(), "Empty value"); ghoul_assert(v != nullptr, "Empty value pointer"); SpiceInt n; bodvrd_c(body.c_str(), value.c_str(), size, &n, v); if (failed_c()) { throwSpiceError( fmt::format("Error getting value '{}' for body '{}'", value, body) ); } } void SpiceManager::getValue(const std::string& body, const std::string& value, double& v) const { getValueInternal(body, value, 1, &v); } void SpiceManager::getValue(const std::string& body, const std::string& value, glm::dvec2& v) const { getValueInternal(body, value, 2, glm::value_ptr(v)); } void SpiceManager::getValue(const std::string& body, const std::string& value, glm::dvec3& v) const { getValueInternal(body, value, 3, glm::value_ptr(v)); } void SpiceManager::getValue(const std::string& body, const std::string& value, glm::dvec4& v) const { getValueInternal(body, value, 4, glm::value_ptr(v)); } void SpiceManager::getValue(const std::string& body, const std::string& value, std::vector& v) const { ghoul_assert(!v.empty(), "Array for values has to be preallocaed"); getValueInternal(body, value, static_cast(v.size()), v.data()); } double SpiceManager::spacecraftClockToET(const std::string& craft, double craftTicks) { ghoul_assert(!craft.empty(), "Empty craft"); int craftId = naifId(craft); double et; sct2e_c(craftId, craftTicks, &et); if (failed_c()) { throwSpiceError(fmt::format( "Error transforming spacecraft clock of '{}' at time {}", craft, craftTicks )); } return et; } double SpiceManager::ephemerisTimeFromDate(const std::string& timeString) const { ghoul_assert(!timeString.empty(), "Empty timeString"); return ephemerisTimeFromDate(timeString.c_str()); } double SpiceManager::ephemerisTimeFromDate(const char* timeString) const { double et; str2et_c(timeString, &et); if (failed_c()) { throwSpiceError(fmt::format("Error converting date '{}'", timeString)); } return et; } std::string SpiceManager::dateFromEphemerisTime(double ephemerisTime, const char* format) { char Buffer[128]; std::memset(Buffer, char(0), 128); timout_c(ephemerisTime, format, 128, Buffer); if (failed_c()) { throwSpiceError(fmt::format( "Error converting ephemeris time '{}' to date with format '{}'", ephemerisTime, format )); } return std::string(Buffer); } glm::dvec3 SpiceManager::targetPosition(const std::string& target, const std::string& observer, const std::string& referenceFrame, AberrationCorrection aberrationCorrection, double ephemerisTime, double& lightTime) const { ghoul_assert(!target.empty(), "Target is not empty"); ghoul_assert(!observer.empty(), "Observer is not empty"); ghoul_assert(!referenceFrame.empty(), "Reference frame is not empty"); bool targetHasCoverage = hasSpkCoverage(target, ephemerisTime); bool observerHasCoverage = hasSpkCoverage(observer, ephemerisTime); if (!targetHasCoverage && !observerHasCoverage) { if (_useExceptions) { throw SpiceException( fmt::format( "Neither target '{}' nor observer '{}' has SPK coverage at time {}", target, observer, ephemerisTime ) ); } else { return glm::dvec3(0.0); } } else if (targetHasCoverage && observerHasCoverage) { glm::dvec3 position = glm::dvec3(0.0); spkpos_c( target.c_str(), ephemerisTime, referenceFrame.c_str(), aberrationCorrection, observer.c_str(), glm::value_ptr(position), &lightTime ); if (failed_c()) { throwSpiceError(fmt::format( "Error getting position from '{}' to '{}' in frame '{}' at time {}", target, observer, referenceFrame, ephemerisTime )); } return position; } else if (targetHasCoverage) { // observer has no coverage return getEstimatedPosition( observer, target, referenceFrame, aberrationCorrection, ephemerisTime, lightTime ) * -1.0; } else { // target has no coverage return getEstimatedPosition( target, observer, referenceFrame, aberrationCorrection, ephemerisTime, lightTime ); } } glm::dvec3 SpiceManager::targetPosition(const std::string& target, const std::string& observer, const std::string& referenceFrame, AberrationCorrection aberrationCorrection, double ephemerisTime) const { double unused = 0.0; return targetPosition( target, observer, referenceFrame, aberrationCorrection, ephemerisTime, unused ); } glm::dmat3 SpiceManager::frameTransformationMatrix(const std::string& from, const std::string& to, double ephemerisTime) const { ghoul_assert(!from.empty(), "From must not be empty"); ghoul_assert(!to.empty(), "To must not be empty"); // get rotation matrix from frame A - frame B glm::dmat3 transform = glm::dmat3(1.0); pxform_c( from.c_str(), to.c_str(), ephemerisTime, reinterpret_cast(glm::value_ptr(transform)) ); if (failed_c()) { throwSpiceError( fmt::format("Error converting from frame '{}' to frame '{}' at time '{}'", from, to, ephemerisTime ) ); } // The rox-major, column-major order are switched in GLM and SPICE, so we have to // transpose the matrix before we can return it return glm::transpose(transform); } SpiceManager::SurfaceInterceptResult SpiceManager::surfaceIntercept( const std::string& target, const std::string& observer, const std::string& fovFrame, const std::string& referenceFrame, AberrationCorrection aberrationCorrection, double ephemerisTime, const glm::dvec3& directionVector) const { ghoul_assert(!target.empty(), "Target must not be empty"); ghoul_assert(!observer.empty(), "Observer must not be empty"); ghoul_assert(target != observer, "Target and observer must be different"); ghoul_assert(!fovFrame.empty(), "FOV frame must not be empty"); ghoul_assert(!referenceFrame.empty(), "Reference frame must not be empty"); ghoul_assert(directionVector != glm::dvec3(0.0), "Direction vector must not be zero"); const std::string ComputationMethod = "ELLIPSOID"; SurfaceInterceptResult result; SpiceBoolean found; sincpt_c(ComputationMethod.c_str(), target.c_str(), ephemerisTime, referenceFrame.c_str(), aberrationCorrection, observer.c_str(), fovFrame.c_str(), glm::value_ptr(directionVector), glm::value_ptr(result.surfaceIntercept), &result.interceptEpoch, glm::value_ptr(result.surfaceVector), &found ); result.interceptFound = (found == SPICETRUE); if (failed_c()) { throwSpiceError(fmt::format( "Error retrieving surface intercept on target '{}' viewed from observer '{}' " "in reference frame '{}' at time '{}'", target, observer, referenceFrame, ephemerisTime )); } return result; } bool SpiceManager::isTargetInFieldOfView(const std::string& target, const std::string& observer, const std::string& referenceFrame, const std::string& instrument, FieldOfViewMethod method, AberrationCorrection aberrationCorrection, double& ephemerisTime) const { ghoul_assert(!target.empty(), "Target must not be empty"); ghoul_assert(!observer.empty(), "Observer must not be empty"); ghoul_assert(target != observer, "Target and observer must be different"); ghoul_assert(!referenceFrame.empty(), "Reference frame must not be empty"); ghoul_assert(!instrument.empty(), "Instrument must not be empty"); int visible; fovtrg_c(instrument.c_str(), target.c_str(), toString(method), referenceFrame.c_str(), aberrationCorrection, observer.c_str(), &ephemerisTime, &visible ); if (failed_c()) { throwSpiceError(fmt::format( "Checking if target '{}' is in view of instrument '{}' failed", target, instrument )); } return visible == SPICETRUE; } SpiceManager::TargetStateResult SpiceManager::targetState(const std::string& target, const std::string& observer, const std::string& referenceFrame, AberrationCorrection aberrationCorrection, double ephemerisTime) const { ghoul_assert(!target.empty(), "Target must not be empty"); ghoul_assert(!observer.empty(), "Observer must not be empty"); ghoul_assert(!referenceFrame.empty(), "Reference frame must not be empty"); TargetStateResult result; result.lightTime = 0.0; double buffer[6]; spkezr_c( target.c_str(), ephemerisTime, referenceFrame.c_str(), aberrationCorrection, observer.c_str(), buffer, &result.lightTime ); if (failed_c()) { throwSpiceError(fmt::format( "Error retrieving state of target '{}' viewed from observer '{}' in " "reference frame '{}' at time '{}'", target, observer, referenceFrame, ephemerisTime )); } memmove(glm::value_ptr(result.position), buffer, sizeof(double) * 3); memmove(glm::value_ptr(result.velocity), buffer + 3, sizeof(double) * 3); return result; } SpiceManager::TransformMatrix SpiceManager::stateTransformMatrix( const std::string& sourceFrame, const std::string& destinationFrame, double ephemerisTime) const { ghoul_assert(!sourceFrame.empty(), "sourceFrame must not be empty"); ghoul_assert(!destinationFrame.empty(), "toFrame must not be empty"); TransformMatrix m; sxform_c( sourceFrame.c_str(), destinationFrame.c_str(), ephemerisTime, reinterpret_cast(m.data()) ); if (failed_c()) { throwSpiceError(fmt::format( "Error retrieved state transform matrix from frame '{}' to frame '{}' at " "time '{}'", sourceFrame, destinationFrame, ephemerisTime )); } return m; } glm::dmat3 SpiceManager::positionTransformMatrix(const std::string& sourceFrame, const std::string& destinationFrame, double ephemerisTime) const { ghoul_assert(!sourceFrame.empty(), "sourceFrame must not be empty"); ghoul_assert(!destinationFrame.empty(), "destinationFrame must not be empty"); glm::dmat3 result = glm::dmat3(1.0); pxform_c( sourceFrame.c_str(), destinationFrame.c_str(), ephemerisTime, reinterpret_cast(glm::value_ptr(result)) ); if (failed_c()) { throwSpiceError(""); } SpiceBoolean success = !(failed_c()); reset_c(); if (!success) { result = getEstimatedTransformMatrix( sourceFrame, destinationFrame, ephemerisTime ); } return glm::transpose(result); } glm::dmat3 SpiceManager::positionTransformMatrix(const std::string& sourceFrame, const std::string& destinationFrame, double ephemerisTimeFrom, double ephemerisTimeTo) const { ghoul_assert(!sourceFrame.empty(), "sourceFrame must not be empty"); ghoul_assert(!destinationFrame.empty(), "destinationFrame must not be empty"); glm::dmat3 result = glm::dmat3(1.0); pxfrm2_c( sourceFrame.c_str(), destinationFrame.c_str(), ephemerisTimeFrom, ephemerisTimeTo, reinterpret_cast(glm::value_ptr(result)) ); if (failed_c()) { throwSpiceError(fmt::format( "Error retrieving position transform matrix from '{}' at time '{}' to frame " "'{}' at time '{}'", sourceFrame, ephemerisTimeFrom, destinationFrame, ephemerisTimeTo )); } return glm::transpose(result); } SpiceManager::FieldOfViewResult SpiceManager::fieldOfView(const std::string& instrument) const { ghoul_assert(!instrument.empty(), "Instrument must not be empty"); return fieldOfView(naifId(instrument)); } SpiceManager::FieldOfViewResult SpiceManager::fieldOfView(int instrument) const { constexpr int MaxBoundsSize = 64; constexpr int BufferSize = 128; FieldOfViewResult res; SpiceInt nrReturned; double boundsArr[MaxBoundsSize][3]; char fovShapeBuffer[BufferSize]; char frameNameBuffer[BufferSize]; getfov_c(instrument, // instrument id MaxBoundsSize, // maximum size for the bounds vector BufferSize, // maximum size for the fov shape buffer BufferSize, // maximum size for the frame name buffer fovShapeBuffer, // the fov shape buffer frameNameBuffer, // the frame name buffer glm::value_ptr(res.boresightVector), // the boresight vector &nrReturned, // the number of returned array values boundsArr // the bounds ); if (failed_c()) { throwSpiceError(fmt::format( "Error getting field-of-view parameters for instrument '{}'", instrument )); return res; } res.bounds.reserve(nrReturned); for (int i = 0; i < nrReturned; ++i) { res.bounds.emplace_back(boundsArr[i][0], boundsArr[i][1], boundsArr[i][2]); } std::string shape = std::string(fovShapeBuffer); static const std::map Map = { { "POLYGON", FieldOfViewResult::Shape::Polygon }, { "RECTANGLE" , FieldOfViewResult::Shape::Rectangle }, { "CIRCLE", FieldOfViewResult::Shape::Circle }, { "ELLIPSE", FieldOfViewResult::Shape::Ellipse } }; res.shape = Map.at(shape); res.frameName = std::string(frameNameBuffer); return res; } SpiceManager::TerminatorEllipseResult SpiceManager::terminatorEllipse( const std::string& target, const std::string& observer, const std::string& frame, const std::string& lightSource, TerminatorType terminatorType, AberrationCorrection aberrationCorrection, double ephemerisTime, int numberOfTerminatorPoints) { ghoul_assert(!target.empty(), "Target must not be empty"); ghoul_assert(!observer.empty(), "Observer must not be empty"); ghoul_assert(!frame.empty(), "Frame must not be empty"); ghoul_assert(!lightSource.empty(), "Light source must not be empty"); ghoul_assert(numberOfTerminatorPoints >= 1, "Terminator points must be >= 1"); TerminatorEllipseResult res; // Warning: This assumes std::vector to have all values memory contiguous res.terminatorPoints.resize(numberOfTerminatorPoints); edterm_c( toString(terminatorType), lightSource.c_str(), target.c_str(), ephemerisTime, frame.c_str(), aberrationCorrection, observer.c_str(), numberOfTerminatorPoints, &res.targetEphemerisTime, glm::value_ptr(res.observerPosition), reinterpret_cast(res.terminatorPoints.data()) ); if (failed_c()) { throwSpiceError(fmt::format( "Error getting terminator ellipse for target '{}' from observer '{}' in " "frame '{}' with light source '{}' at time '{}'", target, observer, frame, lightSource, ephemerisTime )); } return res; } void SpiceManager::findCkCoverage(const std::string& path) { ghoul_assert(!path.empty(), "Empty file path"); ghoul_assert( std::filesystem::is_regular_file(path), fmt::format("File '{}' does not exist", path) ); constexpr unsigned int MaxObj = 1024; constexpr unsigned int WinSiz = 16384; #pragma clang diagnostic push #pragma clang diagnostic ignored "-Wold-style-cast" #pragma GCC diagnostic push #pragma GCC diagnostic ignored "-Wold-style-cast" SPICEINT_CELL(ids, MaxObj); SPICEDOUBLE_CELL(cover, WinSiz); ckobj_c(path.c_str(), &ids); if (failed_c()) { throwSpiceError("Error finding Ck Coverage"); } for (SpiceInt i = 0; i < card_c(&ids); ++i) { const SpiceInt frame = SPICE_CELL_ELEM_I(&ids, i); // NOLINT #pragma clang diagnostic pop #pragma GCC diagnostic pop scard_c(0, &cover); ckcov_c(path.c_str(), frame, SPICEFALSE, "SEGMENT", 0.0, "TDB", &cover); if (failed_c()) { throwSpiceError("Error finding Ck Coverage"); } // Get the number of intervals in the coverage window. const SpiceInt numberOfIntervals = wncard_c(&cover); for (SpiceInt j = 0; j < numberOfIntervals; ++j) { // Get the endpoints of the jth interval. SpiceDouble b, e; wnfetd_c(&cover, j, &b, &e); if (failed_c()) { throwSpiceError("Error finding Ck Coverage"); } _ckCoverageTimes[frame].insert(e); _ckCoverageTimes[frame].insert(b); _ckIntervals[frame].emplace_back(b, e); } } } void SpiceManager::findSpkCoverage(const std::string& path) { ghoul_assert(!path.empty(), "Empty file path"); ghoul_assert( std::filesystem::is_regular_file(path), fmt::format("File '{}' does not exist", path) ); constexpr unsigned int MaxObj = 1024; constexpr unsigned int WinSiz = 16384; #pragma clang diagnostic push #pragma clang diagnostic ignored "-Wold-style-cast" #pragma GCC diagnostic push #pragma GCC diagnostic ignored "-Wold-style-cast" SPICEINT_CELL(ids, MaxObj); SPICEDOUBLE_CELL(cover, WinSiz); spkobj_c(path.c_str(), &ids); if (failed_c()) { throwSpiceError("Error finding Spk ID for coverage"); } for (SpiceInt i = 0; i < card_c(&ids); ++i) { const SpiceInt obj = SPICE_CELL_ELEM_I(&ids, i); // NOLINT #pragma clang diagnostic pop #pragma GCC diagnostic pop scard_c(0, &cover); spkcov_c(path.c_str(), obj, &cover); if (failed_c()) { throwSpiceError("Error finding Spk coverage"); } // Get the number of intervals in the coverage window. const SpiceInt numberOfIntervals = wncard_c(&cover); for (SpiceInt j = 0; j < numberOfIntervals; ++j) { //Get the endpoints of the jth interval. SpiceDouble b, e; wnfetd_c(&cover, j, &b, &e); if (failed_c()) { throwSpiceError("Error finding Spk coverage"); } // insert all into coverage time set, the windows could be merged @AA _spkCoverageTimes[obj].insert(e); _spkCoverageTimes[obj].insert(b); _spkIntervals[obj].emplace_back(b, e); } } } glm::dvec3 SpiceManager::getEstimatedPosition(const std::string& target, const std::string& observer, const std::string& referenceFrame, AberrationCorrection aberrationCorrection, double ephemerisTime, double& lightTime) const { ZoneScoped ghoul_assert(!target.empty(), "Target must not be empty"); ghoul_assert(!observer.empty(), "Observer must not be empty"); ghoul_assert(!referenceFrame.empty(), "Reference frame must not be empty"); ghoul_assert(target != observer, "Target and observer must be different"); int targetId = naifId(target); if (targetId == 0) { // SOLAR SYSTEM BARYCENTER special case, no definition in kernels return glm::dvec3(0.0); } if (_spkCoverageTimes.find(targetId) == _spkCoverageTimes.end()) { if (_useExceptions) { // no coverage throw SpiceException(fmt::format("No position for '{}' at any time", target)); } else { return glm::dvec3(0.0); } } const std::set& coveredTimes = _spkCoverageTimes.find(targetId)->second; glm::dvec3 pos = glm::dvec3(0.0); if (coveredTimes.lower_bound(ephemerisTime) == coveredTimes.begin()) { // coverage later, fetch first position spkpos_c( target.c_str(), *(coveredTimes.begin()), referenceFrame.c_str(), aberrationCorrection, observer.c_str(), glm::value_ptr(pos), &lightTime ); if (failed_c()) { throwSpiceError(fmt::format( "Error estimating position for '{}' with observer '{}' in frame '{}'", target, observer, referenceFrame )); } } else if (coveredTimes.upper_bound(ephemerisTime) == coveredTimes.end()) { // coverage earlier, fetch last position spkpos_c( target.c_str(), *(coveredTimes.rbegin()), referenceFrame.c_str(), aberrationCorrection, observer.c_str(), glm::value_ptr(pos), &lightTime ); if (failed_c()) { throwSpiceError(fmt::format( "Error estimating position for '{}' with observer '{}' in frame '{}'", target, observer, referenceFrame )); } } else { // coverage both earlier and later, interpolate these positions glm::dvec3 posEarlier = glm::dvec3(0.0); double ltEarlier; double timeEarlier = *std::prev((coveredTimes.lower_bound(ephemerisTime))); spkpos_c( target.c_str(), timeEarlier, referenceFrame.c_str(), aberrationCorrection, observer.c_str(), glm::value_ptr(posEarlier), <Earlier ); glm::dvec3 posLater = glm::dvec3(0.0); double ltLater; double timeLater = *(coveredTimes.upper_bound(ephemerisTime)); spkpos_c( target.c_str(), timeLater, referenceFrame.c_str(), aberrationCorrection, observer.c_str(), glm::value_ptr(posLater), <Later ); if (failed_c()) { throwSpiceError(fmt::format( "Error estimating position for '{}' with observer '{}' in frame '{}'", target, observer, referenceFrame )); } // linear interpolation const double t = (ephemerisTime - timeEarlier) / (timeLater - timeEarlier); pos = posEarlier * (1.0 - t) + posLater * t; lightTime = ltEarlier * (1.0 - t) + ltLater * t; } return pos; } glm::dmat3 SpiceManager::getEstimatedTransformMatrix(const std::string& fromFrame, const std::string& toFrame, double time) const { glm::dmat3 result = glm::dmat3(1.0); const int idFrame = frameId(fromFrame); if (_ckCoverageTimes.find(idFrame) == _ckCoverageTimes.end()) { if (_useExceptions) { // no coverage throw SpiceException(fmt::format( "No data available for transform matrix from '{}' to '{}' at any time", fromFrame, toFrame )); } else { return glm::dmat3(1.0); } } std::set coveredTimes = _ckCoverageTimes.find(idFrame)->second; if (coveredTimes.lower_bound(time) == coveredTimes.begin()) { // coverage later, fetch first transform pxform_c( fromFrame.c_str(), toFrame.c_str(), *(coveredTimes.begin()), reinterpret_cast(glm::value_ptr(result)) ); if (failed_c()) { throwSpiceError(fmt::format( "Error estimating transform matrix from '{}' to from '{}' at time '{}'", fromFrame, toFrame, time )); } } else if (coveredTimes.upper_bound(time) == coveredTimes.end()) { // coverage earlier, fetch last transform pxform_c( fromFrame.c_str(), toFrame.c_str(), *(coveredTimes.rbegin()), reinterpret_cast(glm::value_ptr(result)) ); if (failed_c()) { throwSpiceError(fmt::format( "Error estimating transform matrix from frame '{}' to '{}' at time '{}'", fromFrame, toFrame, time )); } } else { // coverage both earlier and later, interpolate these transformations double earlier = *std::prev((coveredTimes.lower_bound(time))); double later = *(coveredTimes.upper_bound(time)); glm::dmat3 earlierTransform = glm::dmat3(1.0); pxform_c( fromFrame.c_str(), toFrame.c_str(), earlier, reinterpret_cast(glm::value_ptr(earlierTransform)) ); if (failed_c()) { throwSpiceError(fmt::format( "Error estimating transform matrix from frame '{}' to '{}' at time '{}'", fromFrame, toFrame, time )); } glm::dmat3 laterTransform = glm::dmat3(1.0); pxform_c( fromFrame.c_str(), toFrame.c_str(), later, reinterpret_cast(glm::value_ptr(laterTransform)) ); if (failed_c()) { throwSpiceError(fmt::format( "Error estimating transform matrix from frame '{}' to '{}' at time '{}'", fromFrame, toFrame, time )); } const double t = (time - earlier) / (later - earlier); result = earlierTransform * (1.0 - t) + laterTransform * t; } return result; } void SpiceManager::setExceptionHandling(UseException useException) { _useExceptions = useException; } SpiceManager::UseException SpiceManager::exceptionHandling() const { return _useExceptions; } scripting::LuaLibrary SpiceManager::luaLibrary() { return { "spice", { codegen::lua::LoadKernel, codegen::lua::UnloadKernel, codegen::lua::SpiceBodies, codegen::lua::RotationMatrix, codegen::lua::Position } }; } } // namespace openspace