diff --git a/src/animation/ease.hpp b/src/animation/ease.hpp index 674900a..5f40c33 100644 --- a/src/animation/ease.hpp +++ b/src/animation/ease.hpp @@ -339,7 +339,6 @@ T ease::in_out_elastic(const T& x, const T& y, S a) if (a < S(0.5)) { return math::lerp(x, y, std::pow(S(2), S(20) * a - S(11)) * std::sin(S(15.5334) - S(27.5293) * a)); - } else { diff --git a/src/astro/illuminance.hpp b/src/astro/illuminance.hpp index 61dfdf6..257eef1 100644 --- a/src/astro/illuminance.hpp +++ b/src/astro/illuminance.hpp @@ -22,8 +22,7 @@ namespace astro { - - + /** * Converts apparent (visual) magnitude to a brightness factor relative to a 0th magnitude star. * diff --git a/src/entity/components/celestial-body.hpp b/src/entity/components/celestial-body.hpp index 691d9d6..025a1f0 100644 --- a/src/entity/components/celestial-body.hpp +++ b/src/entity/components/celestial-body.hpp @@ -26,17 +26,23 @@ namespace component { /// A simple celestial body. struct celestial_body { - /// Body radius, in meters. + /// Mean radius of the body, in meters. double radius; - /// Angle between the body's rotational axis and its orbital axis, in radians. - double axial_tilt; + /// Mass of the body, in kilograms. + double mass; - /// Angle of rotation about the body's rotational axis at epoch, in radians. - double axial_rotation; + /// Right ascension of the body's north pole at epoch, in radians. + double pole_ra; - /// Angular frequency, in radians per day. - double angular_frequency; + /// Declination of the body's north pole at epoch, in radians. + double pole_dec; + + /// Location of the prime meridian at epoch, as a rotation about the north pole, in radians. + double prime_meridian; + + /// Sidereal rotation period, in rotations per day. + double rotation_period; }; } // namespace component diff --git a/src/entity/components/orbit.hpp b/src/entity/components/orbit.hpp index 1d8992d..e19af70 100644 --- a/src/entity/components/orbit.hpp +++ b/src/entity/components/orbit.hpp @@ -20,16 +20,36 @@ #ifndef ANTKEEPER_ENTITY_COMPONENT_ORBIT_HPP #define ANTKEEPER_ENTITY_COMPONENT_ORBIT_HPP +#include "entity/id.hpp" #include "physics/orbit/elements.hpp" #include "physics/orbit/state.hpp" +#include "math/se3.hpp" namespace entity { namespace component { struct orbit { + /// Parent body. + entity::id parent; + + /// Keplerian orbital elements at epoch. physics::orbit::elements elements; - physics::orbit::state state; + + /// Orbital semi-minor axis (b) at epoch. + double semiminor_axis; + + /// Orbital mean motion (n) at epoch. + double mean_motion; + + /// Transformation from the PQW frame to the BCI frame of the parent body. + math::transformation::se3 pqw_to_bci; + + /// Orbital Cartesian position of the body in the BCI frame of the parent body. + math::vector3 bci_position; + + /// Orbital Cartesian position of the body in the ICRF frame. + math::vector3 icrf_position; }; } // namespace component diff --git a/src/entity/systems/astronomy.cpp b/src/entity/systems/astronomy.cpp index cf2fd05..e64a830 100644 --- a/src/entity/systems/astronomy.cpp +++ b/src/entity/systems/astronomy.cpp @@ -22,6 +22,7 @@ #include "entity/components/blackbody.hpp" #include "entity/components/transform.hpp" #include "geom/intersection.hpp" +#include "geom/cartesian.hpp" #include "color/color.hpp" #include "physics/orbit/orbit.hpp" #include "physics/time/ut1.hpp" @@ -57,22 +58,16 @@ astronomy::astronomy(entity::registry& registry): reference_entity(entt::null), observer_location{0, 0, 0}, sun_light(nullptr), + sky_light(nullptr), sky_pass(nullptr) { - // Construct reference frame which transforms coordinates from SEZ to EZS - sez_to_ezs = physics::frame + // Construct transformation which transforms coordinates from ENU to EUS + enu_to_eus = math::transformation::se3 { {0, 0, 0}, - math::normalize - ( - math::quaternion::rotate_x(-math::half_pi) * - math::quaternion::rotate_z(-math::half_pi) - ) + math::quaternion::rotate_x(-math::half_pi) }; - // Construct reference frame which transforms coordinates from EZS to SEZ - ezs_to_sez = sez_to_ezs.inverse(); - registry.on_construct().connect<&astronomy::on_celestial_body_construct>(this); registry.on_replace().connect<&astronomy::on_celestial_body_replace>(this); } @@ -92,47 +87,91 @@ void astronomy::update(double t, double dt) const entity::component::orbit& reference_orbit = registry.get(reference_entity); const entity::component::celestial_body& reference_body = registry.get(reference_entity); - // Determine axial rotation at current time - const double reference_axial_rotation = reference_body.axial_rotation + reference_body.angular_frequency * universal_time; + math::transformation::se3 icrf_to_bci{{-reference_orbit.icrf_position}, math::identity_quaternion}; - // Construct reference frame which transforms coordinates from inertial space to reference body BCBF space - inertial_to_bcbf = physics::orbit::inertial::to_bcbf + // Construct transformation from the ICRF frame to the reference body BCBF frame + icrf_to_bcbf = physics::orbit::frame::bci::to_bcbf ( - reference_orbit.state.r, - reference_orbit.elements.i, - reference_body.axial_tilt, - reference_axial_rotation + reference_body.pole_ra, + reference_body.pole_dec, + reference_body.prime_meridian + (math::two_pi / reference_body.rotation_period) * universal_time ); + icrf_to_bcbf.t = icrf_to_bcbf.r * -reference_orbit.icrf_position; - // Construct reference frame which transforms coordinates from inertial space to reference body topocentric space - inertial_to_topocentric = inertial_to_bcbf * bcbf_to_topocentric; + icrf_to_enu = icrf_to_bcbf * bcbf_to_enu; + icrf_to_eus = icrf_to_enu * enu_to_eus; // Set the transform component translations of orbiting bodies to their topocentric positions registry.view().each( - [&](entity::id entity_id, const auto& celestial_body, const auto& orbit, auto& transform) + [&](entity::id entity_id, const auto& body, const auto& orbit, auto& transform) { - // Transform Cartesian position vector (r) from inertial space to topocentric space - const math::vector3 r_topocentric = inertial_to_topocentric * orbit.state.r; + // Skip reference entity + if (entity_id == this->reference_entity) + return; + + // Transform orbital Cartesian position (r) from the ICRF frame to the EUS frame + const double3 r_eus = icrf_to_eus * orbit.icrf_position; + + // Determine body orientation in the ICRF frame + math::quaternion rotation_icrf = physics::orbit::frame::bcbf::to_bci + ( + body.pole_ra, + body.pole_dec, + body.prime_meridian + (math::two_pi / body.rotation_period) * this->universal_time + ).r; + + // Transform body orientation from the ICRF frame to the EUS frame. + math::quaternion rotation_eus = math::normalize(icrf_to_eus.r * rotation_icrf); // Update local transform - transform.local.translation = math::type_cast(r_topocentric); + if (orbit.parent != entt::null) + { + transform.local.translation = math::normalize(math::type_cast(r_eus)) * 1000.0f; + transform.local.rotation = math::type_cast(rotation_eus); + transform.local.scale = {50.0f, 50.0f, 50.0f}; + } + + /* + if (orbit.parent == entt::null) + { + // RA-DEC + const double3 r_bci = icrf_to_bci * orbit.icrf_position; + double3 r_bci_spherical = physics::orbit::frame::bci::spherical(r_bci); + if (r_bci_spherical.z < 0.0) + r_bci_spherical.z += math::two_pi; + const double dec = math::degrees(r_bci_spherical.y); + const double ra = math::degrees(r_bci_spherical.z); + + // AZ-EL + const double3 r_enu = icrf_to_enu * orbit.icrf_position; + double3 r_enu_spherical = physics::orbit::frame::enu::spherical(r_enu); + if (r_enu_spherical.z < 0.0) + r_enu_spherical.z += math::two_pi; + const double el = math::degrees(r_enu_spherical.y); + const double az = math::degrees(r_enu_spherical.z); + + std::cout << "t: " << this->universal_time << "; ra: " << ra << "; dec: " << dec << std::endl; + std::cout << "t: " << this->universal_time << "; az: " << az << "; el: " << el << std::endl; + } + */ }); // Update blackbody lighting registry.view().each( - [&](entity::id entity_id, const auto& celestial_body, const auto& orbit, const auto& blackbody) + [&](entity::id entity_id, const auto& body, const auto& orbit, const auto& blackbody) { // Calculate blackbody inertial basis - double3 blackbody_forward_inertial = math::normalize(reference_orbit.state.r - orbit.state.r); - double3 blackbody_up_inertial = {0, 0, 1}; + //double3 blackbody_forward_icrf = math::normalize(reference_orbit.icrf_position - orbit.icrf_position); + double3 blackbody_up_icrf = {0, 0, 1}; - // Transform blackbody inertial position and basis into topocentric space - double3 blackbody_position_topocentric = inertial_to_topocentric * orbit.state.r; - double3 blackbody_forward_topocentric = inertial_to_topocentric.rotation * blackbody_forward_inertial; - double3 blackbody_up_topocentric = inertial_to_topocentric.rotation * blackbody_up_inertial; + // Transform blackbody ICRF position and basis into EUS frame + double3 blackbody_position_eus = icrf_to_eus * orbit.icrf_position; + double3 blackbody_position_enu = icrf_to_enu * orbit.icrf_position; + double3 blackbody_forward_eus = math::normalize(-blackbody_position_eus); + double3 blackbody_up_eus = icrf_to_eus.r * blackbody_up_icrf; // Calculate distance from observer to blackbody - double blackbody_distance = math::length(blackbody_position_topocentric) - celestial_body.radius; + double blackbody_distance = math::length(blackbody_position_eus) - body.radius; // Calculate blackbody distance attenuation double distance_attenuation = 1.0 / (blackbody_distance * blackbody_distance); @@ -148,7 +187,7 @@ void astronomy::update(double t, double dt) // Altitude of observer in meters geom::ray sample_ray; sample_ray.origin = {0, reference_body.radius + observer_location[0], 0}; - sample_ray.direction = math::normalize(blackbody_position_topocentric); + sample_ray.direction = math::normalize(blackbody_position_eus); geom::sphere exosphere; exosphere.center = {0, 0, 0}; @@ -172,49 +211,58 @@ void astronomy::update(double t, double dt) if (sun_light != nullptr) { // Update blackbody light transform - sun_light->set_translation(math::normalize(math::type_cast(blackbody_position_topocentric))); + sun_light->set_translation(math::normalize(math::type_cast(blackbody_position_eus))); sun_light->set_rotation ( math::look_rotation ( - math::type_cast(blackbody_forward_topocentric), - math::type_cast(blackbody_up_topocentric) + math::type_cast(blackbody_forward_eus), + math::type_cast(blackbody_up_eus) ) ); + // Sun illuminance at the outer atmosphere + float3 sun_illuminance_outer = math::type_cast(blackbody.luminous_intensity * distance_attenuation); - // Sun color at the outer atmosphere - float3 sun_color_outer = math::type_cast(blackbody.luminous_intensity * distance_attenuation); - - // Sun color at sea level - float3 sun_color_inner = math::type_cast(blackbody.luminous_intensity * distance_attenuation * atmospheric_transmittance); + // Sun illuminance at sea level + float3 sun_illuminance_inner = math::type_cast(blackbody.luminous_intensity * distance_attenuation * atmospheric_transmittance); // Update blackbody light color and intensity - sun_light->set_color(sun_color_inner); + sun_light->set_color(sun_illuminance_inner); sun_light->set_intensity(1.0f); // Upload blackbody params to sky pass if (this->sky_pass) { - this->sky_pass->set_sun_position(math::type_cast(blackbody_position_topocentric)); - this->sky_pass->set_sun_color(sun_color_outer, sun_color_inner); + this->sky_pass->set_sun_position(math::type_cast(blackbody_position_eus)); + this->sky_pass->set_sun_illuminance(sun_illuminance_outer, sun_illuminance_inner); - double blackbody_angular_radius = std::asin((celestial_body.radius * 2.0) / (blackbody_distance * 2.0)); + double blackbody_angular_radius = std::asin((body.radius * 2.0) / (blackbody_distance * 2.0)); this->sky_pass->set_sun_angular_radius(static_cast(blackbody_angular_radius)); } } + + if (sky_light != nullptr) + { + double3 blackbody_position_enu_spherical = physics::orbit::frame::enu::spherical(icrf_to_enu * orbit.icrf_position); + + double illuminance = 25000.0 * std::max(0.0, std::sin(blackbody_position_enu_spherical.y)); + + sky_light->set_color({1.0f, 1.0f, 1.0f}); + sky_light->set_intensity(static_cast(illuminance)); + } }); // Update sky pass topocentric frame if (sky_pass != nullptr) { // Upload topocentric frame to sky pass - sky_pass->set_topocentric_frame + sky_pass->set_icrf_to_eus ( - physics::frame + math::transformation::se3 { - math::type_cast(inertial_to_topocentric.translation), - math::type_cast(inertial_to_topocentric.rotation) + math::type_cast(icrf_to_eus.t), + math::type_cast(icrf_to_eus.r) } ); @@ -247,13 +295,13 @@ void astronomy::set_time_scale(double scale) void astronomy::set_reference_body(entity::id entity_id) { reference_entity = entity_id; - update_bcbf_to_topocentric(); + update_bcbf_to_enu(); } void astronomy::set_observer_location(const double3& location) { observer_location = location; - update_bcbf_to_topocentric(); + update_bcbf_to_enu(); } void astronomy::set_sun_light(scene::directional_light* light) @@ -261,6 +309,11 @@ void astronomy::set_sun_light(scene::directional_light* light) sun_light = light; } +void astronomy::set_sky_light(scene::ambient_light* light) +{ + sky_light = light; +} + void astronomy::set_sky_pass(::render::sky_pass* pass) { this->sky_pass = pass; @@ -269,16 +322,16 @@ void astronomy::set_sky_pass(::render::sky_pass* pass) void astronomy::on_celestial_body_construct(entity::registry& registry, entity::id entity_id, entity::component::celestial_body& celestial_body) { if (entity_id == reference_entity) - update_bcbf_to_topocentric(); + update_bcbf_to_enu(); } void astronomy::on_celestial_body_replace(entity::registry& registry, entity::id entity_id, entity::component::celestial_body& celestial_body) { if (entity_id == reference_entity) - update_bcbf_to_topocentric(); + update_bcbf_to_enu(); } -void astronomy::update_bcbf_to_topocentric() +void astronomy::update_bcbf_to_enu() { double radial_distance = observer_location[0]; @@ -288,13 +341,13 @@ void astronomy::update_bcbf_to_topocentric() radial_distance += registry.get(reference_entity).radius; } - // Construct reference frame which transforms coordinates from BCBF space to topocentric space - bcbf_to_topocentric = physics::orbit::bcbf::to_topocentric + // Construct reference frame which transforms coordinates from a BCBF frame to a horizontal frame + bcbf_to_enu = physics::orbit::frame::bcbf::to_enu ( radial_distance, observer_location[1], observer_location[2] - ) * sez_to_ezs; + ); } } // namespace system diff --git a/src/entity/systems/astronomy.hpp b/src/entity/systems/astronomy.hpp index 445137f..5cf7322 100644 --- a/src/entity/systems/astronomy.hpp +++ b/src/entity/systems/astronomy.hpp @@ -23,8 +23,9 @@ #include "entity/systems/updatable.hpp" #include "entity/id.hpp" #include "scene/directional-light.hpp" +#include "scene/ambient-light.hpp" #include "utility/fundamental-types.hpp" -#include "physics/frame.hpp" +#include "math/se3.hpp" #include "render/passes/sky-pass.hpp" #include "entity/components/atmosphere.hpp" #include "entity/components/celestial-body.hpp" @@ -79,6 +80,7 @@ public: void set_observer_location(const double3& location); void set_sun_light(scene::directional_light* light); + void set_sky_light(scene::ambient_light* light); void set_sky_pass(::render::sky_pass* pass); @@ -86,7 +88,7 @@ private: void on_celestial_body_construct(entity::registry& registry, entity::id entity_id, entity::component::celestial_body& celestial_body); void on_celestial_body_replace(entity::registry& registry, entity::id entity_id, entity::component::celestial_body& celestial_body); - void update_bcbf_to_topocentric(); + void update_bcbf_to_enu(); double universal_time; double time_scale; @@ -95,13 +97,14 @@ private: double3 observer_location; - physics::frame inertial_to_bcbf; - physics::frame bcbf_to_topocentric; - physics::frame inertial_to_topocentric; - physics::frame sez_to_ezs; - physics::frame ezs_to_sez; + math::transformation::se3 icrf_to_bcbf; + math::transformation::se3 bcbf_to_enu; + math::transformation::se3 icrf_to_enu; + math::transformation::se3 enu_to_eus; + math::transformation::se3 icrf_to_eus; scene::directional_light* sun_light; + scene::ambient_light* sky_light; ::render::sky_pass* sky_pass; }; diff --git a/src/entity/systems/orbit.cpp b/src/entity/systems/orbit.cpp index 7835f4a..003733d 100644 --- a/src/entity/systems/orbit.cpp +++ b/src/entity/systems/orbit.cpp @@ -18,9 +18,8 @@ */ #include "entity/systems/orbit.hpp" -#include "entity/components/orbit.hpp" -#include "entity/id.hpp" #include "physics/orbit/orbit.hpp" +#include namespace entity { namespace system { @@ -31,51 +30,46 @@ orbit::orbit(entity::registry& registry): time_scale(1.0), ke_iterations(10), ke_tolerance(1e-6) -{} +{ + registry.on_construct().connect<&orbit::on_orbit_construct>(this); + registry.on_replace().connect<&orbit::on_orbit_replace>(this); +} void orbit::update(double t, double dt) { // Add scaled timestep to current time set_universal_time(universal_time + dt * time_scale); - // Update the orbital state of orbiting bodies + // Propagate orbits registry.view().each( [&](entity::id entity_id, auto& orbit) { - // Calculate semi-minor axis (b) - const double b = physics::orbit::derive_semiminor_axis(orbit.elements.a, orbit.elements.e); + // Determine mean anomaly at current time + const double ma = orbit.elements.ma + orbit.mean_motion * this->universal_time; // Solve Kepler's equation for eccentric anomaly (E) - const double ea = physics::orbit::kepler_ea(orbit.elements.e, orbit.elements.ta, ke_iterations, ke_tolerance); - - // Calculate radial distance and true anomaly (nu) - const double xv = orbit.elements.a * (std::cos(ea) - orbit.elements.e); - const double yv = b * std::sin(ea); - const double distance = std::sqrt(xv * xv + yv * yv); - const double ta = std::atan2(yv, xv); - - // Calculate Cartesian position (r) in perifocal space - const math::vector3 r_perifocal = math::quaternion::rotate_z(ta) * math::vector3{distance, 0, 0}; + const double ea = physics::orbit::anomaly::mean_to_eccentric(orbit.elements.ec, ma, ke_iterations, ke_tolerance); - /// @TODO Calculate Cartesian velocity (v) in perifocal space - //const math::vector3 v_perifocal = ... + // Calculate Cartesian orbital position in the PQW frame + math::vector3 pqw_position = physics::orbit::frame::pqw::cartesian(orbit.elements.ec, orbit.elements.a, ea, orbit.semiminor_axis); - // Construct perifocal to inertial reference frame - const physics::frame perifocal_to_inertial = physics::orbit::inertial::to_perifocal - ( - {0, 0, 0}, - orbit.elements.raan, - orbit.elements.i, - orbit.elements.w - ).inverse(); - - // Transform orbital state vectors from perifocal space to the parent inertial space - const math::vector3 r_inertial = perifocal_to_inertial.transform(r_perifocal); - //const math::vector3 v_inertial = perifocal_frame.transform(v_perifocal); + // Transform orbital position from PQW frame to BCI frame + orbit.bci_position = orbit.pqw_to_bci.transform(pqw_position); + }); + + // Update orbital positions in the ICRF frame + registry.view().each( + [&](entity::id entity_id, auto& orbit) + { + orbit.icrf_position = orbit.bci_position; - // Update orbital state of component - orbit.state.r = r_inertial; - //orbit.state.v = v_inertial; + entity::id parent = orbit.parent; + while (parent != entt::null) + { + const component::orbit& parent_orbit = registry.get(parent); + orbit.icrf_position += parent_orbit.bci_position; + parent = parent_orbit.parent; + } }); } @@ -89,5 +83,27 @@ void orbit::set_time_scale(double scale) time_scale = scale; } +void orbit::on_orbit_construct(entity::registry& registry, entity::id entity_id, entity::component::orbit& component) +{ + component.semiminor_axis = physics::orbit::semiminor_axis(component.elements.a, component.elements.ec); + component.pqw_to_bci = physics::orbit::frame::pqw::to_bci + ( + component.elements.om, + component.elements.in, + component.elements.w + ); +} + +void orbit::on_orbit_replace(entity::registry& registry, entity::id entity_id, entity::component::orbit& component) +{ + component.semiminor_axis = physics::orbit::semiminor_axis(component.elements.a, component.elements.ec); + component.pqw_to_bci = physics::orbit::frame::pqw::to_bci + ( + component.elements.om, + component.elements.in, + component.elements.w + ); +} + } // namespace system } // namespace entity diff --git a/src/entity/systems/orbit.hpp b/src/entity/systems/orbit.hpp index da10d14..7b793b1 100644 --- a/src/entity/systems/orbit.hpp +++ b/src/entity/systems/orbit.hpp @@ -22,6 +22,8 @@ #include "entity/systems/updatable.hpp" #include "utility/fundamental-types.hpp" +#include "entity/id.hpp" +#include "entity/components/orbit.hpp" namespace entity { namespace system { @@ -58,6 +60,9 @@ public: void set_time_scale(double scale); private: + void on_orbit_construct(entity::registry& registry, entity::id entity_id, entity::component::orbit& component); + void on_orbit_replace(entity::registry& registry, entity::id entity_id, entity::component::orbit& component); + double universal_time; double time_scale; std::size_t ke_iterations; diff --git a/src/game/state/nuptial-flight.cpp b/src/game/state/nuptial-flight.cpp index 319c972..c17551c 100644 --- a/src/game/state/nuptial-flight.cpp +++ b/src/game/state/nuptial-flight.cpp @@ -97,14 +97,14 @@ nuptial_flight::nuptial_flight(game::context& ctx): game::world::create_moon(ctx); // Set time to solar noon - game::world::set_time(ctx, 0.0); + game::world::set_time(ctx, 50.0);//107.3); // Freeze time game::world::set_time_scale(ctx, 0.0); } // Load biome - game::load::biome(ctx, "desert-scrub.bio"); + game::load::biome(ctx, "lab.bio"); // Setup and enable sky and ground passes ctx.sky_pass->set_enabled(true); @@ -124,6 +124,20 @@ nuptial_flight::nuptial_flight(game::context& ctx): entity::command::warp_to(*ctx.entity_registry, ruler_10cm_eid, {0, 0, 10}); } + // Create cocoon + { + entity::archetype* cocoon_archetype = ctx.resource_manager->load("cocoon.ent"); + auto cocoon_eid = cocoon_archetype->create(*ctx.entity_registry); + entity::command::warp_to(*ctx.entity_registry, cocoon_eid, {-10, 0, 10}); + } + + // Create larva + { + entity::archetype* larva_archetype = ctx.resource_manager->load("long-neck-larva.ent"); + auto larva_eid = larva_archetype->create(*ctx.entity_registry); + entity::command::warp_to(*ctx.entity_registry, larva_eid, {-13, 0, 10}); + } + // Create keeper if not yet created if (ctx.entities.find("keeper") == ctx.entities.end()) { @@ -152,7 +166,7 @@ nuptial_flight::nuptial_flight(game::context& ctx): locomotion.yaw = 0.0f; ctx.entity_registry->assign(boid_eid, locomotion); - entity::command::warp_to(*ctx.entity_registry, boid_eid, {0, 2, 0}); + entity::command::warp_to(*ctx.entity_registry, boid_eid, {0, 1, 0}); // Set target ant ctx.entities["ant"] = boid_eid; @@ -264,7 +278,7 @@ void nuptial_flight::setup_camera() ctx.entity_registry->assign(camera_eid, constraint_stack); } - float ev100 = 14.5f; + float ev100 = 15.0f; ctx.surface_camera->set_exposure(ev100); } diff --git a/src/game/tools.cpp b/src/game/tools.cpp index ae92534..ad80635 100644 --- a/src/game/tools.cpp +++ b/src/game/tools.cpp @@ -95,7 +95,7 @@ entity::id build_time_tool(game::context& ctx) entity::id planet_eid = ctx.entities["planet"]; entity::component::celestial_body body = ctx.entity_registry->get(planet_eid); - body.axial_rotation = math::radians(360.0f) * ((float)mouse_x / (float)window_w); + //body.axial_rotation = math::radians(360.0f) * ((float)mouse_x / (float)window_w); ctx.entity_registry->replace(planet_eid, body); }; diff --git a/src/game/world.cpp b/src/game/world.cpp index 6023d33..0fc1610 100644 --- a/src/game/world.cpp +++ b/src/game/world.cpp @@ -95,12 +95,8 @@ void create_stars(game::context& ctx) ra = math::wrap_radians(math::radians(ra)); dec = math::wrap_radians(math::radians(dec)); - // Transform spherical equatorial coordinates to rectangular equatorial coordinates - double3 position_bci = geom::spherical::to_cartesian(double3{1.0, dec, ra}); - - // Transform coordinates from equatorial space to inertial space - physics::frame bci_to_inertial = physics::orbit::inertial::to_bci({0, 0, 0}, 0.0, math::radians(23.4393)).inverse(); - double3 position_inertial = bci_to_inertial * position_bci; + // Convert ICRF coordinates from spherical to Cartesian + double3 position = physics::orbit::frame::bci::cartesian(double3{1.0, dec, ra}); // Convert color index to color temperature double cct = color::index::bv_to_cct(bv_color); @@ -118,9 +114,9 @@ void create_stars(game::context& ctx) double3 scaled_color = color_acescg * brightness; // Build vertex - *(star_vertex++) = static_cast(position_inertial.x); - *(star_vertex++) = static_cast(position_inertial.y); - *(star_vertex++) = static_cast(position_inertial.z); + *(star_vertex++) = static_cast(position.x); + *(star_vertex++) = static_cast(position.y); + *(star_vertex++) = static_cast(position.z); *(star_vertex++) = static_cast(scaled_color.x); *(star_vertex++) = static_cast(scaled_color.y); *(star_vertex++) = static_cast(scaled_color.z); @@ -188,22 +184,23 @@ void create_sun(game::context& ctx) entity::id sun_eid = sun_archetype->create(*ctx.entity_registry); ctx.entities["sun"] = sun_eid; - // Create direct sun light scene object - scene::directional_light* sun_direct = new scene::directional_light(); + // Create sun directional light scene object + scene::directional_light* sun_light = new scene::directional_light(); - // Create ambient sun light scene object - scene::ambient_light* sun_ambient = new scene::ambient_light(); - sun_ambient->set_color({1, 1, 1}); - sun_ambient->set_intensity(30000.0f * 0.0f); - sun_ambient->update_tweens(); + // Create sky ambient light scene object + scene::ambient_light* sky_light = new scene::ambient_light(); + sky_light->set_color({1, 1, 1}); + sky_light->set_intensity(0.0f); + sky_light->update_tweens(); // Add sun light scene objects to surface scene - ctx.surface_scene->add_object(sun_direct); - ctx.surface_scene->add_object(sun_ambient); + ctx.surface_scene->add_object(sun_light); + ctx.surface_scene->add_object(sky_light); // Pass direct sun light scene object to shadow map pass and astronomy system - ctx.surface_shadow_map_pass->set_light(sun_direct); - ctx.astronomy_system->set_sun_light(sun_direct); + ctx.surface_shadow_map_pass->set_light(sun_light); + ctx.astronomy_system->set_sun_light(sun_light); + ctx.astronomy_system->set_sky_light(sky_light); } void create_planet(game::context& ctx) @@ -213,6 +210,9 @@ void create_planet(game::context& ctx) entity::id planet_eid = planet_archetype->create(*ctx.entity_registry); ctx.entities["planet"] = planet_eid; + // Assign orbital parent + ctx.entity_registry->get(planet_eid).parent = ctx.entities["sun"]; + // Assign planetary terrain component entity::component::terrain terrain; terrain.elevation = [](double, double) -> double @@ -231,9 +231,13 @@ void create_planet(game::context& ctx) void create_moon(game::context& ctx) { // Create lunar entity - entity::id moon_eid = ctx.entity_registry->create(); + entity::archetype* moon_archetype = ctx.resource_manager->load("moon.ent"); + entity::id moon_eid = moon_archetype->create(*ctx.entity_registry); ctx.entities["moon"] = moon_eid; + // Assign orbital parent + ctx.entity_registry->get(moon_eid).parent = ctx.entities["planet"]; + // Pass moon model to sky pass ctx.sky_pass->set_moon_model(ctx.resource_manager->load("moon.mdl")); } diff --git a/src/math/math.hpp b/src/math/math.hpp index 158f9b6..2c7a45c 100644 --- a/src/math/math.hpp +++ b/src/math/math.hpp @@ -35,6 +35,7 @@ namespace math {} #include "math/quaternion-functions.hpp" #include "math/quaternion-operators.hpp" +#include "math/se3.hpp" #include "math/transform-type.hpp" #include "math/transform-functions.hpp" #include "math/transform-operators.hpp" diff --git a/src/math/se3.hpp b/src/math/se3.hpp new file mode 100644 index 0000000..4141719 --- /dev/null +++ b/src/math/se3.hpp @@ -0,0 +1,139 @@ +/* + * Copyright (C) 2021 Christopher J. Howard + * + * This file is part of Antkeeper source code. + * + * Antkeeper source code is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * Antkeeper source code is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with Antkeeper source code. If not, see . + */ + +#ifndef ANTKEEPER_MATH_TRANSFORMATION_SE3_HPP +#define ANTKEEPER_MATH_TRANSFORMATION_SE3_HPP + +#include "math/vector-type.hpp" +#include "math/vector-operators.hpp" +#include "math/quaternion-type.hpp" +#include "math/quaternion-operators.hpp" +#include "math/quaternion-functions.hpp" + +namespace math { +namespace transformation { + +/** + * 3-dimensional Euclidean proper rigid transformation in SE(3). + * + * @tparam T Scalar type. + */ +template +struct se3 +{ +public: + /// Scalar type. + typedef T scalar_type; + + /// Vector type. + typedef math::vector3 vector_type; + + /// Quaternion type. + typedef math::quaternion quaternion_type; + + /// Transformation matrix type. + typedef math::matrix matrix_type; + + /// Vector representing the translation component of this SE(3) transformation. + vector_type t; + + /// Quaternion representing the rotation component of this SE(3) transformation. + quaternion_type r; + + /// Returns the inverse of this SE(3) transformation. + se3 inverse() const; + + /// Returns a matrix representation of the SE(3) transformation. + matrix_type matrix() const; + + /** + * Transforms a vector by this SE(3) transformation. + * + * @param x Untransformed vector. + * @return Transformed vector. + */ + vector_type transform(const vector_type& x) const; + + /** + * Transforms an SE(3) transformation by this SE(3) transformation. + * + * @param x Other SE(3) transformation. + * @return Frame in this se3's space. + */ + se3 transform(const se3& x) const; + + /// @copydoc se3::transform(const vector_type&) const + vector_type operator*(const vector_type& x) const; + + /// @copydoc se3::transform(const se3&) const + se3 operator*(const se3& x) const; +}; + +template +se3 se3::inverse() const +{ + const quaternion_type inverse_r = math::conjugate(r); + const vector_type inverse_t = -(inverse_r * t); + return se3{inverse_t, inverse_r}; +} + +template +typename se3::matrix_type se3::matrix() const +{ + matrix_type m = math::resize<4, 4>(math::matrix_cast(r)); + + m[3].x = t.x; + m[3].y = t.y; + m[3].z = t.z; + + return m; +} + +template +typename se3::vector_type se3::transform(const vector_type& x) const +{ + return r * x + t; +} + +template +se3 se3::transform(const se3& x) const +{ + return se3 + { + x.transform(t), + math::normalize(x.r * r) + }; +} + +template +typename se3::vector_type se3::operator*(const vector_type& x) const +{ + return transform(x); +} + +template +se3 se3::operator*(const se3& x) const +{ + return transform(x); +} + +} // namespace transformation +} // namespace math + +#endif // ANTKEEPER_MATH_TRANSFORMATION_SE3_HPP diff --git a/src/physics/atmosphere.hpp b/src/physics/atmosphere.hpp index c713b81..2e26cd9 100644 --- a/src/physics/atmosphere.hpp +++ b/src/physics/atmosphere.hpp @@ -124,7 +124,7 @@ T extinction(T ec, T s) template T albedo(T s, T e) { - return s / t; + return s / e; } /** diff --git a/src/physics/frame.hpp b/src/physics/frame.hpp deleted file mode 100644 index f10803e..0000000 --- a/src/physics/frame.hpp +++ /dev/null @@ -1,134 +0,0 @@ -/* - * Copyright (C) 2021 Christopher J. Howard - * - * This file is part of Antkeeper source code. - * - * Antkeeper source code is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * Antkeeper source code is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with Antkeeper source code. If not, see . - */ - -#ifndef ANTKEEPER_PHYSICS_FRAME_HPP -#define ANTKEEPER_PHYSICS_FRAME_HPP - -#include "math/math.hpp" -#include - -namespace physics { - -/** - * 3-dimensional frame of reference. - * - * @tparam T Scalar type. - */ -template -struct frame -{ -public: - /// Scalar type. - typedef T scalar_type; - - /// Vector type. - typedef math::vector3 vector_type; - - /// Quaternion type. - typedef math::quaternion quaternion_type; - - /// Transformation matrix type. - typedef math::matrix matrix_type; - - /// Vector representing a translation from the origin of this frame to the origin of the parent frame. - vector_type translation; - - /// Quaternion representing a rotation from the orientation of this frame to the orientation of the parent frame. - quaternion_type rotation; - - /// Returns the inverse of this frame. - frame inverse() const; - - /// Returns a transformation matrix representation of this frame. - matrix_type matrix() const; - - /** - * Transforms a vector from the parent frame space to this frame's space. - * - * @param v Vector in parent frame space. - * @return Vector in this frame's space. - */ - vector_type transform(const vector_type& v) const; - - /** - * Transforms a frame from the parent frame space to this frame's space. - * - * @param f Frame in parent frame space. - * @return Frame in this frame's space. - */ - frame transform(const frame& f) const; - - /// @copydoc frame::transform(const vector_type&) const - vector_type operator*(const vector_type& v) const; - - /// @copydoc frame::transform(const frame&) const - frame operator*(const frame& f) const; -}; - -template -frame frame::inverse() const -{ - const quaternion_type inverse_rotation = math::conjugate(rotation); - const vector_type inverse_translation = -(inverse_rotation * translation); - return frame{inverse_translation, inverse_rotation}; -} - -template -typename frame::matrix_type frame::matrix() const -{ - matrix_type m = math::resize<4, 4>(math::matrix_cast(rotation)); - - m[3].x = translation.x; - m[3].y = translation.y; - m[3].z = translation.z; - - return m; -} - -template -typename frame::vector_type frame::transform(const vector_type& v) const -{ - return translation + rotation * v; -} - -template -frame frame::transform(const frame& f) const -{ - return frame - { - f.transform(translation), - math::normalize(f.rotation * rotation) - }; -} - -template -typename frame::vector_type frame::operator*(const vector_type& v) const -{ - return transform(v); -} - -template -frame frame::operator*(const frame& f) const -{ - return transform(f); -} - -} // namespace physics - -#endif // ANTKEEPER_PHYSICS_FRAME_HPP diff --git a/src/physics/orbit/anomaly.hpp b/src/physics/orbit/anomaly.hpp new file mode 100644 index 0000000..9533232 --- /dev/null +++ b/src/physics/orbit/anomaly.hpp @@ -0,0 +1,199 @@ +/* + * Copyright (C) 2021 Christopher J. Howard + * + * This file is part of Antkeeper source code. + * + * Antkeeper source code is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * Antkeeper source code is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with Antkeeper source code. If not, see . + */ + +#ifndef ANTKEEPER_PHYSICS_ORBIT_ANOMALY_HPP +#define ANTKEEPER_PHYSICS_ORBIT_ANOMALY_HPP + +#include "math/constants.hpp" +#include + +namespace physics { +namespace orbit { + +/** + * Orbital anomaly functions. + */ +namespace anomaly { + +/** + * Derives the eccentric anomaly given eccentricity and true anomaly. + * + * @param ec Eccentricity (e). + * @param ta True anomaly (nu). + * @return Eccentric anomaly (E). + */ +template +T true_to_eccentric(T ec, T ta); + +/** + * Derives the mean anomaly given eccentricity and true anomaly. + * + * @param ec Eccentricity (e). + * @param ta True anomaly (nu). + * @return Mean anomaly (M). + */ +template +T true_to_mean(T ec, T ta); + +/** + * Derives the true anomaly given eccentricity and eccentric anomaly. + * + * @param ec Eccentricity (e). + * @param ea Eccentric anomaly (E). + * @return True anomaly (nu). + */ +template +T eccentric_to_true(T ec, T ea); + +/** + * Derives the mean anomaly given eccentricity and eccentric anomaly. + * + * @param ec Eccentricity (e). + * @param ea Eccentric anomaly (E). + * @return Mean anomaly (M). + */ +template +T eccentric_to_mean(T ec, T ea); + +/** + * Iteratively derives the eccentric anomaly given eccentricity and mean anomaly. + * + * @param ec Eccentricity (e). + * @param ma Mean anomaly (M). + * @param iterations Maximum number of iterations. + * @param tolerance Solution error tolerance. + * @return Eccentric anomaly (E). + * + * @see Murison, Marc. (2006). A Practical Method for Solving the Kepler Equation. 10.13140/2.1.5019.6808. + */ +template +T mean_to_eccentric(T ec, T ma, std::size_t iterations, T tolerance); + +/** + * Iteratively derives the true anomaly given eccentricity and mean anomaly. + * + * @param ec Eccentricity (e). + * @param ma Mean anomaly (M). + * @param iterations Maximum number of iterations. + * @param tolerance Solution error tolerance. + * @return True anomaly (nu). + */ +template +T mean_to_true(T ec, T ma, std::size_t iterations, T tolerance); + +template +T true_to_eccentric(T ec, T ta) +{ + // Parabolic orbit + if (ec == T(1)) + return std::tan(ta * T(0.5)); + + // Hyperbolic orbit + if (ec > T(1)) + return std::acosh((ec + std::cos(ta)) / (T(1) + ec * std::cos(ta))) * ((ta < T(0)) ? T(-1) : T(1)); + + // Elliptic orbit + return std::atan2(std::sqrt(T(1) - ec * ec) * std::sin(ta), std::cos(ta) + ec); +} + +template +T true_to_mean(T ec, T ta) +{ + return eccentric_to_mean(ec, true_to_eccentric(ec, ta)); +} + +template +T eccentric_to_true(T ec, T ea) +{ + // Parabolic orbit + if (ec == T(1)) + return std::atan(ea) * T(2); + + // Hyperbolic orbit + if (ec > T(1)) + return std::atan(std::sqrt((ec + T(1)) / (ec - T(1))) * std::tanh(ea * T(0.5))) * T(2); + + // Elliptic orbit + return std::atan2(sqrt(T(1) - ec * ec) * std::sin(ea), std::cos(ea) - ec); +} + +template +T eccentric_to_mean(T ec, T ea) +{ + // Parabolic orbit + if (ec == T(1)) + return (ea * ea * ea) / T(6) + ea * T(0.5); + + // Hyperbolic orbit + if (ec > T(1)) + return ec * std::sinh(ea) - ea; + + // Elliptic orbit + return ea - ec * std::sin(ea); +} + +template +T mean_to_eccentric(T ec, T ma, std::size_t iterations, T tolerance) +{ + // Wrap mean anomaly to `[-pi, pi]` + ma = std::remainder(ma, math::two_pi); + + // Third-order approximation of eccentric anomaly starting value, E0 + const T t33 = std::cos(ma); + const T t34 = ec * ec; + const T t35 = t34 * ec; + T ea0 = ma + (T(-0.5) * t35 + ec + (t34 + T(1.5) * t33 * t35) * t33) * std::sin(ma); + + // Iteratively converge E0 and E1 + for (std::size_t i = 0; i < iterations; ++i) + { + // Third-order approximation of eccentric anomaly, E1 + const T t1 = std::cos(ea0); + const T t2 = T(-1) + ec * t1; + const T t3 = std::sin(ea0); + const T t4 = ec * t3; + const T t5 = -ea0 + t4 + ma; + const T t6 = t5 / (T(0.5) * t5 * t4 / t2 + t2); + const T ea1 = ea0 - (t5 / ((T(0.5) * t3 - (T(1) / T(6)) * t1 * t6) * ec * t6 + t2)); + + // Determine solution error + const T error = std::abs(ea1 - ea0); + + // Set E0 to E1 + ea0 = ea1; + + // Break if solution is within error tolerance + if (error < tolerance) + break; + } + + return ea0; +} + +template +T mean_to_true(T ec, T ma, std::size_t iterations, T tolerance) +{ + eccentric_to_true(ec, mean_to_eccentric(ec, ma, iterations, tolerance)); +} + +} // namespace anomaly +} // namespace orbit +} // namespace physics + +#endif // ANTKEEPER_PHYSICS_ORBIT_ANOMALY_HPP diff --git a/src/physics/orbit/elements.hpp b/src/physics/orbit/elements.hpp index 2e26b0e..a16f37b 100644 --- a/src/physics/orbit/elements.hpp +++ b/src/physics/orbit/elements.hpp @@ -21,7 +21,7 @@ #define ANTKEEPER_PHYSICS_ORBIT_ELEMENTS_HPP #include "utility/fundamental-types.hpp" -#include "physics/orbit/state.hpp" +#include "math/constants.hpp" #include namespace physics { @@ -39,54 +39,133 @@ struct elements typedef T scalar_type; /// Eccentricity (e). - scalar_type e; + scalar_type ec; /// Semimajor axis (a). scalar_type a; /// Inclination (i), in radians. - scalar_type i; + scalar_type in; /// Right ascension of the ascending node (OMEGA), in radians. - scalar_type raan; + scalar_type om; /// Argument of periapsis (omega), in radians. scalar_type w; - /// True anomaly (nu) at epoch, in radians. - scalar_type ta; + /// Mean anomaly (M) at epoch, in radians. + scalar_type ma; }; +/** + * Calculates the period of an elliptical orbit according to Kepler's third law. + * + * @param a Semimajor axis (a). + * @param gm Standard gravitational parameter (GM). + * @return Orbital period (T). + */ +template +T period(T a, T gm); + +/** + * Calculates the mean motion (n) of an orbit. + * + * @param a Semimajor axis (a). + * @param gm Standard gravitational parameter (GM). + * @return Mean motion (n), in radians per unit time. + */ +template +T mean_motion(T a, T gm); + +/** + * Calculates the mean motion (n) of an orbit. + * + * @param t Orbital period (T). + * @return Mean motion (n), in radians per unit time. + */ +template +T mean_motion(T t); + +/** + * Derives the argument of the periapsis (omega) of an orbit, given the longitude of periapsis (pomega) and longitude of the ascending node (OMEGA). + * + * @param lp Longitude of the periapsis (pomega), in radians. + * @param om Right ascension of the ascending node (OMEGA), in radians. + * @return Argument of the periapsis (omega), in radians. + */ +template +T argument_periapsis(T om, T lp); + /** * Derives the longitude of the periapsis (pomega) of an orbit, given the argument of periapsis (omega) and longitude of the ascending node (OMEGA). * - * @param w Argument of periapsis (omega), in radians. - * @param raan Right ascension of the ascending node (OMEGA), in radians. + * @param w Argument of the periapsis (omega), in radians. + * @param om Right ascension of the ascending node (OMEGA), in radians. * @return Longitude of the periapsis (pomega), in radians. */ template -T derive_longitude_periapsis(T w, T raan); +T longitude_periapsis(T om, T w); /** * Derives the semiminor axis (b) of an orbit, given the semimajor axis (a) and eccentricity (e). * * @param a Semimajor axis (a). - * @param e Eccentricity (e). + * @param ec Eccentricity (e). * @return Semiminor axis (b). */ template -T derive_semiminor_axis(T a, T e); +T semiminor_axis(T a, T ec); + +/** + * Derives the semi-latus rectum (l) of an orbit, given the semimajor axis (a) and eccentricity (e). + * + * @param a Semimajor axis (a). + * @param ec Eccentricity (e). + * @return Semi-latus rectum (l). + */ +template +T semilatus_rectum(T a, T ec); + +template +T period(T a, T gm) +{ + return math::two_pi * std::sqrt((a * a * a) / gm); +} + +template +T mean_motion(T a, T gm) +{ + return std::sqrt((a * a * a) / gm); +} + +template +T mean_motion(T t) +{ + return math::two_pi / t; +} + +template +T argument_periapsis(T om, T lp) +{ + return lp - om; +} + +template +T longitude_periapsis(T om, T w) +{ + return w + om; +} template -T derive_longitude_periapsis(T w, T raan) +T semiminor_axis(T a, T ec) { - return w + raan; + return a * std::sqrt(T(1) - ec * ec); } template -T derive_semiminor_axis(T a, T e) +T semilatus_rectum(T a, T ec) { - return a * std::sqrt(T(1) - e * e); + return a * (T(1) - ec * ec); } } // namespace orbit diff --git a/src/physics/orbit/frame.hpp b/src/physics/orbit/frame.hpp new file mode 100644 index 0000000..8293694 --- /dev/null +++ b/src/physics/orbit/frame.hpp @@ -0,0 +1,399 @@ +/* + * Copyright (C) 2021 Christopher J. Howard + * + * This file is part of Antkeeper source code. + * + * Antkeeper source code is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * Antkeeper source code is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with Antkeeper source code. If not, see . + */ + +#ifndef ANTKEEPER_PHYSICS_ORBIT_FRAME_HPP +#define ANTKEEPER_PHYSICS_ORBIT_FRAME_HPP + +#include "math/se3.hpp" +#include + +namespace physics { +namespace orbit { + +/// Orbital reference frames. +namespace frame { + +/// Perifocal (PQW) frame. +namespace pqw { + + /** + * Converts PQW coordinates from Cartesian to spherical. + * + * @param v PQW Cartesian coordinates. + * @return PQW spherical coordinates, in the ISO order of radial distance, inclination (radians), and true anomaly (radians). + */ + template + math::vector3 spherical(const math::vector3& v) + { + const T xx_yy = v.x * v.x + v.y * v.y; + + return math::vector3 + { + std::sqrt(xx_yy + v.z * v.z), + std::atan2(v.z, std::sqrt(xx_yy)), + std::atan2(v.y, v.x) + }; + } + + /** + * Constructs spherical PQW coordinates from Keplerian orbital elements. + * + * @param ec Eccentricity (e). + * @param a Semimajor axis (a). + * @param ea Eccentric anomaly (E), in radians. + * @param b Semiminor axis (b). + * @return PQW spherical coordinates, in the ISO order of radial distance, inclination (radians), and true anomaly (radians). + */ + template + math::vector3 spherical(T ec, T a, T ea, T b) + { + const T x = a * (std::cos(ea) - ec); + const T y = b * std::sin(ea); + const T d = std::sqrt(x * x + y * y); + const T ta = std::atan2(y, x); + + return math::vector3{d, T(0), ta}; + } + + /** + * Constructs spherical PQW coordinates from Keplerian orbital elements. + * + * @param ec Eccentricity (e). + * @param a Semimajor axis (a). + * @param ea Eccentric anomaly (E), in radians. + * @return PQW spherical coordinates, in the ISO order of radial distance, inclination (radians), and true anomaly (radians). + */ + template + math::vector3 spherical(T ec, T a, T ea) + { + const T b = a * std::sqrt(T(1) - ec * ec); + return spherical(ec, a, ea, b); + } + + /** + * Converts PQW coordinates from spherical to Cartesian. + * + * @param PQW spherical coordinates, in the ISO order of radial distance, inclination (radians), and true anomaly (radians). + * @return PQW Cartesian coordinates. + */ + template + math::vector3 cartesian(const math::vector3& v) + { + const T x = v[0] * std::cos(v[1]); + + return math::vector3 + { + x * std::cos(v[2]), + x * std::sin(v[2]), + v[0] * std::sin(v[1]), + }; + } + + /** + * Constructs Cartesian PQW coordinates from Keplerian orbital elements. + * + * @param ec Eccentricity (e). + * @param a Semimajor axis (a). + * @param ea Eccentric anomaly (E), in radians. + * @param b Semiminor axis (b). + * @return PQW Cartesian coordinates. + */ + template + math::vector3 cartesian(T ec, T a, T ea, T b) + { + return cartesian(spherical(ec, a, ea, b)); + } + + /** + * Constructs Cartesian PQW coordinates from Keplerian orbital elements. + * + * @param ec Eccentricity (e). + * @param a Semimajor axis (a). + * @param ea Eccentric anomaly (E), in radians. + * @return PQW Cartesian coordinates. + */ + template + math::vector3 cartesian(T ec, T a, T ea) + { + return cartesian(spherical(ec, a, ea)); + } + + /** + * Constructs an SE(3) transformation from a PQW frame to a BCI frame. + * + * @param om Right ascension of the ascending node (OMEGA), in radians. + * @param in Orbital inclination (i), in radians. + * @param w Argument of periapsis (omega), in radians. + * @return PQW to BCI transformation. + */ + template + math::transformation::se3 to_bci(T om, T in, T w) + { + const math::quaternion r = math::normalize + ( + math::quaternion::rotate_z(om) * + math::quaternion::rotate_x(in) * + math::quaternion::rotate_z(w) + ); + + return math::transformation::se3{{T(0), T(0), T(0)}, r}; + } + +} // namespace pqw + +/// Body-centered inertial (BCI) frame. +namespace bci { + + /** + * Converts BCI coordinates from spherical to Cartesian. + * + * @param v BCI spherical coordinates, in the ISO order of radial distance, declination (radians), and right ascension (radians). + * @return BCI Cartesian coordinates. + */ + template + math::vector3 cartesian(const math::vector3& v) + { + const T x = v[0] * std::cos(v[1]); + + return math::vector3 + { + x * std::cos(v[2]), + x * std::sin(v[2]), + v[0] * std::sin(v[1]), + }; + } + + /** + * Converts BCI coordinates from Cartesian to spherical. + * + * @param v BCI Cartesian coordinates. + * @return BCI spherical coordinates, in the ISO order of radial distance, declination (radians), and right ascension (radians). + */ + template + math::vector3 spherical(const math::vector3& v) + { + const T xx_yy = v.x * v.x + v.y * v.y; + + return math::vector3 + { + std::sqrt(xx_yy + v.z * v.z), + std::atan2(v.z, std::sqrt(xx_yy)), + std::atan2(v.y, v.x) + }; + } + + /** + * Constructs an SE(3) transformation from a BCI frame to a BCBF frame. + * + * @param ra Right ascension of the north pole, in radians. + * @param dec Declination of the north pole, in radians. + * @param w Location of the prime meridian, as a rotation about the north pole, in radians. + * @return BCI to BCBF transformation. + */ + template + math::transformation::se3 to_bcbf(T ra, T dec, T w) + { + const math::quaternion r = math::normalize + ( + math::quaternion::rotate_z(-w) * + math::quaternion::rotate_x(-math::half_pi + dec) * + math::quaternion::rotate_z(-math::half_pi - ra) + ); + + return math::transformation::se3{{T(0), T(0), T(0)}, r}; + } + + /** + * Constructs an SE(3) transformation from a BCI frame to a PQW frame. + * + * @param om Right ascension of the ascending node (OMEGA), in radians. + * @param in Orbital inclination (i), in radians. + * @param w Argument of periapsis (omega), in radians. + * @return BCI to PQW transformation. + */ + template + math::transformation::se3 to_pqw(T om, T in, T w) + { + const math::quaternion r = math::normalize + ( + math::quaternion::rotate_z(-w) * + math::quaternion::rotate_x(-in) * + math::quaternion::rotate_z(-om) + ); + + return math::transformation::se3{{T(0), T(0), T(0)}, r}; + } + +} // namespace bci + +/// Body-centered, body-fixed (BCBF) frame. +namespace bcbf { + + /** + * Converts BCBF coordinates from spherical to Cartesian. + * + * @param v BCBF spherical coordinates, in the ISO order of radial distance, latitude (radians), and longitude (radians). + * @return BCBF Cartesian coordinates. + */ + template + math::vector3 cartesian(const math::vector3& v) + { + const T x = v[0] * std::cos(v[1]); + + return math::vector3 + { + x * std::cos(v[2]), + x * std::sin(v[2]), + v[0] * std::sin(v[1]), + }; + } + + /** + * Converts BCBF coordinates from Cartesian to spherical. + * + * @param v BCBF Cartesian coordinates. + * @return BCBF spherical coordinates, in the ISO order of radial distance, latitude (radians), and longitude (radians). + */ + template + math::vector3 spherical(const math::vector3& v) + { + const T xx_yy = v.x * v.x + v.y * v.y; + + return math::vector3 + { + std::sqrt(xx_yy + v.z * v.z), + std::atan2(v.z, std::sqrt(xx_yy)), + std::atan2(v.y, v.x) + }; + } + + /** + * Constructs an SE(3) transformation from a BCBF frame to a BCI frame. + * + * @param ra Right ascension of the north pole, in radians. + * @param dec Declination of the north pole, in radians. + * @param w Location of the prime meridian, as a rotation about the north pole, in radians. + * @return BCBF to BCI transformation. + */ + template + math::transformation::se3 to_bci(T ra, T dec, T w) + { + const math::quaternion r = math::normalize + ( + math::quaternion::rotate_z(math::half_pi + ra) * + math::quaternion::rotate_x(math::half_pi - dec) * + math::quaternion::rotate_z(w) + ); + + return math::transformation::se3{{T(0), T(0), T(0)}, r}; + } + + /** + * Constructs an SE(3) transformation from a BCBF frame to an ENU frame. + * + * @param distance Radial distance of the observer from the center of the body. + * @param latitude Latitude of the observer, in radians. + * @param longitude Longitude of the observer, in radians. + * @return BCBF to ENU transformation. + */ + template + math::transformation::se3 to_enu(T distance, T latitude, T longitude) + { + const math::vector3 t = {T(0), T(0), -distance}; + const math::quaternion r = math::normalize + ( + math::quaternion::rotate_x(-math::half_pi + latitude) * + math::quaternion::rotate_z(-longitude - math::half_pi) + ); + + return math::transformation::se3{t, r}; + } + +} // namespace bcbf + +/// East, North, Up (ENU) horizontal frame. +namespace enu { + + /** + * Converts ENU coordinates from spherical to Cartesian. + * + * @param v ENU spherical coordinates, in the ISO order of radial distance, elevation (radians), and azimuth (radians). + * @return ENU Cartesian coordinates. + */ + template + math::vector3 cartesian(const math::vector3& v) + { + const T x = v[0] * std::cos(v[1]); + const T y = math::half_pi - v[2]; + + return math::vector3 + { + x * std::cos(y), + x * std::sin(y), + v[0] * std::sin(v[1]), + }; + } + + /** + * Converts ENU coordinates from Cartesian to spherical. + * + * @param v ENU Cartesian coordinates. + * @return ENU spherical coordinates, in the ISO order of radial distance, elevation (radians), and azimuth (radians). + */ + template + math::vector3 spherical(const math::vector3& v) + { + const T xx_yy = v.x * v.x + v.y * v.y; + + return math::vector3 + { + std::sqrt(xx_yy + v.z * v.z), + std::atan2(v.z, std::sqrt(xx_yy)), + math::half_pi - std::atan2(v.y, v.x) + }; + } + + /** + * Constructs an SE(3) transformation from an ENU frame to a BCBF frame. + * + * @param distance Radial distance of the observer from the center of the body. + * @param latitude Latitude of the observer, in radians. + * @param longitude Longitude of the observer, in radians. + * @return ENU to BCBF transformation. + */ + template + math::transformation::se3 to_bcbf(T distance, T latitude, T longitude) + { + const math::vector3 t = {T(0), T(0), distance}; + const math::quaternion r = math::normalize + ( + math::quaternion::rotate_z(longitude + math::half_pi) * + math::quaternion::rotate_x(math::half_pi - latitude) + ); + + return math::transformation::se3{r * t, r}; + } + +} // namespace enu + +} // namespace frame +} // namespace orbit +} // namespace physics + +#endif // ANTKEEPER_PHYSICS_ORBIT_FRAME_HPP diff --git a/src/physics/orbit/frames.hpp b/src/physics/orbit/frames.hpp deleted file mode 100644 index 74f6f57..0000000 --- a/src/physics/orbit/frames.hpp +++ /dev/null @@ -1,131 +0,0 @@ -/* - * Copyright (C) 2021 Christopher J. Howard - * - * This file is part of Antkeeper source code. - * - * Antkeeper source code is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * Antkeeper source code is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with Antkeeper source code. If not, see . - */ - -#ifndef ANTKEEPER_PHYSICS_ORBIT_FRAMES_HPP -#define ANTKEEPER_PHYSICS_ORBIT_FRAMES_HPP - -#include "physics/frame.hpp" - -namespace physics { -namespace orbit { - -/// Inertial right-handed coordinate system -namespace inertial { - - /** - * Constucts a reference frame which transforms coordinates from inertial space into perifocal space. - * - * @param focus Cartesian coordinates of the focus of the orbit, in the parent space. - * @param raan Right ascension of the ascending node (OMEGA), in radians. - * @param i Orbital inclination (i), in radians. - * @param w Argument of periapsis (omega), in radians. - * @return Perifocal frame. - */ - template - physics::frame to_perifocal(const math::vector3& focus, T raan, T i, T w) - { - const math::quaternion rotation = math::normalize - ( - math::quaternion::rotate_z(raan) * - math::quaternion::rotate_x(i) * - math::quaternion::rotate_z(w) - ); - - return physics::frame{focus, rotation}.inverse(); - } - - /** - * Constructs a reference frame which transforms coordinates from inertial space to body-centered inertial space. - * - * @param r Cartesian position vector (r) of the center of the body, in inertial space. - * @param i Orbital inclination (i), in radians. - * @param axial_tilt Angle between the body's rotational axis and its orbital axis, in radians. - * @return Body-centered inertial frame. - */ - template - physics::frame to_bci(const math::vector3& r, T i, T axial_tilt) - { - return physics::frame{r, math::quaternion::rotate_x(-axial_tilt - i)}.inverse(); - } - - /** - * Constructs a reference frame which transforms coordinates from inertial space to body-centered, body-fixed space. - * - * @param r Cartesian position vector (r) of the center of the body, in inertial space. - * @param i Orbital inclination (i), in radians. - * @param axial_tilt Angle between the orbiting body's rotational axis and its orbital axis, in radians. - * @param axial_rotation Angle of rotation about the orbiting body's rotational axis, in radians. - */ - template - frame to_bcbf(const math::vector3& r, T i, T axial_tilt, T axial_rotation) - { - const math::quaternion rotation = math::normalize - ( - math::quaternion::rotate_x(-axial_tilt - i) * - math::quaternion::rotate_z(axial_rotation) - ); - - return physics::frame{r, rotation}.inverse(); - } - -} // namespace inertial - -/// Perifocal right-handed coordinate system in which the x-axis points toward the periapsis of the orbit, the y-axis has a true anomaly of 90 degrees past the periapsis, and the z-axis is the angular momentum vector, which is orthogonal to the orbital plane. -namespace perifocal { - - - -} // namespace perifocal - -/// Non-inertial body-centered, body-fixed right-handed coordinate system. The x-axis is orthogonal to the intersection of the prime meridian and the equator. The z-axis points toward the positive pole. The y-axis is right-hand orthogonal to the xz-plane. -namespace bcbf { - - /** - * Constructs a reference frame which transforms coordinates from BCBF space to topocentric space. - * - * @param distance Radial distance of the observer from the center of the body. - * @param latitude Latitude of the observer, in radians. - * @param longitude Longitude of the obserer, in radians. - * @return Topocentric frame. - */ - template - physics::frame to_topocentric(T distance, T latitude, T longitude) - { - const math::vector3 translation = {T(0), T(0), distance}; - - const math::quaternion rotation = math::normalize - ( - math::quaternion::rotate_z(longitude) * - math::quaternion::rotate_y(math::half_pi - latitude) - ); - - return physics::frame{rotation * translation, rotation}.inverse(); - } - -} // namespace bcbf - -/// Non-inertial topocentric right-handed coordinate system. Topocentric frames are constructed as south-east-zenith (SEZ) frames; the x-axis points south, the y-axis points east, and the z-axis points toward the zenith (orthogonal to reference spheroid). -namespace topocentric { - -} // namespace topocentric - -} // namespace orbit -} // namespace physics - -#endif // ANTKEEPER_PHYSICS_ORBIT_FRAMES_HPP diff --git a/src/physics/orbit/kepler.hpp b/src/physics/orbit/kepler.hpp deleted file mode 100644 index 637648f..0000000 --- a/src/physics/orbit/kepler.hpp +++ /dev/null @@ -1,79 +0,0 @@ -/* - * Copyright (C) 2021 Christopher J. Howard - * - * This file is part of Antkeeper source code. - * - * Antkeeper source code is free software: you can redistribute it and/or modify - * it under the terms of the GNU General Public License as published by - * the Free Software Foundation, either version 3 of the License, or - * (at your option) any later version. - * - * Antkeeper source code is distributed in the hope that it will be useful, - * but WITHOUT ANY WARRANTY; without even the implied warranty of - * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the - * GNU General Public License for more details. - * - * You should have received a copy of the GNU General Public License - * along with Antkeeper source code. If not, see . - */ - -#ifndef ANTKEEPER_PHYSICS_ORBIT_KEPLER_HPP -#define ANTKEEPER_PHYSICS_ORBIT_KEPLER_HPP - -#include - -namespace physics { -namespace orbit { - -/** - * Iteratively solves Kepler's equation for eccentric anomaly (E). - * - * @param ec Eccentricity (e). - * @param ma Mean anomaly (M). - * @param iterations Maximum number of iterations. - * @param tolerance Solution error tolerance. - * @return Eccentric anomaly (E). - */ -template -T kepler_ea(T ec, T ma, std::size_t iterations, T tolerance = T(0)); - -/** - * Solves Kepler's equation for mean anomaly (M). - * - * @param ec Eccentricity (e). - * @param ea Eccentric anomaly (E). - * @return Mean anomaly (M). - */ -template -T kepler_ma(T ec, T ea); - -template -T kepler_ea(T ec, T ma, std::size_t iterations, T tolerance) -{ - // Approximate eccentric anomaly, E - T ea0 = ma + ec * std::sin(ma) * (T(1.0) + ec * std::cos(ma)); - - // Iteratively converge E0 and E1 - for (std::size_t i = 0; i < iterations; ++i) - { - const T ea1 = ea0 - (ea0 - ec * std::sin(ea0) - ma) / (T(1.0) - ec * std::cos(ea0)); - const T error = std::abs(ea1 - ea0); - ea0 = ea1; - - if (error < tolerance) - break; - } - - return ea0; -} - -template -T kepler_ma(T ec, T ea) -{ - return ea - ec * std::sin(ea); -} - -} // namespace orbit -} // namespace physics - -#endif // ANTKEEPER_PHYSICS_ORBIT_KEPLER_HPP diff --git a/src/physics/orbit/orbit.hpp b/src/physics/orbit/orbit.hpp index 3871e40..c926c0c 100644 --- a/src/physics/orbit/orbit.hpp +++ b/src/physics/orbit/orbit.hpp @@ -31,9 +31,9 @@ namespace orbit {} } // namespace physics +#include "physics/orbit/anomaly.hpp" #include "physics/orbit/elements.hpp" -#include "physics/orbit/frames.hpp" -#include "physics/orbit/kepler.hpp" +#include "physics/orbit/frame.hpp" #include "physics/orbit/state.hpp" #endif // ANTKEEPER_PHYSICS_ORBIT_HPP diff --git a/src/render/passes/sky-pass.cpp b/src/render/passes/sky-pass.cpp index 74cbbf8..0333573 100644 --- a/src/render/passes/sky-pass.cpp +++ b/src/render/passes/sky-pass.cpp @@ -71,10 +71,10 @@ sky_pass::sky_pass(gl::rasterizer* rasterizer, const gl::framebuffer* framebuffe cloud_shader_program(nullptr), observer_altitude_tween(0.0f, math::lerp), sun_position_tween(float3{1.0f, 0.0f, 0.0f}, math::lerp), - sun_color_outer_tween(float3{1.0f, 1.0f, 1.0f}, math::lerp), - sun_color_inner_tween(float3{1.0f, 1.0f, 1.0f}, math::lerp), - topocentric_frame_translation({0, 0, 0}, math::lerp), - topocentric_frame_rotation(math::quaternion::identity(), math::nlerp) + sun_illuminance_outer_tween(float3{1.0f, 1.0f, 1.0f}, math::lerp), + sun_illuminance_inner_tween(float3{1.0f, 1.0f, 1.0f}, math::lerp), + icrf_to_eus_translation({0, 0, 0}, math::lerp), + icrf_to_eus_rotation(math::quaternion::identity(), math::nlerp) {} sky_pass::~sky_pass() @@ -109,20 +109,20 @@ void sky_pass::render(const render::context& ctx, render::queue& queue) const // Interpolate observer altitude float observer_altitude = observer_altitude_tween.interpolate(ctx.alpha); - // Construct tweened inertial to topocentric frame - physics::frame topocentric_frame = + // Construct tweened ICRF to EUS transformation + math::transformation::se3 icrf_to_eus = { - topocentric_frame_translation.interpolate(ctx.alpha), - topocentric_frame_rotation.interpolate(ctx.alpha) + icrf_to_eus_translation.interpolate(ctx.alpha), + icrf_to_eus_rotation.interpolate(ctx.alpha) }; - // Get topocentric space direction to sun + // Get EUS direction to sun float3 sun_position = sun_position_tween.interpolate(ctx.alpha); float3 sun_direction = math::normalize(sun_position); // Interpolate sun color - float3 sun_color_outer = sun_color_outer_tween.interpolate(ctx.alpha); - float3 sun_color_inner = sun_color_inner_tween.interpolate(ctx.alpha); + float3 sun_illuminance_outer = sun_illuminance_outer_tween.interpolate(ctx.alpha); + float3 sun_illuminance_inner = sun_illuminance_inner_tween.interpolate(ctx.alpha); // Draw atmosphere if (sky_model) @@ -149,8 +149,8 @@ void sky_pass::render(const render::context& ctx, render::queue& queue) const sun_angular_radius_input->upload(sun_angular_radius); // Pre-exposure sun color - if (sun_color_input) - sun_color_input->upload(sun_color_outer * ctx.exposure); + if (sun_illuminance_input) + sun_illuminance_input->upload(sun_illuminance_outer * ctx.exposure); if (scale_height_rm_input) scale_height_rm_input->upload(scale_height_rm); @@ -177,8 +177,8 @@ void sky_pass::render(const render::context& ctx, render::queue& queue) const cloud_model_view_projection_input->upload(model_view_projection); if (cloud_sun_direction_input) cloud_sun_direction_input->upload(sun_direction); - if (cloud_sun_color_input) - cloud_sun_color_input->upload(sun_color_inner); + if (cloud_sun_illuminance_input) + cloud_sun_illuminance_input->upload(sun_illuminance_inner); if (cloud_camera_position_input) cloud_camera_position_input->upload(ctx.camera_transform.translation); if (cloud_camera_exposure_input) @@ -198,7 +198,7 @@ void sky_pass::render(const render::context& ctx, render::queue& queue) const { float star_distance = (clip_near + clip_far) * 0.5f; - model = math::resize<4, 4>(math::matrix_cast(topocentric_frame.rotation)); + model = math::resize<4, 4>(math::matrix_cast(icrf_to_eus.r)); model = math::scale(model, {star_distance, star_distance, star_distance}); model_view = view * model; @@ -218,13 +218,12 @@ void sky_pass::render(const render::context& ctx, render::queue& queue) const rasterizer->draw_arrays(*stars_model_vao, stars_model_drawing_mode, stars_model_start_index, stars_model_index_count); } - // Draw moon model /* - float3 moon_position = {0, 0, 0}; + // Draw moon model + float3 moon_position = moon_position_tween.interpolate(ctx.alpha); + float moon_angular_radius = math::radians(2.0f); if (moon_position.y >= -moon_angular_radius) { - - float moon_distance = (clip_near + clip_far) * 0.5f; float moon_radius = moon_angular_radius * moon_distance; @@ -284,7 +283,7 @@ void sky_pass::set_sky_model(const model* model) observer_altitude_input = sky_shader_program->get_input("observer_altitude"); sun_direction_input = sky_shader_program->get_input("sun_direction"); - sun_color_input = sky_shader_program->get_input("sun_color"); + sun_illuminance_input = sky_shader_program->get_input("sun_illuminance"); sun_angular_radius_input = sky_shader_program->get_input("sun_angular_radius"); scale_height_rm_input = sky_shader_program->get_input("scale_height_rm"); rayleigh_scattering_input = sky_shader_program->get_input("rayleigh_scattering"); @@ -397,7 +396,7 @@ void sky_pass::set_clouds_model(const model* model) { cloud_model_view_projection_input = cloud_shader_program->get_input("model_view_projection"); cloud_sun_direction_input = cloud_shader_program->get_input("sun_direction"); - cloud_sun_color_input = cloud_shader_program->get_input("sun_color"); + cloud_sun_illuminance_input = cloud_shader_program->get_input("sun_illuminance"); cloud_camera_position_input = cloud_shader_program->get_input("camera.position"); cloud_camera_exposure_input = cloud_shader_program->get_input("camera.exposure"); } @@ -413,16 +412,17 @@ void sky_pass::update_tweens() { observer_altitude_tween.update(); sun_position_tween.update(); - sun_color_outer_tween.update(); - sun_color_inner_tween.update(); - topocentric_frame_translation.update(); - topocentric_frame_rotation.update(); + sun_illuminance_outer_tween.update(); + sun_illuminance_inner_tween.update(); + icrf_to_eus_translation.update(); + icrf_to_eus_rotation.update(); + moon_position_tween.update(); } -void sky_pass::set_topocentric_frame(const physics::frame& frame) +void sky_pass::set_icrf_to_eus(const math::transformation::se3& transformation) { - topocentric_frame_translation[1] = frame.translation; - topocentric_frame_rotation[1] = frame.rotation; + icrf_to_eus_translation[1] = transformation.t; + icrf_to_eus_rotation[1] = transformation.r; } void sky_pass::set_sun_position(const float3& position) @@ -430,10 +430,10 @@ void sky_pass::set_sun_position(const float3& position) sun_position_tween[1] = position; } -void sky_pass::set_sun_color(const float3& color_outer, const float3& color_inner) +void sky_pass::set_sun_illuminance(const float3& illuminance_outer, const float3& illuminance_inner) { - sun_color_outer_tween[1] = color_outer; - sun_color_inner_tween[1] = color_inner; + sun_illuminance_outer_tween[1] = illuminance_outer; + sun_illuminance_inner_tween[1] = illuminance_inner; } void sky_pass::set_sun_angular_radius(float radius) @@ -469,6 +469,11 @@ void sky_pass::set_atmosphere_radii(float inner, float outer) atmosphere_radii.z = outer * outer; } +void sky_pass::set_moon_position(const float3& position) +{ + moon_position_tween[1] = position; +} + void sky_pass::handle_event(const mouse_moved_event& event) { mouse_position = {static_cast(event.x), static_cast(event.y)}; diff --git a/src/render/passes/sky-pass.hpp b/src/render/passes/sky-pass.hpp index 4c4593e..723c549 100644 --- a/src/render/passes/sky-pass.hpp +++ b/src/render/passes/sky-pass.hpp @@ -32,7 +32,7 @@ #include "gl/vertex-array.hpp" #include "gl/texture-2d.hpp" #include "gl/drawing-mode.hpp" -#include "physics/frame.hpp" +#include "math/se3.hpp" #include "scene/object.hpp" class resource_manager; @@ -60,16 +60,18 @@ public: void set_stars_model(const model* model); void set_clouds_model(const model* model); - void set_topocentric_frame(const physics::frame& frame); + void set_icrf_to_eus(const math::transformation::se3& transformation); void set_sun_position(const float3& position); - void set_sun_color(const float3& color_outer, const float3& color_inner); + void set_sun_illuminance(const float3& illuminance_outer, const float3& illuminance_inner); void set_sun_angular_radius(float radius); void set_observer_altitude(float altitude); void set_scale_heights(float rayleigh, float mie); void set_scattering_coefficients(const float3& r, const float3& m); void set_mie_anisotropy(float g); void set_atmosphere_radii(float inner, float outer); + + void set_moon_position(const float3& position); private: virtual void handle_event(const mouse_moved_event& event); @@ -82,7 +84,7 @@ private: const gl::shader_input* exposure_input; const gl::shader_input* observer_altitude_input; const gl::shader_input* sun_direction_input; - const gl::shader_input* sun_color_input; + const gl::shader_input* sun_illuminance_input; const gl::shader_input* sun_angular_radius_input; const gl::shader_input* scale_height_rm_input; const gl::shader_input* rayleigh_scattering_input; @@ -132,7 +134,7 @@ private: gl::shader_program* cloud_shader_program; const gl::shader_input* cloud_model_view_projection_input; const gl::shader_input* cloud_sun_direction_input; - const gl::shader_input* cloud_sun_color_input; + const gl::shader_input* cloud_sun_illuminance_input; const gl::shader_input* cloud_camera_position_input; const gl::shader_input* cloud_camera_exposure_input; @@ -143,10 +145,13 @@ private: tween observer_altitude_tween; tween sun_position_tween; - tween sun_color_outer_tween; - tween sun_color_inner_tween; - tween topocentric_frame_translation; - tween> topocentric_frame_rotation; + tween sun_illuminance_outer_tween; + tween sun_illuminance_inner_tween; + tween icrf_to_eus_translation; + tween> icrf_to_eus_rotation; + + tween moon_position_tween; + float sun_angular_radius; float2 scale_height_rm; diff --git a/src/resources/entity-archetype-loader.cpp b/src/resources/entity-archetype-loader.cpp index 97ff38b..55dfa80 100644 --- a/src/resources/entity-archetype-loader.cpp +++ b/src/resources/entity-archetype-loader.cpp @@ -31,6 +31,7 @@ #include "entity/components/celestial-body.hpp" #include "entity/archetype.hpp" #include "entity/ebt.hpp" +#include "physics/orbit/elements.hpp" #include "resources/json.hpp" #include @@ -97,18 +98,24 @@ static bool load_component_celestial_body(entity::archetype& archetype, const js { entity::component::celestial_body component; component.radius = 0.0; - component.axial_tilt = 0.0; - component.axial_rotation = 0.0; - component.angular_frequency = 0.0; + component.mass = 0.0; + component.pole_ra = 0.0; + component.pole_dec = 0.0; + component.prime_meridian = 0.0; + component.rotation_period = 0.0; if (element.contains("radius")) component.radius = element["radius"].get(); - if (element.contains("axial_tilt")) - component.axial_tilt = math::radians(element["axial_tilt"].get()); - if (element.contains("axial_rotation")) - component.axial_rotation = math::radians(element["axial_rotation"].get()); - if (element.contains("angular_frequency")) - component.angular_frequency = math::radians(element["angular_frequency"].get()); + if (element.contains("mass")) + component.mass = element["mass"].get(); + if (element.contains("pole_ra")) + component.pole_ra = math::radians(element["pole_ra"].get()); + if (element.contains("pole_dec")) + component.pole_dec = math::radians(element["pole_dec"].get()); + if (element.contains("prime_meridian")) + component.prime_meridian = math::radians(element["prime_meridian"].get()); + if (element.contains("rotation_period")) + component.rotation_period = element["rotation_period"].get(); archetype.set(component); @@ -151,26 +158,30 @@ static bool load_component_orbit(entity::archetype& archetype, const json& eleme { entity::component::orbit component; - component.elements.e = 0.0; + component.parent = entt::null; + component.elements.ec = 0.0; component.elements.a = 0.0; - component.elements.i = 0.0; - component.elements.raan = 0.0; + component.elements.in = 0.0; + component.elements.om = 0.0; component.elements.w = 0.0; - component.elements.ta = 0.0; + component.elements.ma = 0.0; + component.mean_motion = 0.0; - if (element.contains("e")) - component.elements.e = element["e"].get(); + if (element.contains("ec")) + component.elements.ec = element["ec"].get(); if (element.contains("a")) component.elements.a = element["a"].get(); - if (element.contains("i")) - component.elements.i = math::radians(element["i"].get()); - if (element.contains("raan")) - component.elements.raan = math::radians(element["raan"].get()); + if (element.contains("in")) + component.elements.in = math::radians(element["in"].get()); + if (element.contains("om")) + component.elements.om = math::radians(element["om"].get()); if (element.contains("w")) component.elements.w = math::radians(element["w"].get()); - if (element.contains("ta")) - component.elements.ta = math::radians(element["ta"].get()); - + if (element.contains("ma")) + component.elements.ma = math::radians(element["ma"].get()); + if (element.contains("n")) + component.mean_motion = math::radians(element["n"].get()); + archetype.set(component); return true;