From 94c6e100324b3824d67e11693a857898156cb108 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Gr=C3=A9goire=20Henry?= Date: Thu, 30 Nov 2023 18:05:12 +0100 Subject: [PATCH] module / macro / coords / spice / routines - renamed modules for more explicit name, more will be added and refactored soon - initialized pro macros sub folder for derive trait of body / routines - more work on equatorial coordinates, angle, and conversion to xyz - example with spice for DART impact - reworked and rearranged routines --- Cargo.lock | 10 + Cargo.toml | 1 + .../cfg/cfg.yaml | 26 + examples/viewer-picker/main.rs | 105 +--- macros/Cargo.lock | 47 ++ macros/Cargo.toml | 14 + macros/src/lib.rs | 0 preferences.yaml | 0 src/{ast => body}/asteroid.rs | 4 +- src/{ast => body}/element.rs | 0 src/{ast => body}/interior.rs | 0 src/{ast => body}/material.rs | 0 src/{ast => body}/mesh.rs | 4 +- src/{ast => body}/mod.rs | 0 src/{ast => body}/orbit.rs | 35 +- src/{ast => body}/ray.rs | 12 +- src/{ast => body}/surface.rs | 0 src/cfg/body.rs | 16 + src/cfg/config.rs | 7 +- src/cfg/preferences.rs | 8 + src/cfg/scene.rs | 17 +- src/cfg/simulation.rs | 72 ++- src/lib.rs | 6 +- src/simu/body.rs | 12 +- src/simu/converge.rs | 6 +- src/simu/export.rs | 57 +- src/simu/routines/core.rs | 262 ++++---- src/simu/routines/thermal.rs | 560 ++++++++---------- src/simu/routines/viewer.rs | 132 +---- src/simu/scenario.rs | 121 ++-- src/simu/util.rs | 93 ++- src/win/scene.rs | 4 + src/win/window.rs | 23 +- src/win/window_settings.rs | 5 + 34 files changed, 872 insertions(+), 787 deletions(-) create mode 100644 examples/ground-based-lightcurve-didymos-spice/cfg/cfg.yaml create mode 100644 macros/Cargo.lock create mode 100644 macros/Cargo.toml create mode 100644 macros/src/lib.rs create mode 100644 preferences.yaml rename src/{ast => body}/asteroid.rs (98%) rename src/{ast => body}/element.rs (100%) rename src/{ast => body}/interior.rs (100%) rename src/{ast => body}/material.rs (100%) rename src/{ast => body}/mesh.rs (98%) rename src/{ast => body}/mod.rs (100%) rename src/{ast => body}/orbit.rs (95%) rename src/{ast => body}/ray.rs (95%) rename src/{ast => body}/surface.rs (100%) diff --git a/Cargo.lock b/Cargo.lock index bbd0699..c720d8a 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -1534,6 +1534,7 @@ dependencies = [ "image", "indicatif", "itertools 0.12.0", + "kalast_macros", "lazy_static", "log", "nalgebra 0.32.3", @@ -1558,6 +1559,15 @@ dependencies = [ "uom", ] +[[package]] +name = "kalast_macros" +version = "0.4.0-beta" +dependencies = [ + "proc-macro2", + "quote", + "syn 2.0.31", +] + [[package]] name = "khronos_api" version = "3.1.0" diff --git a/Cargo.toml b/Cargo.toml index 2427a3c..969df3b 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -52,6 +52,7 @@ lazy_static = "1.0" dunce = "1.0" uom = "0.35" regex = "1.10" +kalast_macros = { path = "macros", version = "0.4.0-beta" } [dev-dependencies] rstest = "0.18" diff --git a/examples/ground-based-lightcurve-didymos-spice/cfg/cfg.yaml b/examples/ground-based-lightcurve-didymos-spice/cfg/cfg.yaml new file mode 100644 index 0000000..dfb2527 --- /dev/null +++ b/examples/ground-based-lightcurve-didymos-spice/cfg/cfg.yaml @@ -0,0 +1,26 @@ +window: + export_frames: true + +simulation: + start: 2022-09-26 + step: 60 + duration: 8136 + export: + step: 60 + duration: 8136 + +scene: + spice: true + camera: + !earth 2.0 + +bodies: + - id: Didymos + mesh: + shape: sphere + factor: [0.4095, 0.4005, 0.3035] + smooth: true + material: + albedo: 0.1 + spin: + period: 8136 \ No newline at end of file diff --git a/examples/viewer-picker/main.rs b/examples/viewer-picker/main.rs index d7ae53a..55acf0d 100644 --- a/examples/viewer-picker/main.rs +++ b/examples/viewer-picker/main.rs @@ -1,6 +1,6 @@ use kalast::{ - routines_viewer_default, simu::Scene, util::*, Asteroid, Body, CfgBody, CfgColormap, ColorMode, - FoldersRun, Result, Routines, RoutinesViewerDefault, Scenario, Time, Window, + simu::Scene, util::*, AirlessBody, ColorMode, Result, Routines, RoutinesViewer, Scenario, Time, + Window, }; fn main() -> Result<()> { @@ -29,12 +29,12 @@ fn main() -> Result<()> { fn update_surf_face( routines: &mut RoutinesViewerCustom, - bodies: &mut [Body], + bodies: &mut [AirlessBody], win: &Window, faces: &[usize], ii_body: usize, ) { - let surf = &mut bodies[ii_body].asteroid.surface; + let surf = &mut bodies[ii_body].surface; for &ii_face in faces { if !routines.faces.contains(&ii_face) { @@ -55,95 +55,19 @@ fn update_surf_face( } pub struct RoutinesViewerCustom { - pub default: RoutinesViewerDefault, pub faces: Vec, } -impl Routines for RoutinesViewerCustom { - fn fn_setup_body(&mut self, asteroid: Asteroid, cb: &CfgBody, scene: &Scene) -> Body { - self.default.fn_setup_body(asteroid, cb, scene) - } - fn fn_update_matrix_model( - &self, - ii_body: usize, - ii_other_bodies: &[usize], - cb: &CfgBody, - other_cbs: &[&CfgBody], - bodies: &mut [Body], - time: &Time, - scene: &mut Scene, - ) { - self.default.fn_update_matrix_model( - ii_body, - ii_other_bodies, - cb, - other_cbs, - bodies, - time, - scene, - ) - } - - fn fn_iteration_body( - &mut self, - ii_body: usize, - ii_other_bodies: &[usize], - cb: &CfgBody, - other_cbs: &[&CfgBody], - bodies: &mut [Body], - scene: &Scene, - time: &Time, - ) { - self.default - .fn_iteration_body(ii_body, ii_other_bodies, cb, other_cbs, bodies, scene, time) - } - - fn fn_update_colormap( - &self, - win: &Window, - cmap: &CfgColormap, - ii_body: usize, - body: &mut Body, - scene: &Scene, - ) { - self.default - .fn_update_colormap(win, cmap, ii_body, body, scene) - } - - fn fn_export_iteration( - &self, - cb: &CfgBody, - ii_body: usize, - time: &Time, - folders: &FoldersRun, - is_first_it: bool, - ) { - self.default - .fn_export_iteration(cb, ii_body, time, folders, is_first_it) - } - - fn fn_export_iteration_period( - &self, - cb: &CfgBody, - body: &Body, - ii_body: usize, - folders: &FoldersRun, - exporting_started_elapsed: i64, - is_first_it_export: bool, - ) { - self.default.fn_export_iteration_period( - cb, - body, - ii_body, - folders, - exporting_started_elapsed, - is_first_it_export, - ) +impl RoutinesViewerCustom { + pub fn new() -> Self { + Self { faces: vec![] } } +} +impl Routines for RoutinesViewerCustom { fn fn_end_of_iteration( &mut self, - bodies: &mut [Body], + bodies: &mut [AirlessBody], _time: &Time, _scene: &Scene, win: &Window, @@ -154,11 +78,4 @@ impl Routines for RoutinesViewerCustom { } } -impl RoutinesViewerCustom { - pub fn new() -> Self { - Self { - default: routines_viewer_default(), - faces: vec![], - } - } -} +impl RoutinesViewer for RoutinesViewerCustom {} diff --git a/macros/Cargo.lock b/macros/Cargo.lock new file mode 100644 index 0000000..35206e6 --- /dev/null +++ b/macros/Cargo.lock @@ -0,0 +1,47 @@ +# This file is automatically @generated by Cargo. +# It is not intended for manual editing. +version = 3 + +[[package]] +name = "kalast_macros" +version = "0.4.0-beta" +dependencies = [ + "proc-macro2", + "quote", + "syn", +] + +[[package]] +name = "proc-macro2" +version = "1.0.70" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "39278fbbf5fb4f646ce651690877f89d1c5811a3d4acb27700c1cb3cdb78fd3b" +dependencies = [ + "unicode-ident", +] + +[[package]] +name = "quote" +version = "1.0.33" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "5267fca4496028628a95160fc423a33e8b2e6af8a5302579e322e4b520293cae" +dependencies = [ + "proc-macro2", +] + +[[package]] +name = "syn" +version = "2.0.39" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "23e78b90f2fcf45d3e842032ce32e3f2d1545ba6636271dcbf24fa306d87be7a" +dependencies = [ + "proc-macro2", + "quote", + "unicode-ident", +] + +[[package]] +name = "unicode-ident" +version = "1.0.12" +source = "registry+https://github.com/rust-lang/crates.io-index" +checksum = "3354b9ac3fae1ff6755cb6db53683adb661634f67557942dea4facebec0fee4b" diff --git a/macros/Cargo.toml b/macros/Cargo.toml new file mode 100644 index 0000000..b6b70bc --- /dev/null +++ b/macros/Cargo.toml @@ -0,0 +1,14 @@ +[package] +name = "kalast_macros" +version = "0.4.0-beta" +description = "Kalast Macros" +edition = "2021" +license = "Apache-2.0" + +[lib] +proc-macro = true + +[dependencies] +syn = "2.0" +quote = "1.0" +proc-macro2 = "1.0" \ No newline at end of file diff --git a/macros/src/lib.rs b/macros/src/lib.rs new file mode 100644 index 0000000..e69de29 diff --git a/preferences.yaml b/preferences.yaml new file mode 100644 index 0000000..e69de29 diff --git a/src/ast/asteroid.rs b/src/body/asteroid.rs similarity index 98% rename from src/ast/asteroid.rs rename to src/body/asteroid.rs index 9a3af9c..af501d3 100644 --- a/src/ast/asteroid.rs +++ b/src/body/asteroid.rs @@ -5,13 +5,13 @@ use polars::prelude::{CsvReader, SerReader}; use std::path::Path; #[derive(Debug, Clone)] -pub struct Asteroid { +pub struct AirlessBody { pub surface: Surface, pub interior: Option, pub matrix_model: Mat4, } -impl Asteroid { +impl AirlessBody { pub fn new(surface: Surface) -> Self { Self { surface, diff --git a/src/ast/element.rs b/src/body/element.rs similarity index 100% rename from src/ast/element.rs rename to src/body/element.rs diff --git a/src/ast/interior.rs b/src/body/interior.rs similarity index 100% rename from src/ast/interior.rs rename to src/body/interior.rs diff --git a/src/ast/material.rs b/src/body/material.rs similarity index 100% rename from src/ast/material.rs rename to src/body/material.rs diff --git a/src/ast/mesh.rs b/src/body/mesh.rs similarity index 98% rename from src/ast/mesh.rs rename to src/body/mesh.rs index 8b2a1ba..36f1531 100644 --- a/src/ast/mesh.rs +++ b/src/body/mesh.rs @@ -1,4 +1,4 @@ -use crate::{intersect_surface, util::*, vec3_to_4_one, Asteroid, FaceData}; +use crate::{intersect_surface, util::*, vec3_to_4_one, AirlessBody, FaceData}; use itertools::{iproduct, izip}; @@ -49,7 +49,7 @@ where $a_i$ is the area of the facet $i$, $\Theta_x$ the angle between the cente facets $i$ and $j$ and the normal of the facet $x$, and $r$ is the distance between the facets. */ -pub fn view_factor(b1: &Asteroid, b2: &Asteroid, shadows: bool) -> DMatrix { +pub fn view_factor(b1: &AirlessBody, b2: &AirlessBody, shadows: bool) -> DMatrix { assert!(!b1.surface.is_smooth()); assert!(!b2.surface.is_smooth()); diff --git a/src/ast/mod.rs b/src/body/mod.rs similarity index 100% rename from src/ast/mod.rs rename to src/body/mod.rs diff --git a/src/ast/orbit.rs b/src/body/orbit.rs similarity index 95% rename from src/ast/orbit.rs rename to src/body/orbit.rs index b0ba311..48f87b9 100644 --- a/src/ast/orbit.rs +++ b/src/body/orbit.rs @@ -3,7 +3,10 @@ use crate::{util::*, NEWTON_METHOD_THRESHOLD, NUMBER_ITERATION_FAIL}; use regex::Regex; use serde::{Deserialize, Serialize}; use snafu::{prelude::*, Location}; -use std::str::FromStr; +use std::{ + ops::{Deref, DerefMut}, + str::FromStr, +}; use uom::{ si::{angle::radian, f64::Angle}, str::ParseQuantityError, @@ -64,13 +67,11 @@ pub enum AstronomicalAngleOptions { #[derive(Debug, Clone, Default, Serialize, Deserialize)] #[serde(try_from = "AstronomicalAngleStr")] #[serde(into = "AstronomicalAngleStr")] -pub struct AstronomicalAngle { - pub angle: Angle, -} +pub struct AstronomicalAngle(Angle); impl AstronomicalAngle { pub fn new(angle: Angle) -> Self { - Self { angle } + Self(angle) } pub fn parse(s: &str) -> OrbitResult { @@ -179,6 +180,19 @@ impl FromStr for AstronomicalAngle { } } +impl Deref for AstronomicalAngle { + type Target = Angle; + fn deref(&self) -> &Self::Target { + &self.0 + } +} + +impl DerefMut for AstronomicalAngle { + fn deref_mut(&mut self) -> &mut Self::Target { + &mut self.0 + } +} + #[derive(Serialize, Deserialize)] struct AstronomicalAngleStr(String); @@ -192,7 +206,7 @@ impl TryFrom for AstronomicalAngle { impl From for AstronomicalAngleStr { fn from(value: AstronomicalAngle) -> Self { - Self(value.angle.value.to_string()) + Self(value.value.to_string()) } } @@ -221,6 +235,15 @@ impl Equatorial { pub fn new(ra: AstronomicalAngle, dec: AstronomicalAngle) -> Self { Self { ra, dec } } + + pub fn xyz(&self, distance: Float) -> Vec3 { + distance + * vec3( + (self.ra.cos() * self.dec.cos()).value, + (self.ra.sin() * self.dec.cos()).value, + self.dec.sin().value, + ) + } } pub fn mean_angular_motion(gravitational_parameter: Float, semi_major_axis: Float) -> Float { diff --git a/src/ast/ray.rs b/src/body/ray.rs similarity index 95% rename from src/ast/ray.rs rename to src/body/ray.rs index 126bc3e..c38fbd7 100644 --- a/src/ast/ray.rs +++ b/src/body/ray.rs @@ -1,4 +1,4 @@ -use crate::{util::*, vec3_to_4_one, Asteroid}; +use crate::{util::*, vec3_to_4_one, AirlessBody}; use itertools::izip; @@ -71,11 +71,11 @@ pub fn intersect_surface( pub fn intersect_asteroids( raystart_world: &Vec3, raydir_world: &Vec3, - asteroids: &[&Asteroid], + asteroids: &[AirlessBody], ) -> Option<(Vec3, usize, usize)> { let mut best_intersect: Option<(Vec3, usize, usize)> = None; - for (surface_index, &asteroid) in asteroids.iter().enumerate() { + for (surface_index, asteroid) in asteroids.iter().enumerate() { let raystart = (glm::inverse(&asteroid.matrix_model) * vec3_to_4_one(&raystart_world)).xyz(); @@ -114,7 +114,11 @@ pub fn intersect_asteroids( } /// asteroid1 is shadowed by asteroid2. -pub fn shadows(lightpos_world: &Vec3, asteroid1: &Asteroid, asteroid2: &Asteroid) -> Vec { +pub fn shadows( + lightpos_world: &Vec3, + asteroid1: &AirlessBody, + asteroid2: &AirlessBody, +) -> Vec { if asteroid1.surface.is_smooth() || asteroid2.surface.is_smooth() { unimplemented!("Shadowing by ray-casting is only implemented for flat surface."); } diff --git a/src/ast/surface.rs b/src/body/surface.rs similarity index 100% rename from src/ast/surface.rs rename to src/body/surface.rs diff --git a/src/cfg/body.rs b/src/cfg/body.rs index 2c92694..1ff98d3 100644 --- a/src/cfg/body.rs +++ b/src/cfg/body.rs @@ -576,6 +576,22 @@ pub struct CfgStateCartesian { pub orientation: Mat3, } +impl CfgStateCartesian { + pub fn position_only(position: Vec3) -> Self { + Self { + position, + orientation: Mat3::identity(), + } + } + + pub fn orientation_only(orientation: Mat3) -> Self { + Self { + position: Vec3::zeros(), + orientation, + } + } +} + impl Default for CfgStateCartesian { fn default() -> Self { Self { diff --git a/src/cfg/config.rs b/src/cfg/config.rs index 8cccf1e..c4d9e04 100644 --- a/src/cfg/config.rs +++ b/src/cfg/config.rs @@ -274,8 +274,11 @@ impl Cfg { } } - dbg!(cfg.scene.sun.as_equatorial().unwrap()); - dbg!(cfg.bodies.first().unwrap().state.as_equatorial().unwrap()); + if cfg.preferences.debug_cfg { + println!("\nDebug cfg:"); + dbg!(&cfg); + println!(); + } Ok(cfg) } diff --git a/src/cfg/preferences.rs b/src/cfg/preferences.rs index e5df2af..9784648 100644 --- a/src/cfg/preferences.rs +++ b/src/cfg/preferences.rs @@ -14,6 +14,12 @@ pub struct CfgPreferences { #[serde(default)] pub auto_update: bool, + + #[serde(default)] + pub debug: bool, + + #[serde(default)] + pub debug_cfg: bool, } impl Configuration for CfgPreferences {} @@ -24,6 +30,8 @@ impl Default for CfgPreferences { runs: default_runs(), auto_update: false, do_not_check_latest_version: false, + debug: false, + debug_cfg: false, } } } diff --git a/src/cfg/scene.rs b/src/cfg/scene.rs index 5a848ab..6067c56 100644 --- a/src/cfg/scene.rs +++ b/src/cfg/scene.rs @@ -1,3 +1,5 @@ +use std::path::PathBuf; + use crate::{util::*, CfgBodyError, CfgStateCartesian, Configuration, Equatorial}; use serde::{Deserialize, Serialize}; @@ -16,6 +18,9 @@ pub enum CfgSceneError { #[derive(Clone, Debug, Serialize, Deserialize)] pub struct CfgScene { + #[serde(default)] + pub spice: Option, + #[serde(default)] pub camera: CfgCamera, @@ -26,6 +31,7 @@ pub struct CfgScene { impl Default for CfgScene { fn default() -> Self { Self { + spice: None, camera: CfgCamera::default(), sun: CfgSun::default(), } @@ -48,6 +54,15 @@ pub enum CfgCamera { Earth(Float), } +impl CfgCamera { + pub fn as_earth(&self) -> SceneResult { + match self { + Self::Earth(distance) => Ok(*distance), + _ => panic!("Not in Earth mode."), + } + } +} + impl Default for CfgCamera { fn default() -> Self { Self::Sun(5.0) @@ -77,6 +92,6 @@ impl CfgSun { impl Default for CfgSun { fn default() -> Self { - Self::Cartesian(CfgStateCartesian::default()) + Self::Cartesian(CfgStateCartesian::position_only(Vec3::x())) } } diff --git a/src/cfg/simulation.rs b/src/cfg/simulation.rs index 837a973..484479d 100644 --- a/src/cfg/simulation.rs +++ b/src/cfg/simulation.rs @@ -1,15 +1,24 @@ use crate::{util::*, Configuration}; use serde::{Deserialize, Serialize}; +use snafu::prelude::*; -#[derive(Clone, Debug, Serialize, Deserialize, Default)] +pub type SimulationResult = std::result::Result; + +/// Errors related to Kalast config. +#[derive(Debug, Snafu)] +pub enum CfgSimulationError { + SpiceNotLoaded {}, +} + +#[derive(Clone, Debug, Serialize, Deserialize)] pub struct CfgSimulation { #[serde(default)] pub routines: CfgRoutines, // In seconds. #[serde(default)] - pub start: Float, + pub start: TimeOption, #[serde(default)] pub step: usize, @@ -19,10 +28,56 @@ pub struct CfgSimulation { #[serde(default)] pub export: CfgTimeExport, + + #[serde(default)] + pub pause_after_first_iteration: bool, +} + +impl Default for CfgSimulation { + fn default() -> Self { + Self { + routines: CfgRoutines::default(), + start: TimeOption::default(), + step: 0, + duration: 0, + export: CfgTimeExport::default(), + pause_after_first_iteration: false, + } + } } impl Configuration for CfgSimulation {} +#[derive(Clone, Debug, Serialize, Deserialize, PartialEq)] +#[serde(untagged)] +pub enum TimeOption { + #[serde(rename = "seconds")] + Seconds(Float), + + #[serde(rename = "string")] + String(String), +} + +impl TimeOption { + pub fn seconds(&self) -> SimulationResult { + match self { + Self::Seconds(v) => Ok(*v), + Self::String(_s) => { + #[cfg(feature = "spice")] + return Ok(spice::str2et(_s)); + #[cfg(not(feature = "spice"))] + return Err(CfgSimulationError::SpiceNotLoaded {}); + } + } + } +} + +impl Default for TimeOption { + fn default() -> Self { + Self::Seconds(0.0) + } +} + #[derive(Clone, Debug, Serialize, Deserialize, PartialEq)] pub enum CfgRoutines { #[serde(rename = "viewer")] @@ -38,7 +93,7 @@ impl Default for CfgRoutines { } } -#[derive(Clone, Debug, Serialize, Deserialize, Default)] +#[derive(Clone, Debug, Serialize, Deserialize)] pub struct CfgTimeExport { #[serde(default)] pub step: usize, @@ -52,3 +107,14 @@ pub struct CfgTimeExport { #[serde(default)] pub cooldown_start: Option, } + +impl Default for CfgTimeExport { + fn default() -> Self { + Self { + step: 0, + duration: 0, + period: 0, + cooldown_start: None, + } + } +} diff --git a/src/lib.rs b/src/lib.rs index ddc5db9..213e202 100644 --- a/src/lib.rs +++ b/src/lib.rs @@ -20,21 +20,21 @@ Information on the configuration of your **kalast** scenarios are located at the [releases]: https://github.com/GregoireHENRY/kalast/releases */ -pub mod ast; +pub mod body; pub mod cfg; pub mod simu; pub mod thermal; pub mod util; pub mod win; -pub use ast::*; +pub use body::*; pub use cfg::*; pub use simu::*; pub use thermal::*; pub use util::*; pub use win::*; -#[cfg(feature = "rust_spice")] +#[cfg(feature = "spice")] pub use spice; use snafu::prelude::*; diff --git a/src/simu/body.rs b/src/simu/body.rs index e48edc3..badca39 100644 --- a/src/simu/body.rs +++ b/src/simu/body.rs @@ -1,17 +1,15 @@ -use crate::{matrix_orientation_obliquity, util::*, Asteroid, CfgBody}; +use crate::{matrix_orientation_obliquity, util::*, AirlessBody, CfgBody}; use itertools::Itertools; #[derive(Debug, Clone)] -pub struct Body { - pub id: String, - pub asteroid: Asteroid, +pub struct PreComputedBody { pub mat_orient: Mat4, pub normals: Matrix3xX, } -impl Body { - pub fn new(asteroid: Asteroid, cb: &CfgBody) -> Self { +impl PreComputedBody { + pub fn new(asteroid: &AirlessBody, cb: &CfgBody) -> Self { let mat_orient = matrix_orientation_obliquity(0.0, cb.spin.obliquity * RPD); let normals = Matrix3xX::from_columns( @@ -24,8 +22,6 @@ impl Body { ); Self { - id: cb.id.clone(), - asteroid, mat_orient, normals, } diff --git a/src/simu/converge.rs b/src/simu/converge.rs index 9ab9aea..eebc49e 100644 --- a/src/simu/converge.rs +++ b/src/simu/converge.rs @@ -1,4 +1,4 @@ -use crate::{util::*, Asteroid, CfgBody, CfgTimeExport, ThermalData}; +use crate::{util::*, AirlessBody, CfgBody, CfgTimeExport, ThermalData}; use itertools::izip; use notify_rust::Notification; @@ -9,7 +9,7 @@ use std::path::Path; pub fn check>( id: usize, - asteroid: &mut Asteroid, + asteroid: &mut AirlessBody, cb: &CfgBody, info: &mut ThermalData, path: P, @@ -182,7 +182,7 @@ pub fn check>( } pub fn check_all>( - asteroids: &mut [Asteroid], + asteroids: &mut [AirlessBody], cbs: &[CfgBody], infos: &mut [ThermalData], path: P, diff --git a/src/simu/export.rs b/src/simu/export.rs index cf8f1ed..cdddfac 100644 --- a/src/simu/export.rs +++ b/src/simu/export.rs @@ -1,8 +1,8 @@ use crate::{ - simu::Scene, util::*, Body, Cfg, CfgBody, CfgTimeExport, FoldersRun, Routines, Time, Window, + simu::Scene, util::*, AirlessBody, Cfg, CfgTimeExport, FoldersRun, PreComputedBody, + Routines, Time, Window, }; -use itertools::izip; use polars::prelude::{df, CsvWriter, NamedFrom, SerWriter}; use std::fs; @@ -29,13 +29,14 @@ impl Export { pub fn iteration( &mut self, - time: &mut Time, - folders: &FoldersRun, cfg: &Cfg, - bodies: &[Body], - routines: &dyn Routines, + bodies: &mut [AirlessBody], + precomputed: &[PreComputedBody], + time: &mut Time, scene: &Scene, win: &Window, + folders: &FoldersRun, + routines: &dyn Routines, ) { let dt = time.used_time_step(); let elapsed = time.elapsed_seconds(); @@ -44,14 +45,14 @@ impl Export { fs::create_dir_all(&folders.path).unwrap(); } - for (ii_body, cb) in cfg.bodies.iter().enumerate() { - if ii_body > 0 { + for body in 0..cfg.bodies.len() { + if body > 0 { // print!(" "); } // let np_elapsed = elapsed as Float / cb.spin.period; // print!("{:8.4}", np_elapsed); - routines.fn_export_iteration(cb, ii_body, time, folders, self.is_first_it); + routines.fn_export_iteration(body, cfg, time, folders, self.is_first_it); } // print!(">"); @@ -88,14 +89,22 @@ impl Export { win.export_frame(path); } - for (ii_body, (body, cb)) in izip!(bodies, &cfg.bodies).enumerate() { + for body in 0..cfg.bodies.len() { if self.is_first_it_export { - self.iteration_body_export_start_generic(cb, body, time, folders, scene); + self.iteration_body_export_start_generic( + cfg, + body, + bodies, + precomputed, + time, + scene, + folders, + ); } routines.fn_export_iteration_period( - cb, body, - ii_body, + bodies, + cfg, folders, self.exporting_started_elapsed, self.is_first_it_export, @@ -125,20 +134,24 @@ impl Export { pub fn iteration_body_export_start_generic( &self, - cb: &CfgBody, - body: &Body, + cfg: &Cfg, + body: usize, + _bodies: &mut [AirlessBody], + precomputed: &[PreComputedBody], time: &Time, - folders: &FoldersRun, scene: &Scene, + folders: &FoldersRun, ) { let elapsed = time.elapsed_seconds(); - let np_elapsed = time.elapsed_seconds() as Float / cb.spin.period; + let np_elapsed = time.elapsed_seconds() as Float / cfg.bodies[body].spin.period; let jd = time.jd(); - let folder_state = - folders.simu_rec_time_body_state(self.exporting_started_elapsed as _, &cb.id); - let folder_tpm = - folders.simu_rec_time_body_temperatures(self.exporting_started_elapsed as _, &cb.id); + let folder_state = folders + .simu_rec_time_body_state(self.exporting_started_elapsed as _, &cfg.bodies[body].id); + let folder_tpm = folders.simu_rec_time_body_temperatures( + self.exporting_started_elapsed as _, + &cfg.bodies[body].id, + ); let folder_img = folders.simu_rec_time_frames(self.exporting_started_elapsed as _); fs::create_dir_all(&folder_state).unwrap(); fs::create_dir_all(&folder_tpm).unwrap(); @@ -170,7 +183,7 @@ impl Export { CsvWriter::new(&mut file).finish(&mut df).unwrap(); let mut df = df!( - "spinrot" => body.mat_orient.as_slice(), + "spinrot" => precomputed[body].mat_orient.as_slice(), ) .unwrap(); let mut file = std::fs::File::options() diff --git a/src/simu/routines/core.rs b/src/simu/routines/core.rs index b6b8644..28941b1 100644 --- a/src/simu/routines/core.rs +++ b/src/simu/routines/core.rs @@ -1,150 +1,186 @@ use crate::{ - find_ref_orbit, find_reference_matrix_orientation, matrix_spin, position_in_inertial_frame, - simu::Scene, util::*, Asteroid, Body, CfgBody, CfgColormap, CfgState, CfgStateCartesian, - Equatorial, FoldersRun, Time, Window, + compute_cosine_emission_angle, compute_cosine_incidence_angle, compute_cosine_phase_angle, + find_matrix_spin, find_matrix_translation, find_reference_matrix_orientation, simu::Scene, + update_colormap_scalar, util::*, AirlessBody, Cfg, CfgBody, CfgScalar, FoldersRun, + PreComputedBody, Time, Window, }; +#[cfg(feature = "spice")] +use crate::CfgCamera; + +#[cfg(not(feature = "spice"))] +use crate::{find_ref_orbit, position_in_inertial_frame, CfgState}; + use downcast_rs::{impl_downcast, DowncastSync}; + +#[cfg(not(feature = "spice"))] use itertools::Itertools; pub trait RoutinesData { - fn new(asteroid: &Asteroid, _cb: &CfgBody, _scene: &Scene) -> Self; + fn new(asteroid: &AirlessBody, _cb: &CfgBody, _scene: &Scene) -> Self; } pub trait Routines: DowncastSync { - fn fn_setup_body(&mut self, asteroid: Asteroid, cb: &CfgBody, _scene: &Scene) -> Body { - fn_setup_body_default(asteroid, cb) + fn fn_update_scene(&self, cfg: &Cfg, time: &Time, scene: &mut Scene) { + let elapsed_from_start = time.elapsed_seconds_from_start(); + + #[cfg(not(feature = "spice"))] + if let Some(body) = cfg.bodies.first() { + let other_bodies = cfg.bodies.iter().skip(1).collect_vec(); + + match &body.state { + CfgState::Equatorial(coords) => { + let distance = cfg.scene.camera.as_earth().unwrap(); + let position = coords.xyz(distance); + dbg!(position); + scene.cam_pos = -position; + + let coords_sun = cfg.scene.sun.as_equatorial().unwrap(); + let position_sun = coords_sun.xyz(distance); + dbg!(position_sun); + scene.sun_pos = scene.cam_pos + position_sun; + } + CfgState::Orbit(orb) => { + let (mu_ref, factor) = find_ref_orbit(&orb, &other_bodies); + let pos = position_in_inertial_frame( + orb.a * factor, + orb.e, + orb.i * RPD, + orb.node * RPD, + orb.peri * RPD, + elapsed_from_start as Float, + orb.tp, + mu_ref, + ); + if mu_ref == MU_SUN { + scene.sun_pos = -pos; + } + } + _ => {} + }; + } + + #[cfg(feature = "spice")] + if let Some(body) = cfg.bodies.first() { + let context = match cfg.scene.camera { + CfgCamera::Earth(distance) => ("Earth", "J2000", distance), + _ => panic!("Expecting Earth with distance for camera settings"), + }; + let (position, _lt) = spice::spkpos( + context.0, + elapsed_from_start as f64, + context.1, + "none", + &body.id, + ); + let position = Vec3::from_row_slice(&position); + + let (position_sun, _lt_sun) = spice::spkpos( + "Sun", + elapsed_from_start as f64, + context.1, + "none", + &body.id, + ); + let position_sun = Vec3::from_row_slice(&position_sun); + + scene.cam_pos = position.normalize() * context.2; + scene.sun_pos = position_sun.normalize() * context.2; + } } fn fn_update_matrix_model( &self, - ii_body: usize, - ii_other_bodies: &[usize], - cb: &CfgBody, - other_cbs: &[&CfgBody], - bodies: &mut [Body], + cfg: &Cfg, + body: usize, + bodies: &mut [AirlessBody], + pre_computed_bodies: &mut [PreComputedBody], time: &Time, - scene: &mut Scene, + _scene: &Scene, ) { - fn_update_matrix_model_default( - ii_body, - ii_other_bodies, - cb, - other_cbs, - bodies, - time, - scene, - ); + let mat_orient_ref = find_reference_matrix_orientation(&cfg, body, pre_computed_bodies); + let mat_spin = find_matrix_spin(cfg, body, time); + let mat_translation = find_matrix_translation(cfg, body, time); + + bodies[body].matrix_model = + mat_orient_ref * mat_translation * pre_computed_bodies[body].mat_orient * mat_spin; } fn fn_iteration_body( &mut self, - ii_body: usize, - ii_other_bodies: &[usize], - _cb: &CfgBody, - _other_cbs: &[&CfgBody], - bodies: &mut [Body], - scene: &Scene, - time: &Time, - ); + _body: usize, + _bodies: &mut [AirlessBody], + _pre_computed_bodies: &mut [PreComputedBody], + _time: &Time, + _scene: &Scene, + ) { + } fn fn_update_colormap( &self, - win: &Window, - cmap: &CfgColormap, - ii_body: usize, - body: &mut Body, + body: usize, + bodies: &mut [AirlessBody], + pre_computed_bodies: &[PreComputedBody], + cfg: &Cfg, scene: &Scene, - ); + win: &Window, + ) { + let scalars = match &cfg.window.colormap.scalar { + Some(CfgScalar::AngleIncidence) => compute_cosine_incidence_angle( + &bodies[body], + &pre_computed_bodies[body].normals, + scene, + ) + .map(|a| a.acos() * DPR), + Some(CfgScalar::AngleEmission) => compute_cosine_emission_angle( + &bodies[body], + &pre_computed_bodies[body].normals, + scene, + ) + .map(|a| a.acos() * DPR), + Some(CfgScalar::AnglePhase) => { + compute_cosine_phase_angle(&bodies[body], scene).map(|a| a.acos() * DPR) + } + None => return, + _ => unreachable!(), + }; + + update_colormap_scalar(win, cfg, scalars.as_slice(), &mut bodies[body], body); + } fn fn_export_iteration( &self, - cb: &CfgBody, - ii_body: usize, - time: &Time, - folders: &FoldersRun, - is_first_it: bool, - ); + _body: usize, + _cfg: &Cfg, + _time: &Time, + _folders: &FoldersRun, + _is_first_it: bool, + ) { + } fn fn_export_iteration_period( &self, - cb: &CfgBody, - body: &Body, - ii_body: usize, - folders: &FoldersRun, - exporting_started_elapsed: i64, - is_first_it_export: bool, - ); + _body: usize, + _bodies: &mut [AirlessBody], + _cfg: &Cfg, + _folders: &FoldersRun, + _exporting_started_elapsed: i64, + _is_first_it_export: bool, + ) { + } fn fn_end_of_iteration( &mut self, - bodies: &mut [Body], - time: &Time, - scene: &Scene, + cfg: &Cfg, + _bodies: &mut [AirlessBody], + _time: &Time, + _scene: &Scene, win: &Window, - ); + ) { + if cfg.simulation.pause_after_first_iteration { + win.toggle_pause(); + } + } } impl_downcast!(sync Routines); - -pub fn fn_setup_body_default(asteroid: Asteroid, cb: &CfgBody) -> Body { - Body::new(asteroid, cb) -} - -pub fn fn_update_matrix_model_default( - ii_body: usize, - ii_other_bodies: &[usize], - cb: &CfgBody, - other_cbs: &[&CfgBody], - bodies: &mut [Body], - time: &Time, - scene: &mut Scene, -) { - let elapsed = time.elapsed_seconds(); - let elapsed_from_start = time.elapsed_seconds_from_start(); - - let other_bodies = ii_other_bodies.iter().map(|&ii| &bodies[ii]).collect_vec(); - - let mat_orient_ref = find_reference_matrix_orientation(cb, &other_bodies); - - let mat_spin = { - if cb.spin.period == 0.0 { - Mat4::identity() - } else { - let np_elapsed = elapsed as Float / cb.spin.period; - let spin = (TAU * np_elapsed + cb.spin.spin0 * RPD) % TAU; - matrix_spin(spin) - } - }; - - let mat_translation = match &cb.state { - CfgState::Cartesian(CfgStateCartesian { - position, - orientation: _orientation, - }) => Mat4::new_translation(position), - CfgState::Equatorial(Equatorial { ra, dec }) => Mat4::new_translation(&Vec3::zeros()), - CfgState::File(_p) => Mat4::identity(), - CfgState::Orbit(orb) => { - let (mu_ref, factor) = find_ref_orbit(&orb, &other_cbs); - let pos = position_in_inertial_frame( - orb.a * factor, - orb.e, - orb.i * RPD, - orb.node * RPD, - orb.peri * RPD, - elapsed_from_start as Float, - orb.tp, - mu_ref, - ); - if mu_ref == MU_SUN { - scene.sun_pos = -pos; - Mat4::identity() - } else { - Mat4::new_translation(&(pos * 1e-3)) - } - } - }; - - bodies[ii_body].asteroid.matrix_model = - mat_orient_ref * mat_translation * bodies[ii_body].mat_orient * mat_spin; -} diff --git a/src/simu/routines/thermal.rs b/src/simu/routines/thermal.rs index c7348e6..03bfa87 100644 --- a/src/simu/routines/thermal.rs +++ b/src/simu/routines/thermal.rs @@ -1,9 +1,9 @@ use crate::{ compute_cosine_emission_angle, compute_cosine_incidence_angle, compute_cosine_phase_angle, - effective_temperature, flux_solar_radiation, fn_setup_body_default, newton_method_temperature, - read_surface_low, shadows, simu::Scene, update_colormap_scalar, util::*, Asteroid, Body, - CfgBody, CfgColormap, CfgScalar, CfgTemperatureInit, FoldersRun, Routines, RoutinesData, - Surface, Time, Window, + effective_temperature, flux_solar_radiation, newton_method_temperature, read_surface_low, + shadows, simu::Scene, update_colormap_scalar, util::*, AirlessBody, Cfg, CfgBody, CfgScalar, + CfgTemperatureInit, FoldersRun, PreComputedBody, Routines, RoutinesData, Surface, + Time, Window, }; use itertools::Itertools; @@ -24,7 +24,7 @@ pub struct ThermalData { } impl RoutinesData for ThermalData { - fn new(asteroid: &Asteroid, cb: &CfgBody, scene: &Scene) -> Self { + fn new(asteroid: &AirlessBody, cb: &CfgBody, scene: &Scene) -> Self { let depth_grid = &asteroid.interior.as_ref().unwrap().as_grid().depth; let depth_size = depth_grid.len(); let surf_size = asteroid.surface.faces.len(); @@ -124,28 +124,29 @@ impl RoutinesData for ThermalData { pub trait RoutinesThermal: Routines { fn fn_compute_solar_flux( &self, - body: &Body, + body: &AirlessBody, + precomputed: &PreComputedBody, body_info: &ThermalData, scene: &Scene, ) -> DRVector; fn fn_compute_surface_temperatures( &self, - body: &Body, + body: &AirlessBody, body_info: &ThermalData, surface_fluxes: &DRVector, ) -> DRVector; fn fn_compute_heat_conduction( &self, - body: &Body, + body: &AirlessBody, body_info: &ThermalData, delta_time: Float, ) -> DMatrix; fn fn_compute_bottom_depth_temperatures( &self, - body: &Body, + body: &AirlessBody, body_info: &ThermalData, ) -> DRVector; } @@ -155,38 +156,41 @@ pub struct RoutinesThermalDefault { pub shadows_mutual: bool, } -impl Routines for RoutinesThermalDefault { - fn fn_setup_body(&mut self, asteroid: Asteroid, cb: &CfgBody, scene: &Scene) -> Body { - self.data.push(ThermalData::new(&asteroid, cb, scene)); - fn_setup_body_default(asteroid, cb) +impl RoutinesThermalDefault { + pub fn new() -> Self { + Self { + data: vec![], + shadows_mutual: false, + } } +} +impl Routines for RoutinesThermalDefault { fn fn_iteration_body( &mut self, - ii_body: usize, - ii_other_bodies: &[usize], - _cb: &CfgBody, - _other_cbs: &[&CfgBody], - bodies: &mut [Body], - scene: &Scene, + body: usize, + bodies: &mut [AirlessBody], + pre_computed_bodies: &mut [PreComputedBody], time: &Time, + scene: &Scene, ) { let dt = time.used_time_step(); - let other_bodies = ii_other_bodies.iter().map(|&ii| &bodies[ii]).collect_vec(); - let mut fluxes_solar = - self.fn_compute_solar_flux(&bodies[ii_body], &self.data[ii_body], &scene); + let other_bodies = (0..bodies.len()).filter(|ii| *ii != body).collect_vec(); + + let mut fluxes_solar = self.fn_compute_solar_flux( + &bodies[body], + &pre_computed_bodies[body], + &self.data[body], + &scene, + ); if self.shadows_mutual { let shadows_self: Vec = vec![]; let mut shadows_mutual: Vec = vec![]; for other_body in other_bodies { - shadows_mutual = shadows( - &scene.sun_pos, - &bodies[ii_body].asteroid, - &other_body.asteroid, - ); + shadows_mutual = shadows(&scene.sun_pos, &bodies[body], &bodies[other_body]); } for &index in shadows_mutual.iter().chain(&shadows_self).unique() { @@ -197,355 +201,283 @@ impl Routines for RoutinesThermalDefault { let fluxes = fluxes_solar.clone(); let temperatures_surface = - self.fn_compute_surface_temperatures(&bodies[ii_body], &self.data[ii_body], &fluxes); - self.data[ii_body].tmp.set_row(0, &temperatures_surface); + self.fn_compute_surface_temperatures(&bodies[body], &self.data[body], &fluxes); + self.data[body].tmp.set_row(0, &temperatures_surface); let temperatures_inside = - self.fn_compute_heat_conduction(&bodies[ii_body], &self.data[ii_body], dt as Float); - let curr_size = self.data[ii_body].depth_size - 2; + self.fn_compute_heat_conduction(&bodies[body], &self.data[body], dt as Float); + let curr_size = self.data[body].depth_size - 2; for index in 0..curr_size { - self.data[ii_body] + self.data[body] .tmp .set_row(index + 1, &temperatures_inside.row(index)); } let temperatures_bottom = - self.fn_compute_bottom_depth_temperatures(&bodies[ii_body], &self.data[ii_body]); - self.data[ii_body] + self.fn_compute_bottom_depth_temperatures(&bodies[body], &self.data[body]); + self.data[body] .tmp .set_row(curr_size + 1, &temperatures_bottom); - self.data[ii_body].fluxes = fluxes; - self.data[ii_body].fluxes_solar = fluxes_solar; + self.data[body].fluxes = fluxes; + self.data[body].fluxes_solar = fluxes_solar; } fn fn_update_colormap( &self, - win: &Window, - cmap: &CfgColormap, - ii_body: usize, - body: &mut Body, + body: usize, + bodies: &mut [AirlessBody], + pre_computed_bodies: &[PreComputedBody], + cfg: &Cfg, scene: &Scene, + win: &Window, ) { - fn_update_colormap_thermal_default(&self.data[ii_body], win, cmap, ii_body, body, scene); + let cmap = &cfg.window.colormap; + let scalars = match &cmap.scalar { + Some(CfgScalar::AngleIncidence) => compute_cosine_incidence_angle( + &bodies[body], + &pre_computed_bodies[body].normals, + scene, + ) + .map(|a| a.acos() * DPR), + Some(CfgScalar::AngleEmission) => compute_cosine_emission_angle( + &bodies[body], + &pre_computed_bodies[body].normals, + scene, + ) + .map(|a| a.acos() * DPR), + Some(CfgScalar::AnglePhase) => { + compute_cosine_phase_angle(&bodies[body], scene).map(|a| a.acos() * DPR) + } + Some(CfgScalar::FluxSolar) => self.data[body].fluxes_solar.clone(), + Some(CfgScalar::FluxSurface) => self.data[body].fluxes.clone(), + Some(CfgScalar::FluxEmitted) => unimplemented!(), + Some(CfgScalar::FluxSelf) => unimplemented!(), + Some(CfgScalar::FluxMutual) => unimplemented!(), + None | Some(CfgScalar::Temperature) => self.data[body].tmp.row(0).into_owned(), + }; + + update_colormap_scalar(win, cfg, scalars.as_slice(), &mut bodies[body], body); } fn fn_export_iteration( &self, - cb: &CfgBody, - ii_body: usize, + body: usize, + cfg: &Cfg, time: &Time, folders: &FoldersRun, is_first_it: bool, ) { - fn_export_iteration_thermal_default(&self.data[ii_body], cb, time, folders, is_first_it); + let np_elapsed = time.elapsed_seconds() as Float / cfg.bodies[body].spin.period; + + let data = &self.data[body]; + let tmp_surf_min = data.tmp.row(0).min(); + let tmp_surf_max = data.tmp.row(0).max(); + let tmp_surf_mean = data.tmp.row(0).mean(); + let tmp_surf_stdev = data.tmp.row(0).variance().sqrt(); + let tmp_bottom_min = data.tmp.row(data.depth_size - 1).min(); + let tmp_bottom_max = data.tmp.row(data.depth_size - 1).max(); + let tmp_bottom_mean = data.tmp.row(data.depth_size - 1).mean(); + let tmp_bottom_stdev = data.tmp.row(data.depth_size - 1).variance().sqrt(); + + let mut df = df!( + "time" => [np_elapsed], + "tmp-surf-min" => [tmp_surf_min], + "tmp-surf-max" => [tmp_surf_max], + "tmp-surf-mean" => [tmp_surf_mean], + "tmp-surf-stdev" => [tmp_surf_stdev], + "tmp-bottom-min" => [tmp_bottom_min], + "tmp-bottom-max" => [tmp_bottom_max], + "tmp-bottom-mean" => [tmp_bottom_mean], + "tmp-bottom-stdev" => [tmp_bottom_stdev], + ) + .unwrap(); + + let folder_simu = folders.simu_body(&cfg.bodies[body].id); + fs::create_dir_all(&folder_simu).unwrap(); + + let mut file = std::fs::File::options() + .append(true) + .create(true) + .open(folder_simu.join("progress.csv")) + .unwrap(); + CsvWriter::new(&mut file) + .include_header(is_first_it) + .finish(&mut df) + .unwrap(); } fn fn_export_iteration_period( &self, - cb: &CfgBody, - body: &Body, - ii_body: usize, + body: usize, + _bodies: &mut [AirlessBody], + cfg: &Cfg, folders: &FoldersRun, exporting_started_elapsed: i64, is_first_it_export: bool, ) { - fn_export_iteration_period_thermal_default( - &self.data[ii_body], - cb, - body, - folders, - exporting_started_elapsed, - is_first_it_export, - ); - } + let folder_tpm = folders + .simu_rec_time_body_temperatures(exporting_started_elapsed as _, &cfg.bodies[body].id); + fs::create_dir_all(&folder_tpm).unwrap(); - fn fn_end_of_iteration( - &mut self, - _bodies: &mut [Body], - _time: &Time, - _scene: &Scene, - _win: &Window, - ) { + let data = &self.data[body]; + + if is_first_it_export { + let mut df = df!( + "tmp" => data.tmp.map(|t| t as f32).as_slice(), + ) + .unwrap(); + let mut file = std::fs::File::options() + .append(true) + .create(true) + .open(folder_tpm.join("temperatures-all.csv")) + .unwrap(); + CsvWriter::new(&mut file).finish(&mut df).unwrap(); + } + + if !cfg.bodies[body].record.faces.is_empty() { + let dfcols = cfg.bodies[body] + .record + .faces + .iter() + .map(|&face| Series::new(&format!("{}", face), &vec![data.tmp.row(0)[face]])) + .collect_vec(); + let mut df = DataFrame::new(dfcols).unwrap(); + let p = folder_tpm.join("temperatures-faces.csv"); + let mut file = std::fs::File::options() + .append(true) + .create(true) + .open(&p) + .unwrap(); + CsvWriter::new(&mut file) + .include_header(is_first_it_export) + .finish(&mut df) + .unwrap(); + } + + if !cfg.bodies[body].record.cells.is_empty() { + let dfcols = cfg.bodies[body] + .record + .cells + .iter() + .map(|&cell| Series::new(&format!("{}", cell), &vec![data.tmp[cell]])) + .collect_vec(); + let mut df = DataFrame::new(dfcols).unwrap(); + let p = folder_tpm.join("temperatures-cells.csv"); + let mut file = std::fs::File::options() + .append(true) + .create(true) + .open(&p) + .unwrap(); + CsvWriter::new(&mut file) + .include_header(is_first_it_export) + .finish(&mut df) + .unwrap(); + } + + if !cfg.bodies[body].record.columns.is_empty() { + let dfcols = cfg.bodies[body] + .record + .columns + .iter() + .map(|&column| { + Series::new(&format!("{}", column), data.tmp.column(column).as_slice()) + }) + .collect_vec(); + let mut df = DataFrame::new(dfcols).unwrap(); + let p = folder_tpm.join("temperatures-columns.csv"); + let mut file = std::fs::File::options() + .append(true) + .create(true) + .open(&p) + .unwrap(); + CsvWriter::new(&mut file) + .include_header(is_first_it_export) + .finish(&mut df) + .unwrap(); + } + + if !cfg.bodies[body].record.rows.is_empty() { + let dfcols = cfg.bodies[body] + .record + .rows + .iter() + .map(|&row| { + Series::new( + &format!("{}", row), + data.tmp.row(row).transpose().as_slice(), + ) + }) + .collect_vec(); + let mut df = DataFrame::new(dfcols).unwrap(); + let p = folder_tpm.join("temperatures-rows.csv"); + let mut file = std::fs::File::options() + .append(true) + .create(true) + .open(&p) + .unwrap(); + CsvWriter::new(&mut file) + .include_header(is_first_it_export) + .finish(&mut df) + .unwrap(); + } } } impl RoutinesThermal for RoutinesThermalDefault { fn fn_compute_solar_flux( &self, - body: &Body, + body: &AirlessBody, + precomputed: &PreComputedBody, body_info: &ThermalData, scene: &Scene, ) -> DRVector { - fn_compute_solar_flux_default(body, body_info, scene) + let cos_incidences = compute_cosine_incidence_angle(body, &precomputed.normals, scene); + flux_solar_radiation(&cos_incidences, &body_info.albedos, scene.sun_dist()) } fn fn_compute_surface_temperatures( &self, - body: &Body, + body: &AirlessBody, body_info: &ThermalData, surface_fluxes: &DRVector, ) -> DRVector { - fn_compute_surface_temperatures_default(body, body_info, surface_fluxes) + let depth_grid = &body.interior.as_ref().unwrap().as_grid().depth; + + newton_method_temperature( + body_info.tmp.row(0).as_view(), + &surface_fluxes, + &body_info.emissivities, + &body_info.conductivities, + body_info.tmp.rows(1, 2).as_view(), + depth_grid[2], + ) } fn fn_compute_heat_conduction( &self, - body: &Body, + body: &AirlessBody, body_info: &ThermalData, delta_time: Float, ) -> DMatrix { - fn_compute_heat_conduction_default(body, body_info, delta_time) + let curr_size = body_info.depth_size - 2; + let surf_size = body.surface.faces.len(); + + let prev = body_info.tmp.view((0, 0), (curr_size, surf_size)); + let curr = body_info.tmp.view((1, 0), (curr_size, surf_size)); + let next = body_info.tmp.view((2, 0), (curr_size, surf_size)); + + curr + delta_time + * body_info + .diffu_dx2 + .component_mul(&(prev - 2. * curr + next)) } fn fn_compute_bottom_depth_temperatures( &self, - body: &Body, + _body: &AirlessBody, body_info: &ThermalData, ) -> DRVector { - fn_compute_bottom_depth_temperatures_default(body, body_info) - } -} - -pub fn routines_thermal_default() -> RoutinesThermalDefault { - RoutinesThermalDefault { - data: vec![], - shadows_mutual: false, - } -} - -pub fn fn_compute_solar_flux_default( - body: &Body, - body_info: &ThermalData, - scene: &Scene, -) -> DRVector { - let cos_incidences = compute_cosine_incidence_angle(body, &body.normals, scene); - flux_solar_radiation(&cos_incidences, &body_info.albedos, scene.sun_dist()) -} - -pub fn fn_compute_surface_temperatures_default( - body: &Body, - body_info: &ThermalData, - surface_fluxes: &DRVector, -) -> DRVector { - let depth_grid = &body.asteroid.interior.as_ref().unwrap().as_grid().depth; - - newton_method_temperature( - body_info.tmp.row(0).as_view(), - &surface_fluxes, - &body_info.emissivities, - &body_info.conductivities, - body_info.tmp.rows(1, 2).as_view(), - depth_grid[2], - ) -} - -pub fn fn_compute_heat_conduction_default( - body: &Body, - body_info: &ThermalData, - delta_time: Float, -) -> DMatrix { - let curr_size = body_info.depth_size - 2; - let surf_size = body.asteroid.surface.faces.len(); - - let prev = body_info.tmp.view((0, 0), (curr_size, surf_size)); - let curr = body_info.tmp.view((1, 0), (curr_size, surf_size)); - let next = body_info.tmp.view((2, 0), (curr_size, surf_size)); - - curr + delta_time - * body_info - .diffu_dx2 - .component_mul(&(prev - 2. * curr + next)) -} - -pub fn fn_compute_bottom_depth_temperatures_default( - _body: &Body, - body_info: &ThermalData, -) -> DRVector { - body_info.tmp.row(body_info.depth_size - 2).clone_owned() -} - -pub fn fn_update_colormap_thermal_default( - data: &ThermalData, - win: &Window, - cmap: &CfgColormap, - ii_body: usize, - body: &mut Body, - scene: &Scene, -) { - let scalars = match &cmap.scalar { - Some(CfgScalar::AngleIncidence) => { - compute_cosine_incidence_angle(body, &body.normals, scene).map(|a| a.acos() * DPR) - } - Some(CfgScalar::AngleEmission) => { - compute_cosine_emission_angle(body, &body.normals, scene).map(|a| a.acos() * DPR) - } - Some(CfgScalar::AnglePhase) => { - compute_cosine_phase_angle(body, scene).map(|a| a.acos() * DPR) - } - Some(CfgScalar::FluxSolar) => data.fluxes_solar.clone(), - Some(CfgScalar::FluxSurface) => data.fluxes.clone(), - Some(CfgScalar::FluxEmitted) => unimplemented!(), - Some(CfgScalar::FluxSelf) => unimplemented!(), - Some(CfgScalar::FluxMutual) => unimplemented!(), - None | Some(CfgScalar::Temperature) => data.tmp.row(0).into_owned(), - }; - - update_colormap_scalar(win, cmap, scalars.as_slice(), &mut body.asteroid, ii_body); -} - -pub fn fn_export_iteration_thermal_default( - data: &ThermalData, - cb: &CfgBody, - time: &Time, - folders: &FoldersRun, - is_first_it: bool, -) { - let np_elapsed = time.elapsed_seconds() as Float / cb.spin.period; - - let tmp_surf_min = data.tmp.row(0).min(); - let tmp_surf_max = data.tmp.row(0).max(); - let tmp_surf_mean = data.tmp.row(0).mean(); - let tmp_surf_stdev = data.tmp.row(0).variance().sqrt(); - let tmp_bottom_min = data.tmp.row(data.depth_size - 1).min(); - let tmp_bottom_max = data.tmp.row(data.depth_size - 1).max(); - let tmp_bottom_mean = data.tmp.row(data.depth_size - 1).mean(); - let tmp_bottom_stdev = data.tmp.row(data.depth_size - 1).variance().sqrt(); - - let mut df = df!( - "time" => [np_elapsed], - "tmp-surf-min" => [tmp_surf_min], - "tmp-surf-max" => [tmp_surf_max], - "tmp-surf-mean" => [tmp_surf_mean], - "tmp-surf-stdev" => [tmp_surf_stdev], - "tmp-bottom-min" => [tmp_bottom_min], - "tmp-bottom-max" => [tmp_bottom_max], - "tmp-bottom-mean" => [tmp_bottom_mean], - "tmp-bottom-stdev" => [tmp_bottom_stdev], - ) - .unwrap(); - - let folder_simu = folders.simu_body(&cb.id); - fs::create_dir_all(&folder_simu).unwrap(); - - let mut file = std::fs::File::options() - .append(true) - .create(true) - .open(folder_simu.join("progress.csv")) - .unwrap(); - CsvWriter::new(&mut file) - .include_header(is_first_it) - .finish(&mut df) - .unwrap(); -} - -pub fn fn_export_iteration_period_thermal_default( - data: &ThermalData, - cb: &CfgBody, - _body: &Body, - folders: &FoldersRun, - exporting_started_elapsed: i64, - is_first_it_export: bool, -) { - let folder_tpm = - folders.simu_rec_time_body_temperatures(exporting_started_elapsed as _, &cb.id); - fs::create_dir_all(&folder_tpm).unwrap(); - - if is_first_it_export { - let mut df = df!( - "tmp" => data.tmp.map(|t| t as f32).as_slice(), - ) - .unwrap(); - let mut file = std::fs::File::options() - .append(true) - .create(true) - .open(folder_tpm.join("temperatures-all.csv")) - .unwrap(); - CsvWriter::new(&mut file).finish(&mut df).unwrap(); - } - - if !cb.record.faces.is_empty() { - let dfcols = cb - .record - .faces - .iter() - .map(|&face| Series::new(&format!("{}", face), &vec![data.tmp.row(0)[face]])) - .collect_vec(); - let mut df = DataFrame::new(dfcols).unwrap(); - let p = folder_tpm.join("temperatures-faces.csv"); - let mut file = std::fs::File::options() - .append(true) - .create(true) - .open(&p) - .unwrap(); - CsvWriter::new(&mut file) - .include_header(is_first_it_export) - .finish(&mut df) - .unwrap(); - } - - if !cb.record.cells.is_empty() { - let dfcols = cb - .record - .cells - .iter() - .map(|&cell| Series::new(&format!("{}", cell), &vec![data.tmp[cell]])) - .collect_vec(); - let mut df = DataFrame::new(dfcols).unwrap(); - let p = folder_tpm.join("temperatures-cells.csv"); - let mut file = std::fs::File::options() - .append(true) - .create(true) - .open(&p) - .unwrap(); - CsvWriter::new(&mut file) - .include_header(is_first_it_export) - .finish(&mut df) - .unwrap(); - } - - if !cb.record.columns.is_empty() { - let dfcols = cb - .record - .columns - .iter() - .map(|&column| Series::new(&format!("{}", column), data.tmp.column(column).as_slice())) - .collect_vec(); - let mut df = DataFrame::new(dfcols).unwrap(); - let p = folder_tpm.join("temperatures-columns.csv"); - let mut file = std::fs::File::options() - .append(true) - .create(true) - .open(&p) - .unwrap(); - CsvWriter::new(&mut file) - .include_header(is_first_it_export) - .finish(&mut df) - .unwrap(); - } - - if !cb.record.rows.is_empty() { - let dfcols = cb - .record - .rows - .iter() - .map(|&row| { - Series::new( - &format!("{}", row), - data.tmp.row(row).transpose().as_slice(), - ) - }) - .collect_vec(); - let mut df = DataFrame::new(dfcols).unwrap(); - let p = folder_tpm.join("temperatures-rows.csv"); - let mut file = std::fs::File::options() - .append(true) - .create(true) - .open(&p) - .unwrap(); - CsvWriter::new(&mut file) - .include_header(is_first_it_export) - .finish(&mut df) - .unwrap(); + body_info.tmp.row(body_info.depth_size - 2).clone_owned() } } diff --git a/src/simu/routines/viewer.rs b/src/simu/routines/viewer.rs index 78b8287..829a37a 100644 --- a/src/simu/routines/viewer.rs +++ b/src/simu/routines/viewer.rs @@ -1,15 +1,9 @@ -use crate::{ - compute_cosine_emission_angle, compute_cosine_incidence_angle, compute_cosine_phase_angle, - fn_setup_body_default, simu::Scene, update_colormap_scalar, util::*, Asteroid, Body, CfgBody, - CfgColormap, CfgScalar, FoldersRun, Routines, RoutinesData, Time, Window, -}; - -use itertools::Itertools; +use crate::{simu::Scene, AirlessBody, CfgBody, Routines, RoutinesData}; pub struct ViewerData {} impl RoutinesData for ViewerData { - fn new(_asteroid: &Asteroid, _cb: &CfgBody, _scene: &Scene) -> Self { + fn new(_asteroid: &AirlessBody, _cb: &CfgBody, _scene: &Scene) -> Self { Self {} } } @@ -22,124 +16,12 @@ pub struct RoutinesViewerDefault { pub data: Vec, } -impl Routines for RoutinesViewerDefault { - fn fn_setup_body(&mut self, asteroid: Asteroid, cb: &CfgBody, scene: &Scene) -> Body { - self.data.push(ViewerData::new(&asteroid, cb, scene)); - fn_setup_body_default(asteroid, cb) - } - - fn fn_iteration_body( - &mut self, - _ii_body: usize, - ii_other_bodies: &[usize], - _cb: &CfgBody, - _other_cbs: &[&CfgBody], - bodies: &mut [Body], - _scene: &Scene, - time: &Time, - ) { - let _dt = time.used_time_step(); - let _other_bodies = ii_other_bodies.iter().map(|&ii| &bodies[ii]).collect_vec(); - } - - fn fn_update_colormap( - &self, - win: &Window, - cmap: &CfgColormap, - ii_body: usize, - body: &mut Body, - scene: &Scene, - ) { - fn_update_colormap_viewer_default(&self.data[ii_body], win, cmap, ii_body, body, scene); - } - - fn fn_export_iteration( - &self, - cb: &CfgBody, - ii_body: usize, - time: &Time, - folders: &FoldersRun, - is_first_it: bool, - ) { - fn_export_iteration_viewer_default(&self.data[ii_body], cb, time, folders, is_first_it); - } - - fn fn_export_iteration_period( - &self, - cb: &CfgBody, - body: &Body, - ii_body: usize, - folders: &FoldersRun, - exporting_started_elapsed: i64, - is_first_it_export: bool, - ) { - fn_export_iteration_period_viewer_default( - &self.data[ii_body], - cb, - body, - folders, - exporting_started_elapsed, - is_first_it_export, - ); - } - - fn fn_end_of_iteration( - &mut self, - _bodies: &mut [Body], - _time: &Time, - _scene: &Scene, - _win: &Window, - ) { +impl RoutinesViewerDefault { + pub fn new() -> Self { + RoutinesViewerDefault { data: vec![] } } } -impl RoutinesViewer for RoutinesViewerDefault {} - -pub fn routines_viewer_default() -> RoutinesViewerDefault { - RoutinesViewerDefault { data: vec![] } -} - -pub fn fn_update_colormap_viewer_default( - _data: &D, - win: &Window, - cmap: &CfgColormap, - ii_body: usize, - body: &mut Body, - scene: &Scene, -) { - let scalars = match &cmap.scalar { - Some(CfgScalar::AngleIncidence) => { - compute_cosine_incidence_angle(body, &body.normals, scene).map(|a| a.acos() * DPR) - } - Some(CfgScalar::AngleEmission) => { - compute_cosine_emission_angle(body, &body.normals, scene).map(|a| a.acos() * DPR) - } - Some(CfgScalar::AnglePhase) => { - compute_cosine_phase_angle(body, scene).map(|a| a.acos() * DPR) - } - None => return, - _ => unreachable!(), - }; - - update_colormap_scalar(win, cmap, scalars.as_slice(), &mut body.asteroid, ii_body); -} - -pub fn fn_export_iteration_viewer_default( - _data: &D, - cb: &CfgBody, - time: &Time, - _folders: &FoldersRun, - _is_first_it: bool, -) { - let _np_elapsed = time.elapsed_seconds() as Float / cb.spin.period; -} +impl Routines for RoutinesViewerDefault {} -pub fn fn_export_iteration_period_viewer_default( - _data: &D, - _cb: &CfgBody, - _body: &Body, - _folders: &FoldersRun, - _exporting_started_elapsed: i64, - _is_first_it_export: bool, -) { -} +impl RoutinesViewer for RoutinesViewerDefault {} diff --git a/src/simu/scenario.rs b/src/simu/scenario.rs index 62b30d9..4feec8c 100644 --- a/src/simu/scenario.rs +++ b/src/simu/scenario.rs @@ -1,22 +1,23 @@ use crate::{ - check_if_latest_version, read_surface_main, routines_thermal_default, routines_viewer_default, - simu::Scene, util::*, Asteroid, Body, Cfg, CfgCamera, CfgInterior, CfgInteriorGrid1D, - CfgRoutines, CfgStateCartesian, CfgSun, Export, FoldersRun, FrameEvent, Result, Routines, Time, - Window, Equatorial, + check_if_latest_version, read_surface_main, simu::Scene, util::*, AirlessBody, Cfg, CfgCamera, + CfgInterior, CfgInteriorGrid1D, CfgRoutines, CfgStateCartesian, CfgSun, Equatorial, Export, + FoldersRun, FrameEvent, PreComputedBody, Result, Routines, RoutinesThermalDefault, + RoutinesViewerDefault, Time, Window, }; use chrono::Utc; -use itertools::{izip, Itertools}; +use itertools::Itertools; use std::{env, path::Path}; pub struct Scenario { pub cfg: Cfg, - pub folders: FoldersRun, - pub bodies: Vec, - pub win: Window, + pub bodies: Vec, pub time: Time, pub scene: Scene, + pub win: Window, + pub folders: FoldersRun, pub routines: Box, + pub pre_computed_bodies: Vec, } impl Scenario { @@ -53,16 +54,22 @@ impl Scenario { folders.save_cfgs(&path_cfg); folders.save_src(&path_mainrs); - let bodies: Vec = vec![]; + let bodies = vec![]; + let pre_computed_bodies = vec![]; let routines = match &cfg.simulation.routines { - CfgRoutines::Viewer => Box::new(routines_viewer_default()) as Box, - CfgRoutines::Thermal => Box::new(routines_thermal_default()) as Box, + CfgRoutines::Viewer => Box::new(RoutinesViewerDefault::new()) as Box, + CfgRoutines::Thermal => Box::new(RoutinesThermalDefault::new()) as Box, }; + #[cfg(feature = "spice")] + if let Some(path) = &cfg.scene.spice { + spice::furnsh(path.to_str().unwrap()); + } + let sun_pos = match &cfg.scene.sun { CfgSun::Cartesian(CfgStateCartesian { position, .. }) => position * AU, - CfgSun::Equatorial(Equatorial { ra, dec }) => Vec3::x() * AU, + CfgSun::Equatorial(Equatorial { ra: _ra, dec: _dec }) => Vec3::x() * AU, }; // If camera is from Earth, we cannot give a correct value for position now. @@ -98,21 +105,20 @@ impl Scenario { let time = Time::new() .with_time_step(cfg.simulation.step) - .with_time_start(cfg.simulation.start as _); + .with_time_start(cfg.simulation.start.seconds().unwrap() as _); Ok(Self { cfg, - folders, bodies, - win, time, scene, + win, + folders, routines, + pre_computed_bodies, }) } -} -impl Scenario { pub fn change_routines(&mut self, routines: R) { self.routines = Box::new(routines); } @@ -120,7 +126,7 @@ impl Scenario { pub fn load_bodies(&mut self) -> Result<()> { for (_ii, cb) in self.cfg.bodies.iter().enumerate() { let surface = read_surface_main(cb)?; - let asteroid = Asteroid::new(surface); + let asteroid = AirlessBody::new(surface); let asteroid = match &cb.interior { None => asteroid, @@ -142,13 +148,13 @@ impl Scenario { }, }; - // self.data.push(ThermalData::new(&asteroid, &cb, &self.scene)); - self.bodies - .push(self.routines.fn_setup_body(asteroid, cb, &self.scene)); + self.pre_computed_bodies + .push(PreComputedBody::new(&asteroid, &cb)); + self.bodies.push(asteroid); } self.win - .load_surfaces(self.bodies.iter().map(|b| &b.asteroid.surface)); + .load_surfaces(self.bodies.iter().map(|b| &b.surface)); Ok(()) } @@ -192,8 +198,10 @@ impl Scenario { let mut export = Export::new(&self.cfg.simulation.export); 'main_loop: loop { + // Register keyboard and mouse interactions. let event = self.win.events(); + // Update camera data from keyboard movements. self.scene.cam_pos = self.win.camera_position(); match event { @@ -202,8 +210,7 @@ impl Scenario { }; if self.win.is_paused() { - self.win - .render_asteroids(&self.bodies.iter().map(|b| &b.asteroid).collect_vec()); + self.win.render_asteroids(&self.bodies); self.win.swap_window(); continue; } @@ -216,70 +223,62 @@ impl Scenario { let elapsed = self.time.elapsed_seconds(); let jd = self.time.jd(); - for (indices_bodies, cbs_permut) in izip!( - self.bodies_permutations_indices(), - self.cfg - .bodies - .clone() - .iter() - .permutations(self.bodies.len()) - ) { - if indices_bodies.is_empty() { - break; - } - - let (&ii_body, ii_other_bodies) = indices_bodies.split_first().unwrap(); - let (cb, other_cbs) = cbs_permut.split_first().unwrap(); + self.routines + .fn_update_scene(&self.cfg, &self.time, &mut self.scene); + for body in 0..self.bodies.len() { self.routines.fn_update_matrix_model( - ii_body, - &ii_other_bodies, - cb, - &other_cbs, + &self.cfg, + body, &mut self.bodies, + &mut self.pre_computed_bodies, &self.time, &mut self.scene, ); - self.routines.fn_iteration_body( - ii_body, - ii_other_bodies, - cb, - other_cbs, + body, &mut self.bodies, - &self.scene, + &mut self.pre_computed_bodies, &self.time, + &mut self.scene, ); self.routines.fn_update_colormap( - &self.win, - &self.cfg.window.colormap, - ii_body, - &mut self.bodies[ii_body], + body, + &mut self.bodies, + &self.pre_computed_bodies, + &self.cfg, &self.scene, + &self.win, ); } + self.win.set_camera_position(&self.scene.cam_pos); self.win.set_light_direction(&self.scene.sun_dir()); // self.win.update_vaos(self.bodies.iter_mut().map(|b| &mut b.asteroid_mut().surface)); - self.win - .render_asteroids(&self.bodies.iter().map(|b| &b.asteroid).collect_vec()); + self.win.render_asteroids(&self.bodies); self.win.swap_window(); export.iteration( + &self.cfg, + &mut self.bodies, + &self.pre_computed_bodies, &mut self.time, + &self.scene, + &self.win, &self.folders, - &self.cfg, - &self.bodies, self.routines.as_ref(), + ); + + self.routines.fn_end_of_iteration( + &self.cfg, + &mut self.bodies, + &self.time, &self.scene, &self.win, ); - self.routines - .fn_end_of_iteration(&mut self.bodies, &self.time, &self.scene, &self.win); - if elapsed > self.cfg.simulation.duration { let time_calc = Utc::now().time() - *self.time.real_time(); println!( @@ -303,4 +302,8 @@ impl Scenario { Ok(()) } + + pub fn orientation_reference(&self, body: usize) -> &Mat4 { + &self.pre_computed_bodies[body].mat_orient + } } diff --git a/src/simu/util.rs b/src/simu/util.rs index 41e4b0a..4a3baa2 100644 --- a/src/simu/util.rs +++ b/src/simu/util.rs @@ -1,7 +1,7 @@ use crate::{ - cosine_angle, simu::Scene, util::*, Asteroid, Body, CfgBody, CfgColormap, CfgFrameCenter, - CfgMeshKind, CfgMeshSource, CfgOrbitKepler, CfgOrbitSpeedControl, CfgState, CfgStateCartesian, - Equatorial, Material, Result, Surface, Window, + cosine_angle, matrix_spin, position_in_inertial_frame, simu::Scene, util::*, AirlessBody, Cfg, + CfgBody, CfgFrameCenter, CfgMeshKind, CfgMeshSource, CfgOrbitKepler, CfgOrbitSpeedControl, + CfgState, CfgStateCartesian, Material, PreComputedBody, Result, Surface, Time, Window, }; use itertools::{izip, Itertools}; @@ -48,20 +48,24 @@ pub fn read_surface_low(cb: &CfgBody) -> Result { read_surface(cb, CfgMeshKind::Low) } -pub fn find_reference_matrix_orientation(cb: &CfgBody, other_bodies: &[&Body]) -> Mat4 { - match &cb.state { +pub fn find_reference_matrix_orientation( + cfg: &Cfg, + body: usize, + pre_computed_bodies: &mut [PreComputedBody], +) -> Mat4 { + match &cfg.bodies[body].state { CfgState::Cartesian(CfgStateCartesian { orientation, .. }) => { - glm::mat3_to_mat4(orientation) + glm::mat3_to_mat4(&orientation) } - CfgState::Equatorial(Equatorial { ra, dec }) => Mat4::identity(), + CfgState::Equatorial(_) => Mat4::identity(), CfgState::File(_p) => Mat4::identity(), CfgState::Orbit(orb) => match &orb.frame { CfgFrameCenter::Sun => Mat4::identity(), CfgFrameCenter::Body(id) => { let mut mat = Mat4::identity(); - for other_body in other_bodies { - if other_body.id == *id { - mat = other_body.mat_orient; + for (pre, cb) in izip!(pre_computed_bodies, &cfg.bodies) { + if cb.id == *id { + mat = pre.mat_orient; break; } } @@ -71,6 +75,54 @@ pub fn find_reference_matrix_orientation(cb: &CfgBody, other_bodies: &[&Body]) - } } +pub fn find_matrix_translation(cfg: &Cfg, body: usize, time: &Time) -> Mat4 { + let elapsed_from_start = time.elapsed_seconds_from_start(); + let other_bodies = cfg + .bodies + .iter() + .enumerate() + .filter(|(ii, _)| *ii != body) + .map(|(_, cb)| cb) + .collect_vec(); + + match &cfg.bodies[body].state { + CfgState::Cartesian(CfgStateCartesian { + position, + orientation: _orientation, + }) => Mat4::new_translation(position), + CfgState::Equatorial(_) | CfgState::File(_) => Mat4::identity(), + CfgState::Orbit(orb) => { + let (mu_ref, factor) = find_ref_orbit(&orb, &other_bodies); + if mu_ref == MU_SUN { + Mat4::identity() + } else { + let pos = position_in_inertial_frame( + orb.a * factor, + orb.e, + orb.i * RPD, + orb.node * RPD, + orb.peri * RPD, + elapsed_from_start as Float, + orb.tp, + mu_ref, + ); + Mat4::new_translation(&(pos * 1e-3)) + } + } + } +} + +pub fn find_matrix_spin(cfg: &Cfg, body: usize, time: &Time) -> Mat4 { + let elapsed = time.elapsed_seconds(); + if cfg.bodies[body].spin.period == 0.0 { + Mat4::identity() + } else { + let np_elapsed = elapsed as Float / cfg.bodies[body].spin.period; + let spin = (TAU * np_elapsed + cfg.bodies[body].spin.spin0 * RPD) % TAU; + matrix_spin(spin) + } +} + /// return MU and factor for distances. pub fn find_ref_orbit(orbit: &CfgOrbitKepler, cfgs: &[&CfgBody]) -> (Float, Float) { match &orbit.frame { @@ -99,11 +151,13 @@ pub fn find_ref_orbit(orbit: &CfgOrbitKepler, cfgs: &[&CfgBody]) -> (Float, Floa pub fn update_colormap_scalar( win: &Window, - cmap: &CfgColormap, + cfg: &Cfg, scalars: &[Float], - asteroid: &mut Asteroid, + asteroid: &mut AirlessBody, ii_body: usize, ) { + let cmap = &cfg.window.colormap; + for (element, scalar) in izip!(asteroid.surface.elements_mut(), scalars) { let mut data = (scalar - cmap.vmin) / (cmap.vmax - cmap.vmin); if cmap.reverse { @@ -116,37 +170,34 @@ pub fn update_colormap_scalar( } pub fn compute_cosine_incidence_angle( - body: &Body, + body: &AirlessBody, normals: &Matrix3xX, scene: &Scene, ) -> DRVector { - let matrix_normal: Mat3 = - glm::mat4_to_mat3(&glm::inverse_transpose(body.asteroid.matrix_model)); + let matrix_normal: Mat3 = glm::mat4_to_mat3(&glm::inverse_transpose(body.matrix_model)); let normals_oriented = matrix_normal * normals; cosine_angle(&scene.sun_dir(), &normals_oriented) } pub fn compute_cosine_emission_angle( - body: &Body, + body: &AirlessBody, normals: &Matrix3xX, scene: &Scene, ) -> DRVector { - let matrix_normal: Mat3 = - glm::mat4_to_mat3(&glm::inverse_transpose(body.asteroid.matrix_model)); + let matrix_normal: Mat3 = glm::mat4_to_mat3(&glm::inverse_transpose(body.matrix_model)); let normals_oriented = matrix_normal * normals; cosine_angle(&scene.cam_dir(), &normals_oriented) } -pub fn compute_cosine_phase_angle(body: &Body, scene: &Scene) -> DRVector { +pub fn compute_cosine_phase_angle(body: &AirlessBody, scene: &Scene) -> DRVector { DRVector::from_row_slice( &body - .asteroid .surface .faces .iter() .map(|f| { let v4 = glm::vec3_to_vec4(&f.vertex.position); - let v3_oriented = glm::vec4_to_vec3(&(body.asteroid.matrix_model * v4)); + let v3_oriented = glm::vec4_to_vec3(&(body.matrix_model * v4)); (scene.sun_pos - v3_oriented) .normalize() .dot(&(scene.cam_pos - v3_oriented).normalize()) diff --git a/src/win/scene.rs b/src/win/scene.rs index e831910..a3fd597 100644 --- a/src/win/scene.rs +++ b/src/win/scene.rs @@ -48,6 +48,10 @@ impl WindowScene { self.camera.position = pos; } + pub fn light_direction<>(&self) -> Vec3 { + self.light.direction() + } + pub fn set_light_offset(&mut self, offset: Float) { self.light.set_offset(offset); diff --git a/src/win/window.rs b/src/win/window.rs index 3e6329b..2407aff 100644 --- a/src/win/window.rs +++ b/src/win/window.rs @@ -1,6 +1,6 @@ use crate::{ - intersect_asteroids, util::*, win::WindowScene, Asteroid, Direction, GraphicalPipeline, Shader, - Surface, WindowSettings, WindowState, + intersect_asteroids, util::*, win::WindowScene, AirlessBody, Direction, GraphicalPipeline, + Shader, Surface, WindowSettings, WindowState, }; // use beryllium::events::*; @@ -275,6 +275,14 @@ impl Window { self.scene.borrow().camera_position().clone() } + pub fn set_camera_position(&self, pos: &Vec3) { + self.scene.borrow_mut().set_camera_position(&pos); + } + + pub fn light_direction(&self) -> Vec3 { + self.scene.borrow().light_direction().clone() + } + pub fn set_light_position(&self, pos: &Vec3) { let mut scene = self.scene.borrow_mut(); scene.set_light_position(pos); @@ -513,7 +521,7 @@ impl Window { self.win.win.gl_swap_window(); } - pub fn render_asteroids(&mut self, asteroids: &[&Asteroid]) { + pub fn render_asteroids(&mut self, asteroids: &[AirlessBody]) { let width_viewport = self.settings.borrow().width_viewport(); let height_viewport = self.settings.borrow().height_viewport(); let aspect_ratio = self.settings.borrow().aspect_ratio(); @@ -601,7 +609,12 @@ impl Window { self.render_trajectories(&matrix_projection, &matrix_view, &scene) } - pub fn render_trajectories(&self, matrix_projection: &Mat4, matrix_view: &Mat4, scene: &WindowScene) { + pub fn render_trajectories( + &self, + matrix_projection: &Mat4, + matrix_view: &Mat4, + scene: &WindowScene, + ) { let shader = self.graphical_pipeline.shaders.trajectory(); shader.set_mat4("matrix_projection", &matrix_projection); shader.set_mat4("matrix_view", &matrix_view); @@ -827,7 +840,7 @@ impl Window { fn mouse_click( &self, - asteroids: &[&Asteroid], + asteroids: &[AirlessBody], scene: &WindowScene, matrix_projection: &Mat4, matrix_view: &Mat4, diff --git a/src/win/window_settings.rs b/src/win/window_settings.rs index efe0a05..b87773c 100644 --- a/src/win/window_settings.rs +++ b/src/win/window_settings.rs @@ -71,6 +71,11 @@ impl Default for WindowState { impl WindowState { pub fn toggle_pause(&mut self) -> bool { self.pause = !self.pause; + match self.pause { + true => println!("Paused."), + false => println!("Unpaused."), + } + self.pause }