asgardius/virtualx-engine
1
0
Fork 0
Code Issues Pull requests Projects Releases Packages Wiki Activity Actions

Merge pull request #90701 from permelin/improve-delaunay3d-triangulation

Browse source 
Delaunay3D/LightmapGI: Improve triangulation
This commit is contained in:
Rémi Verschelde 2024-04-22 12:52:12 +02:00
commit 4b6629978e
No known key found for this signature in database
GPG key ID: C3336907360768E1
1 changed files with 38 additions and 39 deletions
Show stats Download patch file Download diff file Expand all files Collapse all files

77
core/math/delaunay_3d.h
View file

@ -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;
@ -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:
@ -215,9 +203,10 @@ public:
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];
@ -226,17 +215,25 @@ public:
} 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;
@ -271,7 +268,7 @@ public:
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;
}
@ -280,7 +277,7 @@ public:
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;) {
@ -300,6 +297,9 @@ public:
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());
@ -307,7 +307,6 @@ public:
}
}
//remove simplex and continue
simplex_list.erase(simplex->SE);
for (const GridPos &gp : simplex->grid_positions) {
@ -334,8 +333,8 @@ public:
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));
@ -370,7 +369,7 @@ public:
break;
}
}
if (invalid || simplex_is_coplanar(points, *simplex)) {
if (invalid || simplex_is_coplanar(src_points, *simplex)) {
memdelete(simplex);
continue;
}