#include "volume_grid_nodes.h"
#include "../core/proc_geometry.h"
#include "../core/proc_mesh.h"
#include "../core/proc_volume.h"
#include "../graph/graph_executor.h"
#include "dc_generator_mesh_sdf.h"
#include "dc_density_grid.h"
#include "dc_mesher.h"
#include <godot_cpp/variant/utility_functions.hpp>
#include <godot_cpp/classes/time.hpp>
#include <thread>
#include <vector>
namespace godot {
// ── Mesh to SDF Grid ────────────────────────────────────────────────────────
void MeshToSDFGridNode::_bind_methods() {}
MeshToSDFGridNode::MeshToSDFGridNode() {
// ⛔ THE DEFAULTS ARE BLENDER'S, and that is not cosmetic. The bridge sends
// only MODIFIED defaults, so a socket left alone in Blender arrives carrying
// OURS. Voxel Size was 0.1 against Blender's 0.30 -- three times the
// resolution on every axis, TWENTY-SEVEN times the samples, and every one of
// them a BVH closest-point query. That alone was most of the 75 s the
// bridged mountain spent in this node. Band Width is an INT on Blender's
// node, not a float.
input_defs.push_back({ "Mesh", SocketType::Geometry, Variant() }); // 0
input_defs.push_back({ "Voxel Size", SocketType::Float, 0.3f }); // 1
input_defs.push_back({ "Band Width", SocketType::Int, 3 }); // 2
output_defs.push_back({ "Volume", SocketType::Geometry, Variant() });
}
void MeshToSDFGridNode::execute(ProcEvalContext &p_ctx, int p_node_idx) {
const uint64_t t_node0 = Time::get_singleton()->get_ticks_usec();
Ref<ProcGeometry> in_geo = p_ctx.read_geometry(p_node_idx, 0);
Ref<ProcGeometry> out_geo;
out_geo.instantiate();
Ref<ProcMesh> mesh = in_geo.is_valid() ? in_geo->get_mesh() : Ref<ProcMesh>();
// `indices` is the CORNER array whenever the mesh carries ngons, so reading
// it three at a time is a triangle list only by accident. Blender's grid,
// circle and cylinder fills all arrive as ngons.
PackedInt32Array tris = mesh.is_valid() ? mesh->triangulated_indices() : PackedInt32Array();
if (mesh.is_null() || mesh->positions.size() == 0 || tris.size() < 3) {
UtilityFunctions::push_warning(
"proc-nodes: Mesh to SDF Grid got no triangles, so there is no "
"surface to measure distance to. Passing an empty volume through.");
out_geo->set_volume(Ref<ProcVolume>(memnew(ProcVolume)));
p_ctx.write_geometry(p_node_idx, 0, out_geo);
return;
}
float voxel = p_ctx.read_float(p_node_idx, 1);
if (voxel <= 0.0f) voxel = 0.1f;
float band = (float)p_ctx.read_int(p_node_idx, 2);
if (band < 0.0f) band = 0.0f;
// libdc's evaluator, not a second one: the pseudonormal method with a BVH
// (Baerentzen & Aanaes), already used by dc-terrain's mesh-SDF generator.
const int vcount = mesh->positions.size();
std::vector<float> pos((size_t)vcount * 3);
for (int i = 0; i < vcount; i++) {
pos[i * 3 + 0] = (float)mesh->positions[i].x;
pos[i * 3 + 1] = (float)mesh->positions[i].y;
pos[i * 3 + 2] = (float)mesh->positions[i].z;
}
// ⛔ DROP THE SLIVER TRIANGLES, AND NOT ONLY THE EXACTLY DEGENERATE ONES.
// Blender's uv-sphere cap is a QUAD with the pole vertex written twice, so
// fanning it about corner 0 leaves one zero-area triangle per segment per
// pole. libdc takes a face normal as normalize((v1-v0)x(v2-v0)) and its
// normalize() only zeroes out below 1e-12 -- so what it does with the sliver
// depends on whether the two pole copies are bit-identical:
//
// north pole (0, 1, 0) exactly, twice |cross| = 0 -> (0,0,0)
// south pole (-8.7e-08, -1, -0) and
// (-8.6e-08, -1, -1.7e-08) |cross| = 3.3e-09 -> a UNIT
// vector of pure rounding noise
//
// and the angle that noise normal is weighted by, at a corner of a sliver, is
// nearly pi -- against ~0.038 rad for the real cap triangle beside it. So the
// noise WINS the vertex pseudonormal: measured 2026-09-06, the south pole's
// pseudonormal came out (-0.108, +0.994, +0.001), pointing INTO the sphere,
// and the field below the south pole read +0.75/+0.50/+0.25 where it should
// read -0.75/-0.50/-0.25. The north pole, whose copies are bit-identical, was
// correct. That asymmetry is the whole tell.
//
// The test is the triangle's HEIGHT -- twice the area over its longest edge --
// against the mesh's own size, so it says "this has no thickness, therefore no
// normal" at any scale. The same argument applies to libdc's own loaders,
// which is where this would belong if that library were only ours to change.
Vector3 bb_min = mesh->positions[0];
Vector3 bb_max = mesh->positions[0];
for (int i = 1; i < vcount; i++) {
bb_min = bb_min.min(mesh->positions[i]);
bb_max = bb_max.max(mesh->positions[i]);
}
const double min_height = (double)(bb_max - bb_min).length() * 1e-6;
std::vector<int> idx;
idx.reserve(tris.size());
for (int i = 0; i + 2 < tris.size(); i += 3) {
const Vector3 &a = mesh->positions[tris[i]];
const Vector3 &b = mesh->positions[tris[i + 1]];
const Vector3 &c = mesh->positions[tris[i + 2]];
const double twice_area = (double)(b - a).cross(c - a).length();
const double longest = (double)MAX(MAX((b - a).length(), (c - b).length()),
(a - c).length());
if (longest <= 0.0 || twice_area / longest < min_height) {
continue;
}
idx.push_back(tris[i]);
idx.push_back(tris[i + 1]);
idx.push_back(tris[i + 2]);
}
if (idx.size() < 3) {
UtilityFunctions::push_warning(
"proc-nodes: Mesh to SDF Grid got no triangle with any area. "
"Passing an empty volume through.");
out_geo->set_volume(Ref<ProcVolume>(memnew(ProcVolume)));
p_ctx.write_geometry(p_node_idx, 0, out_geo);
return;
}
// ⛔ THE WINDING DECIDES THE SIGN, and proc-nodes' primitives do not agree on
// one. libdc calls the side its pseudonormals point to the OUTSIDE, so its
// stated convention -- positive = solid -- holds only for a CLOCKWISE mesh,
// which is Godot's. Measured 2026-09-06 as the enclosed volume under that same
// cross product: cube -1.000 and cylinder -6.243 (clockwise), uv sphere +4.122
// and ico sphere +4.047 (counter-clockwise, because they carry Blender's
// faces). The same call therefore returns distance-outside as POSITIVE for one
// and distance-inside as positive for the other, and a unit sphere came out
// inside out -- 3112 of 3375 samples "solid", most solid at the box corners --
// with the extractor dutifully meshing the box.
//
// A closed mesh's own enclosed volume says which it is: sum v0 . (v1 x v2) is
// positive exactly when that cross product points out of it. Deriving the sign
// per mesh is what makes this node winding-independent, which it has to be:
// the graphs that reach it come from Blender.
double signed_volume = 0.0;
for (size_t t = 0; t + 2 < idx.size(); t += 3) {
const Vector3 &a = mesh->positions[idx[t]];
const Vector3 &b = mesh->positions[idx[t + 1]];
const Vector3 &c = mesh->positions[idx[t + 2]];
signed_volume += (double)a.dot(b.cross(c));
}
signed_volume /= 6.0;
const float sdf_sign = (signed_volume > 0.0) ? -1.0f : 1.0f;
// Closed means the flood from the grid boundary is a trustworthy sign, which
// is what lets the pseudonormals go. Decided here because load_arrays has to
// know before it builds anything.
const bool closed_precheck = Math::abs(signed_volume) >= 1e-9;
if (Math::abs(signed_volume) < 1e-9) {
// Flat or open: there is no inside, so a "signed" distance is only the
// side the faces happen to point. Say so rather than pretend.
UtilityFunctions::push_warning(
"proc-nodes: Mesh to SDF Grid got a mesh enclosing no volume (flat "
"or not closed). The sign of the field is then whichever side the "
"faces point, and the extracted surface will follow it.");
}
const uint64_t t_setup = Time::get_singleton()->get_ticks_usec();
DCGeneratorMeshSDF sdf;
// ⛔ NO PSEUDONORMALS. They exist only to sign a query, and they cost 5554 ms
// of a 6663 ms load_arrays here -- five times the BVH build -- because they
// need per-vertex AND per-EDGE angle-weighted normals over 5.16M triangles.
// This node already determines inside from outside globally: it floods the
// outside from the grid boundary, which is outside by construction because
// the bounds are padded. Asking for unsigned distance and applying that sign
// is the same answer without the 5.5 s.
//
// Only for a CLOSED mesh -- the flood is what supplies the sign, and it
// leaks through a hole. An open mesh keeps the pseudonormal path.
sdf.set_unsigned_only(closed_precheck);
if (!sdf.load_arrays(pos.data(), vcount, idx.data(), (int)(idx.size() / 3))) {
UtilityFunctions::push_error(
"proc-nodes: Mesh to SDF Grid could not build a BVH for this mesh.");
p_ctx.write_geometry(p_node_idx, 0, out_geo);
return;
}
float mn[3], mx[3];
if (!sdf.get_bounds(mn, mx)) {
UtilityFunctions::push_error("proc-nodes: Mesh to SDF Grid got no bounds.");
p_ctx.write_geometry(p_node_idx, 0, out_geo);
return;
}
// The band is the margin the field is sampled OUTSIDE the mesh, in voxels.
// Without it the surface sits exactly on the grid edge and the extractor has
// no outside to find a crossing against.
const float pad = voxel * band;
for (int a = 0; a < 3; a++) { mn[a] -= pad; mx[a] += pad; }
int res[3];
for (int a = 0; a < 3; a++) {
res[a] = (int)((mx[a] - mn[a]) / voxel) + 1;
if (res[a] < 2) res[a] = 2;
// A grid is memory cubed. 256^3 floats is 64 MB; past that a stray
// voxel size turns into an allocation nobody meant to ask for.
if (res[a] > 256) res[a] = 256;
}
Ref<ProcVolume> vol;
vol.instantiate();
vol->grid_origin = Vector3(mn[0], mn[1], mn[2]);
vol->grid_resolution = Vector3i(res[0], res[1], res[2]);
vol->grid_size = Vector3(res[0] * voxel, res[1] * voxel, res[2] * voxel);
vol->density_grid.resize(res[0] * res[1] * res[2]);
float *w = vol->density_grid.ptrw();
const int64_t nx = res[0], ny = res[1], nz = res[2];
const int64_t total = nx * ny * nz;
auto at = [nx, ny](int64_t x, int64_t y, int64_t z) { return (z * ny + y) * nx + x; };
// ── NARROW BAND ─────────────────────────────────────────────────────────
//
// Every voxel used to get an exact closest-point BVH query. For the bridged
// voronoi mountain that is 251 x 251 x 162 = 10.2 MILLION queries against a
// 5.16M triangle BVH -- 53.7 s, threaded, 73% of the whole evaluation.
//
// Dual contouring only reads values NEAR the zero crossing: it finds sign
// changes on edges and interpolates where the magnitude is small. Every sign
// change lies inside a shell around the triangles, so the exact field is
// only needed there. Blender does the same thing and has always done --
// openvdb::tools::meshToLevelSet scan-converts into a sparse tree, and the
// node's "Band Width" IS that shell's width.
//
// ⛔ OURS TOOK THE SAME PARAMETER AND USED IT FOR THE OPPOSITE THING: `band`
// only padded the bounds (`pad` above) and then every voxel inside those
// padded bounds was computed anyway. Same name, inverted meaning.
//
// Requires a closed mesh, because the outside is found by flooding from the
// grid boundary. `signed_volume` already says whether it is closed, and an
// open one falls back to computing every voxel, which is what this has
// always done.
const uint64_t t_bvh = Time::get_singleton()->get_ticks_usec();
const uint64_t t_start = t_bvh;
uint64_t t_stamp = t_start, t_dilate = t_start, t_query = t_start;
const bool closed = closed_precheck;
const int bandv = MAX(1, (int)band);
std::vector<uint8_t> in_band;
if (closed) {
in_band.assign((size_t)total, 0);
// The marked voxels, kept as a WORKLIST. Dilating from a list touches
// only the band; the separable three-pass version this replaced copied
// the whole 10.2M mask twice and scanned it three times to find the
// ~62k voxels that were actually set.
std::vector<int64_t> marked;
marked.reserve(1 << 16);
// 1. Stamp each triangle's voxel AABB. No expansion here: dilating
// afterwards is O(grid) once instead of O(triangles x band^3), and
// these triangles are far smaller than a voxel -- expanding each one
// by 3 would be 5.16M x 343 stamps.
for (size_t t = 0; t + 2 < idx.size(); t += 3) {
const Vector3 &a = mesh->positions[idx[t]];
const Vector3 &b = mesh->positions[idx[t + 1]];
const Vector3 &c = mesh->positions[idx[t + 2]];
const Vector3 lo = a.min(b).min(c);
const Vector3 hi = a.max(b).max(c);
const int x0 = CLAMP((int)Math::floor((lo.x - mn[0]) / voxel), 0, (int)nx - 1);
const int x1 = CLAMP((int)Math::floor((hi.x - mn[0]) / voxel), 0, (int)nx - 1);
const int y0 = CLAMP((int)Math::floor((lo.y - mn[1]) / voxel), 0, (int)ny - 1);
const int y1 = CLAMP((int)Math::floor((hi.y - mn[1]) / voxel), 0, (int)ny - 1);
const int z0 = CLAMP((int)Math::floor((lo.z - mn[2]) / voxel), 0, (int)nz - 1);
const int z1 = CLAMP((int)Math::floor((hi.z - mn[2]) / voxel), 0, (int)nz - 1);
for (int z = z0; z <= z1; z++) {
for (int y = y0; y <= y1; y++) {
const int64_t base = at(0, y, z);
for (int x = x0; x <= x1; x++) {
if (!in_band[(size_t)(base + x)]) {
in_band[(size_t)(base + x)] = 1;
marked.push_back(base + x);
}
}
}
}
}
t_stamp = Time::get_singleton()->get_ticks_usec();
// 2. Dilate by `band` voxels, separably, over the WORKLIST. Three 1D
// passes rather than one 3D neighbourhood, and each pass walks only
// the voxels currently in the band -- appending as it goes and
// iterating to the size the pass started at, which is exactly
// "dilate the current set along this axis". The mask dedupes, so a
// voxel is visited once however many neighbours reach it.
for (int axis = 0; axis < 3; axis++) {
const size_t cur = marked.size();
for (size_t k = 0; k < cur; k++) {
const int64_t o = marked[k];
const int64_t x = o % nx;
const int64_t y = (o / nx) % ny;
const int64_t z = o / (nx * ny);
for (int d = -bandv; d <= bandv; d++) {
if (d == 0) {
continue;
}
int64_t xx = x, yy = y, zz = z;
if (axis == 0) xx += d; else if (axis == 1) yy += d; else zz += d;
if (xx < 0 || xx >= nx || yy < 0 || yy >= ny || zz < 0 || zz >= nz) {
continue;
}
const int64_t nb = at(xx, yy, zz);
if (!in_band[(size_t)nb]) {
in_band[(size_t)nb] = 1;
marked.push_back(nb);
}
}
}
}
}
t_dilate = Time::get_singleton()->get_ticks_usec();
// ⛔ THREADED, because this is the whole cost of the node. libdc parallelises
// its equivalent loops with `#pragma omp parallel for`, but proc-nodes builds
// without -fopenmp so those pragmas are silently ignored and every libdc loop
// this extension compiles is serial. std::thread rather than adding the flag:
// it changes one loop instead of every libdc loop's behaviour in this
// library, and pulls in no runtime.
//
// Safe to share `sdf`: evaluate() is non-const by declaration only -- it
// calls closest_distance_sq/query_bvh, both const, and the BVH is immutable
// after load_arrays.
{
const int nzi = res[2];
unsigned hw = std::thread::hardware_concurrency();
int nthreads = (int)(hw > 0 ? hw : 1u);
if (nthreads > nzi) nthreads = nzi;
if (nthreads < 1) nthreads = 1;
const bool banded = closed && !in_band.empty();
auto slab = [&](int z0, int z1) {
for (int z = z0; z < z1; z++) {
for (int y = 0; y < res[1]; y++) {
for (int x = 0; x < res[0]; x++) {
const int64_t o = at(x, y, z);
if (banded && !in_band[(size_t)o]) {
continue; // filled by the flood below
}
const float d = sdf.evaluate(mn[0] + x * voxel,
mn[1] + y * voxel, mn[2] + z * voxel, voxel);
// banded => unsigned; the flood signs it. Otherwise
// libdc signed it with pseudonormals as before.
w[o] = banded ? d : (sdf_sign * d);
}
}
}
};
if (nthreads == 1) {
slab(0, nzi);
} else {
std::vector<std::thread> pool;
pool.reserve(nthreads);
const int per = (nzi + nthreads - 1) / nthreads;
for (int t = 0; t < nthreads; t++) {
const int z0 = t * per;
const int z1 = MIN(z0 + per, nzi);
if (z0 >= z1) break;
pool.emplace_back(slab, z0, z1);
}
for (std::thread &th : pool) th.join();
}
}
t_query = Time::get_singleton()->get_ticks_usec();
// ── OUTSIDE BY FLOOD, INSIDE BY WHAT IS LEFT ────────────────────────────
//
// A voxel outside the band cannot contain a sign change, so it needs only
// the right SIGN and a magnitude bigger than anything in the band, so no
// edge between the two reads as a crossing. The bounds were padded by
// `pad` above, so the whole grid boundary is outside the mesh and a flood
// from it reaches exactly the outside.
if (closed && !in_band.empty()) {
// positive = solid after this pass (see the winding note above)
const float far_v = voxel * (float)(bandv + 1);
// THE BARRIER is the shell the surface actually lies in: voxels within
// half a voxel of it. Everything else is unambiguously one side or the
// other at this resolution, so a flood from the grid boundary -- which
// is outside, because the bounds are padded -- labels it exactly.
//
// Only the barrier itself is decided by a neighbour, and a shell one
// voxel thick is precisely the width at which a grid cannot tell the
// two sides apart anyway.
const float half = voxel * 0.5f;
std::vector<uint8_t> outside((size_t)total, 0);
std::vector<int64_t> stack;
stack.reserve(4096);
auto passable = [&](int64_t o) {
return !in_band[(size_t)o] || w[o] >= half;
};
auto push_if_open = [&](int64_t x, int64_t y, int64_t z) {
const int64_t o = at(x, y, z);
if (outside[(size_t)o] || !passable(o)) {
return;
}
outside[(size_t)o] = 1;
stack.push_back(o);
};
for (int64_t z = 0; z < nz; z++) {
for (int64_t y = 0; y < ny; y++) {
for (int64_t x = 0; x < nx; x++) {
if (x == 0 || y == 0 || z == 0 || x == nx - 1 || y == ny - 1 || z == nz - 1) {
push_if_open(x, y, z);
}
}
}
}
while (!stack.empty()) {
const int64_t o = stack.back();
stack.pop_back();
const int64_t x = o % nx;
const int64_t y = (o / nx) % ny;
const int64_t z = o / (nx * ny);
if (x > 0) push_if_open(x - 1, y, z);
if (x < nx - 1) push_if_open(x + 1, y, z);
if (y > 0) push_if_open(x, y - 1, z);
if (y < ny - 1) push_if_open(x, y + 1, z);
if (z > 0) push_if_open(x, y, z - 1);
if (z < nz - 1) push_if_open(x, y, z + 1);
}
for (int64_t z = 0; z < nz; z++) {
for (int64_t y = 0; y < ny; y++) {
for (int64_t x = 0; x < nx; x++) {
const int64_t o = at(x, y, z);
if (!in_band[(size_t)o]) {
w[o] = outside[(size_t)o] ? -far_v : far_v;
continue;
}
const float d = w[o];
if (d >= half) {
// Outside the barrier: the flood knows the side.
w[o] = outside[(size_t)o] ? -d : d;
continue;
}
// In the barrier. Take the side from a labelled neighbour;
// prefer an OUTSIDE one, so a voxel touching air is air --
// the alternative closes the surface over thin gaps.
bool any_out = false, any_in = false;
const int64_t nbr[6][3] = {{x - 1, y, z}, {x + 1, y, z}, {x, y - 1, z},
{x, y + 1, z}, {x, y, z - 1}, {x, y, z + 1}};
for (int k = 0; k < 6; k++) {
const int64_t bx = nbr[k][0], by = nbr[k][1], bz = nbr[k][2];
if (bx < 0 || bx >= nx || by < 0 || by >= ny || bz < 0 || bz >= nz) {
continue;
}
const int64_t no = at(bx, by, bz);
if (!in_band[(size_t)no] || w[no] >= half) {
if (outside[(size_t)no]) any_out = true; else any_in = true;
}
}
w[o] = (any_out && !any_in) ? -d : (any_in && !any_out ? d : -d);
}
}
}
}
{
// WHERE THE TIME GOES, per stage. "MeshToSDFGrid: 12.4 s" says nothing
// about whether that is the band being too wide, the queries being slow,
// or the bookkeeping around them.
const uint64_t t_end = Time::get_singleton()->get_ticks_usec();
int64_t band_n = 0;
for (size_t i = 0; i < in_band.size(); i++) {
band_n += in_band[i] ? 1 : 0;
}
UtilityFunctions::print(" [SDF] prep ", (t_setup - t_node0) / 1000, "ms bvh ", (t_bvh - t_setup) / 1000, "ms grid ", (int64_t)nx, "x", (int64_t)ny, "x", (int64_t)nz,
" = ", total, " voxels, band ", band_n,
band_n > 0 ? String(" (") + String::num((double)band_n * 100.0 / (double)total, 1) + "%)" : String(""),
" stamp ", (t_stamp - t_start) / 1000, "ms dilate ", (t_dilate - t_stamp) / 1000,
"ms query ", (t_query - t_dilate) / 1000, "ms flood ", (t_end - t_query) / 1000, "ms");
}
out_geo->set_volume(vol);
p_ctx.write_geometry(p_node_idx, 0, out_geo);
}
// ── SDF Grid Fillet ─────────────────────────────────────────────────────────
void SDFGridFilletNode::_bind_methods() {}
SDFGridFilletNode::SDFGridFilletNode() {
// ⛔ ITERATIONS, not a radius. Blender's node is `Grid` + `Iterations` (INT,
// default 1); this had `Volume` + `Radius` (float, default 1.0), so the
// bridge could only resolve socket 1 by INDEX, and even then it meant a
// different thing -- a Blender `Iterations = 5` landed as radius 5.0 and
// became 5/voxel passes. The socket now says what it does.
input_defs.push_back({ "Grid", SocketType::Geometry, Variant() }); // 0
input_defs.push_back({ "Iterations", SocketType::Int, 1 }); // 1
output_defs.push_back({ "Grid", SocketType::Geometry, Variant() });
}
void SDFGridFilletNode::execute(ProcEvalContext &p_ctx, int p_node_idx) {
Ref<ProcGeometry> in_geo = p_ctx.read_geometry(p_node_idx, 0);
Ref<ProcGeometry> out_geo;
out_geo.instantiate();
Ref<ProcVolume> in_vol = in_geo.is_valid() ? in_geo->get_volume() : Ref<ProcVolume>();
if (in_vol.is_null() || in_vol->density_grid.size() == 0) {
UtilityFunctions::push_warning(
"proc-nodes: SDF Grid Fillet got no grid. Mesh to SDF Grid has to "
"come first -- this filters a field, it does not make one.");
p_ctx.write_geometry(p_node_idx, 0, out_geo);
return;
}
// ⚠ AN APPROXIMATION, AND SAYING SO IS THE POINT. Blender's fillet is
// OpenVDB's level-set filter, which advances the surface by curvature. This
// is iterated Laplacian smoothing of the field, restricted to CONCAVE
// samples -- where the field is below the average of its neighbours, which
// is where a crease is. That rounds creases and leaves convex edges roughly
// alone, which is the visible half of what a fillet does; it is NOT the same
// operator and will not match Blender value for value.
Ref<ProcVolume> vol = in_vol->duplicate_volume();
const Vector3i r = vol->grid_resolution;
int iters = p_ctx.read_int(p_node_idx, 1);
if (iters <= 0) {
p_ctx.write_geometry(p_node_idx, 0, out_geo);
out_geo->set_volume(vol);
return;
}
// Each pass smooths about one voxel, which is what Blender's iteration
// count buys too -- the operator differs, the unit does not.
if (iters > 32) iters = 32;
PackedFloat32Array next = vol->density_grid;
float *dst = next.ptrw();
const int sx = r.x, sy = r.y, sz = r.z;
for (int it = 0; it < iters; it++) {
const float *src = vol->density_grid.ptr();
for (int z = 1; z < sz - 1; z++) {
for (int y = 1; y < sy - 1; y++) {
for (int x = 1; x < sx - 1; x++) {
const int i = (z * sy + y) * sx + x;
const float c = src[i];
const float avg = (src[i - 1] + src[i + 1] + src[i - sx] + src[i + sx] +
src[i - sx * sy] + src[i + sx * sy]) /
6.0f;
// Concave only: the field dips below its neighbourhood.
dst[i] = (avg > c) ? avg : c;
}
}
}
vol->density_grid = next;
}
out_geo->set_volume(vol);
p_ctx.write_geometry(p_node_idx, 0, out_geo);
}
// ── Grid to Mesh ────────────────────────────────────────────────────────────
void GridToMeshNode::_bind_methods() {}
GridToMeshNode::GridToMeshNode() {
// Blender's node is `Grid` + `Threshold` (0.10) + `Adaptivity`. Adaptivity
// has no socket here on purpose -- Dual Contouring's simplifier is
// dc_simplify, which proc-nodes deliberately does not compile (it needs
// meshoptimizer). A Blender graph carrying one lands on nothing, which is
// the same outcome as ignoring it and is visible in the import log.
input_defs.push_back({ "Grid", SocketType::Geometry, Variant() }); // 0
input_defs.push_back({ "Threshold", SocketType::Float, 0.1f }); // 1
output_defs.push_back({ "Mesh", SocketType::Geometry, Variant() });
}
void GridToMeshNode::execute(ProcEvalContext &p_ctx, int p_node_idx) {
Ref<ProcGeometry> in_geo = p_ctx.read_geometry(p_node_idx, 0);
Ref<ProcGeometry> out_geo;
out_geo.instantiate();
Ref<ProcVolume> vol = in_geo.is_valid() ? in_geo->get_volume() : Ref<ProcVolume>();
if (vol.is_null() || vol->density_grid.size() == 0) {
UtilityFunctions::push_warning(
"proc-nodes: Grid to Mesh got no grid, so there is nothing to "
"extract. Mesh to SDF Grid has to come first.");
p_ctx.write_geometry(p_node_idx, 0, out_geo);
return;
}
// DensityGrid is CUBIC -- one `stride` for all three axes -- so a
// non-cubic ProcVolume grid is resampled onto the largest axis rather than
// silently read with the wrong strides, which is a mesh made of garbage.
const Vector3i r = vol->grid_resolution;
const float threshold = p_ctx.read_float(p_node_idx, 1);
const float voxel = (r.x > 0) ? (vol->grid_size.x / (float)r.x) : 1.0f;
// ⛔ THE CELL COUNT MUST BE A POWER OF TWO. construct_octree_grid_pooled
// halves the span at every level, and dc-terrain only ever hands it 32
// (VOXELS_PER_CELL) so nothing had ever exercised another value. A grid of
// 14 cells does not bottom out: the gate HUNG rather than failing, at 15^3,
// with no error -- which is why this reads as "no mesh" from the outside.
int cells = 1;
const int longest = MAX(MAX(r.x, r.y), r.z);
while (cells < longest - 1 && cells < 128) {
cells <<= 1;
}
const int n = cells + 1; // samples per axis
DensityGrid grid;
const float org[3] = { (float)vol->grid_origin.x, (float)vol->grid_origin.y,
(float)vol->grid_origin.z };
grid.allocate(cells, voxel, org);
const float *src = vol->density_grid.ptr();
for (int z = 0; z < n; z++) {
for (int y = 0; y < n; y++) {
for (int x = 0; x < n; x++) {
float v = -1.0f; // outside, for samples past the source grid
if (x < r.x && y < r.y && z < r.z) {
v = src[(z * r.y + y) * r.x + x] - threshold;
}
grid.data[(z * n + y) * n + x] = v;
}
}
}
grid.compute_normals();
// ⛔ IT NEVER HUNG. `demo/sdf_grid_gate.gd` called `vertex_count()` and
// `triangle_count()`, whose BOUND names are `get_vertex_count` /
// `get_triangle_count`, and `m.positions`, which is a C++ member and not a
// property at all. A script error aborts `_init` before it reaches `quit()`,
// and a headless SceneTree with nothing to quit it runs for ever. Disabling
// this node "fixed" the hang only because the empty geometry it wrote made
// `m` null, so the bad call was short-circuited and never evaluated -- which
// is also why 32^3 "hung identically" and why bypassing the fillet changed
// nothing. Extraction takes 0.375 s end to end for the gate's 15^3 grid.
//
// What WAS broken was the field going in, in two ways that only bit the
// extractor: the sign, and the sliver triangles at the sphere's south pole.
// Both are fixed in Mesh to SDF Grid above, and both are documented there.
DCMesherConfig cfg;
DCMeshData md = dc_mesh_from_grid(grid, cfg);
const int vcount = (int)md.positions.size() / 3;
if (vcount == 0 || md.indices.empty()) {
UtilityFunctions::push_warning(
"proc-nodes: Grid to Mesh extracted no surface -- the field does "
"not cross the threshold anywhere inside the grid.");
p_ctx.write_geometry(p_node_idx, 0, out_geo);
return;
}
Ref<ProcMesh> mesh = out_geo->ensure_mesh();
mesh->positions.resize(vcount);
Vector3 *pw = mesh->positions.ptrw();
for (int i = 0; i < vcount; i++) {
pw[i] = Vector3(md.positions[i * 3], md.positions[i * 3 + 1], md.positions[i * 3 + 2]);
}
if ((int)md.normals.size() == vcount * 3) {
mesh->normals.resize(vcount);
Vector3 *nw = mesh->normals.ptrw();
for (int i = 0; i < vcount; i++) {
nw[i] = Vector3(md.normals[i * 3], md.normals[i * 3 + 1], md.normals[i * 3 + 2]);
}
}
mesh->indices.resize((int)md.indices.size());
int32_t *iw = mesh->indices.ptrw();
for (size_t i = 0; i < md.indices.size(); i++) {
iw[i] = md.indices[i];
}
// Dual Contouring emits triangles, so any ngon grouping is gone.
mesh->face_offsets.clear();
p_ctx.write_geometry(p_node_idx, 0, out_geo);
}
} // namespace godot