From 3069b84b8d8c5b8b455452d95e0c7386f577e046 Mon Sep 17 00:00:00 2001 From: Vinh Truong Date: Thu, 9 Feb 2023 17:29:11 +0200 Subject: [PATCH] Replace sky.rs with hw-skymodel --- Cargo.toml | 1 + src/cpupt/mod.rs | 5 +- src/cpupt/sky.rs | 220 +++-------------------------------------------- src/debug/mod.rs | 2 +- src/main.rs | 1 - 5 files changed, 18 insertions(+), 211 deletions(-) diff --git a/Cargo.toml b/Cargo.toml index e525d07..67b0cb1 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -16,6 +16,7 @@ clap = { version = "4.1.4", features = ["derive"] } easer = "0.3.0" env_logger = "0.10.0" gltf = "1.0.0" +hw-skymodel = "0.1.1" image = "0.24.5" imageproc = "0.23.0" imgui = { version = "0.10.0", features = ["tables-api"] } diff --git a/src/cpupt/mod.rs b/src/cpupt/mod.rs index 34b5c36..baddd2a 100644 --- a/src/cpupt/mod.rs +++ b/src/cpupt/mod.rs @@ -150,7 +150,8 @@ impl Raytracer { let mut timer = Instant::now(); let mut ray_stats = intersection::RayBvhHitStats::default(); let mut pixel_buffer = Vec::::new(); - let mut sky_state = sky::ext::StateExt::new(&sky::ext::StateExtParams::default()); + let mut sky_state = + sky::ext::StateExt::new(&sky::ext::StateExtParams::default()).unwrap(); loop { // Check for termination command. @@ -193,7 +194,7 @@ impl Raytracer { camera_position = camera_transform.transform_point(&camera.position()); // Reset sky. - sky_state = sky::ext::StateExt::new(&input.sky_params); + sky_state = sky::ext::StateExt::new(&input.sky_params).unwrap(); // Reset stats. ray_stats = intersection::RayBvhHitStats::default(); diff --git a/src/cpupt/sky.rs b/src/cpupt/sky.rs index 0db4c12..55ba1d5 100644 --- a/src/cpupt/sky.rs +++ b/src/cpupt/sky.rs @@ -7,201 +7,6 @@ use super::*; // License file: sky/hosek-wilkie-license.txt // -// This is a rewrite of the original Hosek & Wilkie sky model. Only the RGB -// model is supported, and the precision has been reduced to f32. Solar radiance -// functions and alien-like worlds are also removed. The state size has also -// been trimmed down from 1088 to only 120 bytes. - -pub enum Channel { - R = 0, - G = 1, - B = 2, -} - -#[derive(Clone, Copy, Debug)] -pub struct Params { - elevation: f32, - turbidity: f32, - albedo: [f32; 3], -} - -impl Default for Params { - fn default() -> Self { - Self { - elevation: 0.0, - turbidity: 1.0, - albedo: [1.0, 1.0, 1.0], - } - } -} - -impl Params { - pub fn elevation(self, elevation: f32) -> Self { - assert!( - (0.0..=FRAC_PI_2).contains(&elevation), - "Solar elevation must be in [0,pi/2], got {elevation} instead" - ); - Self { elevation, ..self } - } - - pub fn turbidity(self, turbidity: f32) -> Self { - assert!( - (1.0..=10.0).contains(&turbidity), - "Turbidity must be in [1,10], got {turbidity} instead" - ); - Self { turbidity, ..self } - } - - pub fn albedo(self, albedo: [f32; 3]) -> Self { - for (channel, albedo) in albedo.iter().enumerate() { - assert!( - (0.0..=1.0).contains(albedo), - "Albedo ({channel}) must be in [0,1], got {albedo} instead", - ); - } - Self { albedo, ..self } - } -} - -#[derive(Clone, Copy, Debug)] -pub struct State { - params: [f32; 27], - radiances: [f32; 3], -} - -impl State { - pub fn new(sky_params: &Params) -> Self { - // Load datasets. - macro_rules! dataset { - ($path:literal) => {{ - // Note: Rust's include_bytes! doesn't guarantee that the byte - // slice is aligned. That's why we use include_bytes_aligned - // crate to guarantee 4-byte alignment, so we can cast into - // f32-slice. - use include_bytes_aligned::include_bytes_aligned; - let bytes = include_bytes_aligned!(4, concat!("sky/", $path)); - bytemuck::cast_slice::<_, f32>(bytes) - }}; - } - let params_r = dataset!("hosek-wilkie-params-r"); - let params_g = dataset!("hosek-wilkie-params-g"); - let params_b = dataset!("hosek-wilkie-params-b"); - let radiances_r = dataset!("hosek-wilkie-radiances-r"); - let radiances_g = dataset!("hosek-wilkie-radiances-g"); - let radiances_b = dataset!("hosek-wilkie-radiances-b"); - - // Init state. - let mut params = [0.0; 27]; - let mut radiances = [0.0; 3]; - let elevation = sky_params.elevation; - let turbidity = sky_params.turbidity; - let albedo = sky_params.albedo; - let t = (elevation / (0.5 * PI)).powf(1.0 / 3.0); - init_params(&mut params[..], params_r, turbidity, albedo[0], t); - init_params(&mut params[9..], params_g, turbidity, albedo[1], t); - init_params(&mut params[(9 * 2)..], params_b, turbidity, albedo[2], t); - init_radiances(&mut radiances[0], radiances_r, turbidity, albedo[0], t); - init_radiances(&mut radiances[1], radiances_g, turbidity, albedo[1], t); - init_radiances(&mut radiances[2], radiances_b, turbidity, albedo[2], t); - - Self { params, radiances } - } - - pub fn radiance(&self, theta: f32, gamma: f32, channel: Channel) -> f32 { - let channel = channel as usize; - let r = self.radiances[channel]; - let p = &self.params[(9 * channel)..]; - let p0 = p[0]; - let p1 = p[1]; - let p2 = p[2]; - let p3 = p[3]; - let p4 = p[4]; - let p5 = p[5]; - let p6 = p[6]; - let p7 = p[7]; - let p8 = p[8]; - - let cos_gamma = gamma.cos(); - let cos_gamma2 = cos_gamma * cos_gamma; - let cos_theta = theta.cos().abs(); - - let exp_m = (p4 * gamma).exp(); - let ray_m = cos_gamma2; - let mie_m_lhs = 1.0 + cos_gamma2; - let mie_m_rhs = (1.0 + p8 * p8 - 2.0 * p8 * cos_gamma).powf(1.5); - let mie_m = mie_m_lhs / mie_m_rhs; - let zenith = cos_theta.sqrt(); - let radiance_lhs = 1.0 + p0 * (p1 / (cos_theta + 0.01)).exp(); - let radiance_rhs = p2 + p3 * exp_m + p5 * ray_m + p6 * mie_m + p7 * zenith; - let radiance_dist = radiance_lhs * radiance_rhs; - r * radiance_dist - } -} - -fn init_params(out_params: &mut [f32], dataset: &[f32], turbidity: f32, albedo: f32, t: f32) { - let turbidity_int = turbidity.trunc() as usize; - let turbidity_rem = turbidity.fract(); - let turbidity_min = turbidity_int.saturating_sub(1); - let turbidity_max = turbidity_int.min(9); - let p0 = &dataset[(9 * 6 * turbidity_min)..]; - let p1 = &dataset[(9 * 6 * turbidity_max)..]; - let p2 = &dataset[(9 * 6 * 10 + 9 * 6 * turbidity_min)..]; - let p3 = &dataset[(9 * 6 * 10 + 9 * 6 * turbidity_max)..]; - let s0 = (1.0 - albedo) * (1.0 - turbidity_rem); - let s1 = (1.0 - albedo) * turbidity_rem; - let s2 = albedo * (1.0 - turbidity_rem); - let s3 = albedo * turbidity_rem; - - for i in 0..9 { - out_params[i] += s0 * quintic::<9>(&p0[i..], t); - out_params[i] += s1 * quintic::<9>(&p1[i..], t); - out_params[i] += s2 * quintic::<9>(&p2[i..], t); - out_params[i] += s3 * quintic::<9>(&p3[i..], t); - } -} - -fn init_radiances(out_radiance: &mut f32, dataset: &[f32], turbidity: f32, albedo: f32, t: f32) { - let turbidity_int = turbidity.trunc() as usize; - let turbidity_rem = turbidity.fract(); - let turbidity_min = turbidity_int.saturating_sub(1); - let turbidity_max = turbidity_int.min(9); - let p0 = &dataset[(6 * turbidity_min)..]; - let p1 = &dataset[(6 * turbidity_max)..]; - let p2 = &dataset[(6 * 10 + 6 * turbidity_min)..]; - let p3 = &dataset[(6 * 10 + 6 * turbidity_max)..]; - let s0 = (1.0 - albedo) * (1.0 - turbidity_rem); - let s1 = (1.0 - albedo) * turbidity_rem; - let s2 = albedo * (1.0 - turbidity_rem); - let s3 = albedo * turbidity_rem; - - *out_radiance += s0 * quintic::<1>(p0, t); - *out_radiance += s1 * quintic::<1>(p1, t); - *out_radiance += s2 * quintic::<1>(p2, t); - *out_radiance += s3 * quintic::<1>(p3, t); -} - -fn quintic(p: &[f32], t: f32) -> f32 { - let t2 = t * t; - let t3 = t2 * t; - let t4 = t2 * t2; - let t5 = t4 * t; - - let inv_t = 1.0 - t; - let inv_t2 = inv_t * inv_t; - let inv_t3 = inv_t2 * inv_t; - let inv_t4 = inv_t2 * inv_t2; - let inv_t5 = inv_t4 * inv_t; - - let m0 = p[0] * inv_t5; - let m1 = p[STRIDE] * 5.0 * inv_t4 * t; - let m2 = p[2 * STRIDE] * 10.0 * inv_t3 * t2; - let m3 = p[3 * STRIDE] * 10.0 * inv_t2 * t3; - let m4 = p[4 * STRIDE] * 5.0 * inv_t * t4; - let m5 = p[5 * STRIDE] * t5; - - m0 + m1 + m2 + m3 + m4 + m5 -} - // // Raydiance-specific extensions // @@ -209,6 +14,8 @@ fn quintic(p: &[f32], t: f32) -> f32 { pub mod ext { use super::*; + use hw_skymodel::rgb::*; + #[derive(Clone, Copy, Debug, PartialEq)] pub struct StateExtParams { pub elevation: f32, @@ -229,12 +36,12 @@ pub mod ext { } pub struct StateExt { - state: State, + state: SkyState, sun_dir: Normal, } impl StateExt { - pub fn new(params: &StateExtParams) -> Self { + pub fn new(params: &StateExtParams) -> Result { // Validate. assert!( (0.0..=TAU).contains(¶ms.azimuth), @@ -243,12 +50,11 @@ pub mod ext { ); // Init state. - let state = State::new( - &super::Params::default() - .elevation(params.elevation) - .turbidity(params.turbidity) - .albedo(params.albedo.into()), - ); + let state = SkyState::new(&SkyParams { + elevation: params.elevation, + turbidity: params.turbidity, + albedo: params.albedo.into(), + })?; // Pre-compute sun direction. let sun_angle_xz = params.azimuth; @@ -259,7 +65,7 @@ pub mod ext { sun_angle_y.sin() * sun_angle_xz.sin() ]; - Self { state, sun_dir } + Ok(Self { state, sun_dir }) } pub fn sun_dir(&self) -> Normal { @@ -271,9 +77,9 @@ pub mod ext { let cos_gamma = ray_dir.dot(&self.sun_dir).clamp(-1.0, 1.0); let gamma = cos_gamma.acos(); ColorRgb::new( - self.state.radiance(theta, gamma, sky::Channel::R), - self.state.radiance(theta, gamma, sky::Channel::G), - self.state.radiance(theta, gamma, sky::Channel::B), + self.state.radiance(theta, gamma, Channel::R), + self.state.radiance(theta, gamma, Channel::G), + self.state.radiance(theta, gamma, Channel::B), ) } } diff --git a/src/debug/mod.rs b/src/debug/mod.rs index d632c69..1daf4ed 100644 --- a/src/debug/mod.rs +++ b/src/debug/mod.rs @@ -477,7 +477,7 @@ fn sky_model_visualizations() -> Result<()> { azimuth, turbidity, albedo, - }); + })?; // Plot. let mut plot = plot::Plot::new(|wo| { diff --git a/src/main.rs b/src/main.rs index 6f601db..32ad229 100644 --- a/src/main.rs +++ b/src/main.rs @@ -76,7 +76,6 @@ use math::*; const PI: f32 = std::f32::consts::PI; const TAU: f32 = std::f32::consts::TAU; const INV_PI: f32 = std::f32::consts::FRAC_1_PI; -const FRAC_PI_2: f32 = std::f32::consts::FRAC_PI_2; // // Input state