//! 3D volumetric forms

use derive_more::{From, TryInto};
#[cfg(feature = "derive_serdes")]
use serde::{Deserialize, Serialize};

use crate::*;

use super::{Aabb3, Capsule3, Cone3, Cylinder3, Hull3, Plane3, Primitive, Sphere3};

////////////////////////////////////////////////////////////////////////////////////////
//  traits                                                                            //
////////////////////////////////////////////////////////////////////////////////////////

/// A trait for bounded and unbounded shapes.
///
/// A 'Shape' implements the following traits:
///
/// - 'Aabb' -- a shape can compute the axis-aligned bounding volume that encloses it
/// - 'Bsphere' -- TODO: a shape can compute its bounding sphere
// TODO: do we want a Stereometric bound here? the shapes in this module all have
// non-zero volume, but if we want to implement Shape for primitives like triangles
// then the volume may be zero
pub trait Shape <S> : Aabb <S> /*+ Bsphere <S>*/ { // TODO: bspheres
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>>;
}

/// Trait for computing the axis-aligned bounding volume of a given shape.
///
/// Note for 'Unbounded' shapes some of the components will contain positive or negative
/// infinity.
pub trait Aabb <S> {
  fn aabb (&self) -> Aabb3 <S>;
}

// TODO: implement bounding spheres
/// A trait for computing bounding spheres
pub trait Bsphere <S> {
  fn sphere (&self) -> Sphere3 <S>;
}

/// Trait for computing volumes of solid figures.
///
/// Note for 'Unbounded' shapes this will be positive infinitiy.
pub trait Stereometric <S> {
  fn volume (&self) -> Positive <S>;
}

////////////////////////////////////////////////////////////////////////////////////////
//  enums                                                                             //
////////////////////////////////////////////////////////////////////////////////////////

/// A shape that is either bounded (encloses a finite volume) or unbounded (delimits an
/// infinite volume)
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Debug, Eq, PartialEq, From, TryInto)]
pub enum Variant <S> {
  Bounded   (Bounded   <S>),
  Unbounded (Unbounded <S>)
}

/// A (totally) bounded shape
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Debug, Eq, PartialEq, From, TryInto)]
pub enum Bounded <S> {
  Sphere   (Sphere   <S>),
  Capsule  (Capsule  <S>),
  Cylinder (Cylinder <S>),
  Cone     (Cone     <S>),
  Cube     (Cube     <S>),
  Cuboid   (Cuboid   <S>),
  Hull     (Hull     <S>)
}

/// A shape that is only partly bounded
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Debug, Eq, PartialEq, From, TryInto)]
pub enum Unbounded <S> {
  Orthant   (Orthant),
  Halfspace (Halfspace <S>)
}

////////////////////////////////////////////////////////////////////////////////////////
//  structs                                                                           //
////////////////////////////////////////////////////////////////////////////////////////

/// An axis-aligned halfspace
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Orthant {
  pub normal_axis : SignedAxis3
}

/// A halfspace defined by a normal vector
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Halfspace <S> {
  pub normal : Unit3 <S>
}

/// A sphere defined by a positive radius
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Sphere <S> {
  pub radius : Positive <S>
}

/// A capsule defined by strictly positive radius and half-height
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Capsule <S> {
  pub radius      : Positive <S>,
  pub half_height : Positive <S>
}

/// A cylinder defined by strictly positive radius and half-height
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Cylinder <S> {
  pub radius      : Positive <S>,
  pub half_height : Positive <S>
}

/// A cone defined by strictly positive radius and half-height
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Cone <S> {
  pub radius      : Positive <S>,
  pub half_height : Positive <S>
}

/// A cube defined by a strictly positive extent (equal to the distance from the cube
/// center to the edges)
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Cube <S> {
  pub extent : Positive <S>
}

/// A box defined by three strictly positive half extents
#[cfg_attr(feature = "derive_serdes", derive(Serialize, Deserialize))]
#[derive(Clone, Copy, Debug, Eq, PartialEq)]
pub struct Cuboid <S> {
  pub extents : Vector3 <Positive <S>>
}

