"""
Renders an IFS fractal with Blender's GPU module. Run inside Blender.
blender -b --factory-startup -noaudio --python blender/ifs_render.py -- --out /tmp/f.png
Pipeline, all on the GPU except the final tonemap:
dispatch the chaos game -> R32UI fixed-point histogram
drain every N dispatches -> RGBA32F accumulator, uint image zeroed
tonemap -> PNG
The drain exists because `imageAtomicAdd` has no float overload, so the histogram
accumulates as scaled integers and uint32 saturates at 4.29e9. Draining folds it
into a float accumulator and zeroes it, which is what lets a render run indefinitely.
The camera convention is the engine's and is NOT the usual one: the view matrix is
built looking along `Position - Forward` (Camera.cs:104), so a visible point has
NEGATIVE w after projection - `project()` discards `p_hom.w >= 0.0`. A matrix built
the conventional way renders an empty image, which looks exactly like a broken
pipeline, so `test_camera.py` asserts the demo's points project inside the frame
with negative w before anything is dispatched.
"""
from __future__ import annotations
import math
import os
import sys
# Importable outside Blender for the camera tests; the GPU work is guarded.
try:
import gpu
import numpy as np
except ImportError: # pragma: no cover - only when imported for the maths
gpu = None
np = None
_HERE = os.path.dirname(os.path.abspath(__file__))
if _HERE not in sys.path:
sys.path.insert(0, _HERE)
from ifs_codegen import apply_portability_fixes, generate_kernel, load_library
from ifs_host import (
build_alias_tables,
flatten_real_params,
flatten_vec3_params,
pack_frame,
pack_iterators,
pack_palette,
pack_real_params,
pack_settings,
pack_vec3_params,
gradient,
sierpinski,
)
import ifs_document
from ifs_layout import DEFAULT_HIST_SCALE, rewrite
REPO = os.path.dirname(_HERE)
STATE_W = 256
WORKGROUPS = 256
LOCAL_SIZE = 64
INVOCATIONS = WORKGROUPS * LOCAL_SIZE
# ------------------------------------------------------------------ camera maths
def _normalize(v: tuple[float, float, float]) -> tuple[float, float, float]:
length = math.sqrt(sum(c * c for c in v))
if length == 0:
raise ValueError("cannot normalize a zero vector")
return (v[0] / length, v[1] / length, v[2] / length)
def _cross(a, b):
return (a[1] * b[2] - a[2] * b[1],
a[2] * b[0] - a[0] * b[2],
a[0] * b[1] - a[1] * b[0])
def _dot(a, b):
return a[0] * b[0] + a[1] * b[1] + a[2] * b[2]
def look_at(position, target, up) -> list[list[float]]:
"""`Matrix4x4.CreateLookAt` — right-handed, row-vector, as .NET defines it."""
zaxis = _normalize((position[0] - target[0], position[1] - target[1],
position[2] - target[2]))
xaxis = _normalize(_cross(up, zaxis))
yaxis = _cross(zaxis, xaxis)
return [
[xaxis[0], yaxis[0], zaxis[0], 0.0],
[xaxis[1], yaxis[1], zaxis[1], 0.0],
[xaxis[2], yaxis[2], zaxis[2], 0.0],
[-_dot(xaxis, position), -_dot(yaxis, position), -_dot(zaxis, position), 1.0],
]
def perspective(fov_degrees: float, aspect: float, near: float,
engine_fov_quirk: bool = True,
flip_x: bool = False) -> list[list[float]]:
"""
`Matrix4x4.CreatePerspectiveFieldOfView` with far = +infinity.
At far = inf the third column degenerates to M33 = -1 and M43 = -near.
`engine_fov_quirk` reproduces Camera.cs:105, which passes
`1 + FieldOfView % 179` rather than the angle itself. Documents are authored
against that behaviour so it is the default and must stay - but it mangles a
real focal length (a 50mm lens is ~39.6 degrees, which the quirk turns into
40.6), so anything driving the camera from OUTSIDE a document, like the
Blender viewport, passes False and gets the angle it asked for.
"""
fov = math.radians(1.0 + (fov_degrees % 179.0)) if engine_fov_quirk \
else math.radians(fov_degrees)
y_scale = 1.0 / math.tan(fov * 0.5)
x_scale = y_scale / aspect
if flip_x:
# Mirrors the HORIZONTAL axis only. IFSEngine's view basis is a horizontal
# mirror of Blender's - its x_ndc = k*x_view/(-z_view) sends camera-right
# to screen-left - so an overlay tracking Blender's camera has to undo
# exactly that, and only that. Negating `up` instead flips both axes and
# puts world +Z below the horizon.
x_scale = -x_scale
return [
[x_scale, 0.0, 0.0, 0.0],
[0.0, y_scale, 0.0, 0.0],
[0.0, 0.0, -1.0, -1.0],
[0.0, 0.0, -near, 0.0],
]
def matmul(a: list[list[float]], b: list[list[float]]) -> list[list[float]]:
return [[sum(a[r][k] * b[k][c] for k in range(4)) for c in range(4)]
for r in range(4)]
def view_projection(position, forward, up, fov_degrees=60.0, near=0.001,
engine_fov_quirk=True, flip_x=False):
"""
The engine's view-projection, flattened in .NET's row-major memory order.
GLSL reads those bytes column-major, i.e. as the transpose, and computes
`M * v` — which equals .NET's row-vector `v * M`. The two conventions cancel,
so the bytes go across unchanged.
"""
target = (position[0] - forward[0], position[1] - forward[1], position[2] - forward[2])
view = look_at(position, target, (-up[0], -up[1], -up[2]))
combined = matmul(view, perspective(fov_degrees, 1.0, near, engine_fov_quirk,
flip_x))
return tuple(value for row in combined for value in row)
def project_point(view_proj, point, width, height):
"""
Mirrors the kernel's `project()` for projection_type 0, for testing.
Returns (x, y, w) in pixels, or None where the kernel would discard.
"""
m = view_proj
p = (point[0], point[1], point[2], 1.0)
# GLSL's M * v against row-major bytes == v * M.
hom = [sum(p[k] * m[k * 4 + c] for k in range(4)) for c in range(4)]
w = hom[3]
if w >= 0.0 or w == 0.0:
return None
ndc = [hom[i] / w for i in range(3)]
ratio = width / float(height)
x = (ndc[0] + 1.0) * 0.5 * width
y = (ndc[1] * ratio + 1.0) * 0.5 * height
if not (0.0 <= x < width and 0.0 <= y < height):
return None
return (x - 0.5, y - 0.5, w)
# ------------------------------------------------------------------ drain shader
_DRAIN = """
void main()
{
ivec2 c = ivec2(gl_GlobalInvocationID.xy);
if (c.x >= out_size.x || c.y >= out_size.y) return;
ivec2 b = ivec2(c.x * 4, c.y);
vec4 v = vec4(
float(imageLoad(histogram_img, b + ivec2(0, 0)).x),
float(imageLoad(histogram_img, b + ivec2(1, 0)).x),
float(imageLoad(histogram_img, b + ivec2(2, 0)).x),
float(imageLoad(histogram_img, b + ivec2(3, 0)).x));
imageStore(accum_img, c, imageLoad(accum_img, c) + v / HIST_SCALE);
imageStore(histogram_img, b + ivec2(0, 0), uvec4(0));
imageStore(histogram_img, b + ivec2(1, 0), uvec4(0));
imageStore(histogram_img, b + ivec2(2, 0), uvec4(0));
imageStore(histogram_img, b + ivec2(3, 0), uvec4(0));
}
"""
def _build_shaders(width, height, hist_scale, depth=False, equirect=False):
"""
`depth` adds the nearest-sample depth image the live overlay needs to be occluded
by scene geometry. It is off for the still path, which has no scene to be behind
— and a declared-but-unused image would still consume one of Blender's eight
slots, so this cannot be left permanently on.
`equirect` swaps the perspective projection for a direction-to-lat/long one, so
the same histogram becomes an environment map. It touches nothing downstream:
the drain, the auto-exposure and the tonemap are shared, which is what makes an
HDRI bake go through the SAME curve as the viewport rather than being a second
renderer that merely resembles it.
"""
library = load_library(os.path.join(REPO, "Transforms"),
os.path.join(REPO, "WpfDisplay", "Library", "Includes"))
template_path = os.path.join(REPO, "IFSEngine", "Rendering", "Shaders",
"ifs_kernel.comp.shader")
with open(template_path, "r", encoding="utf-8-sig") as handle:
kernel = generate_kernel(handle.read(), library)
kernel, _ = apply_portability_fixes(kernel)
body, typedefs, report = rewrite(kernel, hist_scale, STATE_W, depth=depth,
equirect=equirect)
if not report.ok:
raise RuntimeError("kernel rewrite incomplete: " + "; ".join(report.missing))
info = gpu.types.GPUShaderCreateInfo()
info.typedef_source(typedefs)
info.local_group_size(LOCAL_SIZE)
rw = {'NO_RESTRICT', 'READ', 'WRITE'}
info.image(0, 'R32UI', 'UINT_2D', 'histogram_img', qualifiers=rw)
info.image(1, 'RGBA32F', 'FLOAT_2D', 'state_img', qualifiers=rw)
if depth:
info.image(2, 'R32UI', 'UINT_2D', 'depth_img', qualifiers=rw)
# One more push constant on a budget already over the 128-byte portable
# limit (see the add-on's known gaps). Fine on OpenGL; this moves into a UBO
# with the rest when Vulkan or Metal is attempted.
info.push_constant('MAT4', 'blender_view_proj')
info.push_constant('INT', 'depth_slices')
info.image(3, 'R32UI', 'UINT_2D', 'depth_range_img',
qualifiers={'NO_RESTRICT', 'READ'})
if equirect:
info.push_constant('VEC3', 'equirect_centre')
info.sampler(0, 'FLOAT_2D', 'alias_tables_tex')
info.uniform_buf(0, 'FrameParams', 'frame')
info.uniform_buf(1, 'SettingsParams', 'settings')
info.uniform_buf(2, 'IteratorsParams', 'iterators_buf')
info.uniform_buf(3, 'PaletteParams', 'palette_buf')
info.uniform_buf(4, 'RealParams', 'real_params_buf')
info.uniform_buf(5, 'Vec3Params', 'vec3_params_buf')
info.compute_source(body)
kernel_shader = gpu.shader.create_from_info(info)
dinfo = gpu.types.GPUShaderCreateInfo()
dinfo.local_group_size(8, 8)
dinfo.image(0, 'R32UI', 'UINT_2D', 'histogram_img', qualifiers=rw)
dinfo.image(1, 'RGBA32F', 'FLOAT_2D', 'accum_img', qualifiers=rw)
dinfo.push_constant('IVEC2', 'out_size')
dinfo.compute_source(f"#define HIST_SCALE {hist_scale}\n" + _DRAIN)
drain_shader = gpu.shader.create_from_info(dinfo)
return kernel_shader, drain_shader, library
def render(width=512, height=512, dispatches=400, invocation_iters=500,
drain_every=32, hist_scale=DEFAULT_HIST_SCALE, brightness=8.0,
out_path=None, input_path=None, image_name="IFS Render"):
gpu.init()
kernel_shader, drain_shader, library = _build_shaders(width, height, hist_scale)
if input_path:
doc = ifs_document.load(input_path, library)
iterators, camera, settings = doc.iterators, doc.camera, doc.settings
palette = doc.palette(settings.gradient_resolution)
brightness = doc.brightness
up = getattr(camera, "up", (0.0, 1.0, 0.0))
camera.view_proj = view_projection(
camera.position[:3], camera.forward[:3], up,
fov_degrees=getattr(camera, "field_of_view", 60.0))
print(f"loaded {doc.title!r}: {len(iterators)} iterators, "
f"native {doc.width}x{doc.height} level {doc.target_level}")
else:
iterators, camera, settings = sierpinski()
affine_id = library.tf_ids["Affine"]
for it in iterators:
it.tf_id = affine_id
camera.position = (0.0, 0.0, -2.6, 1.0)
camera.forward = (0.0, 0.0, 1.0, 0.0)
camera.view_proj = view_projection((0.0, 0.0, -2.6), (0.0, 0.0, 1.0),
(0.0, 1.0, 0.0))
camera.focus_point = (0.0, 0.0, 0.0, 0.0)
palette = gradient([(1.0, 0.30, 0.10, 1.0),
(1.0, 0.92, 0.35, 1.0),
(0.20, 0.65, 1.0, 1.0)], settings.gradient_resolution)
ubos = {
'settings': gpu.types.GPUUniformBuf(pack_settings(camera, settings)),
'iterators_buf': gpu.types.GPUUniformBuf(pack_iterators(iterators)),
'palette_buf': gpu.types.GPUUniformBuf(pack_palette(palette)),
'real_params_buf': gpu.types.GPUUniformBuf(
pack_real_params(flatten_real_params(iterators))),
'vec3_params_buf': gpu.types.GPUUniformBuf(
pack_vec3_params(flatten_vec3_params(iterators))),
}
alias = build_alias_tables(iterators)
n = len(iterators)
alias_tex = gpu.types.GPUTexture(
(n, n), format='RGBA32F',
data=gpu.types.Buffer('FLOAT', len(alias), alias))
histogram = gpu.types.GPUTexture((width * 4, height), format='R32UI')
histogram.clear(format='UINT', value=(0, 0, 0, 0))
accum = gpu.types.GPUTexture((width, height), format='RGBA32F')
accum.clear(format='FLOAT', value=(0.0, 0.0, 0.0, 0.0))
state_rows = (INVOCATIONS + STATE_W - 1) // STATE_W
state = gpu.types.GPUTexture((STATE_W * 2, state_rows), format='RGBA32F')
state.clear(format='FLOAT', value=(0.0, 0.0, 0.0, 0.0))
def drain():
drain_shader.bind()
drain_shader.image('histogram_img', histogram)
drain_shader.image('accum_img', accum)
drain_shader.uniform_int('out_size', (width, height))
gpu.compute.dispatch(drain_shader, (width + 7) // 8, (height + 7) // 8, 1)
for index in range(dispatches):
frame = gpu.types.GPUUniformBuf(pack_frame(
seed=0x9E3779B9, width=width, height=height, dispatch_cnt=index,
reset_points_state=1 if index == 0 else 0,
invocation_iters=invocation_iters))
kernel_shader.bind()
kernel_shader.image('histogram_img', histogram)
kernel_shader.image('state_img', state)
kernel_shader.uniform_sampler('alias_tables_tex', alias_tex)
kernel_shader.uniform_block('frame', frame)
for name, ubo in ubos.items():
kernel_shader.uniform_block(name, ubo)
gpu.compute.dispatch(kernel_shader, WORKGROUPS, 1, 1)
if (index + 1) % drain_every == 0:
drain()
drain()
pixels = np.array(accum.read(), dtype=np.float32).reshape(height, width, 4)
total = float(pixels[:, :, 3].sum())
peak = float(pixels[:, :, 3].max())
print(f"density: total={total:.4g} peak={peak:.4g} "
f"nonzero={(pixels[:, :, 3] > 0).sum()} of {width * height}")
if peak <= 0:
raise RuntimeError("nothing was plotted - the histogram is empty")
# The engine's tonemap: log density, colour averaged by density.
density = pixels[:, :, 3:4]
safe = np.maximum(density, 1e-12)
scale = np.log10(1.0 + brightness * density / peak) / safe
rgb = np.clip(pixels[:, :, :3] * scale, 0.0, 1.0)
out = np.concatenate([rgb, np.ones_like(density)], axis=2)
_to_image(out, width, height, image_name)
if out_path:
_write_png(out, width, height, out_path)
print(f"wrote {out_path}")
if os.environ.get("IFS_DUMP_DENSITY"):
# Shape without colour: isolates "are the points in the right place" from
# "is the palette being read", which fail identically as a black image.
d = density / peak
d = np.clip(np.log10(1.0 + 500.0 * d), 0.0, 1.0)
_write_png(np.concatenate([d, d, d, np.ones_like(d)], axis=2),
width, height, os.environ["IFS_DUMP_DENSITY"])
return image_name
def _to_image(rgba, width, height, name):
"""
Writes the result into a Blender Image datablock so the Image Editor shows it.
Reused rather than recreated when the size matches, so an editor already
displaying it keeps showing it across re-renders instead of going blank.
"""
import bpy
image = bpy.data.images.get(name)
if image is not None and tuple(image.size) != (width, height):
bpy.data.images.remove(image)
image = None
if image is None:
image = bpy.data.images.new(name, width, height, alpha=True)
image.pixels.foreach_set(rgba.reshape(-1))
image.update()
return image
def _write_png(rgba, width, height, path):
import bpy
image = _to_image(rgba, width, height, "ifs_render_out")
image.filepath_raw = path
image.file_format = 'PNG'
image.save()
def main():
argv = sys.argv[sys.argv.index("--") + 1:] if "--" in sys.argv else []
import argparse
parser = argparse.ArgumentParser()
parser.add_argument("--out", default="/tmp/ifs_blender.png",
help="PNG to write; the result is always left in an Image datablock")
parser.add_argument("--size", type=int, default=512)
parser.add_argument("--dispatches", type=int, default=400)
parser.add_argument("--input", default=None, help="an .ifsjson document")
args = parser.parse_args(argv)
render(width=args.size, height=args.size, dispatches=args.dispatches,
out_path=args.out, input_path=args.input)
if __name__ == "__main__":
main()