DarkRoom/core/dr-pano/src/projection.rs

//! TRACES: FR-MRG-4
//! The surface the composite is drawn on.
//!
//! A panorama is a set of directions; a picture is a plane. The projection
//! is the map between them, and the three offered are the three every
//! stitcher offers because each is right for a different field of view:
//! perspective keeps straight lines straight and cannot reach 180°;
//! cylindrical keeps verticals vertical and stretches nothing horizontally,
//! for the wide single row; spherical for anything that also looks up.
//!
//! Every function here is the *inverse* map — output pixel to direction —
//! because that is what a gather needs (`lens.rs` in `dr-pipeline` says
//! why a warp is written that way), and it is the function the WGSL warp
//! will repeat verbatim. The forward map exists for bounds only.

use crate::linalg::Vec3;

#[derive(Debug, Clone, Copy, PartialEq, Eq)]
pub enum Projection {
    Perspective,
    Cylindrical,
    Spherical,
}

impl Projection {
    /// Which projection a field of view calls for.
    ///
    /// Perspective stretches the edges by `1 / cos` of the angle from the
    /// centre, which is 2× at 60° and unbounded at 90°; the switch is where
    /// that stretch starts to look like a mistake. Spherical is for a set
    /// that spans enough vertically that a cylinder would stretch the top
    /// and bottom the same way.
    pub fn suggest(horizontal_fov: f64, vertical_fov: f64) -> Projection {
        if horizontal_fov < 70f64.to_radians() && vertical_fov < 70f64.to_radians() {
            Projection::Perspective
        } else if vertical_fov < 100f64.to_radians() {
            Projection::Cylindrical
        } else {
            Projection::Spherical
        }
    }

    /// The direction an output point looks along. `scale` is the output's
    /// focal length in pixels: the radius of the cylinder or sphere, or the
    /// plane's distance. Coordinates are centred on the projection's origin
    /// (the direction `+z`).
    pub fn to_direction(self, scale: f64, u: f64, v: f64) -> Vec3 {
        match self {
            Projection::Perspective => Vec3::new(u, v, scale).normalised(),
            Projection::Cylindrical => {
                let theta = u / scale;
                Vec3::new(theta.sin(), v / scale, theta.cos()).normalised()
            }
            Projection::Spherical => {
                let theta = u / scale;
                let phi = v / scale;
                Vec3::new(theta.sin() * phi.cos(), phi.sin(), theta.cos() * phi.cos())
            }
        }
    }

    /// Where a direction lands on the output, or `None` where the
    /// projection cannot show it (behind a perspective plane, at a
    /// cylinder's poles).
    pub fn from_direction(self, scale: f64, d: Vec3) -> Option<(f64, f64)> {
        let (x, y, z) = (d.x(), d.y(), d.z());
        match self {
            Projection::Perspective => (z > 1e-9).then(|| (scale * x / z, scale * y / z)),
            Projection::Cylindrical => {
                let r = (x * x + z * z).sqrt();
                (r > 1e-9).then(|| (scale * x.atan2(z), scale * y / r))
            }
            Projection::Spherical => {
                let r = (x * x + z * z).sqrt();
                Some((scale * x.atan2(z), scale * y.atan2(r)))
            }
        }
    }
}

/// The output rectangle a set of frames covers, in centred output pixels.
#[derive(Debug, Clone, Copy, PartialEq)]
pub struct Bounds {
    pub min_u: f64,
    pub min_v: f64,
    pub max_u: f64,
    pub max_v: f64,
}

impl Bounds {
    pub fn width(&self) -> f64 {
        self.max_u - self.min_u
    }
    pub fn height(&self) -> f64 {
        self.max_v - self.min_v
    }
}

/// Bounds of the frames' footprints under `projection`, by walking each
/// frame's border.
///
/// `frame_size` is the frames' width and height in the same pixels the
/// cameras' focal length is in. The border is sampled rather than only its
/// corners because under a cylinder the widest point of a rolled frame is
/// not a corner.
pub fn bounds(
    projection: Projection,
    scale: f64,
    cameras: &crate::bundle::Cameras,
    frame_size: (f64, f64),
) -> Option<Bounds> {
    let (w, h) = frame_size;
    let mut b: Option<Bounds> = None;
    let steps = 64;
    for k in 0..cameras.rotations.len() {
        for s in 0..steps {
            let t = s as f64 / steps as f64;
            for p in [
                (-w / 2.0 + w * t, -h / 2.0),
                (-w / 2.0 + w * t, h / 2.0),
                (-w / 2.0, -h / 2.0 + h * t),
                (w / 2.0, -h / 2.0 + h * t),
            ] {
                let d = cameras.bearing(k, p);
                let Some((u, v)) = projection.from_direction(scale, d) else {
                    continue;
                };
                b = Some(match b {
                    None => Bounds {
                        min_u: u,
                        min_v: v,
                        max_u: u,
                        max_v: v,
                    },
                    Some(b) => Bounds {
                        min_u: b.min_u.min(u),
                        min_v: b.min_v.min(v),
                        max_u: b.max_u.max(u),
                        max_v: b.max_v.max(v),
                    },
                });
            }
        }
    }
    b
}

#[cfg(test)]
mod tests {
    use super::*;

    #[test]
    fn to_and_from_direction_are_inverses() {
        for proj in [
            Projection::Perspective,
            Projection::Cylindrical,
            Projection::Spherical,
        ] {
            for (u, v) in [(0.0, 0.0), (300.0, -200.0), (-900.0, 450.0)] {
                let d = proj.to_direction(1000.0, u, v);
                let (bu, bv) = proj.from_direction(1000.0, d).expect("in front");
                assert!(
                    (bu - u).abs() < 1e-9 && (bv - v).abs() < 1e-9,
                    "{proj:?} {u} {v}"
                );
            }
        }
    }

    #[test]
    fn the_origin_looks_down_z_in_every_projection() {
        for proj in [
            Projection::Perspective,
            Projection::Cylindrical,
            Projection::Spherical,
        ] {
            let d = proj.to_direction(500.0, 0.0, 0.0);
            assert!((d.z() - 1.0).abs() < 1e-12);
        }
    }

    #[test]
    fn a_cylinder_maps_ninety_degrees_to_a_quarter_turn_of_pixels() {
        let d = Vec3::new(1.0, 0.0, 0.0);
        let (u, v) = Projection::Cylindrical.from_direction(100.0, d).unwrap();
        assert!((u - 100.0 * std::f64::consts::FRAC_PI_2).abs() < 1e-9);
        assert_eq!(v, 0.0);
        assert!(Projection::Perspective.from_direction(100.0, d).is_none());
    }

    #[test]
    fn suggestion_widens_with_the_field() {
        assert_eq!(Projection::suggest(0.5, 0.5), Projection::Perspective);
        assert_eq!(Projection::suggest(2.5, 0.8), Projection::Cylindrical);
        assert_eq!(Projection::suggest(3.0, 2.5), Projection::Spherical);
    }
}