| … | … | @@ -65,6 +65,7 @@ def build_shader(hist_scale: float = HIST_SCALE): |
| 65 | 65 | rw = {'NO_RESTRICT', 'READ', 'WRITE'} |
| 66 | 66 | info.image(0, 'R32UI', 'UINT_3D', 'volume_img', qualifiers=rw) |
| 67 | 67 | info.image(1, 'RGBA32F', 'FLOAT_2D', 'state_img', qualifiers=rw) |
| 68 | info.image(2, 'R32UI', 'UINT_3D', 'colour_img', qualifiers=rw) |
| 68 | 69 | info.sampler(0, 'FLOAT_2D', 'alias_tables_tex') |
| 69 | 70 | info.push_constant('VEC3', 'volume_min') |
| 70 | 71 | info.push_constant('VEC3', 'volume_max') |
| … | … | @@ -96,10 +97,20 @@ def _document_buffers(doc): |
| 96 | 97 | |
| 97 | 98 | def accumulate(shader, doc, ubos, alias_tex, size, bounds, dispatches, |
| 98 | 99 | invocation_iters=500): |
| 99 | | """Runs the chaos game into a fresh 4N x N x N grid and reads it back.""" |
| 100 | """ |
| 101 | Runs the chaos game into a fresh density grid plus a HALF-resolution colour |
| 102 | grid, and reads both back. |
| 103 | |
| 104 | Returns `(density, colour)` — density (z, y, x) float32, colour |
| 105 | (z/2, y/2, x/2, 3) float32 and still premultiplied by opacity. Splitting them |
| 106 | costs 5.5 bytes per voxel against the old interleaved 16. |
| 107 | """ |
| 100 | 108 | lo, hi = bounds |
| 101 | | volume = gpu.types.GPUTexture((size * 4, size, size), format='R32UI') |
| 109 | half = max(1, size // 2) |
| 110 | volume = gpu.types.GPUTexture((size, size, size), format='R32UI') |
| 102 | 111 | volume.clear(format='UINT', value=(0, 0, 0, 0)) |
| 112 | colour_tex = gpu.types.GPUTexture((half * 3, half, half), format='R32UI') |
| 113 | colour_tex.clear(format='UINT', value=(0, 0, 0, 0)) |
| 103 | 114 | rows = (ifs_render.INVOCATIONS + ifs_render.STATE_W - 1) // ifs_render.STATE_W |
| 104 | 115 | state = gpu.types.GPUTexture((ifs_render.STATE_W * 2, rows), format='RGBA32F') |
| 105 | 116 | state.clear(format='FLOAT', value=(0.0, 0.0, 0.0, 0.0)) |
| … | … | @@ -111,6 +122,7 @@ def accumulate(shader, doc, ubos, alias_tex, size, bounds, dispatches, |
| 111 | 122 | invocation_iters=invocation_iters)) |
| 112 | 123 | shader.bind() |
| 113 | 124 | shader.image('volume_img', volume) |
| 125 | shader.image('colour_img', colour_tex) |
| 114 | 126 | shader.image('state_img', state) |
| 115 | 127 | shader.uniform_sampler('alias_tables_tex', alias_tex) |
| 116 | 128 | shader.uniform_float('volume_min', lo) |
| … | … | @@ -121,18 +133,17 @@ def accumulate(shader, doc, ubos, alias_tex, size, bounds, dispatches, |
| 121 | 133 | shader.uniform_block(name, ubo) |
| 122 | 134 | gpu.compute.dispatch(shader, ifs_render.WORKGROUPS, 1, 1) |
| 123 | 135 | |
| 136 | # ONE allocation each, not two. `astype(float32) / HIST_SCALE` materialises the |
| 137 | # float copy AND the quotient; np.divide with an explicit dtype does both in a |
| 138 | # single pass, which is most of a 2.4x peak-RSS overshoot. |
| 124 | 139 | raw = np.array(volume.read(), dtype=np.uint32) |
| 125 | | # (z, y, 4x) as read back; split the interleaved channels out. |
| 126 | | # |
| 127 | | # ONE allocation, not two. `raw.astype(np.float32) / HIST_SCALE` materialises |
| 128 | | # the float copy AND the quotient — at N=384 that is two 906 MB arrays where |
| 129 | | # one will do, and measured peak RSS was 2.4x the grid size because of chains |
| 130 | | # like this. np.divide with an explicit dtype converts and scales in a single |
| 131 | | # pass. |
| 132 | | grid = np.divide(raw.reshape(size, size, size, 4), HIST_SCALE, |
| 133 | | dtype=np.float32) |
| 140 | density = np.divide(raw.reshape(size, size, size), HIST_SCALE, dtype=np.float32) |
| 134 | 141 | del raw |
| 135 | | return grid |
| 142 | raw_colour = np.array(colour_tex.read(), dtype=np.uint32) |
| 143 | colour = np.divide(raw_colour.reshape(half, half, half, 3), HIST_SCALE, |
| 144 | dtype=np.float32) |
| 145 | del raw_colour |
| 146 | return density, colour |
| 136 | 147 | |
| 137 | 148 | |
| 138 | 149 | def auto_bounds(shader, doc, ubos, alias_tex, probe_size=64, dispatches=8, |
| … | … | @@ -156,8 +167,8 @@ def auto_bounds(shader, doc, ubos, alias_tex, probe_size=64, dispatches=8, |
| 156 | 167 | extent = start_extent |
| 157 | 168 | while True: |
| 158 | 169 | box = ((-extent, -extent, -extent), (extent, extent, extent)) |
| 159 | | coarse = accumulate(shader, doc, ubos, alias_tex, probe_size, box, dispatches) |
| 160 | | density = coarse[:, :, :, 3] |
| 170 | density, _coarse_colour = accumulate( |
| 171 | shader, doc, ubos, alias_tex, probe_size, box, dispatches) |
| 161 | 172 | hit = np.argwhere(density > 0) |
| 162 | 173 | if hit.size == 0: |
| 163 | 174 | if extent >= max_extent: |
| … | … | @@ -222,9 +233,8 @@ def auto_bounds(shader, doc, ubos, alias_tex, probe_size=64, dispatches=8, |
| 222 | 233 | probe_lo = tuple(float(lo_w[a] - pad[a]) for a in range(3)) |
| 223 | 234 | probe_hi = tuple(float(hi_w[a] + pad[a]) for a in range(3)) |
| 224 | 235 | |
| 225 | | finer = accumulate(shader, doc, ubos, alias_tex, probe_size, |
| 226 | | (probe_lo, probe_hi), dispatches) |
| 227 | | finer_density = finer[:, :, :, 3] |
| 236 | finer_density, _finer_colour = accumulate( |
| 237 | shader, doc, ubos, alias_tex, probe_size, (probe_lo, probe_hi), dispatches) |
| 228 | 238 | if not (finer_density > 0).any(): |
| 229 | 239 | break # the tighter box found nothing; keep what we had |
| 230 | 240 | |
| … | … | @@ -286,8 +296,8 @@ def _fit(density, box, probe_size, percentile, margin): |
| 286 | 296 | return lo_w - span * margin, hi_w + span * margin |
| 287 | 297 | |
| 288 | 298 | |
| 289 | | def write_vdb(grid: np.ndarray, bounds, path: str, name: str = "density", |
| 290 | | gamma: float = 1.0): |
| 299 | def write_vdb(density_grid: np.ndarray, colour_grid: np.ndarray, bounds, |
| 300 | path: str, name: str = "density", gamma: float = 1.0): |
| 291 | 301 | """ |
| 292 | 302 | Writes density and colour grids to a `.vdb`. |
| 293 | 303 | |
| … | … | @@ -296,7 +306,7 @@ def write_vdb(grid: np.ndarray, bounds, path: str, name: str = "density", |
| 296 | 306 | """ |
| 297 | 307 | import openvdb |
| 298 | 308 | |
| 299 | | size = grid.shape[0] |
| 309 | size = density_grid.shape[0] |
| 300 | 310 | lo, hi = bounds |
| 301 | 311 | |
| 302 | 312 | # Per-axis scale AND a translation to `lo`. Both matter and both were missing: |
| … | … | @@ -323,9 +333,9 @@ def write_vdb(grid: np.ndarray, bounds, path: str, name: str = "density", |
| 323 | 333 | # against the same convention it was written with, so it agreed with itself. |
| 324 | 334 | # test_bake_orientation measures the baked centroid against an analytically |
| 325 | 335 | # known one instead, where it showed up as a clean -0.5 on all three axes. |
| 326 | | def _transform(): |
| 336 | def _transform(scale=1.0): |
| 327 | 337 | """A fresh Transform per grid; sharing one between grids is asking for it.""" |
| 328 | | sx, sy, sz = ((hi[a] - lo[a]) / size for a in range(3)) |
| 338 | sx, sy, sz = ((hi[a] - lo[a]) / size * scale for a in range(3)) |
| 329 | 339 | return openvdb.createLinearTransform(matrix=[ |
| 330 | 340 | [sx, 0.0, 0.0, 0.0], |
| 331 | 341 | [0.0, sy, 0.0, 0.0], |
| … | … | @@ -357,7 +367,7 @@ def write_vdb(grid: np.ndarray, bounds, path: str, name: str = "density", |
| 357 | 367 | # Applied to density ONLY, and only after colour has been divided by the ORIGINAL |
| 358 | 368 | # density below - colour is accumulated weighted by density, so dividing it by a |
| 359 | 369 | # gamma-lifted density collapses it to black. That mistake renders an empty volume. |
| 360 | | density = grid[:, :, :, 3] |
| 370 | density = density_grid |
| 361 | 371 | peak = float(density.max()) |
| 362 | 372 | if peak <= 0: |
| 363 | 373 | raise RuntimeError("nothing was accumulated - the bounds may exclude the fractal") |
| … | … | @@ -377,26 +387,29 @@ def write_vdb(grid: np.ndarray, bounds, path: str, name: str = "density", |
| 377 | 387 | grids.append(dgrid) |
| 378 | 388 | |
| 379 | 389 | # Colour is accumulated weighted by density, so divide it back out. |
| 390 | # Colour is HALF resolution and still premultiplied by opacity, so dividing it |
| 391 | # back out needs the density of the SAME colour voxel — the sum of its 2x2x2 |
| 392 | # density block, not a single voxel. Getting that wrong tints the fractal in a |
| 393 | # way that looks like a palette bug rather than an indexing one. |
| 394 | half = colour_grid.shape[0] |
| 395 | usable = half * 2 |
| 396 | block = density_grid[:usable, :usable, :usable].reshape( |
| 397 | half, 2, half, 2, half, 2).sum(axis=(1, 3, 5)) |
| 398 | safe = np.maximum(block[:, :, :, None], 1e-9) |
| 399 | |
| 380 | 400 | # Same transpose as density - the spatial axes only, leaving the component axis |
| 381 | 401 | # last. A colour grid transposed differently from density would tint the wrong |
| 382 | 402 | # parts of the fractal, which is harder to spot than a whole-volume rotation. |
| 383 | | # |
| 384 | | # Built per SLAB. The whole-array form (divide -> clip -> transpose -> astype) |
| 385 | | # peaked at four simultaneous copies of a 3-channel volume; slabs bound that to |
| 386 | | # one plane at a time while producing byte-identical output. |
| 387 | | colour = np.empty((size, size, size, 3), dtype=np.float32) |
| 388 | | for z in range(size): |
| 389 | | plane = grid[z] # (y, x, 4) |
| 390 | | safe = np.maximum(plane[:, :, 3:4], 1e-9) |
| 391 | | tinted = np.clip(plane[:, :, :3] / safe, 0.0, 1.0) |
| 392 | | # grid is [z][y][x]; the target is [x][y][z], so this plane becomes a |
| 393 | | # column across the first axis. |
| 394 | | colour[:, :, z, :] = tinted.transpose(1, 0, 2) |
| 403 | colour = np.ascontiguousarray( |
| 404 | np.clip(colour_grid / safe, 0.0, 1.0).transpose(2, 1, 0, 3)) |
| 405 | |
| 395 | 406 | try: |
| 396 | 407 | cgrid = openvdb.Vec3SGrid() |
| 397 | 408 | cgrid.copyFromArray(colour) |
| 398 | 409 | cgrid.name = "color" |
| 399 | | cgrid.transform = _transform() |
| 410 | # Twice the voxel size, because it is half the resolution. VDB grids carry |
| 411 | # independent transforms and Blender honours them, so no upsampling. |
| 412 | cgrid.transform = _transform(scale=2.0) |
| 400 | 413 | grids.append(cgrid) |
| 401 | 414 | except Exception: |
| 402 | 415 | # A density-only volume is still useful; do not lose the bake over colour. |
| … | … | @@ -438,10 +451,10 @@ def bake_document(doc, out_path, size=256, dispatches=200, bounds=None, |
| 438 | 451 | lo, hi, fit = auto_bounds(shader, doc, ubos, alias_tex) |
| 439 | 452 | bounds = (lo, hi) |
| 440 | 453 | |
| 441 | | grid = accumulate(shader, doc, ubos, alias_tex, size, bounds, dispatches, |
| 442 | | invocation_iters) |
| 443 | | path, peak = write_vdb(grid, bounds, out_path, gamma=gamma) |
| 444 | | occupied = int((grid[:, :, :, 3] > 0).sum()) |
| 454 | grid, colour_grid = accumulate(shader, doc, ubos, alias_tex, size, bounds, |
| 455 | dispatches, invocation_iters) |
| 456 | path, peak = write_vdb(grid, colour_grid, bounds, out_path, gamma=gamma) |
| 457 | occupied = int((grid > 0).sum()) |
| 445 | 458 | |
| 446 | 459 | # How bright to drive emission, DERIVED rather than guessed. A fixed number cannot |
| 447 | 460 | # work: measured across two shipped templates the required scale differs by ~200x |
| … | … | @@ -458,7 +471,7 @@ def bake_document(doc, out_path, size=256, dispatches=200, bounds=None, |
| 458 | 471 | # probe voxels thick - the fixed axis is the thin direction, so the statistic |
| 459 | 472 | # reports the shortest possible path and derives an exposure far too dark. The |
| 460 | 473 | # camera can look from anywhere, so the metric has to be view-independent. |
| 461 | | normalised = grid[:, :, :, 3] / peak if peak > 0 else grid[:, :, :, 3] |
| 474 | normalised = grid / peak if peak > 0 else grid |
| 462 | 475 | per_axis = [] |
| 463 | 476 | for axis in range(3): |
| 464 | 477 | columns = normalised.sum(axis=axis) |