/// A convex hull
pub type Hull <S> = Hull3 <S>;

////////////////////////////////////////////////////////////////////////////////////////
//  functions                                                                         //
////////////////////////////////////////////////////////////////////////////////////////

pub fn report_sizes() {
  use std::mem::size_of;
  println!("shape report sizes...");

  show!(size_of::<Variant   <f32>>());
  show!(size_of::<Variant   <f64>>());

  show!(size_of::<Bounded   <f32>>());
  show!(size_of::<Bounded   <f64>>());

  show!(size_of::<Unbounded <f32>>());
  show!(size_of::<Unbounded <f64>>());

  show!(size_of::<Sphere    <f32>>());
  show!(size_of::<Capsule   <f32>>());
  show!(size_of::<Cylinder  <f32>>());
  show!(size_of::<Cone      <f32>>());
  show!(size_of::<Cube      <f32>>());
  show!(size_of::<Cuboid    <f32>>());

  show!(size_of::<Orthant>());
  show!(size_of::<Halfspace <f32>>());

  show!(size_of::<Sphere    <f64>>());
  show!(size_of::<Capsule   <f64>>());
  show!(size_of::<Cylinder  <f64>>());
  show!(size_of::<Cone      <f64>>());
  show!(size_of::<Cube      <f64>>());
  show!(size_of::<Cuboid    <f64>>());

  show!(size_of::<Orthant>());
  show!(size_of::<Halfspace <f64>>());

  println!("...shape report sizes");
}

////////////////////////////////////////////////////////////////////////////////////////
//  impls                                                                             //
////////////////////////////////////////////////////////////////////////////////////////

//
//  impl Variant
//
impl <S> Shape <S> for Variant <S> where
  S : Real + num::real::Real + num::float::FloatCore + approx::RelativeEq
    + std::fmt::Debug + MaybeSerDes + 'static
{
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    match self {
      Variant::Bounded   (bounded)   => bounded.to_primitive (position),
      Variant::Unbounded (unbounded) => unbounded.to_primitive (position)
    }
  }
}
impl <S> Aabb <S> for Variant <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  /// Note that unbounded shapes will return a Cuboid with some infinite half extents
  fn aabb (&self) -> Aabb3 <S> {
    match self {
      Variant::Bounded   (bounded)   => bounded.aabb(),
      Variant::Unbounded (unbounded) => unbounded.aabb()
    }
  }
}
impl <S> Stereometric <S> for Variant <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  /// Note that unbounded shapes will return a Cuboid with some infinite half extents
  fn volume (&self) -> Positive <S> {
    match self {
      Variant::Bounded   (bounded)   => bounded.volume(),
      Variant::Unbounded (unbounded) => unbounded.volume()
    }
  }
}

//
//  impl Bounded
//
impl <S> Shape <S> for Bounded <S> where
  S : Real + num::real::Real + approx::RelativeEq + std::fmt::Debug + MaybeSerDes
    + 'static
{
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    match self {
      Bounded::Sphere   (sphere)   => sphere.to_primitive (position),
      Bounded::Capsule  (capsule)  => capsule.to_primitive (position),
      Bounded::Cylinder (cylinder) => cylinder.to_primitive (position),
      Bounded::Cuboid   (cuboid)   => cuboid.to_primitive (position),
      Bounded::Hull     (hull)     => hull.to_primitive (position),
      Bounded::Cone     (cone)     => cone.to_primitive (position),
      Bounded::Cube     (cube)     => cube.to_primitive (position)
    }
  }
}
impl <S> Aabb <S> for Bounded <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    match self {
      Bounded::Sphere   (sphere)   => sphere.aabb(),
      Bounded::Capsule  (capsule)  => capsule.aabb(),
      Bounded::Cylinder (cylinder) => cylinder.aabb(),
      Bounded::Cone     (cone)     => cone.aabb(),
      Bounded::Cube     (cube)     => cube.aabb(),
      Bounded::Cuboid   (cuboid)   => cuboid.aabb(),
      Bounded::Hull     (hull)     => hull.aabb()
    }
  }
}
impl <S : Real> Stereometric <S> for Bounded <S> {
  /// Volume of a bounded solid is always finite
  fn volume (&self) -> Positive <S> {
    match self {
      Bounded::Sphere   (sphere)   => sphere.volume(),
      Bounded::Capsule  (capsule)  => capsule.volume(),
      Bounded::Cylinder (cylinder) => cylinder.volume(),
      Bounded::Cone     (cone)     => cone.volume(),
      Bounded::Cube     (cube)     => cube.volume(),
      Bounded::Cuboid   (cuboid)   => cuboid.volume(),
      Bounded::Hull     (hull)     => hull.volume()
    }
  }
}

