arendur/geometry/
bbox.rs

1// Copyright 2017 Dasein Phaos aka. Luxko
2//
3// Licensed under the Apache License, Version 2.0 <LICENSE-APACHE or
4// http://www.apache.org/licenses/LICENSE-2.0> or the MIT license
5// <LICENSE-MIT or http://opensource.org/licenses/MIT>, at your
6// option. This file may not be copied, modified, or distributed
7// except according to those terms.
8
9//! 2D and 3D bounding box
10use super::float;
11use super::foundamental::*;
12use std::ops;
13use std::mem;
14use super::ray::Ray;
15use num_traits::NumCast;
16
17pub type BBox2f = BBox2<Float>;
18pub type BBox3f = BBox3<Float>;
19
20
21/// A 2D bounding box
22#[derive(PartialEq, Copy, Clone, Serialize, Deserialize, Debug)]
23#[must_use]
24pub struct BBox2<T> {
25    /// min corner of the bounding box
26    pub pmin: Point2<T>,
27    /// max corner of the bounding box
28    pub pmax: Point2<T>,
29}
30
31impl<T: BaseNum> BBox2<T> {
32    /// Construct a new bounding box marked by two corner vertice
33    #[inline]
34    pub fn new(p: Point2<T>, q: Point2<T>) -> BBox2<T> {
35        BBox2{
36            pmin: Point2::new(
37                <T as  PartialOrd>::partial_min(p.x, q.x),
38                <T as  PartialOrd>::partial_min(p.y, q.y),
39            ),
40            pmax: Point2::new(
41                <T as  PartialOrd>::partial_max(p.x, q.x),
42                <T as  PartialOrd>::partial_max(p.y, q.y),
43            )
44        }
45    }
46
47    /// Return the `i`th corner vertex
48    pub fn corner(&self, i: usize) -> Point2<T> {
49        assert!(i < 4, "index out of bound");
50        let x = if (i & 1) == 0 {
51            self.pmin.x
52        } else {
53            self.pmax.x
54        };
55
56        let y = if (i & 2) == 0 {
57            self.pmin.y
58        } else {
59            self.pmax.y
60        };
61
62        Point2::new(x, y)
63    }
64
65    /// Extend the bounding box with `p`, return the resultant new bbox
66    #[inline]
67    pub fn extend(&self, p: Point2<T>) -> Self {
68        BBox2{
69            pmin: Point2::new(
70                <T as  PartialOrd>::partial_min(self.pmin.x, p.x),
71                <T as  PartialOrd>::partial_min(self.pmin.y, p.y),
72            ),
73            pmax: Point2::new(
74                <T as  PartialOrd>::partial_max(self.pmax.x, p.x),
75                <T as  PartialOrd>::partial_max(self.pmax.y, p.y),
76            )        
77        }
78    }
79
80    /// Return the union of two bounding boxes
81    #[inline]
82    pub fn union(&self, other: &Self) -> Self {
83        BBox2{
84            pmin: Point2::new(
85                <T as  PartialOrd>::partial_min(self.pmin.x, other.pmin.x),
86                <T as  PartialOrd>::partial_min(self.pmin.y, other.pmin.y),
87            ),
88            pmax: Point2::new(
89                <T as  PartialOrd>::partial_max(self.pmax.x, other.pmax.x),
90                <T as  PartialOrd>::partial_max(self.pmax.y, other.pmax.y),
91            )        
92        }
93    }
94
95    /// Return the intersection of two bounding boxes
96    #[inline]
97    pub fn intersect(&self, other: &Self) -> Option<Self> {
98        let ret = BBox2 {
99            pmin: Point2::new(
100                <T as  PartialOrd>::partial_max(self.pmin.x, other.pmin.x),
101                <T as  PartialOrd>::partial_max(self.pmin.y, other.pmin.y),
102            ),
103            pmax: Point2::new(
104                <T as  PartialOrd>::partial_min(self.pmax.x, other.pmax.x),
105                <T as  PartialOrd>::partial_min(self.pmax.y, other.pmax.y),
106            )        
107        };
108        if ret.pmin.x > ret.pmax.x || ret.pmin.y > ret.pmax.y {
109            None
110        } else {
111            Some(ret)
112        }
113    }
114
115    /// Return if two bounding boxes overlap
116    #[inline]
117    pub fn overlap(&self, other: &Self) -> bool {
118        (self.pmax.x >= other.pmin.x && self.pmin.x <= other.pmax.x)
119        && (self.pmax.y >= other.pmin.y && self.pmin.y <= other.pmax.y)
120    }
121
122    /// Return if `self` contains `p`
123    #[inline]
124    pub fn contain(&self, p: Point2<T>) -> bool {
125        (p.x <= self.pmax.x && p.x >= self.pmin.x)
126        && (p.y <= self.pmax.y && p.y >= self.pmin.y)
127    }
128
129    /// Return if `self` contains `p`, excluding the case when `p` lies on
130    /// the upper-left boundary
131    #[inline]
132    pub fn contain_lb(&self, p: Point2<T>) -> bool {
133        (p.x < self.pmax.x && p.x >= self.pmin.x)
134        && (p.y < self.pmax.y && p.y >= self.pmin.y)
135    }
136
137    /// Expand each boundary by `delta`. Note that `delta` can be negative
138    #[inline]
139    pub fn expand_by(&self, delta: T) -> Self
140        where T: ops::Neg<Output = T> {
141        BBox2 {
142            pmin: self.pmin + (-Vector2::new(delta, delta)),
143            pmax: self.pmax + Vector2::new(delta, delta),
144        }
145    }
146
147    /// Expand each boundary by `delta`. Note that `delta` can be negative
148    #[inline]
149    pub fn expand_by_vec(&self, delta: Vector2<T>) -> Self
150        where T: ops::Neg<Output = T> {
151        BBox2 {
152            pmin: self.pmin + (-delta),
153            pmax: self.pmax + delta,
154        }
155    }
156
157    /// Return the diagonal vector, from `pmin` to `pmax`
158    #[inline]
159    pub fn diagonal(&self) -> Vector2<T> {
160        self.pmax - self.pmin
161    }
162
163    /// Return the surface area of the bounding box
164    #[inline]
165    pub fn surface_area(&self) -> T {
166        let dx = self.pmax.x - self.pmin.x;
167        let dy = self.pmax.y - self.pmin.y;
168        let zero = <T as Zero>::zero();
169        if dx >= zero && dy >= zero {
170            dx * dy
171        } else {
172            zero
173        }
174    }
175
176    /// Returns the index of the axis along which the bounding
177    /// box has maximum extent
178    pub fn max_extent(&self) -> usize {
179        let d = self.diagonal();
180        if d.x > d.y {
181            0usize
182        } else {
183            1usize
184        }
185    }
186
187    /// Linearly interpolate between the two corners
188    pub fn lerp(&self, t: Vector2<Float>) -> Point2<T>
189        where T: NumCast,
190    {
191        Point2::new(
192            <T as NumCast>::from(
193                t.x * <Float as NumCast>::from(self.pmax.x).unwrap()
194                + (1.0 as Float - t.x) * <Float as NumCast>::from(self.pmin.x).unwrap()
195            ).unwrap(),
196            <T as NumCast>::from(
197                t.y * <Float as NumCast>::from(self.pmax.y).unwrap()
198                + (1.0 as Float - t.y) * <Float as NumCast>::from(self.pmin.y).unwrap()
199            ).unwrap()
200        )
201    }
202
203    /// Returns the bounding circle as `(center, radius)`
204    pub fn bcircle(&self) -> (Point2<T>, T)
205        where T: BaseFloat
206    {
207        let two = <T as One>::one() + <T as One>::one();
208        let zero = <T as Zero>::zero();
209        let pmin: (T, T) = self.pmin.into();
210        let pmin: Vector2<T> = pmin.into();
211        let center = (self.pmax + pmin) / two;
212        let radius = if self.contain(center) {
213            center.distance(self.pmax)
214        } else {
215            zero
216        };
217
218        debug_assert!(radius >= zero);
219        (center, radius)
220    }
221
222    /// Casting to another type of bbox
223    pub fn cast<R: BaseNum>(&self) -> BBox2<R> {
224        BBox2{
225            pmin: self.pmin.cast(),
226            pmax: self.pmax.cast(),
227        }
228    }
229}
230
231#[derive(PartialEq, Eq, Copy, Clone, Hash, Serialize, Deserialize, Debug)]
232#[must_use]
233pub struct BBox3<T> {
234    /// min corner of the bounding box
235    pub pmin: Point3<T>,
236    /// max corner of the bounding box
237    pub pmax: Point3<T>,
238}
239
240/// A 3D bounding box
241impl<T: BaseNum> BBox3<T> {
242    /// Construct a new bounding box marked by two corner vertice
243    #[inline]
244    pub fn new(p: Point3<T>, q: Point3<T>) -> BBox3<T> {
245        BBox3{
246            pmin: Point3::new(
247                <T as  PartialOrd>::partial_min(p.x, q.x),
248                <T as  PartialOrd>::partial_min(p.y, q.y),
249                <T as  PartialOrd>::partial_min(p.z, q.z),
250            ),
251            pmax: Point3::new(
252                <T as  PartialOrd>::partial_max(p.x, q.x),
253                <T as  PartialOrd>::partial_max(p.y, q.y),
254                <T as  PartialOrd>::partial_max(p.z, q.z),
255            )
256        }
257    }
258
259    /// Return the `i`th corner vertex
260    pub fn corner(&self, i: usize) -> Point3<T> {
261        assert!(i < 4, "index out of bound");
262        let x = if (i & 1) == 0 {
263            self.pmin.x
264        } else {
265            self.pmax.x
266        };
267
268        let y = if (i & 2) == 0 {
269            self.pmin.y
270        } else {
271            self.pmax.y
272        };
273
274        let z = if (i & 3) == 0 {
275            self.pmin.z
276        } else {
277            self.pmax.z
278        };
279
280        Point3::new(x, y, z)
281    }
282
283    /// `index==false` for `self.pmin`, else `self.pmax`.
284    #[inline]
285    fn index(&self, index: bool) -> Point3<T> {
286        if index {
287            self.pmax
288        } else {
289            self.pmin
290        }
291    }
292
293    /// Extend the bounding box with `p`, return the resultant new bbox
294    #[inline]
295    pub fn extend(&self, p: Point3<T>) -> Self {
296        BBox3{
297            pmin: Point3::new(
298                <T as  PartialOrd>::partial_min(self.pmin.x, p.x),
299                <T as  PartialOrd>::partial_min(self.pmin.y, p.y),
300                <T as  PartialOrd>::partial_min(self.pmin.z, p.z),
301            ),
302            pmax: Point3::new(
303                <T as  PartialOrd>::partial_max(self.pmax.x, p.x),
304                <T as  PartialOrd>::partial_max(self.pmax.y, p.y),
305                <T as  PartialOrd>::partial_max(self.pmax.z, p.z),
306            )        
307        }
308    }
309
310    /// Return the union of two bounding boxes
311    #[inline]
312    pub fn union(&self, other: &Self) -> Self {
313        BBox3{
314            pmin: Point3::new(
315                <T as  PartialOrd>::partial_min(self.pmin.x, other.pmin.x),
316                <T as  PartialOrd>::partial_min(self.pmin.y, other.pmin.y),
317                <T as  PartialOrd>::partial_min(self.pmin.z, other.pmin.z),
318            ),
319            pmax: Point3::new(
320                <T as  PartialOrd>::partial_max(self.pmax.x, other.pmax.x),
321                <T as  PartialOrd>::partial_max(self.pmax.y, other.pmax.y),
322                <T as  PartialOrd>::partial_max(self.pmax.z, other.pmax.z),
323            )        
324        }
325    }
326
327    /// Return the intersection of two bounding boxes
328    #[inline]
329    pub fn intersect(&self, other: &Self) -> Option<Self> {
330        let ret = BBox3{
331            pmin: Point3::new(
332                <T as  PartialOrd>::partial_max(self.pmin.x, other.pmin.x),
333                <T as  PartialOrd>::partial_max(self.pmin.y, other.pmin.y),
334                <T as  PartialOrd>::partial_max(self.pmin.z, other.pmin.z),
335            ),
336            pmax: Point3::new(
337                <T as  PartialOrd>::partial_min(self.pmax.x, other.pmax.x),
338                <T as  PartialOrd>::partial_min(self.pmax.y, other.pmax.y),
339                <T as  PartialOrd>::partial_min(self.pmax.z, other.pmax.z),
340            )        
341        };
342        if ret.pmin.x > ret.pmax.x || ret.pmin.y > ret.pmax.y || ret.pmin.z > ret.pmax.z {
343            None
344        } else {
345            Some(ret)
346        }
347    }
348
349    /// Return if two bounding boxes overlap
350    #[inline]
351    pub fn overlap(&self, other: &Self) -> bool {
352        (self.pmax.x >= other.pmin.x && self.pmin.x <= other.pmax.x)
353        && (self.pmax.y >= other.pmin.y && self.pmin.y <= other.pmax.y)
354        && (self.pmax.z >= other.pmin.z && self.pmin.z <= other.pmax.z)
355    }
356
357    /// Return if `self` contains `p`
358    #[inline]
359    pub fn contain(&self, p: Point3<T>) -> bool {
360        (p.x <= self.pmax.x && p.x >= self.pmin.x)
361        && (p.y <= self.pmax.y && p.y >= self.pmin.y)
362        && (p.z <= self.pmax.z && p.z >= self.pmin.z)
363    }
364
365    /// Return if `self` contains `p`, excluding the case when `p` lies on
366    /// the greater boundaries
367    #[inline]
368    pub fn contain_lb(&self, p: Point3<T>) -> bool {
369        (p.x < self.pmax.x && p.x >= self.pmin.x)
370        && (p.y < self.pmax.y && p.y >= self.pmin.y)
371        && (p.z < self.pmax.z && p.z >= self.pmin.z)
372    }
373
374    /// Expand each boundary by `delta`. Note that `delta` can be negative
375    #[inline]
376    pub fn expand_by(&self, delta: T) -> Self
377        where T: ops::Neg<Output = T> {
378        BBox3 {
379            pmin: self.pmin + (-Vector3::new(delta, delta, delta)),
380            pmax: self.pmax + Vector3::new(delta, delta, delta),
381        }
382    }
383
384    /// Expand each boundary by `delta`. Note that `delta` can be negative
385    #[inline]
386    pub fn expand_by_vec(&self, delta: Vector3<T>) -> Self
387        where T: ops::Neg<Output = T> {
388        BBox3 {
389            pmin: self.pmin + (-delta),
390            pmax: self.pmax + delta,
391        }
392    }
393
394    /// Return the diagonal vector, from `pmin` to `pmax`
395    #[inline]
396    pub fn diagonal(&self) -> Vector3<T> {
397        self.pmax - self.pmin
398    }
399
400    /// Return the surface area of the bounding box
401    #[inline]
402    pub fn surface_area(&self) -> T {
403        let mut dx = self.pmax.x - self.pmin.x;
404        let mut dy = self.pmax.y - self.pmin.y;
405        let mut dz = self.pmax.z - self.pmin.z;
406        let zero = <T as Zero>::zero();
407        
408        if dx < zero { dx = zero; }
409        if dy < zero { dy = zero; }
410        if dz < zero { dz = zero; }
411        let two = <T as One>::one() + <T as One>::one();
412        two * (dx * dy + dx * dz + dy * dz)
413    }
414
415    /// Renturn the volume of the bounding box
416    #[inline]
417    pub fn volume(&self) -> T {
418        let dx = self.pmax.x - self.pmin.x;
419        let dy = self.pmax.y - self.pmin.y;
420        let dz = self.pmax.z - self.pmin.z;
421        let zero = <T as Zero>::zero();
422        
423        if dx <= zero || dy <= zero || dz <= zero {
424            zero
425        } else {
426            dx * dy * dz
427        }
428    }
429
430    /// Returns the index of the axis along which the bounding
431    /// box has maximum extent
432    pub fn max_extent(&self) -> usize {
433        let d = self.diagonal();
434        if d.x > d.y && d.x > d.z {
435            0usize
436        } else if d.y > d.z {
437            1usize
438        } else {
439            2usize
440        }
441    }
442
443    /// Linearly interpolate between the two corners
444    pub fn lerp(&self, t: Vector3f) -> Point3<T>
445         where T: NumCast,
446    {
447        Point3::new(
448            <T as NumCast>::from(
449                t.x * <Float as NumCast>::from(self.pmax.x).unwrap()
450                + (1.0 as Float - t.x) * <Float as NumCast>::from(self.pmin.x).unwrap()
451            ).unwrap(),
452            <T as NumCast>::from(
453                t.y * <Float as NumCast>::from(self.pmax.y).unwrap()
454                + (1.0 as Float - t.y) * <Float as NumCast>::from(self.pmin.y).unwrap()
455            ).unwrap(),
456            <T as NumCast>::from(
457                t.z * <Float as NumCast>::from(self.pmax.z).unwrap()
458                + (1.0 as Float - t.z) * <Float as NumCast>::from(self.pmin.z).unwrap()
459            ).unwrap(),
460        )
461    }
462
463    pub fn bsphere(&self) -> (Point3<T>, T)
464        where T: BaseFloat
465    {
466        let two = <T as One>::one() + <T as One>::one();
467        let zero = <T as Zero>::zero();
468        let pmin = self.pmin.to_vec();
469        let center = (self.pmax + pmin) / two;
470        let radius = if self.contain(center) {
471            center.distance(self.pmax)
472        } else {
473            zero
474        };
475
476        debug_assert!(radius >= zero);
477        (center, radius)
478    }
479    
480    /// Apply transform `t` on `self`, returning a new bounding box
481    pub fn apply_transform<Tr>(&self, t: &Tr) -> Self
482        where Tr: Transform3<T>
483    {
484        // let bbox = BBox3::new(
485        //     t.transform_point(Point3::new(self.pmin.x, self.pmin.y, self.pmin.z)),
486        //     t.transform_point(Point3::new(self.pmax.x, self.pmin.y, self.pmin.z))
487        // );
488        // bbox.extend(Point3::new(self.pmin.x, self.pmax.y, self.pmin.z))
489        //     .extend(Point3::new(self.pmin.x, self.pmin.y, self.pmax.z))
490        //     .extend(Point3::new(self.pmin.x, self.pmax.y, self.pmax.z))
491        //     .extend(Point3::new(self.pmax.x, self.pmin.y, self.pmax.z))
492        //     .extend(Point3::new(self.pmax.x, self.pmax.y, self.pmin.z))
493        //     .extend(Point3::new(self.pmax.x, self.pmax.y, self.pmax.z))
494        let p = t.transform_point(self.pmin);
495        let diagonal = t.transform_vector(self.diagonal());
496        BBox3::new(
497            p, p + diagonal
498        )
499    }
500
501    /// Casting to another type of bbox
502    pub fn cast<R: BaseNum>(&self) -> BBox3<R> {
503        BBox3{
504            pmin: self.pmin.cast(),
505            pmax: self.pmax.cast(),
506        }
507    }
508}
509
510impl BBox3f {
511    /// Test if the `ray` intersects `self`
512    pub fn intersect_ray<R>(&self, ray: &R) -> Option<(Float, Float)>
513        where R: Ray + ?Sized
514    {
515        let mut t0 = 0.0 as Float;
516        let mut t1 = ray.max_extend();
517
518        let origin = ray.origin();
519        let direction = ray.direction();
520
521        for i in 0..3 {
522            let inv_direction = (1.0 as Float) / direction[i];
523            let mut t_near = (self.pmin[i] - origin[i]) * inv_direction;
524            let mut t_far = (self.pmax[i] - origin[i]) * inv_direction;
525            if t_near > t_far {
526                mem::swap(&mut t_near, &mut t_far);
527            }
528            
529            // Update to ensure robust ray-bounds intersection
530            t_far *= 1. as Float + 2. as Float * float::eb_term(3. as Float);
531
532            if t_near > t0 {
533                t0 = t_near;
534            }
535            if t_far < t1 {
536                t1 = t_far;
537            }
538
539            if t0 > t1 {
540                return None;
541            }
542        }
543
544        Some((t0, t1))
545    }
546
547    /// Test if the `ray` intersects `self`, with cache
548    #[inline]
549    pub fn intersect_ray_cached(&self, cache: &(Point3f, Vector3f, Vector3<bool>, Float)) -> Option<(Float, Float)>
550    {
551        let mut t0 = (self.index(cache.2.x).x - cache.0.x) * cache.1.x;
552        let mut t1 = (self.index(!cache.2.x).x - cache.0.x) * cache.1.x;
553
554        let ty0 = (self.index(cache.2.y).y - cache.0.y) * cache.1.y;
555        let mut ty1 = (self.index(!cache.2.y).y - cache.0.y) * cache.1.y;
556
557        // update for conservative intersection
558        t1 *= 1. as Float + 2. as Float * float::eb_term(3. as Float);
559        ty1 *= 1. as Float + 2. as Float * float::eb_term(3. as Float);
560
561        if t0 > ty1 || ty0 > t1 { return None; }
562        if ty0 > t0 { t0 = ty0; }
563        if ty1 < t1 { t1 = ty1; }
564
565        let tz0 = (self.index(cache.2.z).z - cache.0.z) * cache.1.z;
566        let mut tz1 = (self.index(!cache.2.z).z - cache.0.z) * cache.1.z;
567
568        // update for robustness
569        tz1 *= 1. as Float + 2. as Float * float::eb_term(3. as Float);
570
571        if t0 > tz1 || tz0 > t1 { return None; }
572        if tz0 > t0 { t0 = tz0; }
573        if tz1 < t1 { t1 = tz1; }
574
575        if t0 < cache.3 && t1 > (0.0 as Float) {
576            Some((t0, t1))
577        } else {
578            None
579        }
580    }
581
582    /// Construct cache for chached intersection
583    pub fn construct_ray_cache<R>(ray: &R) -> (Point3f, Vector3f, Vector3<bool>, Float)
584        where R: Ray + ?Sized
585    {
586        let origin = ray.origin();
587        let invert_direction = 1.0 as Float / ray.direction();
588        let zero = 0.0 as Float;
589        let dir_is_neg = Vector3::new(invert_direction.x < zero, invert_direction.y < zero, invert_direction.z < zero);
590        let max_extend = ray.max_extend();
591        (origin, invert_direction, dir_is_neg, max_extend)
592    }
593}
594
595pub struct BBox2iIter {
596    ix: isize,
597    iy: isize,
598    nx: isize,
599    ny: isize,
600    nx_start: isize,
601}
602
603impl Iterator for BBox2iIter {
604    type Item = Point2<isize>;
605
606    #[inline]
607    fn next(&mut self) -> Option<Point2<isize>> {
608        while self.iy < self.ny {
609            if self.ix < self.nx {
610                let ix = self.ix;
611                self.ix += 1;
612                return Some(Point2::new(ix, self.iy))
613            } else {
614                self.iy += 1;
615                self.ix = self.nx_start;
616            }
617        }
618        None
619    }
620}
621
622impl IntoIterator for BBox2<isize> {
623    type Item = Point2<isize>;
624    type IntoIter = BBox2iIter;
625    fn into_iter(self) -> BBox2iIter {
626        BBox2iIter{
627            ix: self.pmin.x,
628            iy: self.pmin.y,
629            nx: self.pmax.x,
630            ny: self.pmax.y,
631            nx_start: self.pmin.x,
632        }
633    }
634}