cadmpeg-asm 0.6.5

// SPDX-License-Identifier: Apache-2.0
//! Decode analytic surfaces and 3D curves, select edge pcurves, reverse
//! curve orientation, and recognize procedural carriers as analytic geometry.

use super::records::TolerantCoedgeExtension;
use crate::ids::{brep_id, IdFormat};
use crate::nurbs::proc_surface::{
    DecodedProceduralSurfaceDefinition, EmbeddedRollingBall, EmbeddedScaledCompoundLoftShape,
};
use crate::nurbs::reader::LEN_TO_MM;
use crate::nurbs::support_chart::{cone_axis_reversed, cone_is_circular_cylinder};
use crate::sab::{Record, Token};
use cadmpeg_ir::features::{FinitePoint3, FiniteVector3};
use cadmpeg_ir::geometry::analytic::{
    CircleCurve, ConeSurface, CylinderSurface, DegenerateCurve, EllipseCurve, LineCurve,
    PlaneSurface, SphereSurface, TorusSurface,
};
use cadmpeg_ir::geometry::{
    Curve, CurveGeometry, SolvedCurveGeometry, SolvedSurfaceGeometry, SurfaceGeometry,
};
use cadmpeg_ir::ids::{CurveId, EdgeId};
use cadmpeg_ir::math::{Point3, Vector3};
use cadmpeg_ir::scalar::{Angle, NonNegativeLength, NonZeroLength, PositiveLength, PositiveReal};
use cadmpeg_ir::topology::{BodyKind, EdgeCarrier, Sense};
use cadmpeg_ir::units::{FiniteVector, OrthonormalFrame3, UnitVector3};
use std::collections::{HashMap, HashSet};

use super::{id, AsmBrep};
const EPS_GEOMETRY_PCURVE_RANGES_ON_DOMAIN_E9: f64 = 1.0e-9;
const EPS_GEOMETRY_ANALYTIC_PROCEDURAL_SURFACE_E10: f64 = 1.0e-10;
const EPS_GEOMETRY_ANALYTIC_ROLLING_BALL_SURFACE_E10: f64 = 1.0e-10;
const EPS_GEOMETRY_LINEAR_NURBS_SPINE_E10: f64 = 1.0e-10;
const EPS_GEOMETRY_RATIONAL_FOUR_ARC_CIRCLE_E12: f64 = 1.0e-12;
const EPS_GEOMETRY_RATIONAL_FOUR_ARC_CIRCLE_E10: f64 = 1.0e-10;
const EPS_GEOMETRY_REDUCE_HOMOGENEOUS_BEZIER_TO_QUADRATIC_E10: f64 = 1.0e-10;
const EPS_GEOMETRY_CLAMP_EDGE_RANGES_TO_CARRIER_DOMAINS_E9: f64 = 1.0e-9;

/// Ordered typed values pulled from a carrier record's payload.
pub(in crate::brep) struct Carrier {
    pub(super) positions: Vec<[f64; 3]>,
    pub(super) vectors: Vec<[f64; 3]>,
    doubles: Vec<f64>,
}

pub(super) fn collect_carrier(rec: &Record) -> Carrier {
    let mut c = Carrier {
        positions: Vec::new(),
        vectors: Vec::new(),
        doubles: Vec::new(),
    };
    for (index, t) in rec.tokens.iter().enumerate() {
        // Cone interval endpoints precede sine and cosine. A present bound
        // is TRUE DOUBLE; it is not a cone carrier scalar.
        if rec.head() == "cone"
            && matches!(t, Token::Double(_))
            && index > 0
            && matches!(rec.tokens[index - 1], Token::True)
        {
            continue;
        }
        match t {
            Token::Position(p) => {
                c.positions.push(*p);
            }
            Token::Vector3(v) => {
                c.vectors.push(*v);
            }
            Token::Double(d) => {
                c.doubles.push(*d);
            }
            _ => {}
        }
    }
    c
}

pub(super) fn scale_point(p: [f64; 3]) -> Point3 {
    Point3::new(p[0] * LEN_TO_MM, p[1] * LEN_TO_MM, p[2] * LEN_TO_MM)
}

pub(super) fn norm3(v: [f64; 3]) -> f64 {
    Vector3::from(v).norm()
}

/// The unit direction of a stored vector, absent when the vector is not
/// finite or its length is within `f64::EPSILON` of zero.
fn direction(v: [f64; 3]) -> Option<UnitVector3> {
    UnitVector3::normalized(Vector3::from(v))
}

/// The frame of two stored directions, absent when either direction is
/// degenerate or the two are not perpendicular.
fn frame(axis: [f64; 3], reference: [f64; 3]) -> Option<OrthonormalFrame3> {
    OrthonormalFrame3::from_units(direction(axis)?, direction(reference)?)
}

/// Whether a record name heads an analytic surface carrier.
pub(super) fn is_analytic_surface(head: &str) -> bool {
    matches!(head, "plane" | "cone" | "sphere" | "torus")
}

/// Whether a record name heads an analytic curve carrier.
pub(super) fn is_analytic_curve(head: &str) -> bool {
    matches!(head, "straight" | "ellipse" | "degenerate_curve")
}

/// Decode an analytic surface carrier. Signed sphere radii and torus minor
/// radii remain in the IR because they are part of the ASM carrier semantics.
/// A negative torus major radius decodes as its magnitude, the same tube
/// circle; its pcurves change chart through `SupportChart::torus`.
pub fn decode_surface(rec: &Record) -> Option<(SolvedSurfaceGeometry, bool)> {
    decode_surface_carrier(rec, &collect_carrier(rec))
}