//
//  impl Unbounded
//
impl <S> Shape <S> for Unbounded <S> where
  S : Real + num::float::FloatCore + approx::RelativeEq + std::fmt::Debug + 'static
{
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    match self {
      Unbounded::Orthant   (orthant)   => orthant.to_primitive (position),
      Unbounded::Halfspace (halfspace) => halfspace.to_primitive (position)
    }
  }
}
impl <S> Aabb <S> for Unbounded <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  fn aabb (&self) -> Aabb3 <S> {
    match self {
      Unbounded::Orthant   (orthant)   => orthant.aabb(),
      Unbounded::Halfspace (halfspace) => halfspace.aabb()
    }
  }
}
impl <S> Stereometric <S> for Unbounded <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  /// Volume of an unbounded solid is always infinite
  fn volume (&self) -> Positive <S> {
    if cfg!(debug_assertions) {
      let volume = match self {
        Unbounded::Orthant   (orthant)   => orthant.volume(),
        Unbounded::Halfspace (halfspace) => halfspace.volume()
      };
      debug_assert_eq!(*volume, S::infinity());
    }
    Positive::infinity()
  }
}

//
//  impl Sphere
//
impl <S : OrderedRing> Sphere <S> {
  #[inline]
  /// Sphere with radius 1.0
  pub fn unit() -> Self where S : Field {
    use num::One;
    Sphere { radius: Positive::one() }
  }
  /// Create a new sphere with the absolute value of the given radius.
  ///
  /// Panics if radius is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Sphere;
  /// let s = Sphere::noisy (0.0);
  /// ```
  #[inline]
  pub fn noisy (radius : S) -> Self where S : std::fmt::Debug {
    assert_ne!(radius, S::zero());
    Sphere { radius: Positive::unchecked (radius.abs()) }
  }
  /// Create a new sphere with the absolute value of the given radius.
  ///
  /// Debug panic if radius is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Sphere;
  /// let s = Sphere::noisy (0.0);
  /// ```
  #[inline]
  pub fn unchecked (radius : S) -> Self where S : std::fmt::Debug {
    debug_assert_ne!(radius, S::zero());
    Sphere { radius: Positive::unchecked (radius.abs()) }
  }
  #[inline]
  pub const fn sphere3 (self, center : Point3 <S>) -> Sphere3 <S> {
    Sphere3 { center, radius: self.radius }
  }
}
impl <S> Shape <S> for Sphere <S> where
  S : Real + num::real::Real + std::fmt::Debug + MaybeSerDes + 'static
{
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    Box::new (self.sphere3 (position))
  }
}
impl <S> Aabb <S> for Sphere <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    Aabb3::with_minmax_unchecked ([-*self.radius; 3].into(), [ *self.radius; 3].into())
  }
}
impl <S : Real> Stereometric <S> for Sphere <S> {
  fn volume (&self) -> Positive <S> {
    let four      = Positive::unchecked (S::four());
    let frac_pi_3 = Positive::unchecked (S::frac_pi_3());
    four * frac_pi_3 * self.radius * self.radius * self.radius
  }
}
impl_numcast_fields!(Sphere, radius);

//
//  impl Capsule
//
impl <S : OrderedRing> Capsule <S> {
  /// Create a new capsule with the absolute values of the given radius and half-height.
  ///
  /// Panics if radius is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Capsule;
  /// let s = Capsule::noisy (0.0, 2.0);
  /// ```
  ///
  /// Panics if half-height is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Capsule;
  /// let s = Capsule::noisy (2.0, 0.0);
  /// ```
  #[inline]
  pub fn noisy (radius : S, half_height : S) -> Self where S : std::fmt::Debug {
    assert_ne!(radius, S::zero());
    assert_ne!(half_height, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Capsule { radius, half_height }
  }
  /// Create a new capsule with the absolute values of the given radius and half-height.
  ///
  /// Debug panic if radius is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Capsule;
  /// let s = Capsule::unchecked (0.0, 2.0);
  /// ```
  ///
  /// Debug panic if half-height is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Capsule;
  /// let s = Capsule::noisy (2.0, 0.0);
  /// ```
  #[inline]
  pub fn unchecked (radius : S, half_height : S) -> Self where
    S : std::fmt::Debug
  {
    debug_assert_ne!(radius, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Capsule { radius, half_height }
  }
  /// Height of the cylinder portion
  #[inline]
  pub fn height (self) -> NonNegative <S> {
    self.half_height * NonNegative::unchecked (S::two())
  }
  #[inline]
  pub const fn capsule3 (self, center : Point3 <S>) -> Capsule3 <S> {
    Capsule3 { center, radius: self.radius, half_height: self.half_height }
  }
}
impl <S> Shape <S> for Capsule <S> where
  S : Real + num::real::Real + std::fmt::Debug + MaybeSerDes + 'static
{
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    Box::new (self.capsule3 (position))
  }
}
impl <S> Aabb <S> for Capsule <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    let r  = *self.radius;
    let hh = *self.half_height;
    Aabb3::with_minmax_unchecked ([-r, -r, -r - hh].into(), [ r,  r,  r + hh].into())
  }
}
impl <S : Real> Stereometric <S> for Capsule <S> {
  fn volume (&self) -> Positive <S> {
    let r               = self.radius;
    let h               = self.height();
    let r2              = r * r;
    let r3              = r2 * r;
    let pi              = Positive::unchecked (S::pi());
    let four            = Positive::unchecked (S::four());
    let frac_pi_3       = Positive::unchecked (S::frac_pi_3());
    let cylinder_volume = pi * r2 * h;
    let sphere_volume   = four * frac_pi_3 * r3;
    sphere_volume + cylinder_volume
  }
}
impl_numcast_fields!(Capsule, radius, half_height);

//
//  impl Cylinder
//
impl <S : OrderedRing> Cylinder <S> {
  #[inline]
  /// Cylinder with radius 1.0 and half-height 1.0
  pub fn unit() -> Self where S : Field {
    use num::One;
    Cylinder {
      radius: Positive::one(), half_height: Positive::one()
    }
  }
  /// Create a new cylinder with the absolute values of the given radius and
  /// half-height.
  ///
  /// Panics if radius or half-height are zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cylinder;
  /// let s = Cylinder::noisy (0.0, 0.0);
  /// ```
  #[inline]
  pub fn noisy (radius : S, half_height : S) -> Self where S : std::fmt::Debug {
    assert_ne!(radius,      S::zero());
    assert_ne!(half_height, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Cylinder { radius, half_height }
  }
  /// Create a new cylinder with the absolute values of the given radius and
  /// half-height.
  ///
  /// Debug panic if radius or half-height are zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cylinder;
  /// let s = Cylinder::unchecked (0.0, 0.0);
  /// ```
  #[inline]
  pub fn unchecked (radius : S, half_height : S) -> Self where S : std::fmt::Debug {
    debug_assert_ne!(radius,      S::zero());
    debug_assert_ne!(half_height, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Cylinder { radius, half_height }
  }
  #[inline]
  pub const fn cylinder3 (self, center : Point3 <S>) -> Cylinder3 <S> {
    Cylinder3 { center, radius: self.radius, half_height: self.half_height }
  }
  #[inline]
  pub fn height (self) -> Positive <S> {
    self.half_height * Positive::unchecked (S::two())
  }
}
impl <S> Shape <S> for Cylinder <S> where
  S : Real + num::real::Real + std::fmt::Debug + MaybeSerDes + 'static
{
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    Box::new (self.cylinder3 (position))
  }
}
impl <S> Aabb <S> for Cylinder <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    let r  = *self.radius;
    let hh = *self.half_height;
    Aabb3::with_minmax_unchecked ([-r, -r, -hh].into(), [ r,  r,  hh].into())
  }
}
impl <S : Real> Stereometric <S> for Cylinder <S> {
  fn volume (&self) -> Positive <S> {
    let pi = Positive::unchecked (S::pi());
    pi * self.radius * self.radius * self.height()
  }
}
impl_numcast_fields!(Cylinder, radius, half_height);

//
//  impl Cone
//
impl <S : OrderedRing> Cone <S> {
  /// Create a new cone with the absolute values of the given radius and half-height.
  ///
  /// Panics if radius or half-height are zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cone;
  /// let s = Cone::noisy (0.0, 0.0);
  /// ```
  #[inline]
  pub fn noisy (radius : S, half_height : S) -> Self where S : std::fmt::Debug {
    assert_ne!(radius,      S::zero());
    assert_ne!(half_height, S::zero());
    let radius      = Positive::unchecked (radius.abs());
    let half_height = Positive::unchecked (half_height.abs());
    Cone { radius, half_height }
  }
}
impl <S> Shape <S> for Cone <S> where
  S : Real + num::real::Real + approx::RelativeEq + std::fmt::Debug + MaybeSerDes
    + 'static
{
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    Box::new (Cone3 {
      center:      position,
      half_height: self.half_height,
      radius:      self.radius
    })
  }
}
impl <S> Aabb <S> for Cone <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    let r  = *self.radius;
    let hh = *self.half_height;
    Aabb3::with_minmax_unchecked ([-r, -r, -hh].into(), [ r,  r,  hh].into())
  }
}
impl <S : Real> Stereometric <S> for Cone <S> {
  fn volume (&self) -> Positive <S> {
    let frac_pi_3 = Positive::unchecked (S::frac_pi_3());
    let two       = Positive::unchecked (S::two());
    frac_pi_3 * self.radius * self.radius * two * self.half_height
  }
}
impl_numcast_fields!(Cone, radius, half_height);

//
//  impl Cube
//
impl <S : OrderedRing> Cube <S> {
  /// Create a new cube with the absolute value of the given extent.
  ///
  /// Panics if extent is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cube;
  /// let s = Cube::noisy (0.0);
  /// ```
  #[inline]
  pub fn noisy (extent : S) -> Self where S : std::fmt::Debug {
    assert_ne!(extent, S::zero());
    Cube { extent: Positive::unchecked (extent.abs()) }
  }
  #[inline]
  pub fn dimensions (self) -> Vector3 <Positive <S>> {
    Vector3::broadcast (self.extent * Positive::unchecked (S::two()))
  }
  #[inline]
  pub fn width (self) -> Positive <S> {
    self.extent * Positive::unchecked (S::two())
  }
  #[inline]
  pub fn height (self) -> Positive <S> {
    self.extent * Positive::unchecked (S::two())
  }
  #[inline]
  pub fn depth (self) -> Positive <S> {
    self.extent * Positive::unchecked (S::two())
  }
}
impl <S> Shape <S> for Cube <S> where S : OrderedField + std::fmt::Debug + 'static {
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    let mut aabb = self.aabb();
    aabb.translate (position.0);
    Box::new (aabb)
  }
}
impl <S> Aabb <S> for Cube <S> where S : OrderedField + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    Aabb3::with_minmax_unchecked (
      [-*self.extent; 3].into(),
      [ *self.extent; 3].into())
  }
}
impl <S : OrderedField> Stereometric <S> for Cube <S> {
  fn volume (&self) -> Positive <S> {
    self.extent * self.extent * self.extent
  }
}
impl_numcast_fields!(Cube, extent);

