From 85896b17a49d3e459d75e6f84056522aed82fa92 Mon Sep 17 00:00:00 2001 From: Connor Johnstone Date: Wed, 14 Oct 2026 21:56:11 -0400 Subject: [PATCH] Earth and moon eclipse --- crates/resolved-cli/Cargo.toml | 4 + crates/resolved-cli/src/main.rs | 4 +- crates/resolved-cli/tests/cli.rs | 74 ++++- crates/resolved/src/fidelity.rs | 7 +- crates/resolved/src/gpu/bodies.rs | 5 +- crates/resolved/src/gpu/earthshine.rs | 4 + crates/resolved/src/render.rs | 3 +- crates/resolved/src/shaders/bodies.wgsl | 14 +- crates/resolved/src/shaders/common.wgsl | 30 ++ crates/resolved/src/shaders/earthshine.wgsl | 12 +- crates/resolved/tests/earthshine.rs | 35 ++- crates/resolved/tests/eclipses.rs | 286 ++++++++++++++++++++ docs/optics.md | 75 ++++- examples/scenes/lunar-eclipse.toml | 19 ++ examples/scenes/solar-eclipse.toml | 16 ++ readme.md | 11 +- 16 files changed, 575 insertions(+), 24 deletions(-) create mode 100644 crates/resolved/tests/eclipses.rs create mode 100644 examples/scenes/lunar-eclipse.toml create mode 100644 examples/scenes/solar-eclipse.toml diff --git a/crates/resolved-cli/Cargo.toml b/crates/resolved-cli/Cargo.toml index fe30369..2ae756b 100644 --- a/crates/resolved-cli/Cargo.toml +++ b/crates/resolved-cli/Cargo.toml @@ -27,3 +27,7 @@ serde = { version = "1.0", features = ["derive"] } toml = "1.1" image = { version = "0.25", default-features = false, features = ["png"] } hifitime = { version = "4.3", optional = true } + +[dev-dependencies] +# JPEG too, so the tests can write a small stand-in for the Blue Marble images. +image = { version = "0.25", default-features = false, features = ["png", "jpeg"] } diff --git a/crates/resolved-cli/src/main.rs b/crates/resolved-cli/src/main.rs index 7425b1f..a3b4c4f 100644 --- a/crates/resolved-cli/src/main.rs +++ b/crates/resolved-cli/src/main.rs @@ -45,8 +45,8 @@ struct RenderArgs { /// 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, shadows on targets and - /// light from Earth and the Moon; needs Blue Marble data). + /// Lighting fidelity tier: L0, or L1 (adds Earth's surface, shadows on targets, + /// light from Earth and the Moon, and eclipses on them; needs Blue Marble data). #[arg(long, default_value = "L0")] lighting: LightingTier, /// Camera fidelity tier: C0. diff --git a/crates/resolved-cli/tests/cli.rs b/crates/resolved-cli/tests/cli.rs index 88e0843..ce71b3b 100644 --- a/crates/resolved-cli/tests/cli.rs +++ b/crates/resolved-cli/tests/cli.rs @@ -27,8 +27,37 @@ fn resolved(args: &[&str]) -> std::process::Output { /// Renders `scene` on `adapter` to a 16-bit PNG and returns its pixels as DN (the /// file holds 12-bit DN shifted up by 4). fn render(scene: &str, adapter: &str) -> image::ImageBuffer, Vec> { + render_at(scene, adapter, "L0") +} + +/// A data directory whose Earth images are small and uniformly gray, for tests at +/// lighting L1 that do not need Blue Marble: loading the real 16384 x 8192 images takes +/// seconds in a debug build, and they are only there once fetched. +fn gray_earth() -> PathBuf { + static DIR: std::sync::OnceLock = std::sync::OnceLock::new(); + DIR.get_or_init(|| { + let dir = out_dir().join("gray-data"); + std::fs::create_dir_all(dir.join("earth")).unwrap(); + let image = image::RgbImage::from_pixel(64, 32, image::Rgb([128; 3])); + for month in 1..=12 { + image + .save(dir.join(format!("earth/bmng-2004-{month:02}.jpg"))) + .unwrap(); + } + dir + }) + .clone() +} + +/// [`render`] at lighting tier `lighting`, with Earth's images from [`gray_earth`]. +fn render_at( + scene: &str, + adapter: &str, + lighting: &str, +) -> image::ImageBuffer, Vec> { + let data = gray_earth(); let out = out_dir().join(format!( - "{}-{adapter}.png", + "{}-{lighting}-{adapter}.png", Path::new(scene).file_stem().unwrap().to_str().unwrap() )); let preview = out.with_extension("preview.png"); @@ -38,6 +67,10 @@ fn render(scene: &str, adapter: &str) -> image::ImageBuffer, Vec scene, "--adapter", adapter, + "--lighting", + lighting, + "--data", + data.to_str().unwrap(), "--out", out.to_str().unwrap(), "--preview", @@ -92,6 +125,45 @@ fn epoch_scene_renders() { assert_eq!(frame.get_pixel(0, 0).0, [64; 3]); } +#[cfg(feature = "astrokit")] +#[test] +fn lunar_eclipse_darkens_the_moon() { + // Mid-totality: at L1 the whole Moon is in Earth's umbra and reads the bias. At L0 + // there are no eclipses on the bodies and it is a full Moon. + let scene = "examples/scenes/lunar-eclipse.toml"; + let l0 = render_at(scene, "gpu", "L0"); + let l1 = render_at(scene, "gpu", "L1"); + assert!( + l0.get_pixel(960, 540).0.iter().all(|&v| v > 150), + "{:?}", + l0.get_pixel(960, 540) + ); + assert!( + l1.pixels().all(|p| p.0 == [64; 3]), + "every pixel at the bias" + ); +} + +#[cfg(feature = "astrokit")] +#[test] +fn solar_eclipse_puts_the_umbra_at_greatest_eclipse() { + // The camera looks straight down at the point of greatest eclipse from NASA's + // bulletin (65.2 N 25.2 W, 17:46 UTC): at L1 it is in the Moon's umbra and reads the + // bias, and the image corner, about 3500 km away in the penumbra, is brighter. At L0 + // the center is lit. + let scene = "examples/scenes/solar-eclipse.toml"; + let l0 = render_at(scene, "gpu", "L0"); + let l1 = render_at(scene, "gpu", "L1"); + assert!( + l0.get_pixel(960, 540).0.iter().all(|&v| v > 100), + "{:?}", + l0.get_pixel(960, 540) + ); + assert_eq!(l1.get_pixel(960, 540).0, [64; 3]); + let corner = l1.get_pixel(40, 1040).0; + assert!(corner.iter().all(|&v| v > 70), "{corner:?}"); +} + #[test] fn errors_are_reported() { let missing = resolved(&["render", "--scene", "no/such.toml", "--preview", "x.png"]); diff --git a/crates/resolved/src/fidelity.rs b/crates/resolved/src/fidelity.rs index b4c0e16..7836166 100644 --- a/crates/resolved/src/fidelity.rs +++ b/crates/resolved/src/fidelity.rs @@ -21,9 +21,10 @@ pub enum LightingTier { /// L0, plus Earth's surface from an image (`RendererOptions::earth_texture`, NASA's /// Blue Marble) with a stand-in for the atmosphere's blue; shadows that targets /// cast on themselves and on each other, with the penumbra the Sun's disk gives - /// them; and sunlight reflected onto targets by Earth and the Moon (not shadowed by - /// the targets). The rest of M2 (the Moon's albedo map, limb darkening, eclipses on - /// the bodies) joins this tier as it lands. + /// them; sunlight reflected onto targets by Earth and the Moon (not shadowed by the + /// targets); and Earth and the Moon shadowing each other, lunar and solar eclipses. + /// The rest of M2 (the Moon's albedo map, limb darkening) joins this tier as it + /// lands. L1, } diff --git a/crates/resolved/src/gpu/bodies.rs b/crates/resolved/src/gpu/bodies.rs index 635aac1..a29757a 100644 --- a/crates/resolved/src/gpu/bodies.rs +++ b/crates/resolved/src/gpu/bodies.rs @@ -35,6 +35,8 @@ pub struct BodyRadiometry { pub earth_haze: [f32; 4], /// Whether Earth is drawn from a texture. pub earth_textured: bool, + /// Whether Earth and the Moon shadow each other: lunar and solar eclipses. + pub bodies_eclipsed: bool, /// Per channel, Earth's albedo from `Appearance::earth_albedo` (no texture). pub earth_albedo: [f32; 4], /// The same as `earth_lambert_1au` for the Moon. @@ -84,6 +86,7 @@ impl BodyRadiometry { earth_rgb: [basis(0), basis(1), basis(2)], earth_haze: albedo(&appearance.earth_haze), earth_textured: false, + bodies_eclipsed: false, earth_albedo: albedo(&appearance.earth_albedo), moon_lambert_1au: lambert(radiometry.effective_albedo(&appearance.moon_albedo)), moon_albedo: albedo(&appearance.moon_albedo), @@ -144,7 +147,7 @@ impl BodiesUniforms { params: [ near as f32, f32::from(u8::from(radiometry.earth_textured)), - 0.0, + f32::from(u8::from(radiometry.bodies_eclipsed)), 0.0, ], earth_center: vec4(rel.earth.translation, EARTH_EQUATORIAL_RADIUS), diff --git a/crates/resolved/src/gpu/earthshine.rs b/crates/resolved/src/gpu/earthshine.rs index dc857da..1d2f413 100644 --- a/crates/resolved/src/gpu/earthshine.rs +++ b/crates/resolved/src/gpu/earthshine.rs @@ -11,6 +11,7 @@ use super::bodies::BodyRadiometry; use super::profile::PassTimestamps; use super::GpuContext; use crate::bodies::{EARTH_EQUATORIAL_RADIUS, EARTH_POLAR_RADIUS, MOON_RADIUS}; +use crate::radiometry::SUN_RADIUS; use crate::scene::{Pose, Scene}; /// Rings of cells over each body's cone, from its center out to its limb. @@ -68,6 +69,7 @@ pub struct TargetLight { moon: Source, // mat3x3 columns, each padded to 16 bytes. earth_axes: [[f32; 4]; 3], + sun: [f32; 4], } impl TargetLight { @@ -81,10 +83,12 @@ impl TargetLight { [c.x as f32, c.y as f32, c.z as f32, 0.0] }; let source = |center, radius| Source::new(point, center, radius, scene.sun, &to_camera); + let sun = to_camera * (scene.sun - point); Self { earth: source(Vector3::zeros(), EARTH_EQUATORIAL_RADIUS), moon: source(scene.moon, MOON_RADIUS), earth_axes: [column(0), column(1), column(2)], + sun: [sun.x as f32, sun.y as f32, sun.z as f32, SUN_RADIUS as f32], } } } diff --git a/crates/resolved/src/render.rs b/crates/resolved/src/render.rs index 16546fc..e7cbf89 100644 --- a/crates/resolved/src/render.rs +++ b/crates/resolved/src/render.rs @@ -111,7 +111,8 @@ impl Renderer { camera.validate()?; let ctx = GpuContext::new(&options.adapter)?; let radiometry = Radiometry::new(&camera); - let body_radiometry = BodyRadiometry::new(&radiometry, &options.appearance)?; + let mut body_radiometry = BodyRadiometry::new(&radiometry, &options.appearance)?; + body_radiometry.bodies_eclipsed = options.fidelity.lighting >= LightingTier::L1; let targets = FrameTargets::new( &ctx, camera.width, diff --git a/crates/resolved/src/shaders/bodies.wgsl b/crates/resolved/src/shaders/bodies.wgsl index bdf949e..758ae9c 100644 --- a/crates/resolved/src/shaders/bodies.wgsl +++ b/crates/resolved/src/shaders/bodies.wgsl @@ -8,7 +8,8 @@ struct Uniforms { // fx, fy, cx, cy in pixels. intrinsics: vec4, - // x: near clip distance. y: 1 if Earth is drawn from the texture. + // x: near clip distance. y: 1 if Earth is drawn from the texture. z: 1 if Earth + // and the Moon shadow each other. params: vec4, // xyz: Earth center. w: equatorial radius. earth_center: vec4, @@ -102,13 +103,22 @@ fn fs_main(@builtin(position) pos: vec4) -> Out { } else { radiance = u.earth_lambert * lit; } + // The Moon's shadow: a solar eclipse. + if u.params.z > 0.5 && lit > 0.0 { + radiance *= sun_visible(te * d, u.sun_center, u.moon_center); + } } let tm = sphere_hit(-u.moon_center.xyz, d, u.moon_center.w); if tm > 0.0 && tm < best { best = tm; let n = normalize(tm * d - u.moon_center.xyz); - radiance = u.moon_lambert * max(dot(n, u.moon_sun.xyz), 0.0); + let lit = max(dot(n, u.moon_sun.xyz), 0.0); + radiance = u.moon_lambert * lit; + // Earth's shadow: a lunar eclipse. + if u.params.z > 0.5 && lit > 0.0 { + radiance *= sun_visible(tm * d, u.sun_center, u.earth_center); + } } let ts = sphere_hit(-u.sun_center.xyz, d, u.sun_center.w); diff --git a/crates/resolved/src/shaders/common.wgsl b/crates/resolved/src/shaders/common.wgsl index 6fd0568..0692f31 100644 --- a/crates/resolved/src/shaders/common.wgsl +++ b/crates/resolved/src/shaders/common.wgsl @@ -57,3 +57,33 @@ fn equirectangular(n: vec3) -> vec2 { fn textured_albedo(rgb: vec3, rgb_to_channels: mat3x4, haze: vec4) -> vec4 { return rgb_to_channels * rgb + haze; } + +// The fraction of the solar disk visible from `p` past one spherical occulter, with +// the Sun and the occulter given as (center, radius), all relative to the same origin. +// +// The overlap of two circles on the sky, as in eclipse.rs (Montenbruck and Gill). The +// lens area uses atan2 where eclipse.rs uses acos: acos loses its precision in f32 +// near 1, which is where a sliver of the Sun is covered. +fn sun_visible(p: vec3, sun: vec4, occulter: vec4) -> f32 { + let to_sun = sun.xyz - p; + let to_body = occulter.xyz - p; + let a = asin(min(sun.w / length(to_sun), 1.0)); + let b = asin(min(occulter.w / length(to_body), 1.0)); + // atan2 keeps precision near zero separation, where acos of a dot product loses it. + let c = atan2(length(cross(to_sun, to_body)), dot(to_sun, to_body)); + if c >= a + b { + return 1.0; + } + if c <= b - a { + return 0.0; + } + if c <= a - b { + return 1.0 - (b / a) * (b / a); + } + // A lens-shaped overlap: the circles cross at distance x along the line between + // their centers from the Sun's, y either side of it. + let x = (c * c + (a - b) * (a + b)) / (2.0 * c); + let y = sqrt(max((a - x) * (a + x), 0.0)); + let area = a * a * atan2(y, x) + b * b * atan2(y, c - x) - c * y; + return clamp(1.0 - area / (PI * a * a), 0.0, 1.0); +} diff --git a/crates/resolved/src/shaders/earthshine.wgsl b/crates/resolved/src/shaders/earthshine.wgsl index 58359fb..a02c5f0 100644 --- a/crates/resolved/src/shaders/earthshine.wgsl +++ b/crates/resolved/src/shaders/earthshine.wgsl @@ -6,8 +6,8 @@ // fixed grid of cells of equal solid angle (rings of equal width in cos(theta) // around the direction to the body, each split into equal sectors). A ray from the // target through each cell's center is shaded where it meets the body, as the -// bodies pass shades what the camera sees. Rays that miss (past an ellipsoid's -// flattened poles) or meet the other body first carry nothing. +// bodies pass shades what the camera sees, eclipses included. Rays that miss (past +// an ellipsoid's flattened poles) or meet the other body first carry nothing. // 2. `irradiance`: for every direction a surface can face, one texel of a small cube // map per target, the sum of what those cells send onto it: // E(n) = sum of L * max(n . d, 0) * dOmega. @@ -38,6 +38,8 @@ struct TargetLight { moon: Source, // Columns: Earth-fixed axes in the camera frame. earth_axes: mat3x3, + // xyz: the Sun's center relative to the target. w: its radius. + sun: vec4, } struct Globals { @@ -130,15 +132,19 @@ fn sources(@builtin(global_invocation_id) id: vec3) { let tm = sphere_hit(-t.moon.center.xyz, d, t.moon.center.w); var albedo = vec4(0.0); var normal = vec3(0.0); + // Each body in the other's shadow (eclipses) as the bodies pass draws it. + var visible = 1.0; if !is_moon && te > 0.0 && !(tm > 0.0 && tm < te) { let surface = earth_surface(t, d, te); albedo = surface[0]; normal = surface[1].xyz; + visible = sun_visible(te * d, t.sun, t.moon.center); } else if is_moon && tm > 0.0 && !(te > 0.0 && te < tm) { albedo = g.moon_albedo; normal = normalize(tm * d - t.moon.center.xyz); + visible = sun_visible(tm * d, t.sun, t.earth.center); } - let lit = max(dot(normal, source.sun.xyz), 0.0); + let lit = max(dot(normal, source.sun.xyz), 0.0) * visible; var out: Sample; out.light = albedo * (lit * source.sun.w * source.cone.y); out.direction = vec4(d, 0.0); diff --git a/crates/resolved/tests/earthshine.rs b/crates/resolved/tests/earthshine.rs index 96cdfff..db4e08d 100644 --- a/crates/resolved/tests/earthshine.rs +++ b/crates/resolved/tests/earthshine.rs @@ -291,7 +291,10 @@ fn quadrature(case: &Case, surfaces: &Surfaces, rings: usize, sectors: usize) -> if other.distance(p, d).is_some_and(|s| s < t) { continue; } - let lit = body.normal(p + t * d).dot(&to_sun).max(0.0); + // The other body's shadow on this one: a lunar or a solar eclipse. + let point = p + t * d; + let shadow = sun_visible_fraction(point, case.sun, &[(other.center, other.radius)]); + let lit = body.normal(point).dot(&to_sun).max(0.0) * shadow; for (s, albedo) in sum.iter_mut().zip(body.albedo) { *s += albedo * lit / PI * ratio * cell * facing; } @@ -534,6 +537,36 @@ fn bodies_hide_each_other() { ); } +#[test] +fn eclipses_dim_the_light_from_the_bodies() { + // The Moon in Earth's umbra sheds nothing; across the umbra's edge, part. Earth + // with the Moon's shadow on it, seen from inside that shadow. + let eclipsed_moon = |name, offset| { + let moon = Vector3::new(-MOON_DISTANCE, offset, 0.0); + Case { + moon, + ..Case::new( + name, + moon + Vector3::new(MOON_RADIUS + 100e3, 0.0, 0.0), + -Vector3::x(), + ) + } + }; + compare( + Surfaces::WHITE_EARTH, + 1e-8, + &[ + eclipsed_moon("moon-in-umbra", 0.0), + eclipsed_moon("moon-in-penumbra", 5.5e6), + ], + ); + let solar = Case { + moon: Vector3::new(3.6e8, 0.0, 0.0), + ..Case::new("earth-under-the-moon", LEO, -Vector3::x()) + }; + compare(Surfaces::WHITE_EARTH, 1e-4, &[solar]); +} + #[test] fn lighting_l0_has_no_earthshine() { let _devices = common::device_lock(); diff --git a/crates/resolved/tests/eclipses.rs b/crates/resolved/tests/eclipses.rs new file mode 100644 index 0000000..4e9edbf --- /dev/null +++ b/crates/resolved/tests/eclipses.rs @@ -0,0 +1,286 @@ +//! Eclipses seen on the bodies at lighting L1, on the GPU and on the CPU (lavapipe): +//! Earth's shadow on the Moon and the Moon's on Earth, pixel by pixel against +//! `eclipse::sun_visible_fraction` at each pixel's surface point. + +mod common; + +use std::sync::Arc; + +use resolved::bodies::{EARTH_EQUATORIAL_RADIUS, EARTH_POLAR_RADIUS, MOON_RADIUS}; +use resolved::eclipse::sun_visible_fraction; +use resolved::gpu::AdapterChoice; +use resolved::nalgebra::{UnitQuaternion, Vector2, Vector3}; +use resolved::radiometry::AU; +use resolved::{ + CameraConfig, EarthTexture, Fidelity, LightingTier, Pose, Renderer, RendererOptions, Scene, +}; + +const SUN: Vector3 = Vector3::new(AU, 0.0, 0.0); +const MOON_DISTANCE: f64 = 3.844e8; + +fn renderer(adapter: &AdapterChoice, lighting: LightingTier) -> Renderer { + let gray = EarthTexture::from_rgba8(1, 1, vec![128, 128, 128, 255]).unwrap(); + let options = RendererOptions { + adapter: adapter.clone(), + fidelity: Fidelity { + lighting, + ..Default::default() + }, + earth_texture: Some(Arc::new(gray)), + read_radiance: true, + ..Default::default() + }; + Renderer::new(CameraConfig::example(), options).unwrap() +} + +fn scene(moon: Vector3) -> Scene { + Scene { + earth_attitude: UnitQuaternion::identity(), + sun: SUN, + moon, + targets: vec![], + } +} + +fn radiance(renderer: &mut Renderer, scene: &Scene, camera: &Pose) -> Vec<[f32; 4]> { + let frame = renderer.submit(scene, camera).unwrap().wait().unwrap(); + frame.radiance.unwrap() +} + +/// Where the line of sight through the center of pixel (u, v) first meets the sphere +/// or (Z-flattened, identity attitude) ellipsoid at `center`, inertial. +fn hit( + camera: &Pose, + u: usize, + v: usize, + center: Vector3, + radius: f64, + flattening: f64, +) -> Option> { + let k = CameraConfig::example().intrinsics(); + let d = camera.attitude * k.ray(&Vector2::new(u as f64 + 0.5, v as f64 + 0.5)); + let stretch = Vector3::new(1.0, 1.0, flattening); + let o = (camera.position - center).component_mul(&stretch); + let ds = d.component_mul(&stretch); + let (a, b, c) = ( + ds.norm_squared(), + o.dot(&ds), + o.norm_squared() - radius * radius, + ); + let disc = b * b - a * c; + let t = (-b - disc.max(0.0).sqrt()) / a; + (disc > 0.0 && t > 0.0).then(|| camera.position + t * d) +} + +/// How the pixels of one body compared: counts of those entirely in shadow, partly and +/// not at all, the least expected fraction among the partly shadowed, and the largest +/// difference between rendered and expected fraction. +#[derive(Debug)] +struct Tally { + dark: usize, + partial: usize, + lit: usize, + least: f64, + worst: f64, +} + +/// The pixels compared: every second one each way, which keeps the f64 work quick in a +/// debug build. +fn sampled() -> impl Iterator { + (0..1080) + .step_by(2) + .flat_map(|v| (0..1920).step_by(2).map(move |u| (u, v))) +} + +/// For each [`sampled`] pixel, the fraction of the Sun `visible` gives at the surface +/// point `point` finds there, if it finds one. Computed once per scene and shared by +/// both devices, on every core. +fn expected( + point: impl Fn(usize, usize) -> Option> + Sync, + visible: impl Fn(Vector3) -> f64 + Sync, +) -> Vec> { + let pixels: Vec<_> = sampled().collect(); + let threads = std::thread::available_parallelism().map_or(1, |n| n.get()); + let chunk = pixels.len().div_ceil(threads); + std::thread::scope(|scope| { + let workers: Vec<_> = pixels + .chunks(chunk) + .map(|part| { + let (point, visible) = (&point, &visible); + scope.spawn(move || { + part.iter() + .map(|&(u, v)| point(u, v).map(visible)) + .collect::>() + }) + }) + .collect(); + workers + .into_iter() + .flat_map(|w| w.join().unwrap()) + .collect() + }) +} + +/// Compares `eclipsed` with `clear` (the same pixels without the shadow) against +/// `expected` at every [`sampled`] pixel that has a value. Pixels at the body's edge, +/// where `clear` is under 1 % of its largest value, are skipped: a rendered edge and +/// an f64 ray need not agree on which pixels they cover. +fn tally(eclipsed: &[[f32; 4]], clear: &[[f32; 4]], expected: &[Option]) -> Tally { + let peak = clear.iter().map(|p| p[1]).fold(0.0, f32::max); + let mut t = Tally { + dark: 0, + partial: 0, + lit: 0, + least: 1.0, + worst: 0.0, + }; + for ((u, v), want) in sampled().zip(expected) { + let (got, base) = (eclipsed[v * 1920 + u], clear[v * 1920 + u]); + let Some(want) = *want else { continue }; + if base[1] < 0.01 * peak { + continue; + } + for k in 0..3 { + let ratio = got[k] as f64 / base[k] as f64; + if want == 0.0 { + assert_eq!(got[k], 0.0, "pixel ({u}, {v}) is in the umbra"); + } else if want == 1.0 { + assert_eq!(got[k], base[k], "pixel ({u}, {v}) is outside the shadow"); + } else { + t.worst = t.worst.max((ratio - want).abs()); + } + } + match want { + 0.0 => t.dark += 1, + 1.0 => t.lit += 1, + _ => { + t.partial += 1; + t.least = t.least.min(want); + } + } + } + t +} + +/// A camera 44,400 km in front of the Moon's sunward face, where the Moon is 712 +/// pixels across, with Earth behind the camera. +fn moon_camera(moon: Vector3) -> Pose { + Pose::looking_at(moon + Vector3::new(4.44e7, 0.0, 0.0), moon, Vector3::z()).unwrap() +} + +#[test] +fn earth_shadows_the_moon() { + let earth = [(Vector3::zeros(), EARTH_EQUATORIAL_RADIUS)]; + // Straight behind Earth: inside the umbra, which reaches 4600 km from the axis at + // the Moon's distance. Then 5500 km off the axis, across the umbra's edge into the + // penumbra (which ends 8200 km out). + let cases = [("total", 0.0), ("partial", 5.5e6)].map(|(label, offset)| { + let moon = Vector3::new(-MOON_DISTANCE, offset, 0.0); + let camera = moon_camera(moon); + let want = expected( + |u, v| hit(&camera, u, v, moon, MOON_RADIUS, 1.0), + |p| sun_visible_fraction(p, SUN, &earth), + ); + (label, moon, camera, want) + }); + let _devices = common::device_lock(); + for adapter in [AdapterChoice::Gpu, AdapterChoice::Cpu] { + let mut l0 = renderer(&adapter, LightingTier::L0); + let mut l1 = renderer(&adapter, LightingTier::L1); + let name = l1.adapter_info().name.clone(); + for (label, moon, camera, want) in &cases { + let clear = radiance(&mut l0, &scene(*moon), camera); + let eclipsed = radiance(&mut l1, &scene(*moon), camera); + common::save( + &format!("eclipse-lunar-{label}"), + l1.adapter_info(), + &eclipsed, + ); + let t = tally(&eclipsed, &clear, want); + match *label { + "total" => assert!( + t.dark > 50_000 && t.partial == 0 && t.lit == 0, + "{name}: {t:?}" + ), + _ => assert!(t.dark > 1000 && t.partial > 10_000, "{name}: {t:?}"), + } + // Measured: 9e-7 on the RTX 3060, 5e-5 on lavapipe. + assert!(t.worst < 1e-4, "{name} {label}: {t:?}"); + } + } +} + +#[test] +fn moon_shadows_earth() { + // Looking down at the subsolar point from 150,000 km, where Earth is 775 pixels + // across. Without the eclipse, the Moon sits behind Earth, out of view. + let camera = Pose::looking_at( + Vector3::new(1.5e8, 0.0, 0.0), + Vector3::zeros(), + Vector3::z(), + ) + .unwrap(); + let hidden = scene(Vector3::new(-MOON_DISTANCE, 0.0, 0.0)); + let flattening = EARTH_EQUATORIAL_RADIUS / EARTH_POLAR_RADIUS; + // Near enough for the umbra to reach the ground (a total eclipse); then at the mean + // distance and off the axis, where it falls short (annular). + let cases = [ + ("total", Vector3::new(3.6e8, 0.0, 0.0)), + ("annular", Vector3::new(MOON_DISTANCE, 2.0e6, 1.0e6)), + ] + .map(|(label, moon)| { + let want = expected( + |u, v| { + hit( + &camera, + u, + v, + Vector3::zeros(), + EARTH_EQUATORIAL_RADIUS, + flattening, + ) + }, + |p| sun_visible_fraction(p, SUN, &[(moon, MOON_RADIUS)]), + ); + (label, moon, want) + }); + let _devices = common::device_lock(); + for adapter in [AdapterChoice::Gpu, AdapterChoice::Cpu] { + let mut l1 = renderer(&adapter, LightingTier::L1); + let name = l1.adapter_info().name.clone(); + let clear = radiance(&mut l1, &hidden, &camera); + for (label, moon, want) in &cases { + let eclipsed = radiance(&mut l1, &scene(*moon), &camera); + common::save( + &format!("eclipse-solar-{label}"), + l1.adapter_info(), + &eclipsed, + ); + let t = tally(&eclipsed, &clear, want); + match *label { + "total" => assert!(t.dark > 10 && t.partial > 10_000, "{name}: {t:?}"), + // In the antumbra the Moon sits inside the Sun's disk and 1 - (b/a)^2, + // about 0.024, of it shows; a partial overlap never gets that low. + _ => assert!( + t.dark == 0 && t.partial > 10_000 && t.least < 0.03, + "{name}: {t:?}" + ), + } + // Measured: 6e-7 on the RTX 3060, 1.6e-5 on lavapipe. + assert!(t.worst < 1e-4, "{name} {label}: {t:?}"); + } + } +} + +#[test] +fn lighting_l0_draws_no_eclipses_on_the_bodies() { + let _devices = common::device_lock(); + let moon = Vector3::new(-MOON_DISTANCE, 0.0, 0.0); + let camera = moon_camera(moon); + for adapter in [AdapterChoice::Gpu, AdapterChoice::Cpu] { + let mut l0 = renderer(&adapter, LightingTier::L0); + let image = radiance(&mut l0, &scene(moon), &camera); + // The Moon in Earth's umbra, drawn as a full Moon. + assert!(image[540 * 1920 + 960][1] > 0.0); + } +} diff --git a/docs/optics.md b/docs/optics.md index c0e7303..bb649bf 100644 --- a/docs/optics.md +++ b/docs/optics.md @@ -690,7 +690,76 @@ The renderer counts both Earth and the Moon as occulters. The fraction is computed once per target, at its center. The whole target dims together, which is accurate unless the target is a sizable fraction of the penumbra's width. Light reaching the target from Earth itself is the -next section. +next section but one. + +### Eclipses on Earth and the Moon + +Code: `shaders/common.wgsl` (`sun_visible`), `shaders/bodies.wgsl`. From +lighting L1. + +Earth and the Moon also shadow each other. When the Moon passes through +Earth's shadow it is a lunar eclipse; when the Moon's shadow falls on Earth, +a solar eclipse. For every pixel that sees Earth or the Moon, the shader +works out the fraction of the Sun's disk visible from that point on the +surface, past the other body, with the same overlap of circles as for the +target, and multiplies the point's sunlight by it. + +The sizes involved: + +- At the Moon's distance, Earth's umbra reaches 4600 km from the shadow's + axis and its penumbra 8200 km. The Moon, 3474 km across, fits inside the + umbra with room to spare, so a total lunar eclipse lasts over an hour. +- The Moon's umbra only just reaches Earth. When the Moon is near its + closest, the umbra touches down as a dark spot up to about 270 km across: + a total solar eclipse. When it is farther, the umbra ends in space, and + beyond its tip the Moon covers only the middle of the Sun's disk, leaving a + ring: an annular eclipse, where a few percent of the Sun still shows. The + penumbra around either is about 7000 km across. + +The overlap is computed in 32-bit floats, which is enough: a lunar surface +point is placed to about 20 m out of 384,000 km, 6e-8 radians against a +solar radius of 4.6e-3 radians seen from there. One formula needed care. The +area where the two circles overlap involves angles whose cosines are close +to 1 when a sliver of the Sun is covered, and the inverse cosine of a number +near 1 loses most of its precision in 32 bits. The shader finds those angles +from their sine and cosine together instead. With the inverse cosine, +pixels at the shadow's edges came out wrong: a point the 64-bit model puts +fully in the umbra read slightly lit, and one just outside the shadow +slightly dimmed. + +Tests (`tests/eclipses.rs`), on the GPU and lavapipe, compare each pixel of +Earth and the Moon against the 64-bit model at that pixel's surface point: +the Moon entirely in Earth's umbra, the Moon across the umbra's edge, a +total solar eclipse and an annular one. Pixels the model puts fully in shadow +read exactly zero, those it puts fully outside are unchanged, and the rest +agree to within 1e-4 of the Sun's light (9e-7 at most on the GPU, 5e-5 on +lavapipe). The command-line tests render two real eclipses from their +epochs: the total lunar eclipse of 2025 March 14 at mid-totality, where the +whole Moon reads the bias, and the total solar eclipse of 2026 August 12, +where the point of greatest eclipse in NASA's bulletin (65.2 N 25.2 W) is in +the umbra (`examples/scenes/lunar-eclipse.toml`, `solar-eclipse.toml`). In +the second, the umbra is 220 by 410 km in the image, stretched along the +direction of the Sun, which is 26 degrees high there; NASA gives the path as +294 km wide. + +The light Earth and the Moon reflect onto targets (next section) includes +these shadows, so an eclipsed Moon sheds less light. + +On the RTX 3060 this adds nothing measurable to the bodies pass (0.202 to +0.205 ms with Earth filling the image). On lavapipe it adds about 1 ms when +Earth or the Moon fills much of the image. + +Limitations: + +- **Earth's atmosphere.** It bends sunlight into the umbra, mostly red, so a + real eclipsed Moon glows a dim copper color instead of going black, and it + makes the umbra about 2 % larger. Both come with the atmosphere (M4). +- **Smooth shapes.** Earth's shadow is cast by a sphere of its equatorial + radius, as for the target's eclipses; its flattening would move the + shadow's edge on the Moon by about 20 km. The Moon's mountains and valleys, + which shape the edge of a solar eclipse's path, are not modeled. +- **Uniform Sun.** The Sun's disk is taken as uniformly bright until limb + darkening (M2) arrives. ### Light from Earth and the Moon @@ -781,9 +850,6 @@ Limitations: here, as they do in the camera's view. The Moon's real reflection is brighter toward the Sun (M4), and Earth has no clouds or atmosphere yet beyond the stand-in haze. -- **No eclipses on the bodies.** Earth's shadow on the Moon and the Moon's - shadow on Earth are not modeled yet (M2), so a Moon in eclipse still lights - the target. ## Where the Sun and Moon appear @@ -815,7 +881,6 @@ velocity as an input, which the scene does not carry yet. | Effect | Where | |---|---| -| Eclipses seen on the bodies: Earth's shadow on the Moon, the Moon's on Earth | M2 | | Earthshine on the Moon's night side | M4 | | The Sun as a disk that is darker toward its edge (limb darkening) | M2 | | Pixel-area sampling (anti-aliasing): graded edges, thin parts | M3 | diff --git a/examples/scenes/lunar-eclipse.toml b/examples/scenes/lunar-eclipse.toml new file mode 100644 index 0000000..2ba60d8 --- /dev/null +++ b/examples/scenes/lunar-eclipse.toml @@ -0,0 +1,19 @@ +# The total lunar eclipse of 2025 March 14 at mid-totality (06:58:44 UTC), seen through +# a 400 mm lens from 7000 km out on the line from Earth's center to the Moon. Needs the +# astrokit feature and its kernels. +# +# resolved render --scene examples/scenes/lunar-eclipse.toml --lighting L1 --preview lunar.png +# +# At lighting L1 the Moon is inside Earth's shadow and dark. At L0 it is drawn without +# Earth's shadow, as a full Moon. Earth's atmosphere bends red light into the shadow +# and a real eclipsed Moon glows a dim red; that comes with the atmosphere (M4). + +epoch = "2025-03-14T06:58:44 UTC" + +[camera] +focal_length = 0.4 + +[observer] +position = [-6956600.88, 697838.87, 344565.04] +look_at = [-399016819.56, 40026652.30, 19763566.70] +up = [0, 0, 1] diff --git a/examples/scenes/solar-eclipse.toml b/examples/scenes/solar-eclipse.toml new file mode 100644 index 0000000..8c1e317 --- /dev/null +++ b/examples/scenes/solar-eclipse.toml @@ -0,0 +1,16 @@ +# The total solar eclipse of 2026 August 12 at greatest eclipse (17:46 UTC), seen from +# geostationary height straight above the point of greatest eclipse, 65.2 N 25.2 W, off +# the west coast of Iceland. Needs the astrokit feature and its kernels. +# +# resolved render --scene examples/scenes/solar-eclipse.toml --lighting L1 --preview solar.png +# +# At lighting L1 the Moon's shadow lies on Earth: the dark umbra, about 290 km across, +# at the image center, inside a penumbra thousands of kilometers wide. At L0 Earth is +# drawn without it. + +epoch = "2026-08-12T17:46:00 UTC" + +[observer] +position = [-16371789.99, -6703593.18, 38273087.50] +look_at = [-2469723.91, -1011253.16, 5773587.33] +up = [0, 0, 1] diff --git a/readme.md b/readme.md index 1dae05e..b6a1fe8 100644 --- a/readme.md +++ b/readme.md @@ -96,9 +96,10 @@ assert!(frame.dn_at(960, 540, 1) > 1000); observer, Sun and Moon positions (or an epoch, with the `astrokit` feature) and targets with their models. `examples/scenes/` has commented examples: `boxwing.toml` and `boxwing-epoch.toml` use a small generated model that is -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`), +in the repository, `lunar-eclipse.toml` and `solar-eclipse.toml` show two +real eclipses at lighting L1, 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 @@ -143,8 +144,8 @@ cargo run --release -p resolved-cli -- render \ M1 is complete: the library renders Earth, Moon, Sun and targets to DN through `Renderer`, the CLI renders scene files, and the benchmarks give the L0 / C0 baseline above. M2 is in progress: lighting L1 so far draws Earth's -surface from Blue Marble, shadows on the target, and light from Earth and the -Moon on it. The work is tracked as +surface from Blue Marble, shadows on the target, light from Earth and the +Moon on it, and eclipses on Earth and the Moon. The work is tracked as Gitea milestones, each split into one issue per change: | Milestone | Scope |