Shared exponent float format

xyz8e5.rs

use glam::{vec3, Vec3};

pub const XYZ8E5_EXPONENT_BITS: i32 = 5;
pub const XYZ8E5_MANTISSA_BITS: i32 = 8;
pub const XYZ8E5_MANTISSA_BITSU: u32 = 8;
pub const XYZ8E5_EXP_BIAS: i32 = 15;
pub const XYZ8E5_MAX_VALID_BIASED_EXP: i32 = 31;

/*
pub const MAX_XYZ8E5_EXP: i32 = 16;
pub const XYZ8E5_MANTISSA_VALUES: i32 = 256;
pub const MAX_XYZ8E5_MANTISSA: i32 = 255;
pub const MAX_XYZ8E5_MANTISSAU: u32 = 255;
pub const MAX_XYZ8E5: f32 = 65280.0;
pub const EPSILON_XYZ8E5: f32 = 1.1920929e-7;
*/

pub const MAX_XYZ8E5_EXP: i32 = XYZ8E5_MAX_VALID_BIASED_EXP - XYZ8E5_EXP_BIAS;
pub const XYZ8E5_MANTISSA_VALUES: i32 = 1 << XYZ8E5_MANTISSA_BITS;
pub const MAX_XYZ8E5_MANTISSA: i32 = XYZ8E5_MANTISSA_VALUES - 1;
pub const MAX_XYZ8E5_MANTISSAU: u32 = (XYZ8E5_MANTISSA_VALUES - 1) as u32;
pub const MAX_XYZ8E5: f32 =
    (MAX_XYZ8E5_MANTISSA as f32) / XYZ8E5_MANTISSA_VALUES as f32 * (1 << MAX_XYZ8E5_EXP) as f32;
pub const EPSILON_XYZ8E5: f32 =
    (1.0 / XYZ8E5_MANTISSA_VALUES as f32) / (1 << XYZ8E5_EXP_BIAS) as f32;

// Similar to https://www.khronos.org/registry/OpenGL/extensions/EXT/EXT_texture_shared_exponent.txt
#[inline]
pub fn vec3_to_xyz8e5(xyz: Vec3) -> u32 {
    let rsign = xyz.x.is_sign_positive() as u32;
    let gsign = xyz.y.is_sign_positive() as u32;
    let bsign = xyz.z.is_sign_positive() as u32;

    let rc = xyz.x.abs().min(MAX_XYZ8E5);
    let gc = xyz.y.abs().min(MAX_XYZ8E5);
    let bc = xyz.z.abs().min(MAX_XYZ8E5);

    let maxxyz = rc.max(gc).max(bc);
    let mut exp_shared =
        (-XYZ8E5_EXP_BIAS - 1).max(maxxyz.log2().floor() as i32) + 1 + XYZ8E5_EXP_BIAS;

    assert!(exp_shared <= XYZ8E5_MAX_VALID_BIASED_EXP);
    assert!(exp_shared >= 0);

    let mut denom = ((exp_shared - XYZ8E5_EXP_BIAS - XYZ8E5_MANTISSA_BITS) as f32).exp2();

    let maxm = (maxxyz / denom + 0.5).floor() as i32;
    if maxm == MAX_XYZ8E5_MANTISSA + 1 {
        denom *= 2.0;
        exp_shared += 1;
        assert!(exp_shared <= XYZ8E5_MAX_VALID_BIASED_EXP);
    } else {
        assert!(maxm <= MAX_XYZ8E5_MANTISSA);
    }

    let rm = (rc / denom + 0.5).floor() as i32;
    let gm = (gc / denom + 0.5).floor() as i32;
    let bm = (bc / denom + 0.5).floor() as i32;

    assert!(rm <= MAX_XYZ8E5_MANTISSA);
    assert!(gm <= MAX_XYZ8E5_MANTISSA);
    assert!(bm <= MAX_XYZ8E5_MANTISSA);
    assert!(rm >= 0);
    assert!(gm >= 0);
    assert!(bm >= 0);

    assert_eq!(rm as u32, rm as u32 & MAX_XYZ8E5_MANTISSAU);
    assert_eq!(gm as u32, gm as u32 & MAX_XYZ8E5_MANTISSAU);
    assert_eq!(bm as u32, bm as u32 & MAX_XYZ8E5_MANTISSAU);

    let rm = rm as u32 | rsign << 8;
    let gm = gm as u32 | gsign << 8;
    let bm = bm as u32 | bsign << 8;
    let exp_shared = exp_shared as u32;

    #[allow(clippy::identity_op)]
    let ret = (exp_shared << 27) | (bm << 18) | (gm << 9) | (rm << 0);

    ret
}

#[inline]
fn bitfield_extract(value: u32, offset: u32, bits: u32) -> u32 {
    let mask = (1u32 << bits) - 1u32;
    (value >> offset) & mask
}

#[inline]
pub fn xyz8e5_to_vec3(v: u32) -> Vec3 {
    let exponent = bitfield_extract(v, 27, XYZ8E5_EXPONENT_BITS as u32) as i32
        - XYZ8E5_EXP_BIAS
        - XYZ8E5_MANTISSA_BITS;
    let scale = (exponent as f32).exp2();

    let rsign = bitfield_extract(v, 8, 1) as f32 * 2.0 - 1.0;
    let gsign = bitfield_extract(v, 17, 1) as f32 * 2.0 - 1.0;
    let bsign = bitfield_extract(v, 26, 1) as f32 * 2.0 - 1.0;

    vec3(
        rsign * bitfield_extract(v, 0, XYZ8E5_MANTISSA_BITS as u32) as f32 * scale,
        gsign * bitfield_extract(v, 9, XYZ8E5_MANTISSA_BITS as u32) as f32 * scale,
        bsign * bitfield_extract(v, 18, XYZ8E5_MANTISSA_BITS as u32) as f32 * scale,
    )
}