//
//  impl Cuboid
//
impl <S : OrderedRing> Cuboid <S> {
  /// Create a new cuboid with the absolute values of the given extents.
  ///
  /// Panics if any extent is zero:
  ///
  /// ```should_panic
  /// # use math_utils::geometry::shape::Cuboid;
  /// let s = Cuboid::noisy ([0.0, 0.0, 0.0]);
  /// ```
  #[inline]
  pub fn noisy (extents : [S; 3]) -> Self where S : std::fmt::Debug {
    assert_ne!(extents[0], S::zero());
    assert_ne!(extents[1], S::zero());
    assert_ne!(extents[2], S::zero());
    let extents = vector3 (
      Positive::unchecked (extents[0].abs()),
      Positive::unchecked (extents[1].abs()),
      Positive::unchecked (extents[2].abs()));
    Cuboid { extents }
  }
  #[inline]
  pub fn dimensions (self) -> Vector3 <Positive <S>> {
    self.extents * Positive::unchecked (S::two())
  }
  /// X dimension
  #[inline]
  pub fn width (self) -> Positive <S> {
    self.extents.x * Positive::unchecked (S::two())
  }
  /// Y dimension
  #[inline]
  pub fn depth (self) -> Positive <S> {
    self.extents.y * Positive::unchecked (S::two())
  }
  /// Z dimension
  #[inline]
  pub fn height (self) -> Positive <S> {
    self.extents.z * Positive::unchecked (S::two())
  }
  #[inline]
  pub fn max (self) -> Point3 <S> {
    self.extents.map (|e| *e).into()
  }
  #[inline]
  pub fn min (self) -> Point3 <S> {
    self.extents.map (|e| -*e).into()
  }
  #[inline]
  pub fn aabb3 (self, center : Point3 <S>) -> Aabb3 <S> where S : std::fmt::Debug {
    Aabb3::with_minmax_unchecked (
      center + self.min().0,
      center + self.max().0)
  }
}
impl <S> TryFrom <Aabb3 <S>> for Cuboid <S> where S : OrderedField {
  type Error = Aabb3 <S>;
  fn try_from (aabb : Aabb3 <S>) -> Result <Self, Self::Error> {
    let extents = {
      let [x, y, z] = aabb.extents().into_array();
      vector3 (
        Positive::new (*x).ok_or (aabb)?,
        Positive::new (*y).ok_or (aabb)?,
        Positive::new (*z).ok_or (aabb)?)
    };
    Ok (Cuboid { extents })
  }
}
impl <S> Shape <S> for Cuboid <S> where S : OrderedField + std::fmt::Debug + 'static {
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    let mut aabb = self.aabb();
    aabb.translate (position.0);
    Box::new (aabb)
  }
}
impl <S> Aabb <S>  for Cuboid <S> where S : OrderedField + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    Aabb3::with_minmax_unchecked (self.min(), self.max())
  }
}
impl <S : OrderedField> Stereometric <S> for Cuboid <S> {
  fn volume (&self) -> Positive <S> {
    self.extents.x * self.extents.y * self.extents.z
  }
}
impl_numcast_fields!(Cuboid, extents);

//
//  impl Hull
//
impl <S> Shape <S> for Hull <S> where
  S : Real + num::real::Real + std::fmt::Debug + MaybeSerDes + 'static
{
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    let mut aabb = self.clone();
    aabb.translate (position.0);
    Box::new (aabb)
  }
}
impl <S> Aabb <S> for Hull <S> where S : Real + std::fmt::Debug {
  fn aabb (&self) -> Aabb3 <S> {
    // TODO: should we re-define a Hull to be a Hull3 + a cached Aabb3 to avoid
    // re-computing every time this is called ?
    Aabb3::containing (self.points()).unwrap()
  }
}
impl <S : Real> Stereometric <S> for Hull <S> {
  fn volume (&self) -> Positive <S> {
    unimplemented!("TODO: hull volume")
  }
}

