Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Delaunay3D/LightmapGI: Improve triangulation #90701

Merged
merged 1 commit into from
Apr 22, 2024
Merged

Delaunay3D/LightmapGI: Improve triangulation #90701

Changes from all commits
Commits
Show all changes
1 commit
Select commit
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
77 changes: 38 additions & 39 deletions core/math/delaunay_3d.h
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -46,7 +46,8 @@ class Delaunay3D {
struct Simplex;

enum {
ACCEL_GRID_SIZE = 16
ACCEL_GRID_SIZE = 16,
QUANTIZATION_MAX = 1 << 16 // A power of two smaller than the 23 bit significand of a float.
};
struct GridPos {
Vector3i pos;
Expand Down Expand Up @@ -173,38 +174,25 @@ class Delaunay3D {

R128 radius2 = rel2_x * rel2_x + rel2_y * rel2_y + rel2_z * rel2_z;

return radius2 < (p_simplex.circum_r2 - R128(0.00001));
return radius2 < (p_simplex.circum_r2 - R128(0.0000000001));
// When this tolerance is too big, it can result in overlapping simplices.
// When it's too small, large amounts of planar simplices are created.
}

static bool simplex_is_coplanar(const Vector3 *p_points, const Simplex &p_simplex) {
Plane p(p_points[p_simplex.points[0]], p_points[p_simplex.points[1]], p_points[p_simplex.points[2]]);
if (ABS(p.distance_to(p_points[p_simplex.points[3]])) < CMP_EPSILON) {
return true;
// Checking every possible distance like this is overkill, but only checking
// one is not enough. If the simplex is almost planar then the vectors p1-p2
// and p1-p3 can be practically collinear, which makes Plane unreliable.
for (uint32_t i = 0; i < 4; i++) {
Plane p(p_points[p_simplex.points[i]], p_points[p_simplex.points[(i + 1) % 4]], p_points[p_simplex.points[(i + 2) % 4]]);
// This tolerance should not be smaller than the one used with
// Plane::has_point() when creating the LightmapGI probe BSP tree.
if (ABS(p.distance_to(p_points[p_simplex.points[(i + 3) % 4]])) < 0.001) {
return true;
}
}

Projection cm;

cm.columns[0][0] = p_points[p_simplex.points[0]].x;
cm.columns[0][1] = p_points[p_simplex.points[1]].x;
cm.columns[0][2] = p_points[p_simplex.points[2]].x;
cm.columns[0][3] = p_points[p_simplex.points[3]].x;

cm.columns[1][0] = p_points[p_simplex.points[0]].y;
cm.columns[1][1] = p_points[p_simplex.points[1]].y;
cm.columns[1][2] = p_points[p_simplex.points[2]].y;
cm.columns[1][3] = p_points[p_simplex.points[3]].y;

cm.columns[2][0] = p_points[p_simplex.points[0]].z;
cm.columns[2][1] = p_points[p_simplex.points[1]].z;
cm.columns[2][2] = p_points[p_simplex.points[2]].z;
cm.columns[2][3] = p_points[p_simplex.points[3]].z;

cm.columns[3][0] = 1.0;
cm.columns[3][1] = 1.0;
cm.columns[3][2] = 1.0;
cm.columns[3][3] = 1.0;

return ABS(cm.determinant()) <= CMP_EPSILON;
return false;
}

public:
Expand All @@ -215,9 +203,10 @@ class Delaunay3D {
static Vector<OutputSimplex> tetrahedralize(const Vector<Vector3> &p_points) {
uint32_t point_count = p_points.size();
Vector3 *points = (Vector3 *)memalloc(sizeof(Vector3) * (point_count + 4));
const Vector3 *src_points = p_points.ptr();
Vector3 proportions;

{
const Vector3 *src_points = p_points.ptr();
AABB rect;
for (uint32_t i = 0; i < point_count; i++) {
Vector3 point = src_points[i];
Expand All @@ -226,17 +215,25 @@ class Delaunay3D {
} else {
rect.expand_to(point);
}
points[i] = point;
}

real_t longest_axis = rect.size[rect.get_longest_axis_index()];
proportions = Vector3(longest_axis, longest_axis, longest_axis) / rect.size;

for (uint32_t i = 0; i < point_count; i++) {
points[i] = (points[i] - rect.position) / rect.size;
// Scale points to the unit cube to better utilize R128 precision
// and quantize to stabilize triangulation over a wide range of
// distances.
points[i] = Vector3(Vector3i((src_points[i] - rect.position) / longest_axis * QUANTIZATION_MAX)) / QUANTIZATION_MAX;
}

const real_t delta_max = Math::sqrt(2.0) * 20.0;
const real_t delta_max = Math::sqrt(2.0) * 100.0;
Vector3 center = Vector3(0.5, 0.5, 0.5);

// any simplex that contains everything is good
// The larger the root simplex is, the more likely it is that the
// triangulation is convex. If it's not absolutely huge, there can
// be missing simplices that are not created for the outermost faces
// of the point cloud if the point density is very low there.
points[point_count + 0] = center + Vector3(0, 1, 0) * delta_max;
points[point_count + 1] = center + Vector3(0, -1, 1) * delta_max;
points[point_count + 2] = center + Vector3(1, -1, -1) * delta_max;
Expand Down Expand Up @@ -271,7 +268,7 @@ class Delaunay3D {
for (uint32_t i = 0; i < point_count; i++) {
bool unique = true;
for (uint32_t j = i + 1; j < point_count; j++) {
if (points[i].is_equal_approx(points[j])) {
if (points[i] == points[j]) {
unique = false;
break;
}
Expand All @@ -280,7 +277,7 @@ class Delaunay3D {
continue;
}

Vector3i grid_pos = Vector3i(points[i] * ACCEL_GRID_SIZE);
Vector3i grid_pos = Vector3i(points[i] * proportions * ACCEL_GRID_SIZE);
grid_pos = grid_pos.clamp(Vector3i(), Vector3i(ACCEL_GRID_SIZE - 1, ACCEL_GRID_SIZE - 1, ACCEL_GRID_SIZE - 1));

for (List<Simplex *>::Element *E = acceleration_grid[grid_pos.x][grid_pos.y][grid_pos.z].front(); E;) {
Expand All @@ -300,14 +297,16 @@ class Delaunay3D {
Triangle t = Triangle(simplex->points[triangle_order[k][0]], simplex->points[triangle_order[k][1]], simplex->points[triangle_order[k][2]]);
uint32_t *p = triangles_inserted.lookup_ptr(t);
if (p) {
// This Delaunay implementation uses the Bowyer-Watson algorithm.
// The rule is that you don't reuse any triangles that were
// shared by any of the retriangulated simplices.
triangles[*p].bad = true;
} else {
triangles_inserted.insert(t, triangles.size());
triangles.push_back(t);
}
}

//remove simplex and continue
simplex_list.erase(simplex->SE);

for (const GridPos &gp : simplex->grid_positions) {
Expand All @@ -334,8 +333,8 @@ class Delaunay3D {

const real_t radius2 = Math::sqrt(double(new_simplex->circum_r2)) + 0.0001;
Vector3 extents = Vector3(radius2, radius2, radius2);
Vector3i from = Vector3i((center - extents) * ACCEL_GRID_SIZE);
Vector3i to = Vector3i((center + extents) * ACCEL_GRID_SIZE);
Vector3i from = Vector3i((center - extents) * proportions * ACCEL_GRID_SIZE);
Vector3i to = Vector3i((center + extents) * proportions * ACCEL_GRID_SIZE);
from = from.clamp(Vector3i(), Vector3i(ACCEL_GRID_SIZE - 1, ACCEL_GRID_SIZE - 1, ACCEL_GRID_SIZE - 1));
to = to.clamp(Vector3i(), Vector3i(ACCEL_GRID_SIZE - 1, ACCEL_GRID_SIZE - 1, ACCEL_GRID_SIZE - 1));

Expand Down Expand Up @@ -370,7 +369,7 @@ class Delaunay3D {
break;
}
}
if (invalid || simplex_is_coplanar(points, *simplex)) {
if (invalid || simplex_is_coplanar(src_points, *simplex)) {
memdelete(simplex);
continue;
}
Expand Down
Loading