fn decode_surface_carrier(rec: &Record, c: &Carrier) -> Option<(SolvedSurfaceGeometry, bool)> {
    let origin = *c.positions.first()?;
    match rec.head() {
        "plane" => Some((
            SolvedSurfaceGeometry::Plane(PlaneSurface::new(
                FinitePoint3::new(scale_point(origin))?,
                frame(*c.vectors.first()?, *c.vectors.get(1)?)?,
            )),
            false,
        )),
        "cone" => {
            let ratio = *c.doubles.first().unwrap_or(&1.0);
            let axis = direction(*c.vectors.first()?)?;
            let major = *c.vectors.get(1)?;
            // Doubles are (ratio, sine, cosine, u_scale). `ratio` is the
            // minor/major radius ratio. `sine` selects cylinder vs cone. The
            // base radius is the major-axis vector's
            // magnitude; the trailing `u_scale` double is the u-parameter
            // scale, which usually coincides with the radius but diverges on
            // offset-derived surfaces. The signed slope `sine / cosine` is the
            // radius change per unit axis distance, and a negative `cosine`
            // points the surface normal toward the axis.
            let sine = *c.doubles.get(1).unwrap_or(&0.0);
            let cosine = *c.doubles.get(2).unwrap_or(&1.0);
            let radius = norm3(major) * LEN_TO_MM;
            (radius > f64::EPSILON).then_some(())?;
            let ref_direction = direction(major)?;
            if cone_is_circular_cylinder(sine, ratio) {
                Some((
                    SolvedSurfaceGeometry::Cylinder(CylinderSurface::new(
                        FinitePoint3::new(scale_point(origin))?,
                        OrthonormalFrame3::from_units(axis, ref_direction)?,
                        PositiveLength::new(radius)?,
                    )),
                    cosine < 0.0,
                ))
            } else {
                // The IR cone's radius grows along `+axis`; a negative native
                // slope shrinks it, so the axis flips to compensate. The
                // outward normal is invariant under the flip; the inward
                // normal of a negative `cosine` folds into the face sense.
                let axis = if cone_axis_reversed(sine, cosine, ratio) {
                    axis.reversed()
                } else {
                    axis
                };
                Some((
                    SolvedSurfaceGeometry::Cone(ConeSurface::new(
                        FinitePoint3::new(scale_point(origin))?,
                        OrthonormalFrame3::from_units(axis, ref_direction)?,
                        NonNegativeLength::new(radius)?,
                        PositiveReal::new(ratio)?,
                        Angle::new(sine.abs().atan2(cosine.abs()))?,
                    )),
                    cosine < 0.0,
                ))
            }
        }
        "sphere" => {
            let signed = *c.doubles.first()?;
            Some((
                SolvedSurfaceGeometry::Sphere(SphereSurface::new(
                    FinitePoint3::new(scale_point(origin))?,
                    frame(*c.vectors.get(1)?, *c.vectors.first()?)?,
                    NonZeroLength::new(signed * LEN_TO_MM)?,
                )),
                false,
            ))
        }
        "torus" => {
            let axis = *c.vectors.first()?;
            let ref_direction = *c.vectors.get(1)?;
            let major = *c.doubles.first()?;
            let minor = *c.doubles.get(1)?;
            Some((
                SolvedSurfaceGeometry::Torus(TorusSurface::new(
                    FinitePoint3::new(scale_point(origin))?,
                    frame(axis, ref_direction)?,
                    PositiveLength::new(major.abs() * LEN_TO_MM)?,
                    NonZeroLength::new(minor * LEN_TO_MM)?,
                )),
                false,
            ))
        }
        _ => None,
    }
}

/// The vertex record's point reference. The modern layout stores the
/// endpoint-index integer at chunk 4 and the point at chunk 5; the
/// save-format 700 layout stores no endpoint index and the point at chunk 4.
/// A modern record always carries the integer, so the legacy branch is
/// unreachable for it.
pub(super) fn vertex_point_ref(record: &Record) -> Option<i64> {
    match record.chunk(4) {
        Some(Token::Long(_)) => record.ref_at(5),
        _ => record.ref_at(4),
    }
}

/// The coedge record's pcurve reference: chunk 10 after the reserved integer
/// in the modern layout, chunk 9 in the save-format 700 layout that stores
/// no reserved integer.
pub(crate) fn coedge_pcurve_ref(record: &Record) -> Option<i64> {
    match record.chunk(9) {
        Some(Token::Long(_)) => record.ref_at(10),
        _ => record.ref_at(9),
    }
}

pub(super) fn is_vertex_record(record: &Record) -> bool {
    matches!(record.head(), "vertex" | "tvertex")
}

pub(super) fn is_edge_record(record: &Record) -> bool {
    matches!(record.head(), "edge" | "tedge")
}

pub(super) fn is_coedge_record(record: &Record) -> bool {
    matches!(record.head(), "coedge" | "tcoedge")
}

pub(super) fn tolerant_coedge_extension(record: &Record) -> Option<TolerantCoedgeExtension> {
    let target = match record.chunk(13)? {
        Token::Ref(target) => (*target >= 0).then_some(*target),
        _ => return None,
    };
    match record.chunk(14)? {
        Token::Long(0) if matches!(record.chunk(15), Some(Token::Long(0))) => {
            Some(TolerantCoedgeExtension::Empty { target })
        }
        Token::Long(1) => {
            let curve_reversed = match record.chunk(15)? {
                Token::True => true,
                Token::False => false,
                _ => return None,
            };
            if !matches!(record.chunk(16), Some(Token::SubtypeOpen)) {
                return None;
            }
            // Raw index of chunk 16: the payload identifiers inside the
            // embedded scope are not chunks, and the serialized token count
            // below is defined over the value tokens alone.
            let open = record
                .tokens
                .iter()
                .enumerate()
                .filter(|(_, token)| !token.is_payload_ident())
                .nth(16)
                .map(|(index, _)| index)?;
            let mut depth = 0usize;
            let mut close = None;
            for (index, token) in record.tokens.iter().enumerate().skip(open) {
                match token {
                    Token::SubtypeOpen => depth += 1,
                    Token::SubtypeClose => {
                        depth = depth.checked_sub(1)?;
                        if depth == 0 {
                            close = Some(index);
                            break;
                        }
                    }
                    _ => {}
                }
            }
            let close = close?;
            // Suffix in chunk space: a record can end with payload
            // identifiers (`null_curve` placeholders for absent curve slots),
            // which are not fields of the extension.
            let mut suffix = record
                .tokens
                .get(close + 1..)?
                .iter()
                .filter(|token| !token.is_payload_ident());
            let parameter_range = match (
                suffix.next(),
                suffix.next(),
                suffix.next(),
                suffix.next(),
                suffix.next(),
                suffix.next(),
            ) {
                (
                    Some(Token::False),
                    Some(Token::False),
                    Some(Token::Long(0)),
                    None,
                    None,
                    None,
                ) => None,
                (
                    Some(Token::True),
                    Some(Token::Double(start)),
                    Some(Token::True),
                    Some(Token::Double(end)),
                    Some(Token::Long(0)),
                    None,
                ) => Some(cadmpeg_ir::units::FiniteVector::new([*start, *end])?),
                _ => return None,
            };
            let payload_token_count = record
                .tokens
                .get(open + 1..close)?
                .iter()
                .filter(|token| !token.is_payload_ident())
                .count();
            Some(TolerantCoedgeExtension::EmbeddedCurve {
                target,
                curve_reversed,
                payload_token_count: u32::try_from(payload_token_count).ok()?,
                parameter_range,
            })
        }
        _ => None,
    }
}

/// Whether a record head belongs to the topology or geometry vocabulary.
///
/// Carrier heads remain known even when no active topology references that
/// particular record. Reachability determines transfer; it does not turn an
/// unreferenced carrier into an application/refinement record.
pub(super) fn is_known_record_head(head: &str) -> bool {
    matches!(
        head,
        "body"
            | "region"
            | "lump"
            | "shell"
            | "subshell"
            | "wire"
            | "face"
            | "loop"
            | "point"
            | "asmheader"
    ) || matches!(
        head,
        "coedge" | "tcoedge" | "edge" | "tedge" | "vertex" | "tvertex"
    ) || is_analytic_surface(head)
        || is_analytic_curve(head)
        || matches!(head, "spline" | "intcurve" | "pcurve")
}

pub(super) fn is_asm_stream_delimiter(name: &str) -> bool {
    matches!(name, "Begin-of-ASM-History-Data" | "End-of-ASM-data")
}

pub(super) fn edge_pcurve_parameter_ranges(edge: &Record) -> Option<[[f64; 2]; 2]> {
    let (Some(Token::Double(start)), Some(Token::Double(end))) = (edge.chunk(4), edge.chunk(6))
    else {
        return None;
    };
    let direct = [*start, *end];
    let negated = [-start, -end];
    Some(if matches!(edge.chunk(9), Some(Token::True)) {
        [negated, direct]
    } else {
        [direct, negated]
    })
}

/// Candidate edge-use intervals whose endpoints lie on this pcurve carrier.
/// Edge sense orders the two signs, but it cannot move a NURBS use outside the
/// carrier's active knot domain. That active domain is the final fallback.
pub(super) fn pcurve_ranges_on_domain(
    candidate: &cadmpeg_ir::geometry::pcurve::PcurveNurbs,
    edge: Option<&Record>,
) -> Option<Vec<[f64; 2]>> {
    let first = *candidate
        .knots()
        .get(usize::try_from(candidate.degree()).ok()?)?;
    let last = *candidate.knots().get(candidate.pole_rows().count())?;
    (first < last).then_some(())?;
    let mut ranges = Vec::new();
    for range in edge
        .and_then(edge_pcurve_parameter_ranges)
        .into_iter()
        .flatten()
    {
        if let Some(range) = (|| {
            range
                .iter()
                .all(|value| {
                    cadmpeg_ir::math::parameter_in_domain(
                        *value,
                        [first, last],
                        EPS_GEOMETRY_PCURVE_RANGES_ON_DOMAIN_E9,
                    )
                })
                .then_some(())?;
            let range = range.map(|value| value.clamp(first, last));
            (range[0] != range[1]).then_some(range)
        })() {
            ranges.push(range);
        }
    }
    if !ranges.contains(&[first, last]) {
        ranges.push([first, last]);
    }
    Some(ranges)
}

/// Decode an analytic curve carrier.
pub fn decode_curve(rec: &Record) -> Option<CurveGeometry> {
    decode_curve_carrier(rec, &collect_carrier(rec))
}

fn decode_curve_carrier(rec: &Record, carrier: &Carrier) -> Option<CurveGeometry> {
    let base = *carrier.positions.first()?;
    match rec.head() {
        "straight" => Some(CurveGeometry::Solved(SolvedCurveGeometry::Line(
            LineCurve::new(
                FinitePoint3::new(scale_point(base))?,
                direction(*carrier.vectors.first()?)?,
            ),
        ))),
        "ellipse" => {
            let axis = *carrier.vectors.first()?;
            let reference = *carrier.vectors.get(1)?;
            let ratio = *carrier.doubles.first()?;
            let major_radius = norm3(reference) * LEN_TO_MM;
            let center = FinitePoint3::new(scale_point(base))?;
            let conic_frame = frame(axis, reference)?;
            if (ratio.abs() - 1.0).abs() <= f64::EPSILON {
                Some(CurveGeometry::Solved(SolvedCurveGeometry::Circle(
                    CircleCurve::new(center, conic_frame, PositiveLength::new(major_radius)?),
                )))
            } else {
                Some(CurveGeometry::Solved(SolvedCurveGeometry::Ellipse(
                    EllipseCurve::try_from_parts(
                        center,
                        conic_frame,
                        PositiveLength::new(major_radius)?,
                        PositiveLength::new(major_radius * ratio.abs())?,
                    )
                    .ok()?,
                )))
            }
        }
        "degenerate_curve" => Some(CurveGeometry::Solved(SolvedCurveGeometry::Degenerate(
            cadmpeg_ir::geometry::analytic::DegenerateCurve::try_new(scale_point(base)).ok()?,
        ))),
        _ => None,
    }
}

pub(super) fn sense_at(rec: &Record, i: usize) -> Sense {
    match rec.chunk(i) {
        Some(Token::True) => Sense::Reversed,
        _ => Sense::Forward,
    }
}

/// The record-level sense bit of an `intcurve` or `spline` carrier: the boolean
/// token immediately before the record's subtype scope ([spec §6.6](https://github.com/cadmpeg/cadmpeg/blob/main/docs/formats/asm.md#66-intcurve-and-spline-subtypes)). `true`
/// marks geometry as the reverse of its cached definition. A reversed intcurve
/// negates the cache parameterization (`C(t) = cache(-t)`), and a reversed
/// spline surface flips the cache normal.
pub(super) fn record_reversed(rec: &Record) -> bool {
    // Adjacency in chunk space: a freestanding payload identifier (e.g. the
    // embedded curve's type name) can sit between value tokens without
    // separating the sense bit from the scope it precedes.
    let mut previous: Option<&Token> = None;
    rec.chunks()
        .find_map(|token| {
            let result = matches!(token, Token::SubtypeOpen)
                .then(|| match previous {
                    Some(Token::True) => Some(true),
                    Some(Token::False) => Some(false),
                    _ => None,
                })
                .flatten();
            previous = Some(token);
            result
        })
        .or_else(|| {
            // A plain `intcurve` companion has no subtype scope after its
            // base header; its sense remains the fourth payload value.
            (rec.head() == "intcurve").then(|| matches!(rec.chunk(3), Some(Token::True)))
        })
        .unwrap_or(false)
}

/// Reverse a curve carrier to its opposite orientation, `C'(t) = C(-t)`.
/// Lines negate their direction, conics negate their plane normal (flipping
/// the angular sweep while keeping the zero-angle direction), and B-splines
/// reverse poles and knots. Carriers without an orientation pass through.
pub(super) fn reverse_curve_geometry(geometry: &mut CurveGeometry) {
    match geometry {
        CurveGeometry::Solved(SolvedCurveGeometry::Line(line_curve)) => {
            line_curve.reverse_parameterization();
        }
        CurveGeometry::Solved(SolvedCurveGeometry::Circle(circle_curve)) => {
            circle_curve.reverse_parameterization();
        }
        CurveGeometry::Solved(SolvedCurveGeometry::Ellipse(ellipse_curve)) => {
            ellipse_curve.reverse_parameterization();
        }
        CurveGeometry::Solved(SolvedCurveGeometry::Nurbs(curve)) => {
            curve.reverse_parameterization();
        }
        _ => {}
    }
}

pub(super) fn reverse_procedural_curve_definition(
    definition: &mut cadmpeg_ir::geometry::ProceduralCurveDefinition,
) -> Result<(), &'static str> {
    if let cadmpeg_ir::geometry::ProceduralCurveDefinition::Helix(helix) = definition {
        helix.try_reverse_parameterization()?;
    }
    Ok(())
}

pub(super) fn double_at(rec: &Record, i: usize) -> Option<f64> {
    match rec.chunk(i) {
        Some(Token::Double(d)) => Some(*d),
        _ => None,
    }
}

