Skip to content

Instantly share code, notes, and snippets.

@DaveH355
Created January 8, 2025 12:37
Show Gist options
  • Select an option

  • Save DaveH355/341a2260f67aabc91a5d51fbae84594d to your computer and use it in GitHub Desktop.

Select an option

Save DaveH355/341a2260f67aabc91a5d51fbae84594d to your computer and use it in GitHub Desktop.
Triangle-AABB intersection as God and Christer Ericson intended
#include <cassert>
#include <cstdlib>
#include <glm/common.hpp>
#include <glm/ext/vector_float3.hpp>
#include <glm/geometric.hpp>
#include <glm/ext/scalar_common.hpp>
#include <limits>
/*
* Credit to the book Real Time Collision Detection
* https://realtimecollisiondetection.net/
*/
struct Plane
{
glm::vec3 n;
float d;
};
// required: min <= max
class AABB
{
public:
glm::vec3 min = glm::vec3(std::numeric_limits<float>::max());
glm::vec3 max = glm::vec3(std::numeric_limits<float>::min());
AABB(glm::vec3 min, glm::vec3 max) : min(min), max(max)
{
assert(min.x <= max.x);
assert(min.y <= max.y);
assert(min.z <= max.z);
}
AABB() = default;
};
bool test_aabb_plane(const AABB &b, const Plane &p)
{
// These two lines not necessary with a (center, extents) AABB representation
glm::vec3 c = (b.max + b.min) * 0.5f; // Compute AABB center
glm::vec3 e = b.max - c; // Compute positive extents
// Compute the projection interval radius of b onto L(t) = b.c + t * p.n
float r = e[0] * std::abs(p.n[0]) + e[1] * std::abs(p.n[1]) +
e[2] * std::abs(p.n[2]);
// Compute distance of box center from plane
float s = dot(p.n, c) - p.d;
// Intersection occurs when distance s falls within [-r,+r] interval
return std::abs(s) <= r;
}
bool test_aabb_triangle(glm::vec3 v0, glm::vec3 v1, glm::vec3 v2, const AABB &b)
{
glm::vec3 a = v0;
// Compute box center and extents (if not already given in that format)
glm::vec3 c = (b.min + b.max) * 0.5f;
float e0 = (b.max.x - b.min.x) * 0.5f;
float e1 = (b.max.y - b.min.y) * 0.5f;
float e2 = (b.max.z - b.min.z) * 0.5f;
// Translate triangle as conceptually moving AABB to origin
v0 = v0 - c;
v1 = v1 - c;
v2 = v2 - c;
// Compute edge vectors for triangle
glm::vec3 f0 = v1 - v0, f1 = v2 - v1, f2 = v0 - v2;
// Test axes a00..a22 (category 3) //
float p0, p1, p2, r;
glm::vec3 a00{0, -f0.z, f0.y};
glm::vec3 a01{0, -f1.z, f1.y};
glm::vec3 a02{0, -f2.z, f2.y};
glm::vec3 a10{f0.z, 0, -f0.x};
glm::vec3 a11{f1.z, 0, -f1.x};
glm::vec3 a12{f2.z, 0, -f2.x};
glm::vec3 a20{-f0.y, f0.x, 0};
glm::vec3 a21{-f1.y, f1.x, 0};
glm::vec3 a22{-f2.y, f2.x, 0};
// Test axis a00 (p0 == p1) delete p1
p0 = dot(v0, a00);
p2 = dot(v2, a00);
r = e1 * std::abs(f0.z) + e2 * std::abs(f0.y);
if (glm::max(-glm::max(p0, p2), glm::min(p0, p2)) > r)
return false; // Axis is a separating axis
// Test axis a01 (p1 == p2) del p2
p0 = dot(v0, a01);
p1 = dot(v1, a01);
r = e1 * std::abs(f1.z) + e2 * std::abs(f1.y);
if (glm::max(-glm::max(p0, p1), glm::min(p0, p1)) > r)
return false;
// Test axis a02 (p0 == p2) del p2
p0 = dot(v0, a02);
p1 = dot(v1, a02);
r = e1 * std::abs(f2.z) + e2 * std::abs(f2.y);
if (glm::max(-glm::max(p0, p1), glm::min(p0, p1)) > r)
return false;
// Test axis a10 (p0 == p1) del p1
p0 = dot(v0, a10);
p2 = dot(v2, a10);
r = e0 * std::abs(f0.z) + e2 * std::abs(f0.x);
if (glm::max(-glm::max(p0, p2), glm::min(p0, p2)) > r)
return false;
// Test axis a11 (p1 == p2) del p2
p0 = dot(v0, a11);
p1 = dot(v1, a11);
r = e0 * std::abs(f1.z) + e2 * std::abs(f1.x);
if (glm::max(-glm::max(p0, p1), glm::min(p0, p1)) > r)
return false;
// Test axis a12 (p0 == p2) del p2
p0 = dot(v0, a12);
p1 = dot(v1, a12);
r = e0 * std::abs(f2.z) + e2 * std::abs(f2.x);
if (glm::max(-glm::max(p0, p1), glm::min(p0, p1)) > r)
return false;
// Test axis 20 (p0 == p1) del p1
p0 = dot(v0, a20);
p2 = dot(v2, a20);
r = e0 * std::abs(f0.y) + e1 * std::abs(f0.x);
if (glm::max(-glm::max(p0, p2), glm::min(p0, p2)) > r)
return false;
// Test axis 21 (p1 == p2) del p2
p0 = dot(v0, a21);
p1 = dot(v1, a21);
r = e0 * std::abs(f1.y) + e1 * std::abs(f1.x);
if (glm::max(-glm::max(p0, p1), glm::min(p0, p1)) > r)
return false;
// Test axis a22 (p0 == p2) del p2
p0 = dot(v0, a22);
p1 = dot(v1, a22);
r = e0 * std::abs(f2.y) + e1 * std::abs(f2.x);
if (glm::max(-glm::max(p0, p1), glm::min(p0, p1)) > r)
return false;
// Test the three axes corresponding to the face normals of AABB b
// (category 1).
// Exit if... [-e0, e0] and [min(v0.x,v1.x,v2.x), max(v0.x,v1.x,v2.x)] do not
// overlap
if (glm::max(v0.x, v1.x, v2.x) < -e0 || glm::min(v0.x, v1.x, v2.x) > e0)
return false;
// ... [-e1, e1] and [min(v0.y,v1.y,v2.y), max(v0.y,v1.y,v2.y)] do not overlap
if (glm::max(v0.y, v1.y, v2.y) < -e1 || glm::min(v0.y, v1.y, v2.y) > e1)
return false;
// ... [-e2, e2] and [min(v0.z,v1.z,v2.z), max(v0.z,v1.z,v2.z)] do not overlap
if (glm::max(v0.z, v1.z, v2.z) < -e2 || glm::min(v0.z, v1.z, v2.z) > e2)
return false;
// Test separating axis corresponding to triangle face normal (category 2)
Plane p{};
p.n = cross(f0, f1);
// using point on original triangle to get correct plane equation.
p.d = dot(p.n, a);
// alternatively, use v0 and move the AABB to origin
// https://realtimecollisiondetection.net/books/rtcd/errata/
return test_aabb_plane(b, p);
}
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment