From fe436748bfb5c4a3241e2e65b779a34b2e21c7fd Mon Sep 17 00:00:00 2001 From: Alexander Bock Date: Fri, 2 Dec 2016 15:05:49 +0100 Subject: [PATCH] Implementing Translation classes for Keplerian elements as well as two-line element formats (#172) --- apps/Launcher/CMakeLists.txt | 2 +- data/scene/earth/earth.mod | 11 + include/openspace/scene/translation.h | 9 + modules/base/CMakeLists.txt | 4 + modules/base/basemodule.cpp | 4 + modules/base/rendering/renderableplanet.cpp | 2 +- modules/base/rendering/renderabletrail.cpp | 1 + .../base/rendering/renderabletrailorbit.cpp | 12 +- .../rendering/renderabletrailtrajectory.cpp | 4 + modules/base/rotation/spicerotation.cpp | 16 +- .../base/translation/keplertranslation.cpp | 378 ++++++++++++++++++ modules/base/translation/keplertranslation.h | 177 ++++++++ modules/base/translation/spicetranslation.cpp | 10 + modules/base/translation/tletranslation.cpp | 374 +++++++++++++++++ modules/base/translation/tletranslation.h | 75 ++++ .../onscreengui/src/guipropertycomponent.cpp | 3 +- src/scene/translation.cpp | 9 + 17 files changed, 1070 insertions(+), 21 deletions(-) create mode 100644 modules/base/translation/keplertranslation.cpp create mode 100644 modules/base/translation/keplertranslation.h create mode 100644 modules/base/translation/tletranslation.cpp create mode 100644 modules/base/translation/tletranslation.h diff --git a/apps/Launcher/CMakeLists.txt b/apps/Launcher/CMakeLists.txt index 34affa1004..a17a8dbe36 100644 --- a/apps/Launcher/CMakeLists.txt +++ b/apps/Launcher/CMakeLists.txt @@ -74,4 +74,4 @@ endif () # Libtorrent include_external_library(${APPLICATION_NAME} libtorrent ${application_path}/ext/libtorrent) -target_include_directories(${APPLICATION_NAME} PUBLIC SYSTEM ${application_path}/ext/libtorrent/include) +target_include_directories(${APPLICATION_NAME} SYSTEM PUBLIC ${application_path}/ext/libtorrent/include) diff --git a/data/scene/earth/earth.mod b/data/scene/earth/earth.mod index 902eefcdfc..f3dbe2bd77 100644 --- a/data/scene/earth/earth.mod +++ b/data/scene/earth/earth.mod @@ -12,6 +12,17 @@ return { } } }, + -- The default reference frame for Earth-orbiting satellites + Name = "EarthInertial", + Parent = "EarthBarycenter", + Transform = { + Rotation = { + Type = "SpiceRotation", + SourceFrame = "J2000", + DestinationFrame = "GALACTIC", + } + }, + }, -- Earth module { Name = "Earth", diff --git a/include/openspace/scene/translation.h b/include/openspace/scene/translation.h index 2fba7f86e8..7d7809e65f 100644 --- a/include/openspace/scene/translation.h +++ b/include/openspace/scene/translation.h @@ -46,7 +46,16 @@ public: glm::dvec3 position(double time); + // Registers a callback that gets called when a significant change has been made that + // invalidates potentially stored points, for example in trails + void onParameterChange(std::function callback); + static openspace::Documentation Documentation(); + +protected: + void notifyObservers(); + + std::function _onParameterChangeCallback; }; } // namespace openspace diff --git a/modules/base/CMakeLists.txt b/modules/base/CMakeLists.txt index 3218c1352d..2a0a2f4f1a 100644 --- a/modules/base/CMakeLists.txt +++ b/modules/base/CMakeLists.txt @@ -43,8 +43,10 @@ set(HEADER_FILES ${CMAKE_CURRENT_SOURCE_DIR}/rendering/simplespheregeometry.h ${CMAKE_CURRENT_SOURCE_DIR}/rendering/screenspaceframebuffer.h ${CMAKE_CURRENT_SOURCE_DIR}/rendering/screenspaceimage.h + ${CMAKE_CURRENT_SOURCE_DIR}/translation/keplertranslation.h ${CMAKE_CURRENT_SOURCE_DIR}/translation/spicetranslation.h ${CMAKE_CURRENT_SOURCE_DIR}/translation/statictranslation.h + ${CMAKE_CURRENT_SOURCE_DIR}/translation/tletranslation.h ${CMAKE_CURRENT_SOURCE_DIR}/rotation/spicerotation.h ${CMAKE_CURRENT_SOURCE_DIR}/rotation/staticrotation.h ${CMAKE_CURRENT_SOURCE_DIR}/scale/staticscale.h @@ -70,8 +72,10 @@ set(SOURCE_FILES ${CMAKE_CURRENT_SOURCE_DIR}/rendering/simplespheregeometry.cpp ${CMAKE_CURRENT_SOURCE_DIR}/rendering/screenspaceframebuffer.cpp ${CMAKE_CURRENT_SOURCE_DIR}/rendering/screenspaceimage.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/translation/keplertranslation.cpp ${CMAKE_CURRENT_SOURCE_DIR}/translation/spicetranslation.cpp ${CMAKE_CURRENT_SOURCE_DIR}/translation/statictranslation.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/translation/tletranslation.cpp ${CMAKE_CURRENT_SOURCE_DIR}/rotation/spicerotation.cpp ${CMAKE_CURRENT_SOURCE_DIR}/rotation/staticrotation.cpp ${CMAKE_CURRENT_SOURCE_DIR}/scale/staticscale.cpp diff --git a/modules/base/basemodule.cpp b/modules/base/basemodule.cpp index 7db855fb09..535a8daa24 100644 --- a/modules/base/basemodule.cpp +++ b/modules/base/basemodule.cpp @@ -48,8 +48,10 @@ #include #include +#include #include #include +#include #include #include @@ -111,8 +113,10 @@ void BaseModule::internalInitialize() { auto fTranslation = FactoryManager::ref().factory(); ghoul_assert(fTranslation, "Ephemeris factory was not created"); + fTranslation->registerClass("KeplerTranslation"); fTranslation->registerClass("StaticTranslation"); fTranslation->registerClass("SpiceTranslation"); + fTranslation->registerClass("TLETranslation"); auto fRotation = FactoryManager::ref().factory(); ghoul_assert(fRotation, "Rotation factory was not created"); diff --git a/modules/base/rendering/renderableplanet.cpp b/modules/base/rendering/renderableplanet.cpp index 1a4ae0d7a9..fdccd9920e 100644 --- a/modules/base/rendering/renderableplanet.cpp +++ b/modules/base/rendering/renderableplanet.cpp @@ -349,7 +349,7 @@ void RenderablePlanet::render(const RenderData& data) { glm::dmat4 rot = glm::rotate(glm::dmat4(1.0), M_PI_2, glm::dvec3(1, 0, 0)); glm::dmat4 roty = glm::rotate(glm::dmat4(1.0), M_PI_2, glm::dvec3(0, -1, 0)); //glm::dmat4 rotProp = glm::rotate(glm::dmat4(1.0), glm::radians(static_cast(_rotation)), glm::dvec3(0, 1, 0)); - modelTransform = modelTransform * rot /** roty*/ /** rotProp*/; + modelTransform = modelTransform * rot * roty /** rotProp*/; glm::dmat4 modelViewTransform = data.camera.combinedViewMatrix() * modelTransform; diff --git a/modules/base/rendering/renderabletrail.cpp b/modules/base/rendering/renderabletrail.cpp index b150cbb871..181a05ce9a 100644 --- a/modules/base/rendering/renderabletrail.cpp +++ b/modules/base/rendering/renderabletrail.cpp @@ -143,6 +143,7 @@ RenderableTrail::RenderableTrail(const ghoul::Dictionary& dictionary) _translation = std::unique_ptr(Translation::createFromDictionary( dictionary.value(KeyTranslation) )); + addPropertySubOwner(_translation.get()); _lineColor = dictionary.value(KeyColor); addProperty(_lineColor); diff --git a/modules/base/rendering/renderabletrailorbit.cpp b/modules/base/rendering/renderabletrailorbit.cpp index ca154b3607..7215c17ab0 100644 --- a/modules/base/rendering/renderabletrailorbit.cpp +++ b/modules/base/rendering/renderabletrailorbit.cpp @@ -131,7 +131,7 @@ openspace::Documentation RenderableTrailOrbit::Documentation() { RenderableTrailOrbit::RenderableTrailOrbit(const ghoul::Dictionary& dictionary) : RenderableTrail(dictionary) , _period("period", "Period in days", 0.0, 0.0, 1e9) - , _resolution("resoluion", "Number of Samples along Orbit", 10000, 1, 1e6) + , _resolution("resolution", "Number of Samples along Orbit", 10000, 1, 1e6) , _needsFullSweep(true) , _indexBufferDirty(true) { @@ -141,6 +141,10 @@ RenderableTrailOrbit::RenderableTrailOrbit(const ghoul::Dictionary& dictionary) "RenderableTrailOrbit" ); + _translation->onParameterChange([this](){ + _needsFullSweep = true; + }); + // Period is in days using namespace std::chrono; int factor = duration_cast(hours(24)).count(); @@ -180,12 +184,6 @@ void RenderableTrailOrbit::update(const UpdateData& data) { // 2. Update floating position // 3. Determine which parts of the array to upload and upload the data - // Early bailout when we don't move in time - if (data.timePaused || data.delta == 0.0) { - return; - } - - // 1 // Update the trails; the report contains whether any of the other values has been // touched and if so, how many diff --git a/modules/base/rendering/renderabletrailtrajectory.cpp b/modules/base/rendering/renderabletrailtrajectory.cpp index 9aa9d1606e..1cbf3e830d 100644 --- a/modules/base/rendering/renderabletrailtrajectory.cpp +++ b/modules/base/rendering/renderabletrailtrajectory.cpp @@ -140,6 +140,10 @@ RenderableTrailTrajectory::RenderableTrailTrajectory(const ghoul::Dictionary& di "RenderableTrailTrajectory" ); + _translation->onParameterChange([this]() { + _needsFullSweep = true; + }); + _startTime = dictionary.value(KeyStartTime); _startTime.onChange([this] { _needsFullSweep = true; }); addProperty(_startTime); diff --git a/modules/base/rotation/spicerotation.cpp b/modules/base/rotation/spicerotation.cpp index 76aaf6dbf3..e1179820d8 100644 --- a/modules/base/rotation/spicerotation.cpp +++ b/modules/base/rotation/spicerotation.cpp @@ -116,17 +116,11 @@ SpiceRotation::SpiceRotation(const ghoul::Dictionary& dictionary) } void SpiceRotation::update(const UpdateData& data) { - try { - _matrix = SpiceManager::ref().positionTransformMatrix( - _sourceFrame, - _destinationFrame, - data.time - ); - } - catch (const ghoul::RuntimeError&) { - // In case of missing coverage - _matrix = glm::dmat3(1); - } + _matrix = SpiceManager::ref().positionTransformMatrix( + _sourceFrame, + _destinationFrame, + data.time + ); } } // namespace openspace diff --git a/modules/base/translation/keplertranslation.cpp b/modules/base/translation/keplertranslation.cpp new file mode 100644 index 0000000000..4e3616f789 --- /dev/null +++ b/modules/base/translation/keplertranslation.cpp @@ -0,0 +1,378 @@ +/***************************************************************************************** + * * + * OpenSpace * + * * + * Copyright (c) 2014-2016 * + * * + * 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 + +namespace { + + const char* KeyEccentricity = "Eccentricity"; + const char* KeySemiMajorAxis = "SemiMajorAxis"; + const char* KeyInclination = "Inclination"; + const char* KeyAscendingNode = "AscendingNode"; + const char* KeyArgumentOfPeriapsis = "ArgumentOfPeriapsis"; + const char* KeyMeanAnomaly = "MeanAnomaly"; + const char* KeyEpoch = "Epoch"; + const char* KeyPeriod = "Period"; + + +template +T solveIteration(Func function, T x0, const T& err = 0.0, int maxIterations = 100) { + T x = 0; + T x2 = x0; + + for (int i = 0; i < maxIterations; ++i) { + x = x2; + x2 = function(x); + if (abs(x2 - x) < err) { + return x2; + } + } + + return x2; +} +} // namespace + +namespace openspace { + +KeplerTranslation::RangeError::RangeError(std::string offender) + : ghoul::RuntimeError("Value '" + offender + "' out of range", "KeplerTranslation") + , offender(std::move(offender)) +{} + +Documentation KeplerTranslation::Documentation() { + using namespace openspace::documentation; + return{ + "Kepler Translation", + "base_transform_kepler", + { + { + "Type", + new StringEqualVerifier("KeplerTranslation"), + "", + Optional::No + }, + { + KeyEccentricity, + new DoubleInRangeVerifier(0.0, 1.0), + "Specifies the eccentricity of the orbit; currently, OpenSpace does not " + "support hyperbolic orbits using Keplerian elements.", + Optional::No + }, + { + KeySemiMajorAxis, + new DoubleVerifier, + "Specifies the semi-major axis of the orbit in kilometers (semi-major " + "axis = average of periapsis and apoapsis).", + Optional::No + }, + { + KeyInclination, + new DoubleInRangeVerifier(0.0, 360.0), + "Specifies the inclination angle (degrees) of the orbit relative to the " + "reference plane (in the case of Earth, the equatorial plane.", + Optional::No + }, + { + KeyAscendingNode, + new DoubleInRangeVerifier(0.0, 360.0), + "Specifies the right ascension of the ascending node (in degrees!) " + "relative to the vernal equinox.", + Optional::No + }, + { + KeyArgumentOfPeriapsis, + new DoubleInRangeVerifier(0.0, 360.0), + "Specifies the argument of periapsis as angle (in degrees) from the " + "ascending.", + Optional::No + }, + { + KeyMeanAnomaly, + new DoubleInRangeVerifier(0.0, 360.0), + "Specifies the position of the orbiting body (in degrees) along the " + "elliptical orbit at epoch time.", + Optional::No + }, + { + KeyEpoch, + new StringVerifier, + "Specifies the epoch time used for position as a string of the form: " + "YYYY MM DD HH:mm:ss", + Optional::No + }, + { + KeyPeriod, + new DoubleGreaterVerifier(0.0), + "Specifies the orbital period (in seconds).", + Optional::No + }, + }, + Exhaustive::Yes + }; +} + +KeplerTranslation::KeplerTranslation() + : Translation() + , _eccentricity("eccentricity", "Eccentricity", 0.0, 0.0, 1.0) + , _semiMajorAxis("semimajorAxis", "Semi-major axis", 0.0, 0.0, 1e6) + , _inclination("inclination", "Inclination", 0.0, 0.0, 360.0) + , _ascendingNode( + "ascendingNode", + "Right ascension of ascending Node", + 0.0, + 0.0, + 360.0 + ) + , _argumentOfPeriapsis( + "argumentOfPeriapsis", + "Argument of Periapsis", + 0.0, + 0.0, + 360.0 + ) + , _meanAnomalyAtEpoch("meanAnomalyAtEpoch", "Mean anomaly at epoch", 0.0, 0.0, 360.0) + , _epoch("epoch", "Epoch", 0.0, 0.0, 1e9) + , _period("period", "Orbit period", 0.0, 0.0, 1e6) + , _orbitPlaneDirty(true) +{ + auto update = [this]() { + _orbitPlaneDirty = true; + }; + + // Only the eccentricity, semimajor axis, inclination, and location of ascending node + // invalidate the shape of the orbit. The other parameters only determine the location + // the spacecraft on that orbit + _eccentricity.onChange(update); + addProperty(_eccentricity); + + _semiMajorAxis.onChange(update); + addProperty(_semiMajorAxis); + + _inclination.onChange(update); + addProperty(_inclination); + + _ascendingNode.onChange(update); + addProperty(_ascendingNode); + + _argumentOfPeriapsis.onChange(update); + addProperty(_argumentOfPeriapsis); + + addProperty(_meanAnomalyAtEpoch); + addProperty(_epoch); + addProperty(_period); +} + +KeplerTranslation::KeplerTranslation(const ghoul::Dictionary& dictionary) + : KeplerTranslation() +{ + documentation::testSpecificationAndThrow( + Documentation(), + dictionary, + "KeplerTranslation" + ); + + setKeplerElements( + dictionary.value(KeyEccentricity), + dictionary.value(KeySemiMajorAxis), + dictionary.value(KeyInclination), + dictionary.value(KeyAscendingNode), + dictionary.value(KeyArgumentOfPeriapsis), + dictionary.value(KeyMeanAnomaly), + dictionary.value(KeyPeriod), + dictionary.value(KeyEpoch) + ); +} + +glm::dvec3 KeplerTranslation::position() const { + return _position; +} + +double KeplerTranslation::eccentricAnomaly(double meanAnomaly) const { + // Compute the eccentric anomaly (the location of the spacecraft taking the + // eccentricity of the orbit into account) using different solves for the regimes in + // which they are most efficient + + if (_eccentricity == 0.0) { + // In a circular orbit, the eccentric anomaly = mean anomaly + return meanAnomaly; + } + else if (_eccentricity < 0.2) { + auto solver = [this, &meanAnomaly](double x) -> double { + // For low eccentricity, using a first order solver sufficient + return meanAnomaly + _eccentricity * sin(x); + }; + return solveIteration(solver, meanAnomaly, 0.0, 5); + } + else if (_eccentricity < 0.9) { + auto solver = [this, &meanAnomaly](double x) -> double { + double e = _eccentricity; + return x + (meanAnomaly + e * sin(x) - x) / (1.0 - e * cos(x)); + }; + return solveIteration(solver, meanAnomaly, 0.0, 6); + } + else if (_eccentricity < 1.0) { + auto sign = [](double val) -> double { + return val > 0.0 ? 1.0 : ((val < 0.0) ? -1.0 : 0.0); + }; + double e = meanAnomaly + 0.85 * _eccentricity * sign(sin(meanAnomaly)); + + auto solver = [this, &meanAnomaly, &sign](double x) -> double { + double e = _eccentricity; + double s = e * sin(x); + double c = e * cos(x); + double f = x - s - meanAnomaly; + double f1 = 1 - c; + double f2 = s; + return x + (-5 * f / (f1 + sign(f1) * sqrt(abs(16 * f1 * f1 - 20 * f * f2)))); + }; + return solveIteration(solver, e, 0.0, 8); + } +} + +void KeplerTranslation::update(const UpdateData& data) { + if (_orbitPlaneDirty) { + computeOrbitPlane(); + _orbitPlaneDirty = false; + } + + double t = data.time - _epoch; + double meanMotion = 2.0 * glm::pi() / _period; + double meanAnomaly = glm::radians(_meanAnomalyAtEpoch.value()) + t * meanMotion; + double e = eccentricAnomaly(meanAnomaly); + + // Use the eccentric anomaly to compute the actual location + double a = _semiMajorAxis / (1.0 - _eccentricity) * 1000.0; + glm::dvec3 p = { + a * (cos(e) - _eccentricity), + a * sqrt(1.0 - _eccentricity * _eccentricity) * sin(e), + 0.0 + }; + _position = _orbitPlaneRotation * p; +} + +void KeplerTranslation::computeOrbitPlane() { + // We assume the following coordinate system: + // z = axis of rotation + // x = pointing towards the first point of Aries + // y completes the righthanded coordinate system + + // Perform three rotations: + // 1. Around the z axis to place the location of the ascending node + // 2. Around the x axis (now aligned with the ascending node) to get the correct + // inclination + // 3. Around the new z axis to place the closest approach to the correct location + + const glm::vec3 ascendingNodeAxisRot = { 0.f, 0.f, 1.f }; + const glm::vec3 inclinationAxisRot = { 1.f, 0.f, 0.f }; + const glm::vec3 argPeriapsisAxisRot = { 0.f, 0.f, 1.f }; + + const double asc = glm::radians(_ascendingNode.value()); + const double inc = glm::radians(_inclination.value()); + const double per = glm::radians(_argumentOfPeriapsis.value()); + + _orbitPlaneRotation = + glm::rotate(asc, glm::dvec3(ascendingNodeAxisRot)) * + glm::rotate(inc, glm::dvec3(inclinationAxisRot)) * + glm::rotate(per, glm::dvec3(argPeriapsisAxisRot)); + + notifyObservers(); + _orbitPlaneDirty = false; +} + +void KeplerTranslation::setKeplerElements(double eccentricity, double semiMajorAxis, + double inclination, double ascendingNode, + double argumentOfPeriapsis, + double meanAnomalyAtEpoch, + double orbitalPeriod, const std::string& epoch) +{ + setKeplerElements( + eccentricity, + semiMajorAxis, + inclination, + ascendingNode, + argumentOfPeriapsis, + meanAnomalyAtEpoch, + orbitalPeriod, + SpiceManager::ref().ephemerisTimeFromDate(epoch) + ); +} + +void KeplerTranslation::setKeplerElements(double eccentricity, double semiMajorAxis, + double inclination, double ascendingNode, + double argumentOfPeriapsis, + double meanAnomalyAtEpoch, + double orbitalPeriod, double epoch) +{ + auto isInRange = [](double val, double min, double max) -> bool { + return val >= min && val <= max; + }; + + if (isInRange(eccentricity, 0.0, 1.0)) { + _eccentricity = eccentricity; + } + else { + throw RangeError("Eccentricity"); + } + + _semiMajorAxis = semiMajorAxis; + + if (isInRange(inclination, 0.0, 360.0)) { + _inclination = inclination; + } + else { + throw RangeError("Inclination"); + } + + if (isInRange(_ascendingNode, 0.0, 360.0)) { + _ascendingNode = ascendingNode; + } + else { + throw RangeError("Ascending Node"); + } + if (isInRange(_argumentOfPeriapsis, 0.0, 360.0)) { + _argumentOfPeriapsis = argumentOfPeriapsis; + } + else { + throw RangeError("Argument of Periapsis"); + } + + if (isInRange(_meanAnomalyAtEpoch, 0.0, 360.0)) { + _meanAnomalyAtEpoch = meanAnomalyAtEpoch; + } + else { + throw RangeError("Mean anomaly at epoch"); + } + + _period = orbitalPeriod; + _epoch = epoch; + + computeOrbitPlane(); +} + +} // namespace openspace diff --git a/modules/base/translation/keplertranslation.h b/modules/base/translation/keplertranslation.h new file mode 100644 index 0000000000..09c84c8dae --- /dev/null +++ b/modules/base/translation/keplertranslation.h @@ -0,0 +1,177 @@ +/***************************************************************************************** + * * + * OpenSpace * + * * + * Copyright (c) 2014-2016 * + * * + * 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. * + ****************************************************************************************/ + +#ifndef __KEPLERTRANSLATION_H__ +#define __KEPLERTRANSLATION_H__ + +#include + +#include + +#include +#include + +namespace openspace { + +/** + * The KeplerTranslation is a concrete Translation implementation that uses the 6 + * Keplerian elements (eccentricity, semi-major axis, inclination, right ascension of the + * ascending node, argument of periapsis, and mean anomaly at epoch) for computing the + * position of a space craft. So far, only eccentricities between [0, 1) are supoorted. + */ +class KeplerTranslation : public Translation { +public: + struct RangeError : public ghoul::RuntimeError { + explicit RangeError(std::string offender); + + std::string offender; + }; + + /** + * The constructor that retrieves the required Keplerian elements from the passed + * \p dictionary. These values are then apssed to the setKeplerElements method for + * further processing. + * The \p dictionary is tested against the Documentation for conformance. + * \param dictionary The ghoul::Dictionary containing all the information about the + * Keplerian elements (see Documentation) + */ + KeplerTranslation(const ghoul::Dictionary& dictionary); + + /// Default destructor + virtual ~KeplerTranslation() = default; + + /** + * This method returns the position of the object at the time that was passed to the + * last call of the #update method. The KeplerTranslation caches its position, thus + * repeated calls to this function are cheap. + * \return The position of the object at the time of last call to the #update method + */ + virtual glm::dvec3 position() const override; + + /** + * Updates the cached position of the object to correspond to the time passed in the + * \p data. + * \param data The UpdateData struct that contains all information necessary to update + * this Translation + */ + void update(const UpdateData& data) override; + + /** + * Method returning the openspace::Documentation that describes the ghoul::Dictinoary + * that can be passed to the constructor. + * \return The openspace::Documentation that describes the ghoul::Dicitonary that can + * be passed to the constructor + */ + static openspace::Documentation Documentation(); + +protected: + /// Default construct that initializes all the properties and member variables + KeplerTranslation(); + + /** + * Sets the internal values for the Keplerian elements and the epoch as a string of + * the form YYYY MM DD HH:mm:ss. + * \param eccentricity The eccentricity of the orbit + * \param semiMajorAxis The semi-major axis of the orbit + * \param inclination The inclination of the orbit relative to the (x-y) reference + * plane (in the case of J2000, the equator) + * \param ascendingNode The right ascension of the ascending node computed relative + * to the x axis (in the case of J2000, the first point of Aries) + * \param argumentOfPeriapsis The location on the orbit with the closes approach + * \param meanAnomalyAtEpoch The location of the space craft on the orbit at the time + * of the \p epoch + * \param orbitalPeriod The period of the orbit in seconds + * \param epoch The epoch to which the orbit is defined as a string of the form: + * YYYY MM DD HH:mm::ss + */ + void setKeplerElements(double eccentricity, double semiMajorAxis, double inclination, + double ascendingNode, double argumentOfPeriapsis, double meanAnomalyAtEpoch, + double orbitalPeriod, const std::string& epoch + ); + + /** + * Sets the internal values for the Keplerian elements and the epoch as seconds past + * J2000 epoch. + * \param eccentricity The eccentricity of the orbit + * \param semiMajorAxis The semi-major axis of the orbit + * \param inclination The inclination of the orbit relative to the (x-y) reference + * plane (in the case of J2000, the equator) + * \param ascendingNode The right ascension of the ascending node computed relative + * to the x axis (in the case of J2000, the first point of Aries) + * \param argumentOfPeriapsis The location on the orbit with the closes approach + * \param meanAnomalyAtEpoch The location of the space craft on the orbit at the time + * of the \p epoch + * \param orbitalPeriod The period of the orbit in seconds + * \param epoch The epoch to which the orbit is defined as number of seconds past the + * J2000 epoch + */ + void setKeplerElements(double eccentricity, double semiMajorAxis, double inclination, + double ascendingNode, double argumentOfPeriapsis, double meanAnomalyAtEpoch, + double orbitalPeriod, double epoch + ); + +private: + /// Recombutes the rotation matrix used in the update method + void computeOrbitPlane(); + + /** + * This method computes the eccentric anomaly (location of the space craft taking the + * eccentricity into acount) based on the mean anomaly (location of the space craft + * assuming an eccentricity of 0.0) + * \param meanAnomaly The mean anomaly for which the eccentric anomaly shall be + * computed + * \return The eccentric anomaly for the provided \p meanAnomaly + */ + double eccentricAnomaly(double meanAnomaly) const; + + /// The eccentricity of the orbit in [0, 1) + properties::DoubleProperty _eccentricity; + /// The semi-major axis in km + properties::DoubleProperty _semiMajorAxis; + /// The inclination of the orbit in [0, 360] + properties::DoubleProperty _inclination; + /// The right ascension of the ascending node in [0, 360] + properties::DoubleProperty _ascendingNode; + /// The argument of periapsis in [0, 360] + properties::DoubleProperty _argumentOfPeriapsis; + /// The mean anomaly at the epoch in [0, 360] + properties::DoubleProperty _meanAnomalyAtEpoch; + + /// The epoch in seconds relative to the J2000 epoch + properties::DoubleProperty _epoch; + /// The period of the orbit in seconds + properties::DoubleProperty _period; + + /// Dirty flag for the _orbitPlaneRotation parameters + bool _orbitPlaneDirty; + /// The rotation matrix that defines the plane of the orbit + glm::dmat3 _orbitPlaneRotation; + + /// The cached position for the last time with which the update method was called + glm::dvec3 _position; +}; + +} // namespace openspace + +#endif // __KEPLERTRANSLATION_H__ diff --git a/modules/base/translation/spicetranslation.cpp b/modules/base/translation/spicetranslation.cpp index 9f3ce6b2c9..5d8a8be05b 100644 --- a/modules/base/translation/spicetranslation.cpp +++ b/modules/base/translation/spicetranslation.cpp @@ -144,8 +144,18 @@ SpiceTranslation::SpiceTranslation(const ghoul::Dictionary& dictionary) } } + auto update = [this](){ + notifyObservers(); + }; + + _target.onChange(update); addProperty(_target); + + _origin.onChange(update); addProperty(_origin); + + _frame.onChange(update); + addProperty(_frame); } glm::dvec3 SpiceTranslation::position() const { diff --git a/modules/base/translation/tletranslation.cpp b/modules/base/translation/tletranslation.cpp new file mode 100644 index 0000000000..13a8f20c0f --- /dev/null +++ b/modules/base/translation/tletranslation.cpp @@ -0,0 +1,374 @@ +/***************************************************************************************** + * * + * OpenSpace * + * * + * Copyright (c) 2014-2016 * + * * + * 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 + +namespace { + const char* KeyFile = "File"; + + // The list of leap years only goes until 2056 as we need to touch this file then + // again anyway ;) + static const std::vector LeapYears = { + 1956, 1960, 1964, 1968, 1972, 1976, 1980, 1984, 1988, 1992, 1996, + 2000, 2004, 2008, 2012, 2016, 2020, 2024, 2028, 2032, 2036, 2040, + 2044, 2048, 2052, 2056 + }; + + // Count the number of full days since the beginning of 2000 to the beginning of + // the parameter 'year' + int countDays(int year) { + // Find the position of the current year in the vector, the difference + // between its position and the position of 2000 (for J2000) gives the + // number of leap years + const int Epoch = 2000; + const int DaysRegularYear = 365; + const int DaysLeapYear = 366; + + if (year == Epoch) { + return 0; + } + + // Get the position of the most recent leap year + auto lb = std::lower_bound(LeapYears.begin(), LeapYears.end(), year); + + // Get the position of the epoch + auto y2000 = std::find(LeapYears.begin(), LeapYears.end(), Epoch); + + // The distance between the two iterators gives us the number of leap years + int nLeapYears = std::abs(std::distance(y2000, lb)); + + int nYears = std::abs(year - Epoch); + int nRegularYears = nYears - nLeapYears; + + // Get the total number of days as the sum of leap years + non leap years + int result = nRegularYears * DaysRegularYear + nLeapYears * DaysLeapYear; + return result; + }; + + // Returns the number of leap seconds that lie between the {year, dayOfYear} + // time point and { 2000, 1 } + int countLeapSeconds(int year, int dayOfYear) { + // Find the position of the current year in the vector; its position in + // the vector gives the number of leap seconds + struct LeapSecond { + int year; + int dayOfYear; + bool operator<(const LeapSecond& rhs) const { + return + std::tie(year, dayOfYear) < std::tie(rhs.year, rhs.dayOfYear); + } + }; + + const LeapSecond Epoch = { 2000, 1}; + + // List taken from: https://www.ietf.org/timezones/data/leap-seconds.list + static const std::vector LeapSeconds = { + { 1972, 1 }, + { 1972, 183 }, + { 1973, 1 }, + { 1974, 1 }, + { 1975, 1 }, + { 1976, 1 }, + { 1977, 1 }, + { 1978, 1 }, + { 1979, 1 }, + { 1980, 1 }, + { 1981, 182 }, + { 1982, 182 }, + { 1983, 182 }, + { 1985, 182 }, + { 1988, 1 }, + { 1990, 1 }, + { 1991, 1 }, + { 1992, 183 }, + { 1993, 182 }, + { 1994, 182 }, + { 1996, 1 }, + { 1997, 182 }, + { 1999, 1 }, + { 2006, 1 }, + { 2009, 1 }, + { 2012, 183 }, + { 2015, 182 }, + { 2017, 1 } + }; + + // Get the position of the last leap second before the desired date + LeapSecond date { year, dayOfYear }; + auto it = std::lower_bound(LeapSeconds.begin(), LeapSeconds.end(), date); + + // Get the position of the Epoch + auto y2000 = std::lower_bound(LeapSeconds.begin(), LeapSeconds.end(), Epoch); + + // The distance between the two iterators gives us the number of leap years + int nLeapSeconds = std::abs(std::distance(y2000, it)); + return nLeapSeconds; + }; + + double epochFromSubstring(const std::string& epochString) { + // The epochString is in the form: + // YYDDD.DDDDDDDD + // With YY being the last two years of the launch epoch, the first DDD the day + // of the year and the remaning a fractional part of the day + + // The main overview of this function: + // 1. Reconstruct the full year from the YY part + // 2. Calculate the number of seconds since the beginning of the year + // 2.a Get the number of full days since the beginning of the year + // 2.b If the year is a leap year, modify the number of days + // 3. Convert the number of days to a number of seconds + // 4. Get the number of leap seconds since January 1st, 2000 and remove them + // 5. Adjust for the fact the epoch starts on 1st Januaray at 12:00:00, not + // midnight + + // According to https://celestrak.com/columns/v04n03/ + // Apparently, US Space Command sees no need to change the two-line element + // set format yet since no artificial earth satellites existed prior to 1957. + // By their reasoning, two-digit years from 57-99 correspond to 1957-1999 and + // those from 00-56 correspond to 2000-2056. We'll see each other again in 2057! + + // 1. Get the full year + std::string yearPrefix = [y = epochString.substr(0, 2)](){ + int year = std::atoi(y.c_str()); + return year >= 57 ? "19" : "20"; + }(); + int year = std::atoi((yearPrefix + epochString.substr(0, 2)).c_str()); + int daysSince2000 = countDays(year); + + // 2. + // 2.a + double daysInYear = std::atof(epochString.substr(2).c_str()); + + // 2.b + bool isInLeapYear = std::find( + LeapYears.begin(), + LeapYears.end(), + year + ) != LeapYears.end(); + if (isInLeapYear && daysInYear >= 60) { + // We are in a leap year, so we have an effective day more if we are + // beyond the end of february (= 31+29 days) + --daysInYear; + } + + // 3 + using namespace std::chrono; + int SecondsPerDay = seconds(hours(24)).count(); + double nSecondsSince2000 = (daysSince2000 + daysInYear) * SecondsPerDay; + + // 4 + // We need to remove additionbal leap seconds past 2000 and add them prior to + // 2000 to sync up the time zones + double nLeapSecondsOffset = -countLeapSeconds(year, std::floor(daysInYear)); + + // 5 + double nSecondsEpochOffset = seconds(hours(12)).count(); + + // Combine all of the values + double epoch = nSecondsSince2000 + nLeapSecondsOffset - nSecondsEpochOffset; + return epoch; + } + + double calculateSemiMajorAxis(double meanMotion) { + using namespace std::chrono; + const double GravitationalConstant = 6.6740831e-11; + const double MassEarth = 5.9721986e24; + const double muEarth = GravitationalConstant * MassEarth; + + // Use Kepler's 3rd law to calculate semimajor axis + // a^3 / P^2 = mu / (2pi)^2 + // <=> a = ((mu * P^2) / (2pi^2))^(1/3) + // with a = semimajor axis + // P = period in seconds + // mu = G*M_earth + using namespace std::chrono; + double period = seconds(hours(24)).count() / meanMotion; + + const double pisq = glm::pi() * glm::pi(); + double semiMajorAxis = pow((muEarth * period*period) / (4 * pisq), 1.0 / 3.0); + + // We need the semi major axis in km instead of m + return semiMajorAxis / 1000.0; + } +} + + +namespace openspace { + +Documentation TLETranslation::Documentation() { + using namespace openspace::documentation; + return { + "TLE Translation", + "base_transform_tle", + { + { + "Type", + new StringEqualVerifier("TLETranslation"), + "", + Optional::No + }, + { + KeyFile, + new StringVerifier, + "Specifies the filename of the Two-Line-Element file", + Optional::No + } + }, + Exhaustive::No + }; +} + + +TLETranslation::TLETranslation(const ghoul::Dictionary& dictionary) + : KeplerTranslation() +{ + documentation::testSpecificationAndThrow( + Documentation(), + dictionary, + "TLETranslation" + ); + + std::string file = dictionary.value(KeyFile); + readTLEFile(file); +} + +void TLETranslation::readTLEFile(const std::string& filename) { + ghoul_assert(FileSys.fileExists(filename), "The filename must exist"); + + std::ifstream file; + file.exceptions(std::ofstream::failbit | std::ofstream::badbit); + file.open(filename); + + // All of the Kepler element information + struct { + double inclination; + double semiMajorAxis; + double ascendingNode; + double eccentricity; + double argumentOfPeriapsis; + double meanAnomaly; + double meanMotion; + double epoch; + } keplerElements; + + std::string line; + while (std::getline(file, line)) { + if (line[0] == '1') { + // First line + // Field Columns Content + // 1 01-01 Line number + // 2 03-07 Satellite number + // 3 08-08 Classification (U = Unclassified) + // 4 10-11 International Designator (Last two digits of launch year) + // 5 12-14 International Designator (Launch number of the year) + // 6 15-17 International Designator(piece of the launch) A + // 7 19-20 Epoch Year(last two digits of year) + // 8 21-32 Epoch(day of the year and fractional portion of the day) + // 9 34-43 First Time Derivative of the Mean Motion divided by two + // 10 45-52 Second Time Derivative of Mean Motion divided by six + // 11 54-61 BSTAR drag term(decimal point assumed)[10] - 11606 - 4 + // 12 63-63 The "Ephemeris type" + // 13 65-68 Element set number.Incremented when a new TLE is generated + // 14 69-69 Checksum (modulo 10) + + keplerElements.epoch = epochFromSubstring(line.substr(18, 14)); + } + else if (line[0] == '2') { + // Second line + //Field Columns Content + // 1 01-01 Line number + // 2 03-07 Satellite number + // 3 09-16 Inclination (degrees) + // 4 18-25 Right ascension of the ascending node (degrees) + // 5 27-33 Eccentricity (decimal point assumed) + // 6 35-42 Argument of perigee (degrees) + // 7 44-51 Mean Anomaly (degrees) + // 8 53-63 Mean Motion (revolutions per day) + // 9 64-68 Revolution number at epoch (revolutions) + // 10 69-69 Checksum (modulo 10) + + std::stringstream stream; + stream.exceptions(std::ios::failbit); + + // Get inclination + stream.str(line.substr(8, 8)); + stream >> keplerElements.inclination; + stream.clear(); + + // Get Right ascension of the ascending node + stream.str(line.substr(17, 8)); + stream >> keplerElements.ascendingNode; + stream.clear(); + + // Get Eccentricity + stream.str("0." + line.substr(26, 7)); + stream >> keplerElements.eccentricity; + stream.clear(); + + // Get argument of periapsis + stream.str(line.substr(34, 8)); + stream >> keplerElements.argumentOfPeriapsis; + stream.clear(); + + // Get mean anomaly + stream.str(line.substr(43, 8)); + stream >> keplerElements.meanAnomaly; + stream.clear(); + + // Get mean motion + stream.str(line.substr(52, 11)); + stream >> keplerElements.meanMotion; + + break; + } + } + + // Calculate the semi major axis based on the mean motion using kepler's laws + keplerElements.semiMajorAxis = calculateSemiMajorAxis(keplerElements.meanMotion); + + // Converting the mean motion (revolutions per day) to period (seconds per revolution) + using namespace std::chrono; + double period = seconds(hours(24)).count() / keplerElements.meanMotion; + + setKeplerElements( + keplerElements.eccentricity, + keplerElements.semiMajorAxis, + keplerElements.inclination, + keplerElements.ascendingNode, + keplerElements.argumentOfPeriapsis, + keplerElements.meanAnomaly, + period, + keplerElements.epoch + ); +} + +} // namespace openspace diff --git a/modules/base/translation/tletranslation.h b/modules/base/translation/tletranslation.h new file mode 100644 index 0000000000..e7df9bb721 --- /dev/null +++ b/modules/base/translation/tletranslation.h @@ -0,0 +1,75 @@ +/***************************************************************************************** + * * + * OpenSpace * + * * + * Copyright (c) 2014-2016 * + * * + * 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. * + ****************************************************************************************/ + +#ifndef __TLETRANSLATION_H__ +#define __TLETRANSLATION_H__ + +#include + +namespace openspace { + +/** + * A specialization of the KeplerTranslation that extracts the Keplerian elements from a + * two-line element as described by the US Space Command + * https://celestrak.com/columns/v04n03 + * The ghoul::Dictionary passed to the constructor must contain the pointer to a file that + * will be read. + */ +class TLETranslation : public KeplerTranslation { +public: + /** + * Constructor for the TLETranslation class. The \p dictionary must contain a key for + * the file that contains the TLE information. The ghoul::Dictionary will be tested + * against the openspace::Documentation returned by Documentation. + * \param The ghoul::Dictionary that contains the information for this TLETranslation + (*/ + TLETranslation(const ghoul::Dictionary& dictionary = ghoul::Dictionary()); + + /** + * Method returning the openspace::Documentation that describes the ghoul::Dictinoary + * that can be passed to the constructor. + * \return The openspace::Documentation that describes the ghoul::Dicitonary that can + * be passed to the constructor + */ + static openspace::Documentation Documentation(); + +private: + /** + * Reads the provided TLE file and calles the KeplerTranslation::setKeplerElments + * method with the correct values. If \p filename is a valid TLE file but contains + * disallowed values (see KeplerTranslation::setKeplerElements), a + * KeplerTranslation::RangeError is thrown. + * \param filename The path to the file that contains the TLE file. + * \throw std::system_error if the TLE file is malformed (does not contain at least + * two lines that start with \c 1 and \c 2. + * \throw KeplerTranslation::RangeError If the Keplerian elements are outside of + * the valid range supported by Kepler::setKeplerElements + * \pre The \p filename must exist + */ + void readTLEFile(const std::string& filename); +}; + +} // namespace openspace + +#endif // __TLETRANSLATION_H__ diff --git a/modules/onscreengui/src/guipropertycomponent.cpp b/modules/onscreengui/src/guipropertycomponent.cpp index 7e247045c7..4093f87573 100644 --- a/modules/onscreengui/src/guipropertycomponent.cpp +++ b/modules/onscreengui/src/guipropertycomponent.cpp @@ -74,6 +74,7 @@ void GuiPropertyComponent::renderPropertyOwner(properties::PropertyOwner* owner) return; } + int nThisProperty = nVisibleProperties(owner->properties(), _visibility); ImGui::PushID(owner->name().c_str()); const auto& subOwners = owner->propertySubOwners(); for (properties::PropertyOwner* subOwner : subOwners) { @@ -82,7 +83,7 @@ void GuiPropertyComponent::renderPropertyOwner(properties::PropertyOwner* owner) if (count == 0) { continue; } - if (subOwners.size() == 1) { + if (subOwners.size() == 1 && (nThisProperty == 0)) { renderPropertyOwner(subOwner); } else { diff --git a/src/scene/translation.cpp b/src/scene/translation.cpp index 8fe9513bce..0f0c19fefd 100644 --- a/src/scene/translation.cpp +++ b/src/scene/translation.cpp @@ -97,5 +97,14 @@ glm::dvec3 Translation::position(double time) { return position(); } +void Translation::notifyObservers() { + if (_onParameterChangeCallback) { + _onParameterChangeCallback(); + } +} + +void Translation::onParameterChange(std::function callback) { + _onParameterChangeCallback = std::move(callback); +} } // namespace openspace