pub(super) fn pcurve_parameter_range(rec: &Record) -> Option<[f64; 2]> {
    // The final two value tokens; a record may end with payload identifiers
    // (e.g. `null_curve` placeholders), which are not fields.
    let mut values = rec.chunks().rev();
    match (values.next(), values.next()) {
        (Some(Token::Double(end)), Some(Token::Double(start))) => Some([*start, *end]),
        _ => None,
    }
}

pub(super) fn pcurve_inline_tail_flags(rec: &Record) -> Option<[bool; 4]> {
    if !matches!(rec.chunk(3), Some(Token::Long(0))) {
        return None;
    }
    // End-relative in chunk space: the four booleans precede the final two
    // value tokens, and trailing payload identifiers are not fields.
    let start = rec.chunks().count().checked_sub(6)?;
    let mut flags = [false; 4];
    for (slot, token) in flags.iter_mut().zip(rec.chunks().skip(start)) {
        *slot = match token {
            Token::True => true,
            Token::False => false,
            _ => return None,
        };
    }
    Some(flags)
}

pub(super) fn procedural_surface_definition_is_exact_carrier(
    definition: &DecodedProceduralSurfaceDefinition,
) -> bool {
    match definition {
        DecodedProceduralSurfaceDefinition::Extrusion { .. }
        | DecodedProceduralSurfaceDefinition::Sum { .. }
        | DecodedProceduralSurfaceDefinition::Helix(_)
        | DecodedProceduralSurfaceDefinition::Ruled { .. }
        | DecodedProceduralSurfaceDefinition::VertexBlend(_)
        | DecodedProceduralSurfaceDefinition::SubSurface { .. } => true,
        DecodedProceduralSurfaceDefinition::Sweep(construction) => matches!(
            &construction.layout,
            crate::nurbs::proc_surface::EmbeddedSweepSurfaceLayout::Revision { form, .. }
                if matches!(form.cache, cadmpeg_ir::geometry::RevisionCacheForm::Parameterization(_))
        ),
        DecodedProceduralSurfaceDefinition::Law(construction) => !matches!(
            construction.tail,
            cadmpeg_ir::geometry::LawSurfaceTail::Full { .. }
        ),
        DecodedProceduralSurfaceDefinition::ScaledCompoundLoft(construction) => matches!(
            construction.shape,
            EmbeddedScaledCompoundLoftShape::None { .. }
        ),
        // Tail form `2` stores no solved cache, so every surface block a blend
        // record holds is a support of its construction.
        DecodedProceduralSurfaceDefinition::Blend {
            native: Some(construction),
            ..
        } => matches!(
            construction.cache,
            cadmpeg_ir::geometry::RevisionCacheForm::Parameterization(_)
        ),
        DecodedProceduralSurfaceDefinition::VariableBlend(construction) => {
            matches!(
                construction.cache,
                cadmpeg_ir::geometry::VariableBlendCache::Parameterization { .. }
                    | cadmpeg_ir::geometry::VariableBlendCache::Stale {}
            )
        }
        _ => false,
    }
}

pub(super) fn analytic_procedural_surface(
    definition: &DecodedProceduralSurfaceDefinition,
) -> Option<SurfaceGeometry> {
    match definition {
        DecodedProceduralSurfaceDefinition::Extrusion {
            directrix,
            direction,
            ..
        } => {
            let (center, normal, ref_direction, radius) = rational_four_arc_circle(directrix)?;
            let axis = UnitVector3::normalized(*direction)?;
            if 1.0 - axis.as_raw().dot(normal).abs() > EPS_GEOMETRY_ANALYTIC_PROCEDURAL_SURFACE_E10
            {
                return None;
            }
            Some(SurfaceGeometry::Solved(SolvedSurfaceGeometry::Cylinder(
                CylinderSurface::new(
                    FinitePoint3::new(center)?,
                    OrthonormalFrame3::from_units(axis, UnitVector3::new(ref_direction)?)?,
                    PositiveLength::new(radius)?,
                ),
            )))
        }
        DecodedProceduralSurfaceDefinition::Blend {
            supports,
            spine: Some(spine),
            radius_offsets: [signed_radius, end_offset],
            cross_section: cadmpeg_ir::geometry::BlendCrossSection::Circular,
            native,
        } if signed_radius == end_offset => {
            analytic_rolling_ball_surface(supports, native.as_deref(), spine, *signed_radius)
        }
        _ => None,
    }
}