#[cfg(test)]
pub mod tests {

    use super::*;
    use rand::Rng;

    #[test]
    fn get_data_for_plot() {
        dbg!(
            MAX_XYZ8E5_EXP,
            XYZ8E5_MANTISSA_VALUES,
            MAX_XYZ8E5_MANTISSA,
            MAX_XYZ8E5_MANTISSAU,
            MAX_XYZ8E5,
            EPSILON_XYZ8E5,
        );
        println!("RANGE   \tMAX      \tAVG");
        for i in 1..65 {
            let mut n = i as f32 * 0.25;
            n = n.exp2() - 1.0;
            let (max, avg) =
                test_conversion(n, 1000000, false, |v| xyz8e5_to_vec3(vec3_to_xyz8e5(v)));
            println!("{:.8}\t{:.8}\t{:.8}", n, max, avg);
        }
    }

    #[test]
    fn get_data_for_plot_f16() {
        for i in 1..57 {
            let mut n = i as f32 * 0.25;
            n = n.exp2() - 1.0;
            let (max, avg) = test_conversion(n, 1000000, false, |v| {
                vec3(
                    half::f16::from_f32(v.x).into(),
                    half::f16::from_f32(v.y).into(),
                    half::f16::from_f32(v.z).into(),
                )
            });
            println!("{} {} {}", n, max, avg);
        }
    }

    #[test]
    fn get_data_for_plot_8unorm_plus_sign() {
        for i in 1..5 {
            let mut n = i as f32 * 0.25;
            n = n.exp2() - 1.0;
            let (max, avg) = test_conversion(n, 1000000, false, |v| {
                ((v.abs() * 255.0 + 0.5).floor() / 255.0)
                    * vec3(v.x.signum(), v.y.signum(), v.z.signum())
            });
            println!("{} {} {}", n, max, avg);
        }
    }

    pub fn test_conversion<F>(distance: f32, iterations: usize, print: bool, proc: F) -> (f32, f32)
    where
        F: Fn(Vec3) -> Vec3,
    {
        let mut rng = rand::thread_rng();

        let mut max_relative_error = 0.0;
        let mut max_error_orig = Vec3::ZERO;
        let mut max_error_xyz_decoded = Vec3::ZERO;

        let mut max_abs_error = 0.0;
        let mut max_abs_error_orig = Vec3::ZERO;
        let mut max_abs_error_decoded = Vec3::ZERO;

        let mut max_nor_dist = 0.0;
        let mut max_nor_dist_orig = Vec3::ZERO;
        let mut max_nor_dist_decoded = Vec3::ZERO;

        let mut max_dist = 0.0;
        let mut max_dist_orig = Vec3::ZERO;
        let mut max_dist_decoded = Vec3::ZERO;

        let mut avg_dist = 0.0;
        for _ in 0..iterations {
            let orig = Vec3::from([
                rng.gen_range(-distance..distance),
                rng.gen_range(-distance..distance),
                rng.gen_range(-distance..distance),
            ]);

            let decoded = proc(orig);

            for i in 0..3 {
                let abs_diff = (orig[i] - decoded[i]).abs();
                avg_dist += abs_diff;
                let relative_error = if orig[i] != 0.0 {
                    abs_diff / orig[i]
                } else {
                    abs_diff
                };

                if relative_error > max_relative_error {
                    max_relative_error = relative_error;
                    max_error_orig = orig;
                    max_error_xyz_decoded = decoded;
                }

                if abs_diff > max_abs_error {
                    max_abs_error = abs_diff;
                    max_abs_error_orig = orig;
                    max_abs_error_decoded = decoded;
                }

                let a = Vec3::from(orig);
                let b = Vec3::from(decoded);
                let nor_dist = a.normalize_or_zero().distance(b.normalize_or_zero());
                if b.normalize_or_zero().length() != 0.0 {
                    if nor_dist > max_nor_dist {
                        max_nor_dist = nor_dist;
                        max_nor_dist_orig = orig;
                        max_nor_dist_decoded = decoded;
                    }
                }

                let dist = a.distance(b);
                if dist > max_dist {
                    max_dist = dist;
                    max_dist_orig = orig;
                    max_dist_decoded = decoded;
                }
            }
        }
        avg_dist /= iterations as f32;
        if print {
            println!("\nMaximum Relative Error:");
            println!("Error:\t {}", max_relative_error);
            println!("Original:\t {:?}", max_error_orig);
            println!("Decoded:\t {:?}", max_error_xyz_decoded);

            println!("\nMaximum Absolute Error:");
            println!("Error:\t {}", max_abs_error);
            println!("Original:\t {:?}", max_abs_error_orig);
            println!("Decoded:\t {:?}", max_abs_error_decoded);

            println!("\nMaximum Normalized Distance:");
            println!("Distance: {}", max_nor_dist);
            println!("Original: {:?}", max_nor_dist_orig);
            println!("Decoded: {:?}", max_nor_dist_decoded);

            println!("\nAt Maximum Distance:");
            println!("Original:\t {:?}", max_dist_orig);
            println!("Decoded:\t {:?}", max_dist_decoded);

            println!("Max Dist:\t {}", max_dist);
            println!("Avg Dist:\t {:?}", avg_dist);
        }
        (max_dist, avg_dist)
    }
}