Effigy.Tests/BvhFaceTests.cs
using System;
using System.Collections.Generic;
using System.Linq;
using Effigy;
using static Effigy.Tests.Report;
namespace Effigy.Tests;
/// <summary>
/// The face query a paint dab depends on, judged against what it has to be: the SAME faces a
/// linear scan over the triangulated surface would find, each once.
///
/// The box test is a pruning step and the triangle test is the verdict. A query that stopped at the
/// box would hand the caller faces the brush never touched — every returned face is rasterised into
/// texels, so a false positive is visible paint, not a loose number. That is why the brute-force
/// agreement below is the check that matters, and everything else is a sharper instance of it.
/// </summary>
public static class BvhFaceTests
{
public static void Run()
{
Section( "bvh faces-in-radius: agrees with a linear scan" );
TestAgreesWithBruteForce();
TestFaceCentreFindsOnlyThatFace();
TestEdgeFindsBothFaces();
TestWholeModelReturnsEveryFaceOnce();
TestOutsideReturnsEmpty();
TestZeroRadiusDoesNotThrow();
}
static void TestAgreesWithBruteForce()
{
// A subdivided sphere is the fixture that would expose a box-only query: its triangles sit
// at every angle, so a sphere centred in a valley must not catch the far wall of that valley
// merely because both bounding boxes overlap. Random points and radii rather than a few
// hand-picked ones, so a pruning bug has to be exactly as lucky as the box to hide.
var mesh = CatmullClark.Subdivide( Primitives.QuadSphere( 1f, 4 ), 2 );
var bvh = MeshBVH.Build( mesh );
var rng = new Random( 20260904 );
var points = new Vec3[100];
for ( var i = 0; i < points.Length; i++ )
{
points[i] = new Vec3(
(float)(rng.NextDouble() * 2.0 - 1.0) * 1.2f,
(float)(rng.NextDouble() * 2.0 - 1.0) * 1.2f,
(float)(rng.NextDouble() * 2.0 - 1.0) * 1.2f );
}
var radii = new[] { 0.05f, 0.2f, 0.5f, 1.5f };
var results = new List<int>();
foreach ( var radius in radii )
{
var mismatches = 0;
foreach ( var point in points )
{
bvh.FacesInRadius( mesh, point, radius, results );
var expected = BruteForceFaces( mesh, point, radius );
results.Sort();
expected.Sort();
if ( !results.SequenceEqual( expected ) )
mismatches++;
}
Check( $"radius {radius:0.##} agrees with a linear scan over {mesh.FaceCount} faces",
mismatches == 0, $"{mismatches} of {points.Length} points disagreed" );
}
}
static void TestFaceCentreFindsOnlyThatFace()
{
// A sphere hugging one face's middle must return that face and nothing else. The opposite
// face's box never overlaps, but the four side faces' boxes reach within the sphere's box —
// it is the triangle test, not the box test, that keeps them out.
var mesh = Primitives.Box( 2, 2, 2 );
var bvh = MeshBVH.Build( mesh );
var centre = mesh.FaceCentroid( mesh.Faces[1] );
var results = new List<int>();
bvh.FacesInRadius( mesh, centre, 0.1f, results );
Check( "a small sphere on a face centre returns exactly that face",
results.Count == 1 && results[0] == 1, $"got [{string.Join( ", ", results )}]" );
}
static void TestEdgeFindsBothFaces()
{
// The edge between the bottom and the -Y side, shared by exactly two faces. A sphere centred
// ON the edge holds a sliver of both triangles, so both faces must come back.
var mesh = Primitives.Box( 2, 2, 2 );
var bvh = MeshBVH.Build( mesh );
var edge = (mesh.Positions[0] + mesh.Positions[1]) * 0.5f;
var results = new List<int>();
bvh.FacesInRadius( mesh, edge, 0.1f, results );
results.Sort();
Check( "a sphere centred on an edge returns both faces sharing it",
results.Count == 2 && results[0] == 0 && results[1] == 2,
$"got [{string.Join( ", ", results )}]" );
}
static void TestWholeModelReturnsEveryFaceOnce()
{
// A sphere big enough to swallow the whole solid must return every face, and none twice —
// each face's two triangles both qualify and the query has to name the face once, not once
// per triangle.
var mesh = Primitives.Box( 2, 2, 2 );
var bvh = MeshBVH.Build( mesh );
var results = new List<int>();
bvh.FacesInRadius( mesh, Vec3.Zero, 2f, results );
Check( "a sphere containing the model returns every face", results.Count == mesh.FaceCount,
$"got {results.Count} of {mesh.FaceCount}" );
Check( "and each face exactly once", results.Distinct().Count() == results.Count,
$"duplicates: [{string.Join( ", ", results )}]" );
}
static void TestOutsideReturnsEmpty()
{
var mesh = Primitives.Box( 2, 2, 2 );
var bvh = MeshBVH.Build( mesh );
var results = new List<int>();
bvh.FacesInRadius( mesh, new Vec3( 10, 10, 10 ), 0.5f, results );
Check( "a sphere entirely outside returns nothing", results.Count == 0,
$"got [{string.Join( ", ", results )}]" );
}
static void TestZeroRadiusDoesNotThrow()
{
// A zero-radius dab is a legal click: it must not divide by anything, and it must still find
// the face the point sits exactly on.
var mesh = Primitives.Box( 2, 2, 2 );
var bvh = MeshBVH.Build( mesh );
var centre = mesh.FaceCentroid( mesh.Faces[1] );
var results = new List<int>();
var ok = true;
try
{
bvh.FacesInRadius( mesh, centre, 0f, results );
}
catch
{
ok = false;
}
Check( "a radius of zero does not throw", ok );
Check( "and still finds the face the point sits on", results.Count == 1 && results[0] == 1,
$"got [{string.Join( ", ", results )}]" );
}
/// <summary>Reference implementation: every face, triangulated, tested by the closest point on
/// each triangle. Independent of the BVH's traversal, so any pruning error shows up as a set
/// that disagrees with this one.</summary>
static List<int> BruteForceFaces( PolyMesh mesh, Vec3 point, float radius )
{
var r2 = radius * radius;
var found = new List<int>();
for ( var fi = 0; fi < mesh.FaceCount; fi++ )
{
var face = mesh.Faces[fi];
if ( face.Count < 3 )
continue;
var corners = new List<Vec3>( face.Count );
for ( var c = 0; c < face.Count; c++ )
corners.Add( mesh.Positions[face.Indices[c]] );
foreach ( var (ia, ib, ic) in Triangulate.Face( corners ) )
{
if ( ClosestPointSq( point, corners[ia], corners[ib], corners[ic] ) <= r2 )
{
found.Add( fi );
break;
}
}
}
return found;
}
/// <summary>Squared distance from a point to a triangle, via the closest point clamped to the
/// vertex, edge and face regions in turn.</summary>
static float ClosestPointSq( Vec3 p, Vec3 a, Vec3 b, Vec3 c )
{
var ab = b - a;
var ac = c - a;
var ap = p - a;
var d1 = Vec3.Dot( ab, ap );
var d2 = Vec3.Dot( ac, ap );
if ( d1 <= 0f && d2 <= 0f )
return ap.LengthSquared;
var bp = p - b;
var d3 = Vec3.Dot( ab, bp );
var d4 = Vec3.Dot( ac, bp );
if ( d3 >= 0f && d4 <= d3 )
return bp.LengthSquared;
var vc = d1 * d4 - d3 * d2;
if ( vc <= 0f && d1 >= 0f && d3 <= 0f )
return (a + ab * (d1 / (d1 - d3)) - p).LengthSquared;
var cp = p - c;
var d5 = Vec3.Dot( ab, cp );
var d6 = Vec3.Dot( ac, cp );
if ( d6 >= 0f && d5 <= d6 )
return cp.LengthSquared;
var vb = d5 * d2 - d1 * d6;
if ( vb <= 0f && d2 >= 0f && d6 <= 0f )
return (a + ac * (d2 / (d2 - d6)) - p).LengthSquared;
var va = d3 * d6 - d5 * d4;
if ( va <= 0f && (d4 - d3) >= 0f && (d5 - d6) >= 0f )
return (b + (c - b) * ((d4 - d3) / ((d4 - d3) + (d5 - d6))) - p).LengthSquared;
var denom = 1f / (va + vb + vc);
return (a + ab * (vb * denom) + ac * (vc * denom) - p).LengthSquared;
}
}