fn analytic_rolling_ball_surface(
    supports: &[Option<SurfaceGeometry>; 2],
    native: Option<&EmbeddedRollingBall>,
    spine: &cadmpeg_ir::geometry::nurbs::NurbsCurve,
    signed_radius: f64,
) -> Option<SurfaceGeometry> {
    let radius = signed_radius.abs();
    if !radius.is_finite() || radius <= f64::EPSILON {
        return None;
    }
    let support = |index: usize| {
        supports[index].as_ref().or_else(|| {
            native.and_then(|native| {
                native.sides[index]
                    .surface
                    .as_ref()
                    .map(|support| &support.surface)
            })
        })
    };
    let first = support(0)?;
    let second = support(1)?;

    if let (
        SurfaceGeometry::Solved(SolvedSurfaceGeometry::Plane(plane_surface)),
        SurfaceGeometry::Solved(SolvedSurfaceGeometry::Plane(plane_surface_2)),
    ) = (first, second)
    {
        let first_origin = plane_surface.origin();
        let first_normal = *plane_surface.frame().axis().as_raw();
        let second_origin = plane_surface_2.origin();
        let second_normal = *plane_surface_2.frame().axis().as_raw();
        let (origin, axis) = linear_nurbs_spine(spine)?;
        let tolerance = EPS_GEOMETRY_ANALYTIC_ROLLING_BALL_SURFACE_E10 * radius;
        let support_intersection = first_normal.cross(second_normal);
        let support_intersection_norm = support_intersection.norm();
        if support_intersection_norm <= EPS_GEOMETRY_ANALYTIC_ROLLING_BALL_SURFACE_E10
            || 1.0
                - axis
                    .as_raw()
                    .dot(support_intersection.scale(1.0 / support_intersection_norm))
                    .abs()
                > EPS_GEOMETRY_ANALYTIC_ROLLING_BALL_SURFACE_E10
        {
            return None;
        }
        for (plane_origin, plane_normal) in [
            (*first_origin, first_normal),
            (*second_origin, second_normal),
        ] {
            if axis.as_raw().dot(plane_normal).abs()
                > EPS_GEOMETRY_ANALYTIC_ROLLING_BALL_SURFACE_E10
                || (point_vector(plane_origin, origin).dot(plane_normal).abs() - radius).abs()
                    > tolerance
            {
                return None;
            }
        }
        return Some(SurfaceGeometry::Solved(SolvedSurfaceGeometry::Cylinder(
            CylinderSurface::new(
                FinitePoint3::new(origin)?,
                OrthonormalFrame3::from_units(
                    axis,
                    UnitVector3::new(cadmpeg_ir::geometry::derive_reference_direction(
                        *axis.as_raw(),
                    ))?,
                )?,
                PositiveLength::new(radius)?,
            ),
        )));
    }

    let (plane_origin, plane_normal, cylinder_origin, cylinder_axis, cylinder_radius) =
        match (first, second) {
            (
                SurfaceGeometry::Solved(SolvedSurfaceGeometry::Plane(plane_surface)),
                SurfaceGeometry::Solved(SolvedSurfaceGeometry::Cylinder(cylinder_surface)),
            ) => {
                let plane_origin = plane_surface.origin().get();
                let plane_normal = *plane_surface.frame().axis().as_raw();
                let cylinder_origin = cylinder_surface.origin().get();
                let cylinder_axis = *cylinder_surface.frame().axis().as_raw();
                let cylinder_radius = cylinder_surface.radius().get();
                (
                    plane_origin,
                    plane_normal,
                    cylinder_origin,
                    cylinder_axis,
                    cylinder_radius,
                )
            }
            (
                SurfaceGeometry::Solved(SolvedSurfaceGeometry::Cylinder(cylinder_surface_2)),
                SurfaceGeometry::Solved(SolvedSurfaceGeometry::Plane(plane_surface_2)),
            ) => {
                let cylinder_origin = cylinder_surface_2.origin().get();
                let cylinder_axis = *cylinder_surface_2.frame().axis().as_raw();
                let cylinder_radius = cylinder_surface_2.radius().get();
                let plane_origin = plane_surface_2.origin().get();
                let plane_normal = *plane_surface_2.frame().axis().as_raw();
                (
                    plane_origin,
                    plane_normal,
                    cylinder_origin,
                    cylinder_axis,
                    cylinder_radius,
                )
            }
            _ => {
                return None;
            }
        };
    let (center, axis, ref_direction, major_radius) = rational_four_arc_circle(spine)?;
    let scale = major_radius.max(radius).max(cylinder_radius);
    let tolerance = EPS_GEOMETRY_ANALYTIC_ROLLING_BALL_SURFACE_E10 * scale;
    let center_offset = point_vector(cylinder_origin, center);
    let axial_offset = center_offset.dot(cylinder_axis);
    let radial_offset = center_offset - cylinder_axis.scale(axial_offset);
    if 1.0 - axis.dot(plane_normal).abs() > EPS_GEOMETRY_ANALYTIC_ROLLING_BALL_SURFACE_E10
        || 1.0 - axis.dot(cylinder_axis).abs() > EPS_GEOMETRY_ANALYTIC_ROLLING_BALL_SURFACE_E10
        || (point_vector(plane_origin, center).dot(plane_normal).abs() - radius).abs() > tolerance
        || radial_offset.norm() > tolerance
        || ((major_radius - cylinder_radius).abs() - radius).abs() > tolerance
    {
        return None;
    }
    Some(SurfaceGeometry::Solved(SolvedSurfaceGeometry::Torus(
        TorusSurface::try_new(center, axis, ref_direction, major_radius, signed_radius).ok()?,
    )))
}

fn linear_nurbs_spine(
    curve: &cadmpeg_ir::geometry::nurbs::NurbsCurve,
) -> Option<(Point3, UnitVector3)> {
    if curve.degree() == 0 || curve.periodic() {
        return None;
    }
    if let Some(weights) = curve.weights() {
        let first_sign = weights.first()?.get().signum();
        if weights
            .iter()
            .any(|weight| weight.get().signum() != first_sign)
        {
            return None;
        }
    }
    let origin = curve.control_points()[0];
    let (_, farthest) = curve
        .control_points()
        .iter()
        .copied()
        .map(|point| (point_vector(origin.get(), point.get()).norm(), point))
        .max_by(|left, right| left.0.total_cmp(&right.0))?;
    let extent = point_vector(origin.get(), farthest.get()).norm();
    if !extent.is_finite() || extent <= f64::EPSILON {
        return None;
    }
    let axis = UnitVector3::normalized(point_vector(origin.get(), farthest.get()))?;
    // This admits an analytic replacement, not a model-length approximation.
    if curve.control_points().iter().any(|point| {
        let relative = point_vector(origin.get(), point.get());
        let relative = Vector3::new(
            relative.x / extent,
            relative.y / extent,
            relative.z / extent,
        );
        axis.as_raw().cross(relative).norm() > EPS_GEOMETRY_LINEAR_NURBS_SPINE_E10
    }) {
        return None;
    }
    Some((origin.get(), axis))
}