//
//  impl Orthant
//
impl <S : Field> TryFrom <Halfspace <S>> for Orthant {
  type Error = Halfspace <S>;
  fn try_from (halfspace : Halfspace <S>) -> Result <Self, Self::Error> {
    SignedAxis3::try_from (*halfspace.normal)
      .map (|normal_axis| Orthant { normal_axis })
      .map_err (|_| halfspace)
  }
}
impl From <SignedAxis3> for Orthant {
  fn from (normal_axis : SignedAxis3) -> Self {
    Orthant { normal_axis }
  }
}
impl <S> Shape <S> for Orthant where
  S : Real + num::float::FloatCore + approx::RelativeEq + std::fmt::Debug + 'static
{
  #[expect(clippy::renamed_function_params)]
  fn to_primitive (&self, base : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    Box::new (Plane3 { base, normal: self.normal_axis.to_unit() })
  }
}
impl <S> Aabb <S> for Orthant where
  S : OrderedField + num::float::FloatCore + std::fmt::Debug
{
  /// Computes the containing Aabb3.
  ///
  /// ```
  /// # use math_utils::geometry::shape::Orthant;
  /// # use math_utils::geometry::Aabb3;
  /// # use math_utils::SignedAxis3;
  /// use math_utils::geometry::shape::Aabb;
  /// let s = Orthant::from (SignedAxis3::PosZ);
  /// assert_eq!(s.aabb(), Aabb3::with_minmax (
  ///   [f32::NEG_INFINITY, f32::NEG_INFINITY, f32::NEG_INFINITY].into(),
  ///   [f32::INFINITY,     f32::INFINITY,     0.0].into()
  /// ).unwrap());
  /// let s = Orthant::from (SignedAxis3::NegZ);
  /// assert_eq!(s.aabb(), Aabb3::with_minmax (
  ///   [f32::NEG_INFINITY, f32::NEG_INFINITY, 0.0].into(),
  ///   [f32::INFINITY,     f32::INFINITY,     f32::INFINITY].into()
  /// ).unwrap());
  /// ```
  fn aabb (&self) -> Aabb3 <S> {
    use num::float::FloatCore;
    let axis_vec = self.normal_axis.to_vec::<S>();
    let surface  = Vector3::broadcast (S::one()) - axis_vec.map (FloatCore::abs);

    let do_min = |i| if surface[i] == S::one() || axis_vec[i] == S::one() {
      S::neg_infinity()
    } else {
      debug_assert_eq!(axis_vec[i], -S::one());
      S::zero()
    };
    let do_max = |i| if surface[i] == S::one() {
      S::infinity()
    } else if axis_vec[i] == S::one() {
      S::zero()
    } else {
      debug_assert_eq!(axis_vec[i], -S::one());
      S::infinity()
    };
    let min = [ do_min (0), do_min (1), do_min (2) ].into();
    let max = [ do_max (0), do_max (1), do_max (2) ].into();

    Aabb3::with_minmax_unchecked (min, max)
  }
}
impl <S> Stereometric <S> for Orthant where S : OrderedRing + num::float::FloatCore {
  /// Volume of an unbounded solid is always infinite
  fn volume (&self) -> Positive <S> {
    Positive::infinity()
  }
}

//
//  impl Halfspace
//
impl <S : Copy> Halfspace <S> {
  #[inline]
  pub const fn plane3 (self, base : Point3 <S>) -> Plane3 <S> {
    Plane3 { base, normal: self.normal }
  }
}
impl <S : Real> From <Unit3 <S>> for Halfspace <S> {
  fn from (normal_vector : Unit3 <S>) -> Self {
    Halfspace { normal: normal_vector }
  }
}
impl <S> Shape <S> for Halfspace <S> where
  S : Real + num::float::FloatCore + approx::RelativeEq + std::fmt::Debug + 'static
{
  fn to_primitive (&self, position : Point3 <S>) -> Box <dyn Primitive <S, Point3 <S>>> {
    Box::new (self.plane3 (position))
  }
}
impl <S> Aabb <S> for Halfspace <S> where
  S : Real + num::float::FloatCore + std::fmt::Debug
{
  fn aabb (&self) -> Aabb3 <S> {
    if let Ok (orthant) = Orthant::try_from (*self) {
      orthant.aabb()
    } else {
      Aabb3::with_minmax_unchecked (
        [S::neg_infinity(); 3].into(),
        [S::infinity(); 3].into())
    }
  }
}
impl <S> Stereometric <S> for Halfspace <S> where S : Real + num::float::FloatCore {
  /// Volume of an unbounded solid is always infinite
  fn volume (&self) -> Positive <S> {
    Positive::infinity()
  }
}
impl_numcast_fields!(Halfspace, normal);
