Self-intersection detection and resolution

What it does

Two TBPs that bracket the self-intersection problem:

• surface%find_self_intersections(pairs, status, include_touching) — detection. Returns the list of every facet pair (a, b) whose triangles cross or overlap in a common plane, decided exactly, with the kind of contact and the 3D segment endpoints of the intersection.
• surface%resolve_self_intersections(status) — resolution. Cleans the surface up by collapsing every self-crossing through a self-boolean union (boolean(self, self, BOOL_UNION)).

Real-world STL exports routinely have self-intersections. They slip through is_watertight() (the topological-edge test is satisfied even when triangles geometrically cross), so any downstream operation that assumes a clean solid — booleans, signed distance, voxelization for cut-cell — silently produces wrong answers. Detection is the first step; resolution is what you reach for when sanitize alone can't fix the input.

Pipeline (detection)

1. Tree-vs-tree broad phase. Walk the AABB tree against itself, collecting every leaf-pair whose AABBs overlap. Reuses the existing enumerate_children dispatcher so SAH BVH and octree both share the traversal code.

2. Pool-aware adjacency filter. The shared vertices of the pair are counted by vertex-pool ID (exact coordinate equality when the pool is not in use). Two facets sharing an edge meet along it, as expected; they are a defect only when exactly coplanar and folded over it (both free vertices on the same side of the edge, exact orient2d): SELFINT_EDGE_FOLD (issue #69).

3. Narrow phase. No shared vertex (issue #69): the exact pair_cut of the boolean (fossil_tritri_exact, exact orient3d signs and implicit section points). A common segment of positive length is a crossing (SELFINT_CROSSING), reported with its endpoints rounded to doubles. A single common point is a touch (SELFINT_TOUCHING), reported only with include_touching=.true.. Two exactly coplanar triangles are compared in the projection that keeps them non-degenerate (exact orient2d, separating-axis test on the six edge lines): overlapping interiors are a defect (SELFINT_COPLANAR, p and q the two centroids), boundaries meeting without area a touch. Before #69 this was the floating Möller test, which never reported a coplanar overlap. One shared vertex (issue #46): the exact pair_classify (orient3d); a common segment of positive length is a crossing, the vertex alone is a touch and is not reported. Pairs whose other vertices lie strictly on one side of the other plane can only touch and are skipped. Exactly coplanar pairs sharing a vertex (issue #69; skipped to 1e-12 rad before, so crossings tilted by less went unseen) are decided in the projection: overlapping interiors are reported (SELFINT_COPLANAR), a common edge segment from the vertex is a touch (include_touching).

4. Non-manifold edges (issue #69). Sheets crossing along a shared edge (e.g. two cubes merged so that one passes through the other along a common edge) meet only on the edge, so no pair test sees them. A pass over the edges with 4 or more facets orders the facets around each edge by the angle of their half-planes, exactly (orient3d on input vertices), and notes whether each traverses the edge forwards or backwards. Embedded oriented sheets alternate around the edge (the winding number of the wedges stays in {w, w + 1}); two adjacent facets with the same direction bound a wedge wound twice, where the sheets cross. Each such pair is reported as SELFINT_SHEETS, with the edge as p–q. Fans with unequal numbers of forward and backward facets (an orientation or boundary defect, not a crossing) and fans with two facets in the same half-plane (folds, reported pairwise) are skipped. The criterion assumes consistently oriented sheets: on a misoriented input a reported pair can be an orientation defect rather than a crossing (sanitize orients the normals first).

Not reported: duplicate facets (three shared vertices), which remove_duplicate_facets handles.

Pipeline (resolution)

resolve_self_intersections is a one-line wrapper over the boolean engine:

call self%boolean(other=self, op=BOOL_UNION, status=status)

The §1.1 arrangement code treats self vs self symmetrically: every self-crossing becomes an arrangement segment, every input triangle gets retriangulated, the per-sub-triangle winding-number tagging keeps only the outer-manifold sub-triangles. Inherits §1.1's axis-aligned coplanar limitation.

API

Detection

call surface%find_self_intersections(pairs=pairs, status=status, include_touching=touching)
   type(intersection_pair_t), allocatable, intent(out) :: pairs(:)
   integer(I4P),                           intent(out), optional :: status
   logical,                                intent(in),  optional :: include_touching ! default .false.

pairs is always allocated (zero-length if no intersections found). Each record carries:

• pairs(i)%a, pairs(i)%b — facet IDs with a < b.
• pairs(i)%kind — SELFINT_CROSSING (a common segment of positive length), SELFINT_COPLANAR (coplanar, overlapping interiors), SELFINT_EDGE_FOLD (sharing an edge, coplanar, folded over it), SELFINT_SHEETS (adjacent around a non-manifold edge with the same direction: the sheets cross along the edge) or, with include_touching, SELFINT_TOUCHING (one common point, or coplanar boundaries meeting without area).
• pairs(i)%p, pairs(i)%q — segment endpoints in world coordinates (the common point twice for a touch, the two centroids for a coplanar overlap).

Resolution

call surface%resolve_self_intersections(status=status)
   integer(I4P), intent(out), optional :: status  ! BOOL_STATUS_*

Mutates self in place. Output is self-intersection-free modulo the §1.1 coplanar limitation — perturb the input by ~1e-6 * bbox_diag along each axis if you suspect coplanar self-overlaps.

Example: detection

The dragon fixture in this repo is known dirty — it ships with 1191 self-intersection pairs: 1167 crossings (759 before issue #46 added the pairs sharing one vertex) and, since issue #69, 24 pairs of sheets crossing at its non-manifold edges. It is the regression workhorse for §1.2:

program ex_find_self_intersections
use fossil
use penf,   only : I4P
implicit none

type(surface_stl_object)               :: dragon
type(intersection_pair_t), allocatable :: pairs(:)
integer(I4P)                           :: status

call dragon%load_from_file(file_name='src/tests/dragon.stl', guess_format=.true.)
call dragon%find_self_intersections(pairs=pairs, status=status)
print '(A,I0,A,I0)', 'dragon: ', dragon%get_facets_number(), &
      ' facets, ', size(pairs), ' self-intersection pairs'
endprogram ex_find_self_intersections

For a clean mesh (cube.stl, bunny.stl), the same call returns size(pairs) == 0 — the no-false-positive guarantee of the adjacency filter (every shared edge and every fan around a vertex would otherwise register as a pseudo-intersection).

Example: resolution

program ex_resolve_self_intersections
use fossil
use penf,   only : I4P
implicit none

type(surface_stl_object)               :: dragon
type(intersection_pair_t), allocatable :: pairs(:)
integer(I4P)                           :: status

call dragon%load_from_file(file_name='src/tests/dragon.stl', guess_format=.true.)
call dragon%find_self_intersections(pairs=pairs, status=status)
print '(A,I0)', 'before: ', size(pairs)

call dragon%resolve_self_intersections(status=status)

call dragon%find_self_intersections(pairs=pairs, status=status)
print '(A,I0,A,I0)', 'after:  ', size(pairs), '  status=', status
endprogram ex_resolve_self_intersections

Visual reference

The dragon fixture (6588 facets), drawn translucent, with the segment p–q of each of the 1191 pairs find_self_intersections reports:

Drawn opaque, the dragon looks fine: the crossings are small folds the eye does not catch. That's the whole point of find_self_intersections — these are the defects that escape eye inspection and that is_watertight() doesn't catch.

Known limitations

• Resolution inherits §1.1's coplanar limitation. Self- intersections that involve coplanar facet pairs (axis-aligned cases especially) may produce a resolve output with the same silent-wrong-volume behaviour documented for the boolean. Workaround is the same: perturb by ~1e-6 * bbox_diag. Most real-world dirty STL has generic-position self-intersections (random angles, no coplanar overlap) and the resolution works cleanly.
• resolve_self_intersections returns the topology of BOOL_UNION applied to the surface against itself. Interior cavities formed by self-intersections are dropped, not re-meshed — if the input had meaningful interior structure, it's gone after resolve.
• Detection is O(n_pairs) in the worst case, where n_pairs is the number of overlapping leaf-pairs found in the broad phase. For the dragon fixture (~6.6k facets) this is ~10× faster than brute- force O(n²) pair enumeration; for cleaner meshes the speedup is much larger because few leaves overlap at all.

See also

• booleans — the resolution path is a one-liner over the boolean engine; same coplanar limitation applies.
• surface%sanitize — local defect repair (degenerate / duplicate facets, normal orientation). Run this before find_self_intersections to filter out intersection records that exist only because of degenerate-facet noise.
• alpha wrapping — the heavyweight alternative for inputs too dirty for resolve to handle (large holes, missing facets that aren't just self-intersection artifacts).

References

• Möller, A Fast Triangle-Triangle Intersection Test, Journal of Graphics Tools 2(2), 1997. The narrow-phase test until issue #69, replaced by the exact pair_cut (see booleans).
• Zhou, Grinspun, Zorin & Jacobson, Mesh Arrangements for Solid Geometry, SIGGRAPH 2016. The reference algorithm for the §1.1 boolean that resolve_self_intersections reuses.
• libigl igl::SelfIntersectMesh — the closest open-source reference; FOSSIL's detector matches its broad-phase + narrow-phase decomposition.