pub(super) fn rational_four_arc_circle(
    curve: &cadmpeg_ir::geometry::nurbs::NurbsCurve,
) -> Option<(Point3, Vector3, Vector3, f64)> {
    let cadmpeg_ir::geometry::nurbs::NurbsPoles3::Rational { points } = curve.pole_rows() else {
        return None;
    };
    let degree = usize::try_from(curve.degree()).ok()?;
    if degree < 2 || curve.periodic() || points.len() != 4 * degree + 1 {
        return None;
    }
    let knot_tolerance = EPS_GEOMETRY_RATIONAL_FOUR_ARC_CIRCLE_E12
        * (curve.knots()[curve.knots().len() - 1] * 0.5 - curve.knots()[0] * 0.5).abs()
        * 2.0;
    let spans = [
        curve.knots()[0],
        curve.knots()[degree + 1],
        curve.knots()[2 * degree + 1],
        curve.knots()[3 * degree + 1],
        curve.knots()[4 * degree + 1],
    ];
    if spans
        .windows(2)
        .any(|pair| pair[1] - pair[0] <= knot_tolerance)
        || (0..5).any(|span| {
            let range = if span == 0 {
                0..degree + 1
            } else if span == 4 {
                4 * degree + 1..curve.knots().len()
            } else {
                span * degree + 1..(span + 1) * degree + 1
            };
            curve.knots()[range]
                .iter()
                .any(|value| (*value - spans[span]).abs() > knot_tolerance)
        })
    {
        return None;
    }
    let weight_scale = points
        .iter()
        .map(|pole| pole.weight.get().abs())
        .fold(0.0, f64::max);
    let mut homogeneous = Vec::new();
    for pole in points {
        let point = pole.point;
        let weight = pole.weight.get() / weight_scale;
        let homogeneous_pole = [point.x * weight, point.y * weight, point.z * weight, weight];
        if !(weight.is_finite()
            && weight != 0.0
            && homogeneous_pole.iter().all(|value| value.is_finite()))
        {
            return None;
        }
        homogeneous.push(homogeneous_pole);
    }
    let mut quadratics = [None; 4];
    for (span, quadratic) in quadratics.iter_mut().enumerate() {
        *quadratic = Some(reduce_homogeneous_bezier_to_quadratic(
            &homogeneous[span * degree..=span * degree + degree],
        )?);
    }
    let quadratics = [
        quadratics[0]?,
        quadratics[1]?,
        quadratics[2]?,
        quadratics[3]?,
    ];
    let base_weight = quadratics[0][0][3];
    let weight_scale = base_weight.abs();
    let weight_tolerance = EPS_GEOMETRY_RATIONAL_FOUR_ARC_CIRCLE_E10 * weight_scale;
    if !base_weight.is_finite()
        || base_weight == 0.0
        || quadratics.iter().any(|span| {
            (span[0][3] - base_weight).abs() > weight_tolerance
                || (span[2][3] - base_weight).abs() > weight_tolerance
                || (span[1][3] - base_weight * std::f64::consts::FRAC_1_SQRT_2).abs()
                    > weight_tolerance
        })
    {
        return None;
    }
    let quadratic_points = quadratics.map(|span| {
        span.map(|point| {
            Point3::new(
                point[0] / point[3],
                point[1] / point[3],
                point[2] / point[3],
            )
        })
    });
    let point_distance = |left: Point3, right: Point3| point_vector(left, right).norm();
    let scale = quadratic_points
        .iter()
        .flat_map(|span| span.windows(2))
        .map(|pair| point_distance(pair[0], pair[1]))
        .fold(0.0_f64, f64::max);
    let tolerance = EPS_GEOMETRY_RATIONAL_FOUR_ARC_CIRCLE_E10 * scale;
    if point_distance(quadratic_points[0][0], quadratic_points[3][2]) > tolerance
        || quadratic_points
            .windows(2)
            .any(|pair| point_distance(pair[0][2], pair[1][0]) > tolerance)
    {
        return None;
    }
    let first_center = point_sum_difference(
        quadratic_points[0][0],
        quadratic_points[0][2],
        quadratic_points[0][1],
    );
    for span in &quadratic_points {
        let [start, control, end] = *span;
        let center = point_sum_difference(start, end, control);
        if point_distance(center, first_center) > tolerance {
            return None;
        }
    }
    let first_radial = point_vector(first_center, quadratic_points[0][0]);
    let radius = first_radial.norm();
    if !radius.is_finite() || radius <= tolerance {
        return None;
    }
    let mut normal = None;
    for span in &quadratic_points {
        let radial = point_vector(first_center, span[0]);
        let next = point_vector(first_center, span[2]);
        let radial_unit = FiniteVector3::new(radial)?.unit_nonzero()?;
        let next_unit = FiniteVector3::new(next)?.unit_nonzero()?;
        if (radial.norm() - radius).abs() > tolerance
            || radial_unit.dot(next_unit).abs() > EPS_GEOMETRY_RATIONAL_FOUR_ARC_CIRCLE_E10
        {
            return None;
        }
        let span_normal = FiniteVector3::new(radial_unit.cross(next_unit))?.unit_nonzero()?;
        if normal.is_some_and(|normal: Vector3| {
            normal.dot(span_normal) < 1.0 - EPS_GEOMETRY_RATIONAL_FOUR_ARC_CIRCLE_E10
        }) {
            return None;
        }
        normal.get_or_insert(span_normal);
    }
    Some((
        first_center,
        normal?,
        FiniteVector3::new(first_radial)?.unit_nonzero()?,
        radius,
    ))
}

fn reduce_homogeneous_bezier_to_quadratic(input: &[[f64; 4]]) -> Option<[[f64; 4]; 3]> {
    let mut control = Vec::new();
    control.extend_from_slice(input);
    while control.len() > 3 {
        let degree = control.len() - 1;
        let mut reduced = Vec::new();
        reduced.push(control[0]);
        for index in 1..degree {
            let alpha = cadmpeg_core::convert::f64_from_index(index)?
                / cadmpeg_core::convert::f64_from_index(degree)?;
            let denominator = 1.0 - alpha;
            reduced.push(std::array::from_fn(|coordinate| {
                (control[index][coordinate] - alpha * reduced[index - 1][coordinate]) / denominator
            }));
        }
        if reduced.iter().flatten().any(|value| !value.is_finite()) {
            return None;
        }
        if (0..4).any(|coordinate| {
            let scale = control
                .iter()
                .map(|point| point[coordinate].abs())
                .fold(0.0, f64::max);
            (reduced[degree - 1][coordinate] - control[degree][coordinate]).abs()
                > EPS_GEOMETRY_REDUCE_HOMOGENEOUS_BEZIER_TO_QUADRATIC_E10 * scale
        }) {
            return None;
        }
        control = reduced;
    }
    control.try_into().ok()
}

pub(super) fn point_vector(origin: Point3, point: Point3) -> Vector3 {
    Vector3::new(point.x - origin.x, point.y - origin.y, point.z - origin.z)
}

fn point_sum_difference(first: Point3, second: Point3, subtract: Point3) -> Point3 {
    Point3::new(
        first.x + second.x - subtract.x,
        first.y + second.y - subtract.y,
        first.z + second.z - subtract.z,
    )
}

/// Snap edge parameter ranges that overshoot their B-spline carrier's knot
/// domain by floating-point noise back onto the domain boundary. Native edge
/// ranges and cache knot vectors are stored independently and can disagree in
/// their last few bits; a genuine domain violation is left for validation.
pub(super) fn clamp_edge_ranges_to_carrier_domains(
    out: &mut AsmBrep,
) -> Result<(), cadmpeg_core::CodecError> {
    let domains: HashMap<&str, [f64; 2]> = out
        .curves
        .iter()
        .filter_map(|curve| match &curve.geometry {
            CurveGeometry::Solved(SolvedCurveGeometry::Nurbs(nurbs)) => {
                let first = nurbs.knots().get(usize::try_from(nurbs.degree()).ok()?)?;
                let last = nurbs.knots().get(nurbs.pole_count())?;
                Some((curve.id.as_str(), [*first, *last]))
            }
            _ => None,
        })
        .collect::<std::collections::HashMap<_, _>>();
    for edge in &mut out.edges {
        let Some([mut start, mut end]) = edge.param_range().map(FiniteVector::get) else {
            continue;
        };
        let Some([first, last]) = edge.curve().and_then(|curve| domains.get(curve.as_str())) else {
            continue;
        };
        let tolerance = (last * 0.5 - first * 0.5).abs()
            * (2.0 * EPS_GEOMETRY_CLAMP_EDGE_RANGES_TO_CARRIER_DOMAINS_E9);
        if start < *first && *first - start <= tolerance {
            start = *first;
        }
        if end > *last && end - *last <= tolerance {
            end = *last;
        }
        edge.carrier =
            cadmpeg_ir::topology::EdgeCarrier::new(edge.curve().cloned(), Some([start, end]))
                .map_err(cadmpeg_core::CodecError::malformed)?;
    }
    Ok(())
}

/// Solid or sheet kind of a faced body without wire topology, from the edge
/// and sense of each coedge use.
///
/// The faces bound a volume when every edge that does not collapse to a point
/// has as many forward as reversed uses. A non-manifold edge closes the faces
/// in the same way as a two-use edge. A collapsed edge bounds no area, so its
/// uses do not open the faces. A body with no counted edge is a sheet.
pub fn faced_body_kind<'a>(
    uses: impl IntoIterator<Item = (&'a EdgeId, Sense)>,
    collapsed: impl Fn(&EdgeId) -> bool,
) -> BodyKind {
    let mut balance = HashMap::<&EdgeId, (usize, usize)>::new();
    for (edge, sense) in uses {
        if collapsed(edge) {
            continue;
        }
        let (forward, reversed) = balance.entry(edge).or_default();
        match sense {
            Sense::Forward => *forward += 1,
            Sense::Reversed => *reversed += 1,
        }
    }
    if !balance.is_empty()
        && balance
            .values()
            .all(|(forward, reversed)| forward == reversed)
    {
        BodyKind::Solid
    } else {
        BodyKind::Sheet
    }
}

/// Kept edges that collapse to one point: both vertex references name the
/// same vertex and the curve reference is null or names a `degenerate_curve`.
///
/// An undecoded curve also leaves the neutral edge without a carrier, so the
/// test reads the native curve reference.
pub(super) fn collapsed_edges(
    records: &[Record],
    by_index: &HashMap<i64, &Record>,
    kept_edges: &HashSet<i64>,
    format: IdFormat,
) -> Result<HashSet<EdgeId>, cadmpeg_core::CodecError> {
    let mut collapsed = HashSet::new();
    for record in records {
        if !is_edge_record(record) {
            continue;
        }
        let index = i64::try_from(record.index).map_err(|_| {
            cadmpeg_core::CodecError::malformed("ASM record index exceeds the i64 range")
        })?;
        if !kept_edges.contains(&index) {
            continue;
        }
        let (Some(start), Some(end)) = (record.ref_at(3), record.ref_at(5)) else {
            continue;
        };
        let point_carrier = match record.chunk(8) {
            Some(Token::Ref(curve)) if *curve < 0 => true,
            Some(Token::Ref(curve)) => by_index
                .get(curve)
                .is_some_and(|curve| curve.head() == "degenerate_curve"),
            _ => false,
        };
        if start == end && point_carrier {
            collapsed.insert(EdgeId::from(id(format, index)));
        }
    }
    Ok(collapsed)
}

/// Give each native point edge a degenerate carrier at its vertex. A plain
/// `edge` record whose curve reference is null and whose two vertex
/// references are identical collapses to that vertex
/// ([spec §6.3](https://github.com/cadmpeg/cadmpeg/blob/main/docs/formats/asm.md#63-analytic-curve-byte-layouts)).
pub(super) fn attach_point_edge_carriers(out: &mut AsmBrep, records: &[Record], format: IdFormat) {
    let point_edges = records
        .iter()
        .filter(|record| {
            record.head() == "edge"
                && matches!(record.tokens.get(8), Some(Token::Ref(curve)) if *curve < 0)
                && record.ref_at(3).is_some()
                && record.ref_at(3) == record.ref_at(5)
        })
        .filter_map(|record| {
            let index = i64::try_from(record.index).ok()?;
            Some((EdgeId::from(id(format, index)), index))
        })
        .collect::<HashMap<_, _>>();
    if point_edges.is_empty() {
        return;
    }
    let positions = out
        .points
        .iter()
        .map(|point| (&point.id, point.position()))
        .collect::<HashMap<_, _>>();
    let vertex_positions = out
        .vertices
        .iter()
        .filter_map(|vertex| Some((&vertex.id, *positions.get(&vertex.point)?)))
        .collect::<HashMap<_, _>>();
    let mut curves = Vec::new();
    for edge in &mut out.edges {
        let Some(index) = point_edges.get(&edge.id) else {
            continue;
        };
        if edge.curve().is_some() || edge.start != edge.end {
            continue;
        }
        let Some(position) = vertex_positions.get(&edge.start) else {
            continue;
        };
        let curve = brep_id!(format, CurveId, "point-edge-curve", *index);
        edge.carrier = EdgeCarrier::Curve(curve.clone());
        curves.push(Curve {
            id: curve,
            geometry: CurveGeometry::Solved(SolvedCurveGeometry::Degenerate(DegenerateCurve::new(
                *position,
            ))),
            parameter_range: None,
            source_object: None,
        });
    }
    out.curves.extend(curves);
}

pub(super) fn classify_body_kinds(out: &mut AsmBrep, collapsed: &HashSet<EdgeId>) {
    let mut shell_bodies = HashMap::new();
    for region in &out.regions {
        for shell in &region.shells {
            shell_bodies.insert(shell, &region.body);
        }
    }
    let mut body_has_faces = HashSet::new();
    let mut body_has_wires = HashSet::new();
    let mut face_bodies = HashMap::new();
    for shell in &out.shells {
        let Some(body) = shell_bodies.get(&shell.id).copied() else {
            continue;
        };
        if !shell.wire_edges().is_empty() || !shell.free_vertices().is_empty() {
            body_has_wires.insert(body);
        }
        if !shell.faces().is_empty() {
            body_has_faces.insert(body);
        }
        for face in shell.faces() {
            face_bodies.insert(face, body);
        }
    }
    let mut loop_bodies = HashMap::new();
    for face in &out.faces {
        let Some(body) = face_bodies.get(&face.id).copied() else {
            continue;
        };
        for loop_id in &face.loops {
            loop_bodies.insert(loop_id, body);
        }
    }
    let mut coedge_bodies = HashMap::new();
    for loop_ in &out.loops {
        let Some(body) = loop_bodies.get(&loop_.id).copied() else {
            continue;
        };
        for coedge in loop_.coedges() {
            coedge_bodies.insert(coedge, body);
        }
    }
    let mut body_uses = HashMap::<_, Vec<(&EdgeId, Sense)>>::new();
    for coedge in &out.coedges {
        if let Some(body) = coedge_bodies.get(&coedge.id).copied() {
            body_uses
                .entry(body)
                .or_default()
                .push((&coedge.edge, coedge.sense));
        }
    }
    for body in &mut out.bodies {
        if !body_has_faces.contains(&body.id) {
            body.kind = BodyKind::Wire;
            continue;
        }
        if body_has_wires.contains(&body.id) {
            body.kind = BodyKind::General;
            continue;
        }
        body.kind = faced_body_kind(
            body_uses.get(&body.id).into_iter().flatten().copied(),
            |edge| collapsed.contains(edge),
        );
    }
}

#[cfg(test)]
mod analytic_surface_tests {
    use super::{collect_carrier, decode_surface};
    use crate::sab::{Record, Token};
    use std::sync::Arc;

    fn surface_record(head: &str, tokens: Vec<Token>) -> Record {
        Record {
            index: 1,
            name: format!("{head}-surface"),

            tokens: Arc::from(tokens),
            offset: 0,
            len: 0,
        }
    }

