diff --git a/.gitignore b/.gitignore index dd4ff23..123d8d7 100644 --- a/.gitignore +++ b/.gitignore @@ -1,2 +1,3 @@ /target /examples/models/nasa/ +/data/ diff --git a/Cargo.lock b/Cargo.lock index db87dbd..1887508 100644 --- a/Cargo.lock +++ b/Cargo.lock @@ -696,6 +696,7 @@ version = "2.7.1" source = "registry+https://github.com/rust-lang/crates.io-index" checksum = "6ea2d84b969582b4b1864a92dc5d27cd2b77b622a8d79306834f1be5ba20d84b" dependencies = [ + "bytemuck", "cfg-if", "crunchy", "num-traits", @@ -1555,6 +1556,7 @@ dependencies = [ "bytemuck", "criterion", "gltf", + "half", "hifitime", "image", "log", diff --git a/crates/resolved-cli/src/main.rs b/crates/resolved-cli/src/main.rs index 5912424..ff56f96 100644 --- a/crates/resolved-cli/src/main.rs +++ b/crates/resolved-cli/src/main.rs @@ -4,11 +4,14 @@ mod scene_file; use std::path::{Path, PathBuf}; use std::process::ExitCode; +use std::sync::Arc; use std::time::Instant; -use clap::{Parser, Subcommand}; +use clap::{Args, Parser, Subcommand}; use resolved::gpu::AdapterChoice; -use resolved::{Frame, Renderer, RendererOptions}; +use resolved::{ + CameraTier, EarthTexture, Fidelity, Frame, LightingTier, Renderer, RendererOptions, +}; use crate::scene_file::SceneFile; @@ -25,32 +28,46 @@ struct Cli { #[derive(Subcommand)] enum Command { /// Render a scene file. - Render { - /// The scene, as TOML (see examples/scenes/). - #[arg(long)] - scene: PathBuf, - /// Raw DN as a 16-bit PNG (RGB for three channels, gray for one), shifted up to - /// fill the 16-bit range: DN = value >> (16 - ADC bits), so 12-bit DN are x 16. - #[arg(long)] - out: Option, - /// An 8-bit sRGB PNG for viewing: bias subtracted, saturation as white. - #[arg(long)] - preview: Option, - /// Device: auto, gpu, cpu, or part of an adapter's name. - #[arg(long, default_value = "auto")] - adapter: String, - /// Render this many frames and report setup and per-frame times. - #[arg(long, value_name = "FRAMES")] - timings: Option, - /// Directory of astrokit's kernels, for scenes given by epoch. Without it: - /// ASTROKIT_KERNELS, then ./data/kernels, then the sibling astrokit checkout. - #[arg(long, value_name = "DIR")] - kernels: Option, - }, + Render(RenderArgs), /// List the adapters this machine offers. Adapters, } +#[derive(Args)] +struct RenderArgs { + /// The scene, as TOML (see examples/scenes/). + #[arg(long)] + scene: PathBuf, + /// Raw DN as a 16-bit PNG (RGB for three channels, gray for one), shifted up to + /// fill the 16-bit range: DN = value >> (16 - ADC bits), so 12-bit DN are x 16. + #[arg(long)] + out: Option, + /// An 8-bit sRGB PNG for viewing: bias subtracted, saturation as white. + #[arg(long)] + preview: Option, + /// Lighting fidelity tier: L0, or L1 (adds Earth's surface; needs Blue Marble data). + #[arg(long, default_value = "L0")] + lighting: LightingTier, + /// Camera fidelity tier: C0. + #[arg(long, default_value = "C0")] + camera: CameraTier, + /// Device: auto, gpu, cpu, or part of an adapter's name. + #[arg(long, default_value = "auto")] + adapter: String, + /// Render this many frames and report setup and per-frame times. + #[arg(long, value_name = "FRAMES")] + timings: Option, + /// Directory of data the tiers need (earth/bmng-2004-MM.jpg from + /// scripts/fetch-earth.sh). Without it: RESOLVED_DATA, then ./data, then the + /// repository's data/. + #[arg(long, value_name = "DIR")] + data: Option, + /// Directory of astrokit's kernels, for scenes given by epoch. Without it: + /// ASTROKIT_KERNELS, then ./data/kernels, then the sibling astrokit checkout. + #[arg(long, value_name = "DIR")] + kernels: Option, +} + fn main() -> ExitCode { let result = match Cli::parse().command { Command::Adapters => { @@ -62,21 +79,7 @@ fn main() -> ExitCode { } Ok(()) } - Command::Render { - scene, - out, - preview, - adapter, - timings, - kernels, - } => render( - &scene, - out.as_deref(), - preview.as_deref(), - &adapter, - timings, - kernels.as_deref(), - ), + Command::Render(args) => render(&args), }; match result { Ok(()) => ExitCode::SUCCESS, @@ -87,37 +90,79 @@ fn main() -> ExitCode { } } -fn render( - path: &Path, - out: Option<&Path>, - preview: Option<&Path>, - adapter: &str, - timings: Option, - kernels: Option<&Path>, -) -> Result<(), String> { +/// Where to look for data: `explicit` alone if given, else `RESOLVED_DATA`, `./data` +/// and the repository's `data/`. +fn data_directories(explicit: Option<&Path>) -> Vec { + if let Some(dir) = explicit { + return vec![dir.to_owned()]; + } + let repository = Path::new(env!("CARGO_MANIFEST_DIR")).join("../../data"); + std::env::var_os("RESOLVED_DATA") + .map(PathBuf::from) + .into_iter() + .chain([PathBuf::from("data"), repository]) + .collect() +} + +/// The Earth image for `month` (1-12), from the first data directory that has it. +fn earth_image(data: Option<&Path>, month: u32) -> Result { + let name = format!("earth/bmng-2004-{month:02}.jpg"); + let tried = data_directories(data); + tried.iter().map(|dir| dir.join(&name)).find(|p| p.is_file()).ok_or_else(|| { + let tried: Vec<_> = tried.iter().map(|d| d.display().to_string()).collect(); + format!( + "lighting L1 needs {name} (run scripts/fetch-earth.sh, or pass --data DIR); looked in {}", + tried.join(", ") + ) + }) +} + +fn render(args: &RenderArgs) -> Result<(), String> { + let (out, preview, timings) = (args.out.as_deref(), args.preview.as_deref(), args.timings); if out.is_none() && preview.is_none() && timings.is_none() { return Err("nothing to do: give --out, --preview or --timings".into()); } + let path = args.scene.as_path(); let file = SceneFile::read(path)?; let camera = file.camera(); + + // What the tiers need, gathered before the renderer exists so it can check them. + let start = Instant::now(); + let earth_texture = if args.lighting >= LightingTier::L1 { + let month = file.earth_month()?.ok_or( + "lighting L1 needs the month for Earth's surface: give the scene an epoch, or [appearance] earth_month", + )?; + let image = earth_image(args.data.as_deref(), month)?; + Some(Arc::new( + EarthTexture::load(&image).map_err(|e| e.to_string())?, + )) + } else { + None + }; + let texturing = start.elapsed(); + let options = RendererOptions { - adapter: match adapter { + adapter: match args.adapter.as_str() { "auto" => AdapterChoice::Auto, "gpu" => AdapterChoice::Gpu, "cpu" => AdapterChoice::Cpu, name => AdapterChoice::Named(name.to_owned()), }, + fidelity: Fidelity { + lighting: args.lighting, + camera: args.camera, + }, + earth_texture, appearance: file.appearance(), read_radiance: false, profile: timings.is_some(), }; - let start = Instant::now(); let mut renderer = Renderer::new(camera, options).map_err(|e| e.to_string())?; let setup = start.elapsed(); let start = Instant::now(); let dir = path.parent().unwrap_or(Path::new(".")); - let scene = file.scene(&mut renderer, dir, kernels)?; + let scene = file.scene(&mut renderer, dir, args.kernels.as_deref())?; let observer = file.observer()?; let loading = start.elapsed(); @@ -138,7 +183,15 @@ fn render( if timings.is_some() { println!("adapter {}", renderer.adapter_info().name); - println!("setup {:.1} ms", setup.as_secs_f64() * 1e3); + println!("tiers {} {}", args.lighting, args.camera); + println!( + "setup {:.1} ms (includes the Earth texture's upload)", + setup.as_secs_f64() * 1e3 + ); + println!( + "earth {:.1} ms (texture decode)", + texturing.as_secs_f64() * 1e3 + ); println!( "scene {:.1} ms (models and ephemeris)", loading.as_secs_f64() * 1e3 diff --git a/crates/resolved-cli/src/scene_file.rs b/crates/resolved-cli/src/scene_file.rs index a9c6679..7086dcb 100644 --- a/crates/resolved-cli/src/scene_file.rs +++ b/crates/resolved-cli/src/scene_file.rs @@ -75,6 +75,9 @@ pub struct ObserverFile { pub struct AppearanceFile { pub earth_albedo: Option, pub moon_albedo: Option, + /// Which month (1-12) Earth's surface shows at lighting L1, for scenes without an + /// epoch; with one, the epoch's month is the default. + pub earth_month: Option, } /// One spacecraft. @@ -154,6 +157,32 @@ impl SceneFile { a } + /// The month (1-12) for Earth's surface: `earth_month`, else the epoch's month. + pub fn earth_month(&self) -> Result> { + let month = match self.appearance.earth_month { + Some(m) => Some(m), + None => self.epoch_month()?, + }; + match month { + Some(m) if !(1..=12).contains(&m) => Err(format!("earth_month {m} is not 1-12")), + m => Ok(m), + } + } + + #[cfg(feature = "astrokit")] + fn epoch_month(&self) -> Result> { + let Some(epoch) = &self.epoch else { + return Ok(None); + }; + let epoch: hifitime::Epoch = epoch.parse().map_err(|e| format!("epoch {epoch:?}: {e}"))?; + Ok(Some(epoch.to_gregorian_utc().1 as u32)) + } + + #[cfg(not(feature = "astrokit"))] + fn epoch_month(&self) -> Result> { + Ok(None) + } + /// The observer's pose. pub fn observer(&self) -> Result { let o = &self.observer; @@ -329,6 +358,28 @@ mod tests { assert_eq!(camera.detector.adc_bits, 10); } + #[test] + fn earth_month_comes_from_the_setting() { + let with = |extra: &str| -> SceneFile { + toml::from_str(&format!("{MINIMAL}\n[appearance]\n{extra}\n")).unwrap() + }; + assert_eq!(with("earth_month = 3").earth_month().unwrap(), Some(3)); + assert_eq!( + with("").earth_month().unwrap(), + None, + "no month and no epoch" + ); + assert!(with("earth_month = 13").earth_month().is_err()); + } + + #[cfg(feature = "astrokit")] + #[test] + fn earth_month_defaults_to_the_epoch() { + let file: SceneFile = + toml::from_str(&format!("epoch = \"2026-10-14T12:00:00 UTC\"\n{MINIMAL}")).unwrap(); + assert_eq!(file.earth_month().unwrap(), Some(10)); + } + #[test] fn mistakes_are_reported() { // A misspelled field. diff --git a/crates/resolved-cli/tests/cli.rs b/crates/resolved-cli/tests/cli.rs index 9fb66d1..88e0843 100644 --- a/crates/resolved-cli/tests/cli.rs +++ b/crates/resolved-cli/tests/cli.rs @@ -123,6 +123,60 @@ fn missing_kernels_say_where_they_looked() { ); } +#[test] +fn lighting_l1_says_what_it_is_missing() { + // boxwing.toml has no epoch and no earth_month: L1 cannot pick a month. + let no_month = resolved(&[ + "render", + "--scene", + "examples/scenes/boxwing.toml", + "--lighting", + "L1", + "--preview", + "x.png", + ]); + assert!(!no_month.status.success()); + let stderr = String::from_utf8_lossy(&no_month.stderr); + assert!(stderr.contains("earth_month"), "{stderr}"); + + // With a month but no data: the error names the file and the fetch script. + let scene = out_dir().join("month.toml"); + let models = repo().join("examples/models").canonicalize().unwrap(); + let text = std::fs::read_to_string(repo().join("examples/scenes/boxwing.toml")) + .unwrap() + .replace("../models", models.to_str().unwrap()) + + "\n[appearance]\nearth_month = 7\n"; + std::fs::write(&scene, text).unwrap(); + let no_data = resolved(&[ + "render", + "--scene", + scene.to_str().unwrap(), + "--lighting", + "L1", + "--data", + "no/such/dir", + "--preview", + "x.png", + ]); + assert!(!no_data.status.success()); + let stderr = String::from_utf8_lossy(&no_data.stderr); + assert!( + stderr.contains("bmng-2004-07.jpg") && stderr.contains("fetch-earth.sh"), + "{stderr}" + ); + + let unknown = resolved(&[ + "render", + "--scene", + "examples/scenes/boxwing.toml", + "--lighting", + "L7", + "--preview", + "x.png", + ]); + assert!(!unknown.status.success()); +} + #[test] fn adapters_lists_something() { let result = resolved(&["adapters"]); diff --git a/crates/resolved/Cargo.toml b/crates/resolved/Cargo.toml index 5bd1492..93e5545 100644 --- a/crates/resolved/Cargo.toml +++ b/crates/resolved/Cargo.toml @@ -30,6 +30,8 @@ pollster = "1.0" bytemuck = { version = "1.0", features = ["derive"] } log = "0.4" gltf = "1.4" +image = { version = "0.25", default-features = false, features = ["png", "jpeg"] } +half = { version = "2.7", features = ["bytemuck"] } stl_io = "0.11" [dev-dependencies] diff --git a/crates/resolved/src/bodies.rs b/crates/resolved/src/bodies.rs index 4c6be67..fc53bf9 100644 --- a/crates/resolved/src/bodies.rs +++ b/crates/resolved/src/bodies.rs @@ -17,8 +17,12 @@ pub const MOON_RADIUS: f64 = 1_737_400.0; /// and better reflectance models in later milestones. #[derive(Clone, Debug, PartialEq)] pub struct Appearance { - /// Earth's reflectance. + /// Earth's reflectance, used when no Earth texture is set. pub earth_albedo: Spectrum, + /// Added to an Earth texture's reflectance: a stand-in for the light Earth's + /// atmosphere scatters toward the camera, until the atmosphere is modeled (M4). + /// Not used with `earth_albedo`, which already includes the atmosphere's effect. + pub earth_haze: Spectrum, /// The Moon's reflectance. pub moon_albedo: Spectrum, } @@ -42,11 +46,24 @@ impl Default for Appearance { Self { earth_albedo: Spectrum::sampled(nm.map(|x| x * 1e-9).to_vec(), earth.to_vec()) .expect("fixed samples"), + earth_haze: rayleigh_haze(), moon_albedo: Spectrum::constant(0.12), } } } +/// Reflectance shaped like Rayleigh scattering, falling as 1/lambda^4: 0.07 at 450 nm, +/// 0.017 at 650 nm. Sampled every 10 nm from 300 to 1100 nm, zero outside. +/// +/// The size is a round figure of the right order: clear air's Rayleigh optical depth is +/// about 0.2 at 450 nm, and the light it scatters toward a camera above it amounts to a +/// reflectance of a few hundredths to a tenth, depending on the angles. +fn rayleigh_haze() -> Spectrum { + let nm: Vec = (30..=110).map(|i| i as f64 * 10.0).collect(); + let value = nm.iter().map(|l| 0.07 * (450.0 / l).powi(4)).collect(); + Spectrum::sampled(nm.iter().map(|l| l * 1e-9).collect(), value).expect("increasing samples") +} + #[cfg(test)] mod tests { use super::*; diff --git a/crates/resolved/src/earth.rs b/crates/resolved/src/earth.rs new file mode 100644 index 0000000..a93227f --- /dev/null +++ b/crates/resolved/src/earth.rs @@ -0,0 +1,192 @@ +//! Earth's surface as an image: albedo by latitude and longitude. +//! +//! The image is equirectangular: longitude runs from -180 degrees at the left edge to +//! +180 at the right, latitude from +90 at the top to -90 at the bottom, in Earth-fixed +//! (ITRF) axes with geodetic latitude. Pixel values are sRGB-encoded red, green and +//! blue reflectance, as in NASA's Blue Marble (`scripts/fetch-earth.sh`). + +use std::path::Path; + +use half::f16; + +use crate::Error; + +/// An equirectangular image of Earth's surface reflectance. +#[derive(Clone)] +pub struct EarthTexture { + width: u32, + height: u32, + /// RGBA8, sRGB-encoded; alpha unused. + pixels: Vec, +} + +impl std::fmt::Debug for EarthTexture { + fn fmt(&self, f: &mut std::fmt::Formatter<'_>) -> std::fmt::Result { + f.debug_struct("EarthTexture") + .field("width", &self.width) + .field("height", &self.height) + .finish_non_exhaustive() + } +} + +impl EarthTexture { + /// Loads an image file (JPEG or PNG). + pub fn load(path: &Path) -> Result { + let image = + image::open(path).map_err(|e| Error::Texture(format!("{}: {e}", path.display())))?; + let rgba = image.to_rgba8(); + Self::from_rgba8(rgba.width(), rgba.height(), rgba.into_raw()) + } + + /// An image from sRGB-encoded RGBA8 pixels, rows from the top. + pub fn from_rgba8(width: u32, height: u32, pixels: Vec) -> Result { + if width == 0 || height == 0 || pixels.len() != 4 * width as usize * height as usize { + return Err(Error::Texture(format!( + "{width} x {height} needs {} bytes of RGBA, got {}", + 4 * width as usize * height as usize, + pixels.len() + ))); + } + Ok(Self { + width, + height, + pixels, + }) + } + + /// Width in pixels. + pub fn width(&self) -> u32 { + self.width + } + + /// Height in pixels. + pub fn height(&self) -> u32 { + self.height + } + + /// Every mip level from full size down to 1 x 1, each half the size of the one + /// before (rounded down, at least 1), as linear reflectance in RGBA f16. + /// + /// Each level averages 2 x 2 blocks of the one above, in linear values, so a far view + /// of Earth averages reflectance. The GPU stores and filters these linear values too: + /// sRGB textures would let the hardware decode, but lavapipe's sRGB decode is off by + /// up to 3 % for dark values like ocean, and filtering encoded values instead darkens + /// the blend between contrasting texels (a coastline) by up to half. + pub(crate) fn mip_chain(&self) -> Vec<(u32, u32, Vec)> { + let decode: Vec = (0..=255u8) + .map(|v| srgb_to_linear(v as f32 / 255.0)) + .collect(); + let one = f16::from_f32(1.0); + let mut first = vec![one; self.pixels.len()]; + parallel_rows(&mut first, 4 * self.width as usize, |y, row| { + let src = &self.pixels[y * row.len()..(y + 1) * row.len()]; + for (dst, src) in row.chunks_exact_mut(4).zip(src.chunks_exact(4)) { + for c in 0..3 { + dst[c] = f16::from_f32(decode[src[c] as usize]); + } + } + }); + let mut levels = vec![(self.width, self.height, first)]; + loop { + let (w, h, above) = levels.last().expect("at least one level"); + let (w, h) = (*w, *h); + if w == 1 && h == 1 { + break; + } + let (nw, nh) = ((w / 2).max(1), (h / 2).max(1)); + let mut below = vec![one; 4 * nw as usize * nh as usize]; + parallel_rows(&mut below, 4 * nw as usize, |y, row| { + let y = y as u32; + for x in 0..nw { + for c in 0..3 { + let mut sum = 0.0; + for (dx, dy) in [(0, 0), (1, 0), (0, 1), (1, 1)] { + let sx = (2 * x + dx).min(w - 1); + let sy = (2 * y + dy).min(h - 1); + sum += above[4 * (sy * w + sx) as usize + c].to_f32(); + } + row[4 * x as usize + c] = f16::from_f32(sum / 4.0); + } + } + }); + levels.push((nw, nh, below)); + } + levels + } +} + +/// Runs `f(row_index, row)` over the rows of `data` (rows of `row_len` values) on all +/// cores. +fn parallel_rows(data: &mut [T], row_len: usize, f: impl Fn(usize, &mut [T]) + Sync) { + let rows = data.len() / row_len; + let threads = std::thread::available_parallelism() + .map_or(1, |n| n.get()) + .min(rows.max(1)); + let per_thread = rows.div_ceil(threads); + std::thread::scope(|scope| { + for (t, block) in data.chunks_mut(per_thread * row_len).enumerate() { + let f = &f; + scope.spawn(move || { + for (i, row) in block.chunks_mut(row_len).enumerate() { + f(t * per_thread + i, row); + } + }); + } + }); +} + +fn srgb_to_linear(v: f32) -> f32 { + if v <= 0.040_45 { + v / 12.92 + } else { + ((v + 0.055) / 1.055).powf(2.4) + } +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn mip_chain_averages_in_linear_light() { + // A 2 x 2 checker of black and white averages to 50 % reflectance, not the 22 % + // that averaging encoded values (sRGB 128) would give. + let px = |v: u8| [v, v, v, 255]; + let pixels = [px(0), px(255), px(255), px(0)].concat(); + let texture = EarthTexture::from_rgba8(2, 2, pixels).unwrap(); + let chain = texture.mip_chain(); + assert_eq!(chain.len(), 2); + assert_eq!((chain[1].0, chain[1].1), (1, 1)); + assert_eq!( + chain[1].2[..3] + .iter() + .map(|v| v.to_f32()) + .collect::>(), + [0.5; 3] + ); + assert_eq!(chain[0].2[4].to_f32(), 1.0, "sRGB 255 is reflectance 1"); + } + + #[test] + fn mip_chain_handles_non_square_sizes() { + let texture = EarthTexture::from_rgba8(8, 4, vec![100; 8 * 4 * 4]).unwrap(); + let sizes: Vec<_> = texture + .mip_chain() + .iter() + .map(|(w, h, _)| (*w, *h)) + .collect(); + assert_eq!(sizes, [(8, 4), (4, 2), (2, 1), (1, 1)]); + // A uniform image stays uniform at every level. + let first = texture.mip_chain()[0].2[0]; + assert!(texture + .mip_chain() + .iter() + .all(|(_, _, p)| p.chunks(4).all(|t| t[..3] == [first; 3]))); + } + + #[test] + fn wrong_sizes_are_rejected() { + assert!(EarthTexture::from_rgba8(4, 2, vec![0; 31]).is_err()); + assert!(EarthTexture::from_rgba8(0, 2, vec![]).is_err()); + } +} diff --git a/crates/resolved/src/error.rs b/crates/resolved/src/error.rs index 7ee3c11..e2a7e4f 100644 --- a/crates/resolved/src/error.rs +++ b/crates/resolved/src/error.rs @@ -21,10 +21,18 @@ pub enum Error { #[error("ephemeris: {0}")] Ephemeris(String), + /// An image that could not be read or has the wrong size. + #[error("texture: {0}")] + Texture(String), + /// A model file that could not be read or describes no usable geometry. #[error("model: {0}")] Model(String), + /// The chosen fidelity tier needs an input that was not supplied. + #[error("missing input: {0}")] + MissingInput(String), + /// A request this version cannot carry out yet. #[error("unsupported: {0}")] Unsupported(String), diff --git a/crates/resolved/src/fidelity.rs b/crates/resolved/src/fidelity.rs new file mode 100644 index 0000000..1d16300 --- /dev/null +++ b/crates/resolved/src/fidelity.rs @@ -0,0 +1,100 @@ +//! Fidelity tiers: how much of the physics a frame models, chosen separately for +//! lighting and for the camera. +//! +//! A tier is a fixed bundle of features, so a tier name says exactly what an image +//! includes, and each tier's cost is measured as a whole (the readme's performance +//! table). A tier that needs data, such as L1's Earth image, makes +//! [`Renderer::new`](crate::Renderer::new) fail if that data is not supplied. + +use std::fmt; +use std::str::FromStr; + +/// Lighting tiers, from cheapest to most complete. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, PartialOrd, Ord, Hash)] +#[non_exhaustive] +pub enum LightingTier { + /// Lambertian surfaces with one reflectance per body, exact Earth and Moon limbs, + /// eclipses of the target by Earth and the Moon. No shadows on the target, no light + /// from Earth or the Moon onto it. + #[default] + L0, + /// L0, plus Earth's surface from an image (`RendererOptions::earth_texture`, NASA's + /// Blue Marble) with a stand-in for the atmosphere's blue. The rest of M2 (shadows, + /// Earthshine, the Moon's albedo map) joins this tier as it lands. + L1, +} + +/// Camera tiers, from cheapest to most complete. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, PartialOrd, Ord, Hash)] +#[non_exhaustive] +pub enum CameraTier { + /// A pinhole camera with absolute radiometry to DN: no distortion, blur or noise. + #[default] + C0, +} + +/// The lighting and camera tiers a renderer draws at. +#[derive(Clone, Copy, Debug, Default, PartialEq, Eq, Hash)] +pub struct Fidelity { + /// How light reaches and leaves surfaces. + pub lighting: LightingTier, + /// How the camera turns light into DN. + pub camera: CameraTier, +} + +impl fmt::Display for LightingTier { + fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result { + f.write_str(match self { + Self::L0 => "L0", + Self::L1 => "L1", + }) + } +} + +impl fmt::Display for CameraTier { + fn fmt(&self, f: &mut fmt::Formatter<'_>) -> fmt::Result { + f.write_str(match self { + Self::C0 => "C0", + }) + } +} + +impl FromStr for LightingTier { + type Err = String; + + /// `L0` or `L1`, either case. + fn from_str(s: &str) -> Result { + match s.to_ascii_uppercase().as_str() { + "L0" => Ok(Self::L0), + "L1" => Ok(Self::L1), + _ => Err(format!("unknown lighting tier {s:?}: L0 or L1")), + } + } +} + +impl FromStr for CameraTier { + type Err = String; + + /// `C0`, either case. + fn from_str(s: &str) -> Result { + match s.to_ascii_uppercase().as_str() { + "C0" => Ok(Self::C0), + _ => Err(format!("unknown camera tier {s:?}: C0")), + } + } +} + +#[cfg(test)] +mod tests { + use super::*; + + #[test] + fn tiers_parse_and_print() { + assert_eq!("l1".parse::().unwrap(), LightingTier::L1); + assert_eq!(LightingTier::L1.to_string(), "L1"); + assert_eq!("C0".parse::().unwrap(), CameraTier::C0); + assert!("L9".parse::().is_err()); + assert!(LightingTier::L0 < LightingTier::L1); + assert_eq!(Fidelity::default().lighting, LightingTier::L0); + } +} diff --git a/crates/resolved/src/gpu/bodies.rs b/crates/resolved/src/gpu/bodies.rs index a1277c8..2ad20d5 100644 --- a/crates/resolved/src/gpu/bodies.rs +++ b/crates/resolved/src/gpu/bodies.rs @@ -11,6 +11,8 @@ use super::targets::{FrameTargets, DEPTH_COMPARE, DEPTH_FORMAT, RADIANCE_FORMAT} use super::GpuContext; use crate::bodies::{Appearance, EARTH_EQUATORIAL_RADIUS, EARTH_POLAR_RADIUS, MOON_RADIUS}; use crate::camera::Intrinsics; +use crate::earth::EarthTexture; +use crate::model::rgb_spectrum; use crate::radiometry::{Radiometry, AU, SUN_RADIUS}; use crate::scene::{Pose, Scene}; use crate::Error; @@ -21,8 +23,18 @@ pub const MAX_CHANNELS: usize = 4; /// Per-channel constants for the bodies that do not change from frame to frame. #[derive(Clone, Copy, Debug, PartialEq)] pub struct BodyRadiometry { - /// Channel radiance of Earth's surface facing the Sun at 1 AU: albedo S_k / pi. + /// Channel radiance of Earth's surface facing the Sun at 1 AU: albedo S_k / pi, + /// with `Appearance::earth_albedo`. Used when Earth has no texture. pub earth_lambert_1au: [f32; 4], + /// The same with albedo 1: what a texture's per-channel albedo multiplies. + pub earth_irradiance_1au: [f32; 4], + /// Per channel, the albedo of each RGB basis spectrum (see [`rgb_spectrum`]): column + /// j maps a texture's channel j (red, green, blue) to the camera's channels. + pub earth_rgb: [[f32; 4]; 3], + /// Per channel, the albedo of `Appearance::earth_haze`, added to a texture's. + pub earth_haze: [f32; 4], + /// Whether Earth is drawn from a texture. + pub earth_textured: bool, /// The same for the Moon. pub moon_lambert_1au: [f32; 4], /// Channel radiance of the solar disk. @@ -46,12 +58,28 @@ impl BodyRadiometry { } out }; + let albedo = |spectrum: &crate::spectrum::Spectrum| { + let mut out = [0.0; 4]; + for (k, a) in radiometry.effective_albedo(spectrum).iter().enumerate() { + out[k] = *a as f32; + } + out + }; + let basis = |j: usize| { + let mut rgb = [0.0; 3]; + rgb[j] = 1.0; + albedo(&rgb_spectrum(rgb)) + }; let mut sun_radiance = [0.0; 4]; for (k, s) in sun_radiance.iter_mut().take(n).enumerate() { *s = radiometry.sun_radiance(k) as f32; } Ok(Self { earth_lambert_1au: lambert(radiometry.effective_albedo(&appearance.earth_albedo)), + earth_irradiance_1au: lambert(vec![1.0; n]), + earth_rgb: [basis(0), basis(1), basis(2)], + earth_haze: albedo(&appearance.earth_haze), + earth_textured: false, moon_lambert_1au: lambert(radiometry.effective_albedo(&appearance.moon_albedo)), sun_radiance, }) @@ -70,6 +98,10 @@ pub struct BodiesUniforms { earth_shape: [f32; 4], earth_sun: [f32; 4], earth_lambert: [f32; 4], + earth_irradiance: [f32; 4], + // mat3x4 columns: texture red, green, blue to channel albedo. + earth_rgb: [[f32; 4]; 3], + earth_haze: [f32; 4], moon_center: [f32; 4], moon_sun: [f32; 4], moon_lambert: [f32; 4], @@ -103,7 +135,12 @@ impl BodiesUniforms { Self { intrinsics: [intrinsics.fx, intrinsics.fy, intrinsics.cx, intrinsics.cy] .map(|x| x as f32), - params: [near as f32, 0.0, 0.0, 0.0], + params: [ + near as f32, + f32::from(u8::from(radiometry.earth_textured)), + 0.0, + 0.0, + ], earth_center: vec4(rel.earth.translation, EARTH_EQUATORIAL_RADIUS), earth_axes: [column(0), column(1), column(2)], earth_shape: [ @@ -114,6 +151,9 @@ impl BodiesUniforms { ], earth_sun: vec4(sun_from(Vector3::zeros()), 0.0), earth_lambert: scaled(radiometry.earth_lambert_1au, scene.sun.norm()), + earth_irradiance: scaled(radiometry.earth_irradiance_1au, scene.sun.norm()), + earth_rgb: radiometry.earth_rgb, + earth_haze: radiometry.earth_haze, moon_center: vec4(rel.moon, MOON_RADIUS), moon_sun: vec4(sun_from(scene.moon), 0.0), moon_lambert: scaled(radiometry.moon_lambert_1au, (scene.sun - scene.moon).norm()), @@ -128,6 +168,7 @@ impl BodiesUniforms { pub struct BodiesPass { pipeline: wgpu::RenderPipeline, uniforms: wgpu::Buffer, + sampler: wgpu::Sampler, bind_group: wgpu::BindGroup, } @@ -172,21 +213,99 @@ impl BodiesPass { usage: wgpu::BufferUsages::UNIFORM | wgpu::BufferUsages::COPY_DST, mapped_at_creation: false, }); - let bind_group = device.create_bind_group(&wgpu::BindGroupDescriptor { - label: Some("bodies"), - layout: &pipeline.get_bind_group_layout(0), - entries: &[wgpu::BindGroupEntry { - binding: 0, - resource: uniforms.as_entire_binding(), - }], + // Repeat across the 180 degree meridian, clamp at the poles. Not anisotropic: + // lavapipe's anisotropic filtering with explicit gradients falls back to coarse + // mips in places (measured as a color 26x off 13 px from the seam). + let sampler = device.create_sampler(&wgpu::SamplerDescriptor { + label: Some("earth"), + address_mode_u: wgpu::AddressMode::Repeat, + address_mode_v: wgpu::AddressMode::ClampToEdge, + mag_filter: wgpu::FilterMode::Linear, + min_filter: wgpu::FilterMode::Linear, + mipmap_filter: wgpu::MipmapFilterMode::Linear, + ..Default::default() }); + // A 1 x 1 placeholder until a texture is set; the shader ignores it then. + let placeholder = EarthTexture::from_rgba8(1, 1, vec![0, 0, 0, 255]).expect("1 x 1"); + let bind_group = Self::bind(ctx, &pipeline, &uniforms, &sampler, &placeholder); Self { pipeline, uniforms, + sampler, bind_group, } } + /// Uploads `texture` with its mip chain and binds it for later frames. + pub fn set_earth_texture(&mut self, ctx: &GpuContext, texture: &EarthTexture) { + self.bind_group = Self::bind(ctx, &self.pipeline, &self.uniforms, &self.sampler, texture); + } + + fn bind( + ctx: &GpuContext, + pipeline: &wgpu::RenderPipeline, + uniforms: &wgpu::Buffer, + sampler: &wgpu::Sampler, + texture: &EarthTexture, + ) -> wgpu::BindGroup { + let levels = texture.mip_chain(); + let gpu = ctx.device.create_texture(&wgpu::TextureDescriptor { + label: Some("earth"), + size: wgpu::Extent3d { + width: texture.width(), + height: texture.height(), + depth_or_array_layers: 1, + }, + mip_level_count: levels.len() as u32, + sample_count: 1, + dimension: wgpu::TextureDimension::D2, + // Linear reflectance, so filtering averages reflectance (see `mip_chain`). + format: wgpu::TextureFormat::Rgba16Float, + usage: wgpu::TextureUsages::TEXTURE_BINDING | wgpu::TextureUsages::COPY_DST, + view_formats: &[], + }); + for (level, (w, h, pixels)) in levels.iter().enumerate() { + ctx.queue.write_texture( + wgpu::TexelCopyTextureInfo { + texture: &gpu, + mip_level: level as u32, + origin: wgpu::Origin3d::ZERO, + aspect: wgpu::TextureAspect::All, + }, + bytemuck::cast_slice(pixels), + wgpu::TexelCopyBufferLayout { + offset: 0, + bytes_per_row: Some(8 * w), + rows_per_image: Some(*h), + }, + wgpu::Extent3d { + width: *w, + height: *h, + depth_or_array_layers: 1, + }, + ); + } + let view = gpu.create_view(&Default::default()); + ctx.device.create_bind_group(&wgpu::BindGroupDescriptor { + label: Some("bodies"), + layout: &pipeline.get_bind_group_layout(0), + entries: &[ + wgpu::BindGroupEntry { + binding: 0, + resource: uniforms.as_entire_binding(), + }, + wgpu::BindGroupEntry { + binding: 1, + resource: wgpu::BindingResource::TextureView(&view), + }, + wgpu::BindGroupEntry { + binding: 2, + resource: wgpu::BindingResource::Sampler(sampler), + }, + ], + }) + } + /// Uploads this frame's uniforms and records the pass into `encoder`, writing /// `timestamps` around it if given. pub fn encode( diff --git a/crates/resolved/src/lib.rs b/crates/resolved/src/lib.rs index 24ad84b..dcd405b 100644 --- a/crates/resolved/src/lib.rs +++ b/crates/resolved/src/lib.rs @@ -19,10 +19,12 @@ pub use wgpu; pub mod bodies; pub mod camera; +pub mod earth; pub mod eclipse; #[cfg(feature = "astrokit")] pub mod ephemeris; pub mod error; +pub mod fidelity; pub mod gpu; pub mod model; pub mod radiometry; @@ -32,7 +34,9 @@ mod solar_data; pub mod spectrum; pub use camera::{CameraConfig, Channel, Detector, Intrinsics}; +pub use earth::EarthTexture; pub use error::Error; +pub use fidelity::{CameraTier, Fidelity, LightingTier}; pub use model::{LoadOptions, Material, Mesh, ModelId}; pub use radiometry::Radiometry; pub use render::{Frame, Pending, Renderer, RendererOptions}; diff --git a/crates/resolved/src/model.rs b/crates/resolved/src/model.rs index e4e106b..8d2fc02 100644 --- a/crates/resolved/src/model.rs +++ b/crates/resolved/src/model.rs @@ -39,12 +39,9 @@ impl Material { /// does not determine a spectrum, so a material whose reflectance matters should be /// given one measured spectrum instead. pub fn from_rgb(name: &str, rgb: [f64; 3]) -> Self { - let [r, g, b] = rgb; - let nm = [200.0, 490.0, 500.0, 580.0, 590.0, 2800.0].map(|x| x * 1e-9); Self { name: name.to_owned(), - albedo: Spectrum::sampled(nm.to_vec(), vec![b, b, g, g, r, r]) - .expect("fixed wavelengths"), + albedo: rgb_spectrum(rgb), } } } @@ -55,6 +52,17 @@ impl Material { #[derive(Clone, Copy, Debug, PartialEq, Eq, Hash)] pub struct ModelId(pub usize); +/// The stepped spectrum an RGB reflectance stands for: blue below 490 nm, green from +/// 500 to 580 nm, red above 590 nm (into the infrared), with linear steps between. +/// +/// It is linear in `rgb`, so per-channel albedos of an RGB image can be folded into a +/// 3 x N matrix once (see [`rgb_spectrum`] with unit vectors) instead of per texel. +pub fn rgb_spectrum(rgb: [f64; 3]) -> Spectrum { + let [r, g, b] = rgb; + let nm = [200.0, 490.0, 500.0, 580.0, 590.0, 2800.0].map(|x| x * 1e-9); + Spectrum::sampled(nm.to_vec(), vec![b, b, g, g, r, r]).expect("fixed wavelengths") +} + /// A named piece of a model, drawn with one material. #[derive(Clone, Debug, PartialEq)] pub struct Part { diff --git a/crates/resolved/src/render.rs b/crates/resolved/src/render.rs index f2a5058..8a0a3ee 100644 --- a/crates/resolved/src/render.rs +++ b/crates/resolved/src/render.rs @@ -12,7 +12,11 @@ use nalgebra::Vector3; use crate::bodies::{Appearance, MOON_RADIUS}; use crate::camera::CameraConfig; +use std::sync::Arc; + +use crate::earth::EarthTexture; use crate::eclipse::sun_visible_fraction; +use crate::fidelity::{Fidelity, LightingTier}; use crate::gpu::bodies::{BodiesPass, BodiesUniforms, BodyRadiometry}; use crate::gpu::mesh::{GpuMesh, MeshDraw, MeshPass}; use crate::gpu::profile::{PassTimings, Profiler}; @@ -25,10 +29,14 @@ use crate::scene::{Pose, Scene}; use crate::Error; /// Choices made once, when the renderer is created. -#[derive(Clone, Debug, Default, PartialEq)] +#[derive(Clone, Debug, Default)] pub struct RendererOptions { /// Which device to render on. pub adapter: AdapterChoice, + /// The lighting and camera tiers to draw at. + pub fidelity: Fidelity, + /// Earth's surface, needed from lighting tier L1 up; ignored below it. + pub earth_texture: Option>, /// How Earth and the Moon reflect light. pub appearance: Appearance, /// Also read back the channel radiance image each frame. Costs another 16 bytes per @@ -73,6 +81,7 @@ impl Frame { #[derive(Debug)] pub struct Renderer { ctx: GpuContext, + fidelity: Fidelity, camera: CameraConfig, radiometry: Radiometry, body_radiometry: BodyRadiometry, @@ -117,17 +126,24 @@ impl Renderer { row, ) }); + if options.fidelity.lighting >= LightingTier::L1 && options.earth_texture.is_none() { + return Err(Error::MissingInput(format!( + "lighting {} draws Earth's surface from an image: set RendererOptions::earth_texture", + options.fidelity.lighting + ))); + } let profiler = if options.profile { Profiler::new(&ctx) } else { None }; - Ok(Self { + let mut renderer = Self { bodies: BodiesPass::new(&ctx), meshes: MeshPass::new(&ctx), sensor: SensorPass::new(&ctx), sensor_uniforms: SensorUniforms::new(&camera, &radiometry), ctx, + fidelity: options.fidelity, camera, radiometry, body_radiometry, @@ -136,7 +152,11 @@ impl Renderer { dn_staging, radiance_staging, profiler, - }) + }; + if let Some(texture) = &options.earth_texture { + renderer.set_earth_texture(texture); + } + Ok(renderer) } /// The device this renders on. @@ -149,6 +169,29 @@ impl Renderer { &self.camera } + /// Draws Earth's surface from `texture` from now on, replacing the one given in + /// [`RendererOptions::earth_texture`]: for switching months. Uploads the texture + /// with its mip chain, which for a 16384 x 8192 image is about 1.4 GB of GPU memory + /// and 0.4 s. + /// + /// Below lighting tier L1 there is no Earth texture, and this does nothing. + pub fn set_earth_texture(&mut self, texture: &EarthTexture) { + if self.fidelity.lighting < LightingTier::L1 { + log::warn!( + "lighting {} draws no Earth texture; ignoring it", + self.fidelity.lighting + ); + return; + } + self.bodies.set_earth_texture(&self.ctx, texture); + self.body_radiometry.earth_textured = true; + } + + /// The tiers this renderer draws at. + pub fn fidelity(&self) -> Fidelity { + self.fidelity + } + /// Uploads a model, to be drawn by any [`Target`](crate::Target) that names the /// returned id. pub fn load_model(&mut self, mesh: &Mesh) -> Result { diff --git a/crates/resolved/src/shaders/bodies.wgsl b/crates/resolved/src/shaders/bodies.wgsl index 9a7e7f6..5e8efe4 100644 --- a/crates/resolved/src/shaders/bodies.wgsl +++ b/crates/resolved/src/shaders/bodies.wgsl @@ -8,7 +8,7 @@ struct Uniforms { // fx, fy, cx, cy in pixels. intrinsics: vec4, - // x: near clip distance. + // x: near clip distance. y: 1 if Earth is drawn from the texture. params: vec4, // xyz: Earth center. w: equatorial radius. earth_center: vec4, @@ -18,8 +18,15 @@ struct Uniforms { earth_shape: vec4, // xyz: unit vector from Earth toward the Sun. earth_sun: vec4, - // Per channel: channel radiance of the surface facing the Sun (albedo S_k / (pi d^2)). + // Per channel: channel radiance of the surface facing the Sun (albedo S_k / (pi d^2)), + // for the untextured Earth. earth_lambert: vec4, + // The same with albedo 1, for the textured Earth. + earth_irradiance: vec4, + // Columns: per-channel albedo of texture red, green and blue reflectance. + earth_rgb: mat3x4, + // Per channel: albedo added to the texture's (stand-in for the atmosphere). + earth_haze: vec4, // xyz: Moon center. w: radius. moon_center: vec4, moon_sun: vec4, @@ -31,6 +38,11 @@ struct Uniforms { } @group(0) @binding(0) var u: Uniforms; +// Equirectangular reflectance: longitude -180..180 left to right, latitude +90 at the top. +@group(0) @binding(1) var earth_texture: texture_2d; +@group(0) @binding(2) var earth_sampler: sampler; + +const PI: f32 = 3.14159265358979; struct Out { @location(0) radiance: vec4, @@ -63,6 +75,27 @@ fn sphere_hit(o: vec3, d: vec3, r: f32) -> f32 { return select(MISS, t, t > 0.0); } +// Linear RGB reflectance at the Earth-fixed unit normal `n`. +// +// Longitude wraps at 180 degrees, where the texture coordinate jumps from 1 to 0 and its +// screen-space derivative with it, which would select the coarsest mip along the seam. +// So the coordinate is computed two ways, wrapping at 180 and at 0 degrees, and the one +// whose derivative is smaller here is used: both name the same texel under the repeat +// sampler (Tarini's method). +fn sample_earth(n: vec3) -> vec3 { + let x = atan2(n.y, n.x) / (2.0 * PI) + 0.5; + let v = 0.5 - asin(clamp(n.z, -1.0, 1.0)) / PI; + let ua = fract(x); + let ub = fract(x + 0.5) - 0.5; + let ga = vec2(dpdx(ua), dpdy(ua)); + let gb = vec2(dpdx(ub), dpdy(ub)); + let gv = vec2(dpdx(v), dpdy(v)); + let use_b = dot(gb, gb) < dot(ga, ga); + let uu = select(ua, ub, use_b); + let gu = select(ga, gb, use_b); + return textureSampleGrad(earth_texture, earth_sampler, vec2(uu, v), vec2(gu.x, gv.x), vec2(gu.y, gv.y)).rgb; +} + @fragment fn fs_main(@builtin(position) pos: vec4) -> Out { let k = u.intrinsics; @@ -80,12 +113,20 @@ fn fs_main(@builtin(position) pos: vec4) -> Out { let ds = dir * stretch; let ls = length(ds); let te = sphere_hit(to_cam * stretch, ds / ls, u.earth_center.w) / ls; + // The ellipsoid normal in Earth-fixed axes; its direction gives geodetic latitude. + // Computed for every pixel (misses use the camera's position and ignore the result) + // because texture gradients may only be taken in uniform control flow. + let n_fixed = normalize((to_cam + max(te, 0.0) * dir) * stretch * stretch); + let earth_rgb = sample_earth(n_fixed); if te > 0.0 && te < best { best = te; - let p = to_cam + te * dir; - let n_fixed = normalize(p * stretch * stretch); let n = u.earth_axes * n_fixed; - radiance = u.earth_lambert * max(dot(n, u.earth_sun.xyz), 0.0); + let lit = max(dot(n, u.earth_sun.xyz), 0.0); + if u.params.y > 0.5 { + radiance = u.earth_irradiance * (u.earth_rgb * earth_rgb + u.earth_haze) * lit; + } else { + radiance = u.earth_lambert * lit; + } } let tm = sphere_hit(-u.moon_center.xyz, d, u.moon_center.w); diff --git a/crates/resolved/tests/earth_texture.rs b/crates/resolved/tests/earth_texture.rs new file mode 100644 index 0000000..b98a7db --- /dev/null +++ b/crates/resolved/tests/earth_texture.rs @@ -0,0 +1,308 @@ +//! Earth drawn from a texture: lookup by latitude and longitude, Earth's rotation, the +//! 180 degree seam and minification, on the GPU and on the CPU (lavapipe). + +mod common; + +use resolved::bodies::Appearance; +use resolved::gpu::AdapterChoice; +use resolved::model::rgb_spectrum; +use resolved::nalgebra::{UnitQuaternion, Vector3}; +use resolved::radiometry::{Radiometry, AU}; +use std::sync::Arc; + +use resolved::{ + CameraConfig, EarthTexture, Fidelity, LightingTier, Pose, Renderer, RendererOptions, Scene, + Spectrum, +}; + +const ADAPTERS: [AdapterChoice; 2] = [AdapterChoice::Gpu, AdapterChoice::Cpu]; + +/// The camera's distance from Earth's center: Earth is about 580 px across. +const RANGE: f64 = 2.0e8; + +/// A lighting tier L1 renderer drawing Earth from `texture`. +fn renderer(adapter: &AdapterChoice, haze: bool, texture: &EarthTexture) -> Renderer { + let mut appearance = Appearance::default(); + if !haze { + appearance.earth_haze = Spectrum::constant(0.0); + } + let options = RendererOptions { + adapter: adapter.clone(), + fidelity: Fidelity { + lighting: LightingTier::L1, + ..Default::default() + }, + earth_texture: Some(Arc::new(texture.clone())), + appearance, + read_radiance: true, + ..Default::default() + }; + Renderer::new(CameraConfig::example(), options).unwrap() +} + +/// A unit vector at geocentric latitude and longitude, in degrees. +fn direction(lat: f64, lon: f64) -> Vector3 { + let (lat, lon) = (lat.to_radians(), lon.to_radians()); + Vector3::new(lat.cos() * lon.cos(), lat.cos() * lon.sin(), lat.sin()) +} + +/// Earth turned by `earth_attitude`, the camera over `lat`, `lon` (in inertial axes), +/// the Sun straight behind the camera. +fn view(lat: f64, lon: f64, earth_attitude: UnitQuaternion) -> (Scene, Pose) { + let up = direction(lat, lon); + let scene = Scene { + earth_attitude, + sun: up * AU, + moon: Vector3::new(0.0, 0.0, -3.844e8), + targets: vec![], + }; + let north = if lat.abs() > 80.0 { + Vector3::x() + } else { + Vector3::z() + }; + let camera = Pose::looking_at(up * RANGE, Vector3::zeros(), north).unwrap(); + (scene, camera) +} + +fn linear(byte: u8) -> f64 { + let v = byte as f64 / 255.0; + if v <= 0.040_45 { + v / 12.92 + } else { + ((v + 0.055) / 1.055).powf(2.4) + } +} + +/// Channel radiance of Earth's surface with sRGB-encoded reflectance `rgb`, facing the +/// Sun at 1 AU, without haze. +fn expected(rgb: [u8; 3]) -> [f64; 3] { + let r = Radiometry::new(&CameraConfig::example()); + let albedo = r.effective_albedo(&rgb_spectrum([rgb[0], rgb[1], rgb[2]].map(linear))); + [0, 1, 2].map(|k| r.lambert_radiance(k, albedo[k], 1.0, AU)) +} + +/// Eight 90-degree regions: columns are longitude -180..-90, -90..0, 0..90, 90..180; +/// rows are the northern then the southern hemisphere. +const REGIONS: [[u8; 3]; 8] = [ + [200, 40, 40], + [40, 200, 40], + [40, 40, 200], + [200, 200, 40], + [40, 200, 200], + [200, 40, 200], + [120, 120, 120], + [250, 250, 250], +]; + +/// The eight regions at 16 x 16 texels each, so a sample near a region's center has +/// only that region's color around it. (With one texel per region, the half-pixel +/// between the boresight and the center pixel's line of sight is enough for bilinear +/// filtering to mix in 0.1 % of a neighbor, a few percent in a dark channel.) +fn regions() -> EarthTexture { + let (w, h) = (64u32, 32u32); + let pixels = (0..h) + .flat_map(|y| (0..w).map(move |x| REGIONS[(y / 16 * 4 + x / 16) as usize])) + .flat_map(|[r, g, b]| [r, g, b, 255]) + .collect(); + EarthTexture::from_rgba8(w, h, pixels).unwrap() +} + +fn center(frame: &resolved::Frame) -> [f32; 4] { + frame.radiance.as_ref().unwrap()[540 * 1920 + 960] +} + +fn assert_close(got: [f32; 4], want: [f64; 3], context: &str) { + for c in 0..3 { + let error = got[c] as f64 / want[c] - 1.0; + assert!( + error.abs() < 1e-3, + "{context} channel {c}: {} vs {} ({:+.4} %)", + got[c], + want[c], + 100.0 * error + ); + } +} + +#[test] +fn texel_under_the_camera_sets_the_radiance() { + // Over 45 N 45 E the camera sees the middle of region 2 (longitude 0..90, north). + let _devices = common::device_lock(); + for adapter in &ADAPTERS { + let mut renderer = renderer(adapter, false, ®ions()); + let (scene, camera) = view(45.0, 45.0, UnitQuaternion::identity()); + let frame = renderer.submit(&scene, &camera).unwrap().wait().unwrap(); + assert_close( + center(&frame), + expected(REGIONS[2]), + &renderer.adapter_info().name, + ); + } +} + +#[test] +fn earth_rotation_moves_the_surface() { + // Turning Earth 90 degrees east about its axis puts longitude -45 under a camera + // fixed in inertial space at 45 E: region 1. And the south: region 5. + let _devices = common::device_lock(); + let turn = UnitQuaternion::from_axis_angle(&Vector3::z_axis(), 90f64.to_radians()); + for adapter in &ADAPTERS { + let mut renderer = renderer(adapter, false, ®ions()); + for (lat, region) in [(45.0, 1), (-45.0, 5)] { + let (scene, camera) = view(lat, 45.0, turn); + let frame = renderer.submit(&scene, &camera).unwrap().wait().unwrap(); + let name = &renderer.adapter_info().name; + assert_close( + center(&frame), + expected(REGIONS[region]), + &format!("{name} lat {lat}"), + ); + } + } +} + +#[test] +fn haze_adds_its_spectrum() { + let _devices = common::device_lock(); + let r = Radiometry::new(&CameraConfig::example()); + let haze = r.effective_albedo(&Appearance::default().earth_haze); + let surface = r.effective_albedo(&rgb_spectrum(REGIONS[2].map(linear))); + let want = [0, 1, 2].map(|k| r.lambert_radiance(k, surface[k] + haze[k], 1.0, AU)); + // Bluer haze: it raises blue the most. + assert!(haze[2] > haze[1] && haze[1] > haze[0], "{haze:?}"); + for adapter in &ADAPTERS { + let mut renderer = renderer(adapter, true, ®ions()); + let (scene, camera) = view(45.0, 45.0, UnitQuaternion::identity()); + let frame = renderer.submit(&scene, &camera).unwrap().wait().unwrap(); + assert_close(center(&frame), want, &renderer.adapter_info().name); + } +} + +#[test] +fn no_seam_at_180_degrees() { + // The four texel columns either side of 180 degrees (22.5 degrees each way) are one + // color and everything else another, so a seam that fell back to a coarse mip would + // show the average. Pixels near the image center, over 180 degrees, all read the + // edge color. + let _devices = common::device_lock(); + let (w, h) = (64u32, 32u32); + let (edge, middle) = ([30u8, 60, 200], [220u8, 220, 220]); + let pixels = (0..h) + .flat_map(|_| { + (0..w).map(|x| { + if !(4..w - 4).contains(&x) { + edge + } else { + middle + } + }) + }) + .flat_map(|[r, g, b]| [r, g, b, 255]) + .collect(); + let texture = EarthTexture::from_rgba8(w, h, pixels).unwrap(); + let want = expected(edge); + for adapter in &ADAPTERS { + let mut renderer = renderer(adapter, false, &texture); + // Over 179.8 E the seam is about one pixel right of center, inside a 2 x 2 pixel + // quad (screen-space derivatives are taken per quad). Over exactly 180 it would + // fall on a quad boundary, where even unhandled it shows nothing. + let (scene, camera) = view(0.0, 179.8, UnitQuaternion::identity()); + let frame = renderer.submit(&scene, &camera).unwrap().wait().unwrap(); + let radiance = frame.radiance.unwrap(); + let name = &renderer.adapter_info().name; + // 20 px from the disk center is 4 degrees of longitude, well inside the edge + // columns' texel centers. + for v in (440..640).step_by(10) { + for u in 940..980 { + let got = radiance[v * 1920 + u]; + // Shading falls off away from the sub-solar point; compare the color + // (channel ratios), which the texture alone sets. + let ratio = |a: usize, b: usize| got[a] as f64 / got[b] as f64; + let want_ratio = |a: usize, b: usize| want[a] / want[b]; + for (a, b) in [(0, 2), (1, 2)] { + let error = ratio(a, b) / want_ratio(a, b) - 1.0; + assert!( + error.abs() < 0.01, + "{name} pixel ({u}, {v}): ratio {a}/{b} off by {:+.2} %", + 100.0 * error + ); + } + } + } + } +} + +#[test] +fn far_views_average_rather_than_alias() { + // A checker of black and white texels, each 0.18 degrees (20 km), seen from 400,000 + // km where a pixel spans 44 km: mipmapping averages the checker to 50 % reflectance. + // Point sampling would make noise of 0 % and 100 %. + let _devices = common::device_lock(); + let (w, h) = (2048u32, 1024u32); + let pixels = (0..h) + .flat_map(|y| (0..w).map(move |x| if (x + y) % 2 == 0 { 0u8 } else { 255 })) + .flat_map(|v| [v, v, v, 255]) + .collect(); + let texture = EarthTexture::from_rgba8(w, h, pixels).unwrap(); + // 50 % linear reflectance is sRGB 188. + let gray = EarthTexture::from_rgba8(1, 1, vec![188, 188, 188, 255]).unwrap(); + for adapter in &ADAPTERS { + let mut renderer = renderer(adapter, false, &gray); + let (scene, mut camera) = view(10.0, 30.0, UnitQuaternion::identity()); + camera.position *= 2.0; + let mut block = |texture: &EarthTexture| { + renderer.set_earth_texture(texture); + let frame = renderer.submit(&scene, &camera).unwrap().wait().unwrap(); + let radiance = frame.radiance.unwrap(); + (520..560) + .flat_map(|v| (940..980).map(move |u| (u, v))) + .map(|(u, v)| radiance[v * 1920 + u][1] as f64) + .collect::>() + }; + let checker = block(&texture); + let reference = block(&gray); + let name = renderer.adapter_info().name.clone(); + for (c, g) in checker.iter().zip(&reference) { + let error = c / g - 1.0; + assert!( + error.abs() < 0.03, + "{name}: {c} vs {g} ({:+.2} %)", + 100.0 * error + ); + } + } +} + +#[test] +fn lighting_l1_needs_an_earth_texture() { + let _devices = common::device_lock(); + let options = RendererOptions { + fidelity: Fidelity { + lighting: LightingTier::L1, + ..Default::default() + }, + ..Default::default() + }; + let error = Renderer::new(CameraConfig::example(), options) + .unwrap_err() + .to_string(); + assert!(error.contains("earth_texture"), "{error}"); +} + +#[test] +fn lighting_l0_ignores_an_earth_texture() { + // L0 draws Earth from Appearance::earth_albedo whether or not a texture is given. + let _devices = common::device_lock(); + let (scene, camera) = view(45.0, 45.0, UnitQuaternion::identity()); + let center_at = |earth_texture: Option>| { + let options = RendererOptions { + earth_texture, + read_radiance: true, + ..Default::default() + }; + let mut renderer = Renderer::new(CameraConfig::example(), options).unwrap(); + center(&renderer.submit(&scene, &camera).unwrap().wait().unwrap()) + }; + assert_eq!(center_at(Some(Arc::new(regions()))), center_at(None)); +} diff --git a/docs/images/l0-c0.png b/docs/images/l0-c0.png index f040e39..441994f 100644 Binary files a/docs/images/l0-c0.png and b/docs/images/l0-c0.png differ diff --git a/docs/optics.md b/docs/optics.md index 1c7f561..2cb72ca 100644 --- a/docs/optics.md +++ b/docs/optics.md @@ -459,6 +459,64 @@ lit) right. It misses some things that are visible in real images: - Earth's atmosphere softens the limb and adds a thin bright band around it (M4). Here Earth's edge is a sharp line. +### Earth's surface + +Code: `earth.rs`, `shaders/bodies.wgsl`; data from `scripts/fetch-earth.sh`. + +At lighting tier L1 (`RendererOptions::fidelity`, with the image in +`RendererOptions::earth_texture`; `--lighting L1` on the command line), Earth's +reflectance comes from an image of the whole surface instead of one spectrum. The image is NASA's +Blue Marble Next Generation: a cloud-free mosaic of the 2004 MODIS satellite +record, one image per month, so snow cover and vegetation follow the season. +The CLI picks the month from the scene's epoch. + +The image is equirectangular: each column is a longitude and each row a +latitude. For every pixel that sees Earth, the shader finds the surface +point's geodetic latitude and longitude from the ellipsoid's normal there, in +Earth-fixed axes, and reads the image at that place. Earth's orientation from +the scene turns the continents with the epoch. + +The full image is 500 m per pixel at the equator, resized to 660 m to fit the +largest texture lavapipe allows (16384 x 8192). Whether that is fine or coarse +depends on the view: + +- From low orbit looking straight down, a camera pixel of the example camera + covers 77 m of ground, so each texture pixel spans about 8 image pixels and + looks soft. +- From geostationary orbit, a camera pixel covers about 4 km, so each image + pixel averages about 35 texture pixels (6 x 6). + +For the second case the texture is stored with mipmaps: copies at half, a +quarter, an eighth of the size, and so on. The GPU reads from the copy whose +pixels best match the image pixel, so a far view shows the true average color +of a region rather than a shimmering sample of individual texture pixels. + +Two details matter for brightness: + +- Averaging has to happen on reflectance, not on the encoded values image files + store (sRGB, which compresses the darks). Averaging encoded black and white + gives a reflectance of 21 %, not 50 %. So the texture is decoded to linear + reflectance once at load and kept in 16-bit floating point, and both the + mipmaps and the GPU's blending between neighboring texture pixels work on + those values. A coastline seen up close blends correctly from sea to land. +- Blue Marble is surface reflectance with no atmosphere, and open ocean is + dark, a few percent. From orbit the ocean looks bright blue because the air + above it scatters sunlight toward the camera: Rayleigh scattering, which + grows as 1 / wavelength^4, which is also why the sky is blue. Until the + atmosphere is modeled (M4), a stand-in adds a reflectance of that shape, + 0.07 at 450 nm falling to 0.017 at 650 nm (`Appearance::earth_haze`). It + brightens and blues the whole lit disk evenly, where the real effect grows + toward the limb, where the line of sight crosses more air. + +The texture's red, green and blue are reflectances in the same three bands a +glTF base color uses (see the target model section), folded into the camera's +channels once at setup as a small matrix, so any channel layout works. + +Blue Marble is a visual product built from calibrated MODIS reflectances, not a +calibrated reflectance product itself, so its values are close but not exact +reflectances. A calibrated alternative (MODIS's albedo product) is a possible +later refinement. + ### Depth from the bodies The pass also writes each hit's depth, near / z, as described in the @@ -591,7 +649,8 @@ velocity as an input, which the scene does not carry yet. | Noise: photon shot noise, read noise, dark current, pixel-to-pixel variation | M3 | | Saturation spilling into neighbors (blooming) | M3 | | Non-Lambertian materials (BRDFs), glints, lunar backscatter | M4 | -| Earth and Moon albedo maps | M2 | +| Moon albedo map | M2 | +| Clouds on Earth | M4 | | Earth's atmosphere: scattering, the bright limb | M4 | | Motion blur, rolling shutter, stray light | M6 | | Apparent positions: light time and aberration | M3 | diff --git a/examples/scenes/ace.toml b/examples/scenes/ace.toml index f455d66..448b7ff 100644 --- a/examples/scenes/ace.toml +++ b/examples/scenes/ace.toml @@ -4,15 +4,22 @@ # comes from scripts/fetch-models.sh. # # resolved render --scene examples/scenes/ace.toml --preview ace.png +# resolved render --scene examples/scenes/ace.toml --lighting L1 --preview ace.png -sun = [111172883513, -100100513931, 0] +sun = [145154174600, -36191000457, 0] moon = [0, 0, -3.844e8] +# Earth-fixed to inertial: 25 N 33 E behind the spacecraft, looking north. +earth_attitude = [0.660019964, -0.691724929, -0.146685455, -0.253719624] [observer] position = [7000036.38, -73.51, 0.00] look_at = [7000000.00, 0.00, 0.00] up = [1, 0, 0] +# At lighting L1, Earth's surface shows July (from scripts/fetch-earth.sh). +[appearance] +earth_month = 7 + [[target]] model = "../models/nasa/ace.glb" scale = 0.00485 diff --git a/examples/scenes/boxwing-epoch.toml b/examples/scenes/boxwing-epoch.toml index 3db8968..fecc3f8 100644 --- a/examples/scenes/boxwing-epoch.toml +++ b/examples/scenes/boxwing-epoch.toml @@ -4,6 +4,8 @@ # this epoch: same limb, same isometric view, same lighting. # # resolved render --scene examples/scenes/boxwing-epoch.toml --preview epoch.png +# +# At lighting L1 (--lighting L1), Earth's surface shows the epoch's month. epoch = "2026-10-14T12:00:00 UTC" diff --git a/examples/scenes/icesat2.toml b/examples/scenes/icesat2.toml index a087b60..0d0b0f5 100644 --- a/examples/scenes/icesat2.toml +++ b/examples/scenes/icesat2.toml @@ -3,15 +3,22 @@ # scripts/make-example-scenes.py; the model comes from scripts/fetch-models.sh. # # resolved render --scene examples/scenes/icesat2.toml --preview icesat2.png +# resolved render --scene examples/scenes/icesat2.toml --lighting L1 --preview icesat2.png -sun = [105781668823, -105781668823, 0] +sun = [143802702898, -41234761516, 0] moon = [0, 0, -3.844e8] +# Earth-fixed to inertial: 25 N 33 E behind the spacecraft, looking north. +earth_attitude = [0.655810260, -0.694017033, -0.135426580, -0.264410481] [observer] position = [6874182.69, -103.11, 0.00] look_at = [6874137.00, 0.00, 0.00] up = [1, 0, 0] +# At lighting L1, Earth's surface shows July (from scripts/fetch-earth.sh). +[appearance] +earth_month = 7 + [[target]] model = "../models/nasa/icesat2.glb" scale = 1.0 diff --git a/examples/scenes/landsat8.toml b/examples/scenes/landsat8.toml index 22fcbd8..e389243 100644 --- a/examples/scenes/landsat8.toml +++ b/examples/scenes/landsat8.toml @@ -3,15 +3,22 @@ # comes from scripts/fetch-models.sh. # # resolved render --scene examples/scenes/landsat8.toml --preview landsat8.png +# resolved render --scene examples/scenes/landsat8.toml --lighting L1 --preview landsat8.png -sun = [111172883513, -100100513931, 0] +sun = [145154174600, -36191000457, 0] moon = [0, 0, -3.844e8] +# Earth-fixed to inertial: 25 N 33 E behind the spacecraft, looking north. +earth_attitude = [0.662448702, -0.690270337, -0.153384684, -0.247308951] [observer] position = [7083216.12, -150.19, 0.00] look_at = [7083137.00, 0.00, 0.00] up = [1, 0, 0] +# At lighting L1, Earth's surface shows July (from scripts/fetch-earth.sh). +[appearance] +earth_month = 7 + [[target]] model = "../models/nasa/landsat8.glb" scale = 1.0 diff --git a/examples/scenes/ostm.toml b/examples/scenes/ostm.toml index d4434b6..8d7303b 100644 --- a/examples/scenes/ostm.toml +++ b/examples/scenes/ostm.toml @@ -3,15 +3,22 @@ # comes from scripts/fetch-models.sh. # # resolved render --scene examples/scenes/ostm.toml --preview ostm.png +# resolved render --scene examples/scenes/ostm.toml --lighting L1 --preview ostm.png -sun = [125463331183, -81476839938, 0] +sun = [148141994506, -20819999562, 0] moon = [0, 0, -3.844e8] +# Earth-fixed to inertial: 25 N 33 E behind the spacecraft, looking north. +earth_attitude = [0.675405688, -0.680462708, -0.192277152, -0.209349365] [observer] position = [7714190.15, -72.55, 0.00] look_at = [7714137.00, 0.00, 0.00] up = [1, 0, 0] +# At lighting L1, Earth's surface shows July (from scripts/fetch-earth.sh). +[appearance] +earth_month = 7 + [[target]] model = "../models/nasa/ostm.glb" scale = 3.5 diff --git a/examples/scenes/ssl1300.toml b/examples/scenes/ssl1300.toml index f7c22ce..c821de1 100644 --- a/examples/scenes/ssl1300.toml +++ b/examples/scenes/ssl1300.toml @@ -3,15 +3,22 @@ # scripts/make-example-scenes.py; the model comes from scripts/fetch-models.sh. # # resolved render --scene examples/scenes/ssl1300.toml --preview ssl1300.png +# resolved render --scene examples/scenes/ssl1300.toml --lighting L1 --preview ssl1300.png -sun = [144500446867, 38718778044, 0] +sun = [146328798300, 31103146240, 0] moon = [0, 0, -3.844e8] +# Earth-fixed to inertial: 25 N 33 E behind the spacecraft, looking north. +earth_attitude = [0.701664364, -0.636154985, -0.308718051, -0.087562094] [observer] position = [42164291.13, -34.20, 0.00] look_at = [42164000.00, 0.00, 0.00] up = [1, 0, 0] +# At lighting L1, Earth's surface shows July (from scripts/fetch-earth.sh). +[appearance] +earth_month = 7 + [[target]] model = "../models/nasa/ssl1300.glb" scale = 0.1 diff --git a/readme.md b/readme.md index 64f356d..ff4a629 100644 --- a/readme.md +++ b/readme.md @@ -21,9 +21,12 @@ on the CPU through Mesa's lavapipe (`vulkan-swrast` on Arch, ## Performance and fidelity Lighting and camera fidelity are chosen separately: lighting tiers L0 to L3, -camera tiers C0 to C3. A row is added here as each tier lands, with its -measured cost and a sample image from the same scene, so the cost of each step -up is visible at a glance. +camera tiers C0 to C3, set with `RendererOptions::fidelity` in the library and +`--lighting` and `--camera` on the command line. A tier is a fixed bundle of +features, and a tier that needs data (L1 needs an image of Earth's surface) +refuses to start without it. A row is added here as each tier is complete, with +its measured cost and a sample image from the same scene, so the cost of each +step up is visible at a glance. | Tier | Lighting | Camera | RTX 3060 | lavapipe (CPU) | Sample | |---|---|---|---|---|---| @@ -95,6 +98,10 @@ in the repository, and the rest use NASA spacecraft models (Landsat 8, ICESat-2, OSTM/Jason-2, an SSL-1300 and ACE) that `scripts/fetch-models.sh` downloads into `examples/models/nasa/`. That script needs Node (for `npx`), because most NASA models are Draco-compressed and are decoded on the way in. +At lighting L1, Earth's surface comes from NASA's Blue Marble, which +`scripts/fetch-earth.sh` downloads and resizes (with ImageMagick) into +`data/earth/`, one image per month. The month comes from the scene's epoch, or +from `[appearance] earth_month` in a scene without one. ```sh cargo run --release -p resolved-cli -- render \ @@ -108,6 +115,12 @@ cargo run --release -p resolved-cli -- render \ - `--preview` writes an 8-bit sRGB image for looking at, with the bias as black and saturation as white. - `--timings N` renders N frames and reports setup and per-frame times. +- `--lighting L0|L1` and `--camera C0` pick the fidelity tiers (default L0 + and C0). A tier the scene or the data cannot support stops with an error + naming what is missing. +- `--data DIR` says where the tiers' data is (`earth/bmng-2004-MM.jpg`). + Without it the CLI tries `RESOLVED_DATA`, then `./data`, then the + repository's `data/`. - `--adapter` takes `auto`, `gpu`, `cpu` or part of an adapter's name; `resolved adapters` lists them. - `--kernels DIR` says where astrokit's kernels are, for scenes given by diff --git a/scripts/fetch-earth.sh b/scripts/fetch-earth.sh new file mode 100755 index 0000000..d0a67a9 --- /dev/null +++ b/scripts/fetch-earth.sh @@ -0,0 +1,53 @@ +#!/bin/sh +# Downloads NASA's Blue Marble Next Generation monthly Earth images and resizes +# them for the renderer. +# +# Source: Blue Marble: Next Generation (2004), base map without topography or +# bathymetry, 21600 x 10800 (500 m per pixel), one image per month, from +# https://science.nasa.gov/earth/earth-observatory/blue-marble-next-generation/. +# Produced by Reto Stöckli, NASA Earth Observatory; NASA asks that users credit +# "NASA Earth Observatory". +# +# The originals go to data/earth/source/ and are checked by sha256. Each is +# resized to 16384 x 8192 (660 m per pixel), the largest texture lavapipe takes, +# into data/earth/bmng-2004-MM.jpg. Resizing uses ImageMagick (`magick`). +# Both directories are gitignored. + +set -eu + +base='https://assets.science.nasa.gov/content/dam/science/esd/eo/images/bmng/bmng-base' +root="$(cd "$(dirname "$0")/.." && pwd)/data/earth" +mkdir -p "$root/source" + +# month number, folder name, sha256 of the original +while read -r n month sum; do + original="$root/source/world.2004$n.3x21600x10800.jpg" + if [ ! -f "$original" ]; then + echo "fetching $month" + curl -fsSL "$base/$month/world.2004$n.3x21600x10800.jpg" -o "$original.part" + mv "$original.part" "$original" + fi + echo "$sum $original" | sha256sum -c --quiet - + resized="$root/bmng-2004-$n.jpg" + if [ ! -f "$resized" ]; then + echo "resizing $month" + # Resample in linear light, then back to sRGB, so averaging is of reflectance. + magick "$original" -colorspace RGB -filter Lanczos -resize 16384x8192! \ + -colorspace sRGB -quality 95 "$resized" + fi +done <<'EOF' +01 january c2f2c80a352518cfea23d462dad6f6e778375ddd74bbec5786bcfc85cc49c61d +02 february 8800081f38b53bfa1e3f0c2fa17be694b7646942fe5add009ea0dcd64c03e70b +03 march af37f61a9715c59dc86763df19f588a20c9d92b3d0a8034a0b00acd6b97ccb7d +04 april 3bd13f173555eb24451c45d3a83ca177baca02e9a8bdd9fa3e4989f7b496ce6b +05 may 0de40c868680ab9e86fef42406d76386b243b32444ff50edfd406b296638c206 +06 june 282608580a290723655bfbe917e69ea2f93f4da068dfa96ff4a152b1c6ec36be +07 july dea8b4dc8a4f93f5f8bce0c8c85a508a178e7901e9ed8e6bf86e6ce7ef6d61e2 +08 august 926be7253be90a1fd57071a33afbfe7a712a2dfe0d5b5b4602bc5dfbcfdc2c89 +09 september 7cf788e13a3a7b4a926b524f8f71d635c56a18c48ba0bfe35695b05a305db3f3 +10 october bf568ceca2ad55e9c8dd7bcedfd0e419201d03bb560206a7bfa2cc57c1e04c64 +11 november 1959bbe72897d3e9526c440dae3ec04cfdbdc4080c8ff0b14632925ae31167b1 +12 december b0ffc9d12205689aeaec66116bbae2c4c263fcd98a0cbadba8d5488567b53278 +EOF + +echo "Earth images in $root" diff --git a/scripts/make-example-scenes.py b/scripts/make-example-scenes.py index c14a69f..a0b2c35 100755 --- a/scripts/make-example-scenes.py +++ b/scripts/make-example-scenes.py @@ -5,9 +5,15 @@ Every scene has the same composition: the spacecraft in an isometric view (the camera looks along its body (+1, +1, +1) diagonal, body +Z toward image up) with Earth's limb crossing the image 2 degrees above center, Earth below and space above. The camera distance makes the model's bounding sphere span about 5.6 -degrees, half the example camera's width. The Sun direction is the one, on a -3-degree grid, that best lights the three visible body faces and Earth near the -limb, without putting the spacecraft in Earth's shadow. +degrees, half the example camera's width. + +Earth is turned so the ground behind the spacecraft, at the image center, is +25 N 33 E (Egypt, the Nile and the Red Sea), with the camera looking north +toward the Mediterranean on the limb. The Sun is the direction, on a 2-degree +grid, that puts it at least 58 degrees above that ground (cos i >= 0.85) and +best lights the spacecraft's three visible faces, without putting the +spacecraft in Earth's shadow. With the Sun that high over the ground, one face +is lit at a slant. Model centers and radii were measured from the decoded files at the listed scales (the bounding box of every vertex); remeasure if a model changes. @@ -18,6 +24,10 @@ import pathlib import textwrap RE = 6378137.0 +# The ground at the image center, geocentric latitude and longitude in degrees. +GROUND = (25.0, 33.0) +# The least cos(incidence) the Sun may have on that ground. +GROUND_SUN = 0.85 AU = 1.495978707e11 LEO = 7.0e6 GEO = 42_164_000.0 @@ -105,20 +115,29 @@ def scene(name, file, scale, center, radius, orbit, description): offset = mul(attitude, center) position = [look[i] - offset[i] for i in range(3)] - to_limb = [-math.cos(limb), math.sin(limb), 0.0] - reach = math.sqrt(orbit * orbit - RE * RE) - limb_normal = norm([look[i] + reach * to_limb[i] for i in range(3)]) + # The ground behind the spacecraft at the image center: where the boresight, + # continued past the spacecraft, meets Earth. + b = dot(look, boresight) + reach = -b - math.sqrt(b * b - (dot(look, look) - RE * RE)) + ground = norm([look[i] + reach * boresight[i] for i in range(3)]) + # Turn Earth so GROUND is there, with north along the camera's view across it. + lat, lon = (math.radians(x) for x in GROUND) + fixed = [math.cos(lat) * math.cos(lon), math.cos(lat) * math.sin(lon), math.sin(lat)] + north = norm([boresight[i] - dot(boresight, ground) * ground[i] for i in range(3)]) + earth = triad(fixed, [0, 0, 1], ground, north) + best = None - for az in range(0, 360, 3): - for el in range(-87, 90, 3): + for az in range(0, 360, 2): + for el in range(-88, 90, 2): a, e = math.radians(az), math.radians(el) sun = [math.cos(e) * math.cos(a), math.cos(e) * math.sin(a), math.sin(e)] if dot(sun, [-1.0, 0.0, 0.0]) > math.cos(limb): continue # behind Earth: the spacecraft would be in shadow + if dot(ground, sun) < GROUND_SUN: + continue faces = [dot(mul(attitude, n), sun) for n in ([1, 0, 0], [0, 1, 0], [0, 0, 1])] - score = min(min(faces), dot(limb_normal, sun)) - if best is None or score > best[0]: - best = (score, sun) + if best is None or min(faces) > best[0]: + best = (min(faces), sun) sun = best[1] fmt = lambda v, d=2: "[" + ", ".join(f"{x:.{d}f}" for x in v) + "]" @@ -131,15 +150,22 @@ def scene(name, file, scale, center, radius, orbit, description): return f"""{header} # # resolved render --scene examples/scenes/{name}.toml --preview {name}.png +# resolved render --scene examples/scenes/{name}.toml --lighting L1 --preview {name}.png sun = {fmt([AU * x for x in sun], 0)} moon = [0, 0, -3.844e8] +# Earth-fixed to inertial: 25 N 33 E behind the spacecraft, looking north. +earth_attitude = {fmt(quaternion(earth), 9)} [observer] position = {fmt(camera)} look_at = {fmt(look)} up = [1, 0, 0] +# At lighting L1, Earth's surface shows July (from scripts/fetch-earth.sh). +[appearance] +earth_month = 7 + [[target]] model = "../models/nasa/{file}" scale = {scale}