diff --git a/lib/common/src/bvh.rs b/lib/common/src/bvh.rs index c5345c3..11e4694 100644 --- a/lib/common/src/bvh.rs +++ b/lib/common/src/bvh.rs @@ -1,4 +1,10 @@ +use std::cmp::Ordering; +use std::collections::BTreeMap; + use crate::aabb::Aabb; +use crate::ray::Ray; +use crate::integer::{Ratio,Planar64}; +use crate::instruction::{InstructionCollector,TimedInstruction}; //da algaritum //lista boxens @@ -10,6 +16,96 @@ use crate::aabb::Aabb; //sort the centerpoints on each axis (3 lists) //bv is put into octant based on whether it is upper or lower in each list + +pub fn intersect_aabb(ray:&Ray,aabb:&Aabb)->Option<Ratio<Planar64,Planar64>>{ + // n.(o+d*t)==n.p + // n.o + n.d * t == n.p + // t == (n.p - n.o)/n.d + let mut hit=None; + match ray.direction.x.cmp(&Planar64::ZERO){ + Ordering::Less=>{ + let rel_min=aabb.min()-ray.origin; + let rel_max=aabb.max()-ray.origin; + let dy=rel_max.x*ray.direction.y; + let dz=rel_max.x*ray.direction.z; + // x is negative, so inequalities are flipped + if rel_min.y*ray.direction.x>dy&&dy>rel_max.y*ray.direction.x + &&rel_min.z*ray.direction.x>dz&&dz>rel_max.z*ray.direction.x{ + let t=rel_max.x/ray.direction.x; + hit=Some(hit.map_or(t,|best_t|t.min(best_t))); + } + }, + Ordering::Equal=>(), + Ordering::Greater=>{ + let rel_min=aabb.min()-ray.origin; + let rel_max=aabb.max()-ray.origin; + let dy=rel_min.x*ray.direction.y; + let dz=rel_min.x*ray.direction.z; + // x is positive, so inequalities are normal + if rel_min.y*ray.direction.x<dy&&dy<rel_max.y*ray.direction.x + &&rel_min.z*ray.direction.x<dz&&dz<rel_max.z*ray.direction.x{ + let t=rel_min.x/ray.direction.x; + hit=Some(hit.map_or(t,|best_t|t.min(best_t))); + } + }, + } + match ray.direction.z.cmp(&Planar64::ZERO){ + Ordering::Less=>{ + let rel_min=aabb.min()-ray.origin; + let rel_max=aabb.max()-ray.origin; + let dx=rel_max.z*ray.direction.x; + let dy=rel_max.z*ray.direction.y; + // z is negative, so inequalities are flipped + if rel_min.x*ray.direction.z>dx&&dx>rel_max.x*ray.direction.z + &&rel_min.y*ray.direction.z>dy&&dy>rel_max.y*ray.direction.z{ + let t=rel_max.z/ray.direction.z; + hit=Some(hit.map_or(t,|best_t|t.min(best_t))); + } + }, + Ordering::Equal=>(), + Ordering::Greater=>{ + let rel_min=aabb.min()-ray.origin; + let rel_max=aabb.max()-ray.origin; + let dx=rel_min.z*ray.direction.x; + let dy=rel_min.z*ray.direction.y; + // z is positive, so inequalities are normal + if rel_min.x*ray.direction.z<dx&&dx<rel_max.x*ray.direction.z + &&rel_min.y*ray.direction.z<dy&&dy<rel_max.y*ray.direction.z{ + let t=rel_min.z/ray.direction.z; + hit=Some(hit.map_or(t,|best_t|t.min(best_t))); + } + }, + } + match ray.direction.y.cmp(&Planar64::ZERO){ + Ordering::Less=>{ + let rel_min=aabb.min()-ray.origin; + let rel_max=aabb.max()-ray.origin; + let dz=rel_max.y*ray.direction.z; + let dx=rel_max.y*ray.direction.x; + // y is negative, so inequalities are flipped + if rel_min.z*ray.direction.y>dz&&dz>rel_max.z*ray.direction.y + &&rel_min.x*ray.direction.y>dx&&dx>rel_max.x*ray.direction.y{ + let t=rel_max.y/ray.direction.y; + hit=Some(hit.map_or(t,|best_t|t.min(best_t))); + } + }, + Ordering::Equal=>(), + Ordering::Greater=>{ + let rel_min=aabb.min()-ray.origin; + let rel_max=aabb.max()-ray.origin; + let dz=rel_min.y*ray.direction.z; + let dx=rel_min.y*ray.direction.x; + // y is positive, so inequalities are normal + if rel_min.z*ray.direction.y<dz&&dz<rel_max.z*ray.direction.y + &&rel_min.x*ray.direction.y<dx&&dx<rel_max.x*ray.direction.y{ + let t=rel_min.y/ray.direction.y; + hit=Some(hit.map_or(t,|best_t|t.min(best_t))); + } + }, + } + hit +} + pub enum RecursiveContent<N,L>{ Branch(Vec<N>), Leaf(L), @@ -44,6 +140,92 @@ impl<L> BvhNode<L>{ }, } } + fn populate_nodes<'a,T,F>( + &'a self, + collector:&mut InstructionCollector<&'a L,Ratio<Planar64,Planar64>>, + nodes:&mut BTreeMap<Ratio<Planar64,Planar64>,&'a BvhNode<L>>, + ray:&Ray, + start_time:Ratio<Planar64,Planar64>, + f:&F, + ) + where + T:Ord+Copy, + Ratio<Planar64,Planar64>:From<T>, + F:Fn(&L,&Ray)->Option<T>, + { + match &self.content{ + RecursiveContent::Leaf(leaf)=>if let Some(time)=f(leaf,ray){ + let ins=TimedInstruction{time:time.into(),instruction:leaf}; + if start_time.lt_ratio(ins.time)&&ins.time.lt_ratio(collector.time()){ + collector.collect(Some(ins)); + } + }, + RecursiveContent::Branch(children)=>for child in children{ + if child.aabb.contains(ray.origin){ + child.populate_nodes(collector,nodes,ray,start_time,f); + }else{ + // Am I an upcoming superstar? + if let Some(t)=intersect_aabb(ray,&child.aabb){ + if start_time.lt_ratio(t)&&t.lt_ratio(collector.time()){ + nodes.insert(t,child); + } + } + } + }, + } + } + pub fn sample_ray<T,F>( + &self, + ray:&Ray, + start_time:T, + time_limit:T, + f:F, + )->Option<(T,&L)> + where + T:Ord+Copy, + T:From<Ratio<Planar64,Planar64>>, + Ratio<Planar64,Planar64>:From<T>, + F:Fn(&L,&Ray)->Option<T>, + { + // source of nondeterminism when Aabb boundaries are coplanar + let mut nodes=BTreeMap::new(); + + let start_time=start_time.into(); + let time_limit=time_limit.into(); + let mut collector=InstructionCollector::new(time_limit); + // break open all nodes that contain ray.origin and populate nodes with future intersection times + self.populate_nodes(&mut collector,&mut nodes,ray,start_time,&f); + + // swim through nodes one at a time + while let Some((t,node))=nodes.pop_first(){ + if collector.time()<t{ + break; + } + match &node.content{ + RecursiveContent::Leaf(leaf)=>if let Some(time)=f(leaf,ray){ + let ins=TimedInstruction{time:time.into(),instruction:leaf}; + // this lower bound can also be omitted + // but it causes type inference errors lol + if start_time.lt_ratio(ins.time)&&ins.time.lt_ratio(collector.time()){ + collector.collect(Some(ins)); + } + }, + // break open the node and predict collisions with the child nodes + RecursiveContent::Branch(children)=>for child in children{ + // Am I an upcoming superstar? + if let Some(t)=intersect_aabb(ray,&child.aabb){ + // we don't need to check the lower bound + // because child aabbs are guaranteed to be within the parent bounds. + if t<collector.time(){ + nodes.insert(t,child); + } + } + }, + } + } + + collector.take().map(|TimedInstruction{time,instruction:leaf}|(time.into(),leaf)) + } pub fn into_inner(self)->(RecursiveContent<BvhNode<L>,L>,Aabb){ (self.content,self.aabb) }