1use half::f16;
2use serde::{Deserialize, Serialize};
3
4pub const AXES: [[f32; 3]; 6] = [
5 [1.0, 0.0, 0.0],
6 [-1.0, 0.0, 0.0],
7 [0.0, 1.0, 0.0],
8 [0.0, -1.0, 0.0],
9 [0.0, 0.0, 1.0],
10 [0.0, 0.0, -1.0],
11];
12const MAGIC: &[u8; 8] = b"PITOPRB1";
13const HEADER: usize = 64;
14const STRIDE_V1: usize = 60;
15const STRIDE_V2: usize = 80;
16
17fn packed(value: f32) -> f16 {
18 f16::from_f32(value.clamp(0.0, 65504.0))
19}
20
21#[derive(Clone, Copy, Debug, PartialEq, Serialize, Deserialize)]
22pub struct GridSpec {
23 pub min: [f32; 3],
24 pub max: [f32; 3],
25 pub spacing: f32,
26}
27
28impl GridSpec {
29 pub fn dimensions(self) -> Result<[u32; 3], String> {
30 if !self.spacing.is_finite() || self.spacing <= 0.0 {
31 return Err("invalid probe spacing".into());
32 }
33 let mut dims = [0; 3];
34 for (axis, dim) in dims.iter_mut().enumerate() {
35 if !self.min[axis].is_finite()
36 || !self.max[axis].is_finite()
37 || self.max[axis] < self.min[axis]
38 {
39 return Err("invalid probe bounds".into());
40 }
41 let ratio = (self.max[axis] as f64 - self.min[axis] as f64) / self.spacing as f64;
42 let steps = ratio.round();
43 if (ratio - steps).abs() > 1e-5 {
44 return Err("probe bounds must be a multiple of spacing".into());
45 }
46 if steps >= u32::MAX as f64 {
47 return Err("probe dimensions overflow".into());
48 }
49 *dim = steps as u32 + 1;
50 }
51 let count = u64::from(dims[0]) * u64::from(dims[1]) * u64::from(dims[2]);
52 if count > u32::MAX as u64 || count > usize::MAX as u64 / STRIDE_V2 as u64 {
53 return Err("probe count overflow".into());
54 }
55 Ok(dims)
56 }
57
58 pub fn index(dims: [u32; 3], xyz: [u32; 3]) -> usize {
59 (xyz[2] as usize * dims[1] as usize + xyz[1] as usize) * dims[0] as usize + xyz[0] as usize
60 }
61
62 pub fn position(self, dims: [u32; 3], index: usize) -> [f32; 3] {
63 let xyz = [
64 index as u32 % dims[0],
65 index as u32 / dims[0] % dims[1],
66 index as u32 / dims[0] / dims[1],
67 ];
68 [0, 1, 2].map(|axis| self.min[axis] + xyz[axis] as f32 * self.spacing)
69 }
70}
71
72#[derive(Clone, Debug, PartialEq)]
73pub struct Probe {
74 pub lobes: [[f32; 3]; 6],
75 pub visibility: [f32; 6],
76}
77
78#[derive(Clone, Debug, PartialEq)]
79pub struct ProbeExtra {
80 pub second: [f32; 6],
81 pub offset: [f32; 3],
82 pub enabled: bool,
83 pub backfaces: u16,
84}
85
86impl Default for Probe {
87 fn default() -> Self {
88 Self {
89 lobes: [[0.0; 3]; 6],
90 visibility: [65504.0; 6],
91 }
92 }
93}
94
95#[derive(Clone, Debug, PartialEq)]
96pub struct Grid {
97 pub spec: GridSpec,
98 pub dims: [u32; 3],
99 pub probes: Vec<Probe>,
100 pub extras: Option<Vec<ProbeExtra>>,
101}
102
103impl Grid {
104 pub fn new(spec: GridSpec, probes: Vec<Probe>) -> Result<Self, String> {
105 let dims = spec.dimensions()?;
106 if probes.len() != dims.iter().map(|&v| v as usize).product::<usize>() {
107 return Err("probe count does not match grid".into());
108 }
109 if probes.iter().any(|probe| {
110 probe
111 .lobes
112 .iter()
113 .flatten()
114 .any(|value| !value.is_finite() || *value < 0.0)
115 || probe
116 .visibility
117 .iter()
118 .any(|value| !value.is_finite() || *value < 0.0)
119 }) {
120 return Err("probe values must be finite and nonnegative".into());
121 }
122 Ok(Self {
123 spec,
124 dims,
125 probes,
126 extras: None,
127 })
128 }
129
130 pub fn with_extras(
131 spec: GridSpec,
132 probes: Vec<Probe>,
133 extras: Vec<ProbeExtra>,
134 ) -> Result<Self, String> {
135 let mut grid = Self::new(spec, probes)?;
136 if extras.len() != grid.probes.len()
137 || extras.iter().any(|extra| {
138 extra
139 .second
140 .iter()
141 .any(|value| !value.is_finite() || *value < 0.0)
142 || extra
143 .offset
144 .iter()
145 .any(|value| !value.is_finite() || value.abs() > spec.spacing * 0.5)
146 || extra.backfaces > 32767
147 })
148 {
149 return Err("invalid probe depth moments or state".into());
150 }
151 grid.extras = Some(extras);
152 Ok(grid)
153 }
154
155 pub fn version(&self) -> u32 {
156 if self.extras.is_some() { 2 } else { 1 }
157 }
158
159 pub fn sample(&self, position: [f32; 3], normal: [f32; 3], fallback: [f32; 3]) -> [f32; 3] {
160 let length = normal.iter().map(|v| v * v).sum::<f32>().sqrt();
161 if !length.is_finite() || length < 1e-6 || position.iter().any(|v| !v.is_finite()) {
162 return fallback;
163 }
164 let n = normal.map(|v| v / length);
165 let point: [f32; 3] = std::array::from_fn(|axis| position[axis] + n[axis] * 0.006);
166 let mut base = [0; 3];
167 let mut t = [0.0; 3];
168 let last: [f32; 3] = std::array::from_fn(|axis| (self.dims[axis] - 1) as f32);
169 let outside = |at: [f32; 3]| {
170 (0..3).any(|axis| {
171 let coord = (at[axis] - self.spec.min[axis]) / self.spec.spacing;
172 coord < 0.0 || coord > last[axis]
173 })
174 };
175 if outside(position) && outside(point) {
176 return fallback;
177 }
178 for axis in 0..3 {
179 let coord =
180 ((point[axis] - self.spec.min[axis]) / self.spec.spacing).clamp(0.0, last[axis]);
181 base[axis] = (coord.floor() as u32).min(self.dims[axis].saturating_sub(2));
182 t[axis] = if self.dims[axis] == 1 {
183 0.0
184 } else {
185 (coord - base[axis] as f32).clamp(0.0, 1.0)
186 };
187 }
188 let mut color = [0.0; 3];
189 let mut total = 0.0;
190 for corner in 0..8 {
191 let mut xyz = base;
192 let mut weight = 1.0;
193 for axis in 0..3 {
194 let high = corner & (1 << axis) != 0;
195 if high && self.dims[axis] > 1 {
196 xyz[axis] += 1;
197 }
198 weight *= if high { t[axis] } else { 1.0 - t[axis] };
199 }
200 if weight == 0.0 {
201 continue;
202 }
203 let index = GridSpec::index(self.dims, xyz);
204 let probe = &self.probes[index];
205 let extra = self.extras.as_ref().map(|extras| &extras[index]);
206 if extra.is_some_and(|extra| !extra.enabled) {
207 continue;
208 }
209 let mut center = self.spec.position(self.dims, index);
210 if let Some(extra) = extra {
211 for (axis, value) in center.iter_mut().enumerate() {
212 *value += packed(extra.offset[axis].abs() / self.spec.spacing).to_f32()
213 * self.spec.spacing
214 * extra.offset[axis].signum();
215 }
216 }
217 let toward_probe: [f32; 3] =
218 std::array::from_fn(|axis| center[axis] - position[axis] + n[axis] * 0.0001);
219 let side_length = toward_probe.iter().map(|v| v * v).sum::<f32>().sqrt();
220 let facing = ((toward_probe
221 .iter()
222 .zip(n)
223 .map(|(v, axis)| v * axis)
224 .sum::<f32>()
225 / side_length.max(1e-6))
226 + 1.0)
227 .clamp(0.0, 2.0)
228 * 0.5;
229 weight *= facing * facing * facing;
230 if weight == 0.0 {
231 continue;
232 }
233 let ray: [f32; 3] = std::array::from_fn(|axis| point[axis] - center[axis]);
234 let distance = ray.iter().map(|v| v * v).sum::<f32>().sqrt();
235 let mut irradiance = [0.0; 3];
236 let mut visibility = 0.0;
237 let mut mean = 0.0;
238 let mut second = 0.0;
239 for (lobe, axis_direction) in AXES.iter().enumerate() {
240 let nw = (0..3)
241 .map(|axis| n[axis] * axis_direction[axis])
242 .sum::<f32>()
243 .max(0.0)
244 .powi(2);
245 for (channel, value) in irradiance.iter_mut().enumerate() {
246 *value += packed(probe.lobes[lobe][channel]).to_f32() * nw;
247 }
248 let direction_weight = (0..3)
249 .map(|axis| ray[axis] * axis_direction[axis])
250 .sum::<f32>()
251 .max(0.0)
252 / distance.max(1e-6);
253 let ray_weight = direction_weight * direction_weight;
254 if let Some(extra) = extra {
255 mean += ray_weight
256 * packed(probe.visibility[lobe] / self.spec.spacing).to_f32()
257 * self.spec.spacing;
258 second += ray_weight
259 * packed(extra.second[lobe] / self.spec.spacing.powi(2)).to_f32()
260 * self.spec.spacing.powi(2);
261 } else {
262 let reach = packed(probe.visibility[lobe]).to_f32() + self.spec.spacing * 0.25;
263 visibility +=
264 ray_weight * ((reach - distance) / self.spec.spacing).clamp(0.0, 1.0);
265 }
266 }
267 if extra.is_some() {
268 let variance = (second - mean * mean).max((self.spec.spacing * 0.1).powi(2));
269 let delta = (distance - mean - self.spec.spacing * 0.1).max(0.0);
270 let chance = variance / (variance + delta * delta);
271 weight *= chance * chance;
272 } else {
273 weight *= if distance > 1e-6 { visibility } else { 1.0 };
274 }
275 for channel in 0..3 {
276 color[channel] += irradiance[channel] * weight;
277 }
278 total += weight;
279 }
280 if total <= 1e-6 {
281 fallback
282 } else {
283 color.map(|v| v / total)
284 }
285 }
286
287 pub fn bytes(&self) -> Vec<u8> {
288 let version = self.version();
289 let stride = if version == 2 { STRIDE_V2 } else { STRIDE_V1 };
290 let mut bytes = Vec::with_capacity(HEADER + self.probes.len() * stride);
291 bytes.extend_from_slice(MAGIC);
292 for value in [version, self.dims[0], self.dims[1], self.dims[2]] {
293 bytes.extend_from_slice(&value.to_le_bytes());
294 }
295 for value in self
296 .spec
297 .min
298 .into_iter()
299 .chain(self.spec.max)
300 .chain([self.spec.spacing])
301 {
302 bytes.extend_from_slice(&value.to_le_bytes());
303 }
304 bytes.resize(HEADER, 0);
305 for (index, probe) in self.probes.iter().enumerate() {
306 for lobe in probe.lobes {
307 for value in lobe.into_iter().chain([0.0]) {
308 bytes.extend_from_slice(&packed(value).to_bits().to_le_bytes());
309 }
310 }
311 for value in probe.visibility {
312 let stored = if version == 2 {
313 value / self.spec.spacing
314 } else {
315 value
316 };
317 bytes.extend_from_slice(&packed(stored).to_bits().to_le_bytes());
318 }
319 if let Some(extras) = &self.extras {
320 let extra = &extras[index];
321 for value in extra.second {
322 bytes.extend_from_slice(
323 &packed(value / self.spec.spacing.powi(2))
324 .to_bits()
325 .to_le_bytes(),
326 );
327 }
328 for value in extra.offset {
329 bytes.extend_from_slice(
330 &f16::from_f32((value / self.spec.spacing).clamp(-0.5, 0.5))
331 .to_bits()
332 .to_le_bytes(),
333 );
334 }
335 let flags = (extra.backfaces << 1) | u16::from(extra.enabled);
336 bytes.extend_from_slice(&flags.to_le_bytes());
337 }
338 }
339 bytes
340 }
341
342 pub fn from_bytes(bytes: &[u8]) -> Result<Self, String> {
343 if bytes.len() < HEADER || &bytes[..8] != MAGIC {
344 return Err("invalid probe header".into());
345 }
346 let word = |at| u32::from_le_bytes(bytes[at..at + 4].try_into().unwrap());
347 let version = word(8);
348 if !matches!(version, 1 | 2) || bytes[52..64].iter().any(|&v| v != 0) {
349 return Err("unsupported probe version".into());
350 }
351 let dims = [word(12), word(16), word(20)];
352 let float = |at| f32::from_bits(word(at));
353 let spec = GridSpec {
354 min: [float(24), float(28), float(32)],
355 max: [float(36), float(40), float(44)],
356 spacing: float(48),
357 };
358 if spec.dimensions()? != dims {
359 return Err("invalid probe dimensions".into());
360 }
361 let count = dims.iter().map(|&v| v as usize).product::<usize>();
362 let stride = if version == 2 { STRIDE_V2 } else { STRIDE_V1 };
363 if bytes.len() != HEADER + count * stride {
364 return Err("invalid probe byte length".into());
365 }
366 let mut probes = Vec::with_capacity(count);
367 let mut extras = Vec::with_capacity(count);
368 for chunk in bytes[HEADER..].chunks_exact(stride) {
369 let half = |at| f16::from_bits(u16::from_le_bytes([chunk[at], chunk[at + 1]])).to_f32();
370 let mut probe = Probe::default();
371 for lobe in 0..6 {
372 for channel in 0..3 {
373 probe.lobes[lobe][channel] = half(lobe * 8 + channel * 2);
374 }
375 let stored = half(48 + lobe * 2);
376 probe.visibility[lobe] = if version == 2 {
377 stored * spec.spacing
378 } else {
379 stored
380 };
381 }
382 probes.push(probe);
383 if version == 2 {
384 let flags = u16::from_le_bytes([chunk[78], chunk[79]]);
385 extras.push(ProbeExtra {
386 second: std::array::from_fn(|lobe| half(60 + lobe * 2) * spec.spacing.powi(2)),
387 offset: std::array::from_fn(|axis| half(72 + axis * 2) * spec.spacing),
388 enabled: flags & 1 != 0,
389 backfaces: flags >> 1,
390 });
391 }
392 }
393 if version == 2 {
394 Grid::with_extras(spec, probes, extras)
395 } else {
396 Grid::new(spec, probes)
397 }
398 }
399}
400
401#[derive(Clone, Copy, Debug, PartialEq)]
402pub struct AnchorBlend {
403 pub lower: usize,
404 pub upper: usize,
405 pub upper_weight: f32,
406}
407
408pub fn blend_anchors(hours: &[f32], hour: f32) -> Result<AnchorBlend, String> {
409 if hours.is_empty()
410 || !hour.is_finite()
411 || hours.iter().any(|v| !v.is_finite())
412 || hours.windows(2).any(|v| v[0] >= v[1])
413 {
414 return Err("anchors must be finite, sorted and distinct".into());
415 }
416 if hour <= hours[0] {
417 return Ok(AnchorBlend {
418 lower: 0,
419 upper: 0,
420 upper_weight: 0.0,
421 });
422 }
423 for upper in 1..hours.len() {
424 if hour < hours[upper] {
425 return Ok(AnchorBlend {
426 lower: upper - 1,
427 upper,
428 upper_weight: (hour - hours[upper - 1]) / (hours[upper] - hours[upper - 1]),
429 });
430 }
431 }
432 let last = hours.len() - 1;
433 Ok(AnchorBlend {
434 lower: last,
435 upper: last,
436 upper_weight: 0.0,
437 })
438}
439
440#[cfg(test)]
441mod tests {
442 use super::*;
443 #[test]
444 fn indexing_and_binary_round_trip() {
445 let spec = GridSpec {
446 min: [-1.0, 0.0, 2.0],
447 max: [1.0, 2.0, 4.0],
448 spacing: 1.0,
449 };
450 let dims = spec.dimensions().unwrap();
451 assert_eq!(dims, [3, 3, 3]);
452 assert_eq!(GridSpec::index(dims, [2, 1, 1]), 14);
453 assert_eq!(spec.position(dims, 14), [1.0, 1.0, 3.0]);
454 let grid = Grid::new(spec, vec![Probe::default(); 27]).unwrap();
455 assert_eq!(Grid::from_bytes(&grid.bytes()).unwrap(), grid);
456 assert_eq!(
457 Grid::from_bytes(&grid.bytes()).unwrap().bytes(),
458 grid.bytes()
459 );
460 assert_eq!(&grid.bytes()[..8], b"PITOPRB1");
461 }
462 #[test]
463 fn anchors_and_sampling() {
464 assert_eq!(
465 blend_anchors(&[8.0, 12.0, 16.0, 20.0], 14.0).unwrap(),
466 AnchorBlend {
467 lower: 1,
468 upper: 2,
469 upper_weight: 0.5
470 }
471 );
472 for (hour, lower) in [(9.0, 0), (20.0, 1), (22.0, 2)] {
473 let blend = blend_anchors(&[9.0, 20.0, 22.0], hour).unwrap();
474 assert_eq!((blend.lower, blend.upper_weight), (lower, 0.0));
475 }
476 let spec = GridSpec {
477 min: [0.0; 3],
478 max: [1.0; 3],
479 spacing: 1.0,
480 };
481 let mut probes = vec![Probe::default(); 8];
482 for (index, probe) in probes.iter_mut().enumerate() {
483 probe.lobes[0] = [index as f32; 3];
484 }
485 let grid = Grid::new(spec, probes).unwrap();
486 let sample = grid.sample([0.5; 3], [1.0, 0.0, 0.0], [9.0; 3]);
487 assert!((3.95..4.0).contains(&sample[0]));
488 assert_eq!(grid.sample([2.0; 3], [1.0, 0.0, 0.0], [9.0; 3]), [9.0; 3]);
489 }
490
491 #[test]
492 fn thin_wall_keeps_the_unlit_side_dark_at_each_spacing() {
493 for spacing in [0.02, 0.03, 0.05] {
494 let spec = GridSpec {
495 min: [0.0, 0.0, -spacing],
496 max: [0.0, 0.0, spacing],
497 spacing,
498 };
499 let probes = [-spacing, 0.0, spacing]
500 .map(|z| Probe {
501 lobes: [[if z > 0.0 { 1.0 } else { 0.0 }; 3]; 6],
502 ..Probe::default()
503 })
504 .to_vec();
505 let grid = Grid::new(spec, probes).unwrap();
506 assert_eq!(
507 grid.sample([0.0, 0.0, -0.001], [0.0, 0.0, -1.0], [0.0; 3]),
508 [0.0; 3]
509 );
510 assert!(grid.sample([0.0, 0.0, 0.001], [0.0, 0.0, 1.0], [0.0; 3])[0] > 0.9);
511 }
512 }
513
514 #[test]
515 fn a_surface_inside_a_thin_grid_keeps_its_probe_after_the_normal_offset() {
516 let spec = GridSpec {
517 min: [0.0, 0.01, 0.0],
518 max: [0.0, 0.01, 0.0],
519 spacing: 1.0,
520 };
521 let mut lit = Probe::default();
522 lit.lobes[2] = [0.5; 3];
523 lit.visibility = [4.0; 6];
524 let grid = Grid::new(spec, vec![lit]).unwrap();
525 assert_eq!(
526 grid.sample([0.0, 0.01, 0.0], [0.0, 1.0, 0.0], [9.0; 3]),
527 [0.5; 3]
528 );
529 assert_eq!(
530 grid.sample([0.0, 0.5, 0.0], [0.0, 1.0, 0.0], [9.0; 3]),
531 [9.0; 3]
532 );
533 }
534
535 #[test]
536 fn visibility_uses_the_probe_to_point_direction() {
537 let spec = GridSpec {
538 min: [0.0, 0.006, 0.0],
539 max: [1.0, 0.006, 0.0],
540 spacing: 1.0,
541 };
542 let mut blocked = Probe::default();
543 blocked.lobes[2] = [1.0; 3];
544 blocked.visibility[0] = 0.1;
545 let grid = Grid::new(spec, vec![blocked, Probe::default()]).unwrap();
546 assert_eq!(
547 grid.sample([0.5, 0.0, 0.0], [0.0, 1.0, 0.0], [9.0; 3]),
548 [0.0; 3]
549 );
550 }
551
552 #[test]
553 fn version_two_round_trip_keeps_probe_state_and_moments() {
554 let spec = GridSpec {
555 min: [0.0; 3],
556 max: [1.0, 0.0, 0.0],
557 spacing: 1.0,
558 };
559 let probe = Probe {
560 visibility: [4.0; 6],
561 ..Probe::default()
562 };
563 let extra = ProbeExtra {
564 second: [16.0; 6],
565 offset: [0.25, 0.0, 0.0],
566 enabled: true,
567 backfaces: 3,
568 };
569 let grid = Grid::with_extras(spec, vec![probe; 2], vec![extra; 2]).unwrap();
570 let bytes = grid.bytes();
571 assert_eq!(bytes.len(), 64 + 2 * 80);
572 assert_eq!(u32::from_le_bytes(bytes[8..12].try_into().unwrap()), 2);
573 let decoded = Grid::from_bytes(&bytes).unwrap();
574 assert_eq!(decoded.version(), 2);
575 assert_eq!(decoded.bytes(), bytes);
576 assert_eq!(decoded.extras.unwrap()[0].backfaces, 3);
577 }
578
579 #[test]
580 fn chebyshev_moments_reject_a_blocked_probe() {
581 let spec = GridSpec {
582 min: [0.0, 0.006, 0.0],
583 max: [1.0, 0.006, 0.0],
584 spacing: 1.0,
585 };
586 let mut blocked = Probe::default();
587 blocked.lobes[2] = [1.0; 3];
588 blocked.visibility = [0.1; 6];
589 let open = Probe {
590 visibility: [4.0; 6],
591 ..Probe::default()
592 };
593 let visible = ProbeExtra {
594 second: [0.01; 6],
595 offset: [0.0; 3],
596 enabled: true,
597 backfaces: 0,
598 };
599 let clear = ProbeExtra {
600 second: [16.0; 6],
601 ..visible.clone()
602 };
603 let grid = Grid::with_extras(
604 spec,
605 vec![blocked.clone(), open.clone()],
606 vec![visible.clone(), clear.clone()],
607 )
608 .unwrap();
609 let sample = grid.sample([0.5, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0; 3]);
610 assert!(sample[0] < 0.01);
611 let disabled = Grid::with_extras(
612 spec,
613 vec![blocked, open],
614 vec![
615 ProbeExtra {
616 enabled: false,
617 ..visible
618 },
619 clear,
620 ],
621 )
622 .unwrap();
623 assert_eq!(
624 disabled.sample([0.5, 0.0, 0.0], [0.0, 1.0, 0.0], [0.0; 3]),
625 [0.0; 3]
626 );
627 }
628}