    #[test]
    fn finite_cone_intervals_do_not_shift_the_carrier_scalars() {
        let record = surface_record(
            "cone",
            vec![
                Token::Position([0.0, 0.0, 0.0]),
                Token::Vector3([0.0, 0.0, 1.0]),
                Token::Vector3([7.0, 0.0, 0.0]),
                Token::Double(1.0),
                Token::True,
                Token::Double(2.0),
                Token::True,
                Token::Double(3.0),
                Token::Double(0.0),
                Token::Double(1.0),
                Token::Double(7.0),
                Token::False,
                Token::True,
                Token::Double(4.0),
                Token::True,
                Token::Double(5.0),
                Token::False,
                Token::False,
            ],
        );
        {
            let carrier = collect_carrier(&record);
            assert_eq!(carrier.doubles, [1.0, 0.0, 1.0, 7.0]);
            let (surface, _) = decode_surface(&record).unwrap();
            assert!(matches!(
                surface,
                cadmpeg_ir::geometry::SolvedSurfaceGeometry::Cylinder(_)
            ));
        };
    }

    #[test]
    fn analytic_surfaces_require_complete_serialized_frames() {
        let origin = Token::Position([0.0, 0.0, 0.0]);
        let axis = Token::Vector3([0.0, 0.0, 1.0]);

        let plane = surface_record("plane", vec![origin.clone(), axis.clone()]);
        let cone_without_major = surface_record(
            "cone",
            vec![
                origin.clone(),
                axis.clone(),
                Token::Double(1.0),
                Token::Double(0.0),
                Token::Double(1.0),
                Token::Double(7.0),
            ],
        );
        let cone_with_zero_major = surface_record(
            "cone",
            vec![
                origin.clone(),
                axis.clone(),
                Token::Vector3([0.0, 0.0, 0.0]),
                Token::Double(1.0),
            ],
        );
        let sphere = surface_record(
            "sphere",
            vec![origin.clone(), Token::Double(2.0), axis.clone()],
        );
        let torus = surface_record(
            "torus",
            vec![origin, axis, Token::Double(3.0), Token::Double(1.0)],
        );

        for record in [
            plane,
            cone_without_major,
            cone_with_zero_major,
            sphere,
            torus,
        ] {
            assert!(decode_surface(&record).is_none());
        }
    }
}

#[cfg(test)]
mod sense_tests {
    use super::record_reversed;
    use crate::sab::{Record, Token};

    #[test]
    fn intcurve_sense_falls_back_only_when_the_scope_has_no_boolean() {
        for (header, before_scope, expected) in [
            (Token::True, Token::Long(7), true),
            (Token::False, Token::Long(7), false),
            (Token::True, Token::False, false),
            (Token::False, Token::True, true),
        ] {
            let record = Record {
                index: 0,
                name: "intcurve".into(),
                tokens: vec![
                    Token::Ref(-1),
                    Token::Long(-1),
                    Token::Ref(-1),
                    header,
                    before_scope,
                    Token::SubtypeOpen,
                    Token::SubtypeClose,
                ]
                .into(),
                offset: 0,
                len: 0,
            };
            assert_eq!(record_reversed(&record), expected);
        }
    }
}

#[cfg(test)]
mod tests {
    mod numerical_ranges;
    use super::Point3;
    const SMALL_CURVED_SPINE_EXTENT: f64 = 1.0e-10;

    #[test]
    fn numerical_seventh_analytic_spine_rejects_relative_curvature_at_small_scale() {
        for scale in [1.0, SMALL_CURVED_SPINE_EXTENT, 1.0e100] {
            let spine = |height| {
                cadmpeg_ir::geometry::nurbs::NurbsCurve::from_lanes(
                    2,
                    vec![0.0, 0.0, 0.0, 1.0, 1.0, 1.0],
                    vec![
                        Point3::new(0.0, 0.0, 0.0),
                        Point3::new(0.5 * scale, height * scale, 0.0),
                        Point3::new(scale, 0.0, 0.0),
                    ],
                    None,
                    false,
                )
                .unwrap()
            };
            assert!(super::linear_nurbs_spine(&spine(0.4)).is_none());
            assert!(super::linear_nurbs_spine(&spine(0.0)).is_some());
        }
    }
    #[test]
    fn audit_regression_circle_recognition_ignores_knot_units() {
        use cadmpeg_ir::geometry::nurbs::NurbsCurve;

        for scale in [1.0, 1e-13] {
            let knots = [0., 0., 0., 1., 1., 2., 2., 3., 3., 4., 4., 4.]
                .map(|knot| knot * scale)
                .to_vec();
            let poles = [
                (1., 0.),
                (1., 1.),
                (0., 1.),
                (-1., 1.),
                (-1., 0.),
                (-1., -1.),
                (0., -1.),
                (1., -1.),
                (1., 0.),
            ]
            .into_iter()
            .map(|(x, y)| Point3::new(x, y, 0.))
            .collect();
            let weights = (0..9)
                .map(|i| {
                    if i % 2 == 0 {
                        1.
                    } else {
                        std::f64::consts::FRAC_1_SQRT_2
                    }
                })
                .collect();
            let circle = NurbsCurve::from_lanes(2, knots, poles, Some(weights), false).unwrap();
            let (_, _, _, radius) = super::rational_four_arc_circle(&circle).unwrap();
            assert!((radius - 1.0).abs() <= 8.0 * f64::EPSILON);
        }
    }

    #[test]
    fn audit_regression_edge_clamping_preserves_real_domain_violation() {
        use cadmpeg_ir::geometry::nurbs::NurbsCurve;
        use cadmpeg_ir::geometry::{Curve, CurveGeometry, SolvedCurveGeometry};
        use cadmpeg_ir::ids::{CurveId, EdgeId, VertexId};
        use cadmpeg_ir::topology::{Edge, EdgeCarrier};

        let id = CurveId::mint("sat:audit:curve#domain").unwrap();
        let curve = NurbsCurve::from_lanes(
            1,
            vec![0., 0., 1e-12, 1e-12],
            vec![Point3::new(0., 0., 0.), Point3::new(1., 0., 0.)],
            None,
            false,
        )
        .unwrap();
        let edge = Edge {
            id: EdgeId::mint("sat:audit:edge#domain").unwrap(),
            carrier: EdgeCarrier::new(Some(id.clone()), Some([-1e-10, 5e-13])).unwrap(),
            start: VertexId::mint("sat:audit:vertex#a").unwrap(),
            end: VertexId::mint("sat:audit:vertex#b").unwrap(),
            tolerance: None,
        };
        let mut out = super::AsmBrep {
            curves: vec![Curve {
                parameter_range: None,
                id,
                geometry: CurveGeometry::Solved(SolvedCurveGeometry::Nurbs(curve)),
                source_object: None,
            }],
            edges: vec![edge],
            ..Default::default()
        };
        super::clamp_edge_ranges_to_carrier_domains(&mut out).unwrap();
        assert_eq!(
            out.edges[0]
                .param_range()
                .map(cadmpeg_ir::units::FiniteVector::get),
            Some([-1e-10, 5e-13])
        );
        out.edges[0].set_param_range(Some(
            cadmpeg_ir::topology::ParameterInterval::new([-1e-23, 5e-13]).unwrap(),
        ));
        super::clamp_edge_ranges_to_carrier_domains(&mut out).unwrap();
        assert_eq!(
            out.edges[0]
                .param_range()
                .map(cadmpeg_ir::units::FiniteVector::get),
            Some([0., 5e-13])
        );
    }
}