skills/shader-dev/techniques/csg-boolean-operations.md

WebGL2 Adaptation Requirements

The code templates in this document use ShaderToy GLSL style. When generating standalone HTML pages, you must adapt for WebGL2:

• Use canvas.getContext("webgl2")
• First line of shaders: #version 300 es, add precision highp float; in fragment shaders
• Vertex shader: attribute -> in, varying -> out
• Fragment shader: varying -> in, gl_FragColor -> custom output variable (must be declared before void main(), e.g. out vec4 outColor;), texture2D() -> texture()
• ShaderToy's void mainImage(out vec4 fragColor, in vec2 fragCoord) must be adapted to the standard void main() entry point

CSG Boolean Operations

Core Principles

CSG boolean operations are per-point value operations on two distance fields:

Operation	Expression	Meaning
Union	min(d1, d2)	Take nearest surface, keeping both shapes
Intersection	max(d1, d2)	Take farthest surface, keeping only the overlap
Subtraction	max(d1, -d2)	Cut d1 using the interior of d2

Smooth booleans (smooth min/max) introduce a blending band in the transition region. The parameter k controls the blend band width (larger = rounder, k=0 degenerates to hard boolean). Multiple variants exist with different mathematical properties.

Implementation Steps

Step 1: Hard Boolean Operations

float opUnion(float d1, float d2) { return min(d1, d2); }
float opIntersection(float d1, float d2) { return max(d1, d2); }
float opSubtraction(float d1, float d2) { return max(d1, -d2); }

Step 2: Smooth Union (Polynomial Version)

// k: blend radius, typical values 0.05~0.5
float opSmoothUnion(float d1, float d2, float k) {
    float h = clamp(0.5 + 0.5 * (d2 - d1) / k, 0.0, 1.0);
    return mix(d2, d1, h) - k * h * (1.0 - h);
}

Step 3: Smooth Subtraction and Intersection (Polynomial Version)

float opSmoothSubtraction(float d1, float d2, float k) {
    float h = clamp(0.5 - 0.5 * (d2 + d1) / k, 0.0, 1.0);
    return mix(d2, -d1, h) + k * h * (1.0 - h);
}

float opSmoothIntersection(float d1, float d2, float k) {
    float h = clamp(0.5 - 0.5 * (d2 - d1) / k, 0.0, 1.0);
    return mix(d2, d1, h) + k * h * (1.0 - h);
}

Step 4: Quadratic Optimized Version (Recommended as Default)

float smin(float a, float b, float k) {
    float h = max(k - abs(a - b), 0.0);
    return min(a, b) - h * h * 0.25 / k;
}

float smax(float a, float b, float k) {
    float h = max(k - abs(a - b), 0.0);
    return max(a, b) + h * h * 0.25 / k;
}

// Subtraction via smax
float sSub(float d1, float d2, float k) {
    return smax(d1, -d2, k);
}

Step 4b: Smooth Minimum Variant Library

Different smin implementations have different mathematical properties. Choose based on your needs:

Variant	Rigid	Associative	Best For
Quadratic (default above)	Yes	No	General use, fastest
Cubic	Yes	No	Smoother C2 transitions
Quartic	Yes	No	Highest quality blending
Exponential	No	Yes	Multi-body blending (order-independent)
Circular Geometric	Yes	Yes	Strict local blending

Rigid: preserves original SDF shape outside the blend region (no under-estimation). Associative: smin(a, smin(b, c)) == smin(smin(a, b), c) — important when blending many objects where evaluation order varies.

// --- Cubic Polynomial smin (C2 continuous, smoother transitions) ---
float sminCubic(float a, float b, float k) {
    k *= 6.0;
    float h = max(k - abs(a - b), 0.0) / k;
    return min(a, b) - h * h * h * k * (1.0 / 6.0);
}

// --- Quartic Polynomial smin (C3 continuous, highest quality) ---
float sminQuartic(float a, float b, float k) {
    k *= 16.0 / 3.0;
    float h = max(k - abs(a - b), 0.0) / k;
    return min(a, b) - h * h * h * (4.0 - h) * k * (1.0 / 16.0);
}

// --- Exponential smin (associative — order independent for multi-body blending) ---
float sminExp(float a, float b, float k) {
    float r = exp2(-a / k) + exp2(-b / k);
    return -k * log2(r);
}

// --- Circular Geometric smin (rigid + local + associative) ---
float sminCircle(float a, float b, float k) {
    k *= 1.0 / (1.0 - sqrt(0.5));
    return max(k, min(a, b)) - length(max(k - vec2(a, b), 0.0));
}

// --- Gradient-aware smin (carries material/color through blending) ---
// x = distance, yzw = material properties or color components
vec4 sminColor(vec4 a, vec4 b, float k) {
    k *= 4.0;
    float h = max(k - abs(a.x - b.x), 0.0) / (2.0 * k);
    return vec4(
        min(a.x, b.x) - h * h * k,
        mix(a.yzw, b.yzw, (a.x < b.x) ? h : 1.0 - h)
    );
}

// --- Smooth maximum from any smin variant ---
// smax(a, b, k) = -smin(-a, -b, k)
// Smooth subtraction: smax(d1, -d2, k)
// Smooth intersection: smax(d1, d2, k)

Step 5: Basic SDF Primitives

float sdSphere(vec3 p, float r) {
    return length(p) - r;
}

float sdBox(vec3 p, vec3 b) {
    vec3 d = abs(p) - b;
    return length(max(d, 0.0)) + min(max(d.x, max(d.y, d.z)), 0.0);
}

float sdCylinder(vec3 p, float h, float r) {
    vec2 d = abs(vec2(length(p.xz), p.y)) - vec2(r, h);
    return min(max(d.x, d.y), 0.0) + length(max(d, 0.0));
}

Step 6: CSG Composition for Scene Building

float mapScene(vec3 p) {
    float cube = sdBox(p, vec3(1.0));
    float sphere = sdSphere(p, 1.2);
    float cylX = sdCylinder(p.yzx, 2.0, 0.4);
    float cylY = sdCylinder(p.xyz, 2.0, 0.4);
    float cylZ = sdCylinder(p.zxy, 2.0, 0.4);

    // (cube intersect sphere) - three cylinders = nut
    float shape = opIntersection(cube, sphere);
    float holes = opUnion(cylX, opUnion(cylY, cylZ));
    return opSubtraction(shape, holes);
}

Step 7: Smooth CSG Modeling for Organic Forms

// Use different k values for different body parts: large k for major joints, small k for fine details
float mapCreature(vec3 p) {
    float body = sdSphere(p, 0.5);
    float head = sdSphere(p - vec3(0.0, 0.6, 0.3), 0.25);
    float d = smin(body, head, 0.15);          // large blend

    float leg = sdCylinder(p - vec3(0.2, -0.5, 0.0), 0.3, 0.08);
    d = smin(d, leg, 0.08);                    // medium blend

    float eye = sdSphere(p - vec3(0.05, 0.75, 0.4), 0.05);
    d = smax(d, -eye, 0.02);                  // small blend for subtraction
    return d;
}

Step 8: Ray Marching Main Loop

float rayMarch(vec3 ro, vec3 rd, float maxDist) {
    float t = 0.0;
    for (int i = 0; i < MAX_STEPS; i++) {
        vec3 p = ro + rd * t;
        float d = mapScene(p);
        if (d < SURF_DIST) return t;
        t += d;
        if (t > maxDist) break;
    }
    return -1.0;
}

Step 9: Normal Calculation (Tetrahedral Sampling, 4 Samples More Efficient Than 6 with Central Differences)

vec3 calcNormal(vec3 pos) {
    vec2 e = vec2(0.001, -0.001);
    return normalize(
        e.xyy * mapScene(pos + e.xyy) +
        e.yyx * mapScene(pos + e.yyx) +
        e.yxy * mapScene(pos + e.yxy) +
        e.xxx * mapScene(pos + e.xxx)
    );
}

Full Code Template

// === CSG Boolean Operations - WebGL2 Full Template ===
// Note: When generating HTML with this template, pass iTime, iResolution, etc. via uniforms

#define MAX_STEPS 128
#define MAX_DIST 50.0
#define SURF_DIST 0.001
#define SMOOTH_K 0.1

// === Hard Boolean Operations ===
float opUnion(float d1, float d2) { return min(d1, d2); }
float opIntersection(float d1, float d2) { return max(d1, d2); }
float opSubtraction(float d1, float d2) { return max(d1, -d2); }

// === Smooth Boolean Operations (Quadratic Optimized) ===
float smin(float a, float b, float k) {
    float h = max(k - abs(a - b), 0.0);
    return min(a, b) - h * h * 0.25 / k;
}

float smax(float a, float b, float k) {
    float h = max(k - abs(a - b), 0.0);
    return max(a, b) + h * h * 0.25 / k;
}

// === SDF Primitives ===
float sdSphere(vec3 p, float r) {
    return length(p) - r;
}

float sdBox(vec3 p, vec3 b) {
    vec3 d = abs(p) - b;
    return length(max(d, 0.0)) + min(max(d.x, max(d.y, d.z)), 0.0);
}

float sdCylinder(vec3 p, float h, float r) {
    vec2 d = abs(vec2(length(p.xz), p.y)) - vec2(r, h);
    return min(max(d.x, d.y), 0.0) + length(max(d, 0.0));
}

float sdEllipsoid(vec3 p, vec3 r) {
    float k0 = length(p / r);
    float k1 = length(p / (r * r));
    return k0 * (k0 - 1.0) / k1;
}

float sdCapsule(vec3 p, vec3 a, vec3 b, float r) {
    vec3 pa = p - a, ba = b - a;
    float h = clamp(dot(pa, ba) / dot(ba, ba), 0.0, 1.0);
    return length(pa - ba * h) - r;
}

// === Scene Definition ===
float mapScene(vec3 p) {
    // Rotation animation
    float angle = iTime * 0.3;
    float c = cos(angle), s = sin(angle);
    p.xz = mat2(c, -s, s, c) * p.xz;

    // Primitives
    float cube = sdBox(p, vec3(1.0));
    float sphere = sdSphere(p, 1.25);
    float cylR = 0.45;
    float cylX = sdCylinder(p.yzx, 2.0, cylR);
    float cylY = sdCylinder(p.xyz, 2.0, cylR);
    float cylZ = sdCylinder(p.zxy, 2.0, cylR);

    // Hard boolean combination: nut = (cube intersect sphere) - three cylinders
    float nut = opSubtraction(
        opIntersection(cube, sphere),
        opUnion(cylX, opUnion(cylY, cylZ))
    );

    // Organic spheres -- smooth union blending
    float blob1 = sdSphere(p - vec3(1.8, 0.0, 0.0), 0.4);
    float blob2 = sdSphere(p - vec3(-1.8, 0.0, 0.0), 0.4);
    float blob3 = sdSphere(p - vec3(0.0, 1.8, 0.0), 0.4);
    float blobs = smin(blob1, smin(blob2, blob3, 0.3), 0.3);

    return smin(nut, blobs, 0.15);
}

// === Normal Calculation (Tetrahedral Sampling) ===
vec3 calcNormal(vec3 pos) {
    vec2 e = vec2(0.001, -0.001);
    return normalize(
        e.xyy * mapScene(pos + e.xyy) +
        e.yyx * mapScene(pos + e.yyx) +
        e.yxy * mapScene(pos + e.yxy) +
        e.xxx * mapScene(pos + e.xxx)
    );
}

// === Ray Marching ===
float rayMarch(vec3 ro, vec3 rd) {
    float t = 0.0;
    for (int i = 0; i < MAX_STEPS; i++) {
        vec3 p = ro + rd * t;
        float d = mapScene(p);
        if (d < SURF_DIST) return t;
        t += d;
        if (t > MAX_DIST) break;
    }
    return -1.0;
}

// === Soft Shadows ===
float calcSoftShadow(vec3 ro, vec3 rd, float k) {
    float res = 1.0;
    float t = 0.02;
    for (int i = 0; i < 64; i++) {
        float h = mapScene(ro + rd * t);
        res = min(res, k * h / t);
        t += clamp(h, 0.01, 0.2);
        if (res < 0.001 || t > 20.0) break;
    }
    return clamp(res, 0.0, 1.0);
}

// === AO (Ambient Occlusion) ===
float calcAO(vec3 pos, vec3 nor) {
    float occ = 0.0;
    float sca = 1.0;
    for (int i = 0; i < 5; i++) {
        float h = 0.01 + 0.12 * float(i);
        float d = mapScene(pos + h * nor);
        occ += (h - d) * sca;
        sca *= 0.95;
    }
    return clamp(1.0 - 3.0 * occ, 0.0, 1.0);
}

// === Main Function (WebGL2 Adapted) ===
out vec4 outColor;
void main() {
    vec2 uv = (gl_FragCoord.xy - 0.5 * iResolution.xy) / iResolution.y;

    // Camera
    float camDist = 4.0;
    float camAngle = 0.3;
    vec3 ro = vec3(
        camDist * cos(iTime * 0.2),
        camDist * sin(camAngle),
        camDist * sin(iTime * 0.2)
    );
    vec3 ta = vec3(0.0, 0.0, 0.0);

    // Camera matrix
    vec3 ww = normalize(ta - ro);
    vec3 uu = normalize(cross(ww, vec3(0.0, 1.0, 0.0)));
    vec3 vv = cross(uu, ww);
    vec3 rd = normalize(uv.x * uu + uv.y * vv + 2.0 * ww);

    // Background color
    vec3 col = vec3(0.4, 0.5, 0.6) - 0.3 * rd.y;

    // Ray marching
    float t = rayMarch(ro, rd);
    if (t > 0.0) {
        vec3 pos = ro + rd * t;
        vec3 nor = calcNormal(pos);

        vec3 lightDir = normalize(vec3(0.8, 0.6, -0.3));
        float dif = clamp(dot(nor, lightDir), 0.0, 1.0);
        float sha = calcSoftShadow(pos + nor * 0.01, lightDir, 16.0);
        float ao = calcAO(pos, nor);
        float amb = 0.5 + 0.5 * nor.y;

        vec3 mate = vec3(0.2, 0.3, 0.4);
        col = vec3(0.0);
        col += mate * 2.0 * dif * sha;
        col += mate * 0.3 * amb * ao;
    }

    col = pow(col, vec3(0.4545));
    outColor = vec4(col, 1.0);
}

Common Variants

Variant 1: Exponential Smooth Union

float sminExp(float a, float b, float k) {
    float res = exp(-k * a) + exp(-k * b);
    return -log(res) / k;
}

Variant 2: Smooth Operations with Color Blending

// Returns blend factor for the caller to blend colors
float sminWithFactor(float a, float b, float k, out float blend) {
    float h = clamp(0.5 + 0.5 * (b - a) / k, 0.0, 1.0);
    blend = h;
    return mix(b, a, h) - k * h * (1.0 - h);
}
// float blend;
// float d = sminWithFactor(d1, d2, 0.1, blend);
// vec3 color = mix(color2, color1, blend);

// vec3 overload of smax
vec3 smax(vec3 a, vec3 b, float k) {
    vec3 h = max(k - abs(a - b), 0.0);
    return max(a, b) + h * h * 0.25 / k;
}

Variant 3: Stepwise CSG Modeling (Architectural/Industrial)

float sdBuilding(vec3 p) {
    float walls = sdBox(p, vec3(1.0, 0.8, 1.0));
    vec3 roofP = p;
    roofP.y -= 0.8;
    float roof = sdBox(roofP, vec3(1.2, 0.3, 1.2));
    float d = opUnion(walls, roof);

    // Cut windows (exploiting symmetry)
    vec3 winP = abs(p);
    winP -= vec3(1.01, 0.3, 0.4);
    float window = sdBox(winP, vec3(0.1, 0.15, 0.12));
    d = opSubtraction(d, window);

    // Hollow out interior
    float hollow = sdBox(p, vec3(0.95, 0.75, 0.95));
    d = opSubtraction(d, hollow);
    return d;
}

Variant 4: Large-Scale Organic Character Modeling

float mapCharacter(vec3 p) {
    float body = sdEllipsoid(p, vec3(0.5, 0.4, 0.6));
    float head = sdEllipsoid(p - vec3(0.0, 0.5, 0.5), vec3(0.25));
    float d = smin(body, head, 0.2);           // large k: wide blend

    float ear = sdEllipsoid(p - vec3(0.3, 0.6, 0.3), vec3(0.15, 0.2, 0.05));
    d = smin(d, ear, 0.08);                    // medium blend

    float nostril = sdSphere(p - vec3(0.0, 0.4, 0.7), 0.03);
    d = smax(d, -nostril, 0.02);               // small k: fine sculpting
    return d;
}

Performance & Composition Tips

Performance:

• Bounding volume acceleration: use AABB/bounding spheres to skip distant sub-scenes, reducing mapScene() calls
• Tetrahedral sampling normals (4 samples) outperform central differences (6 samples)
• Step scaling t += d * 0.9 can reduce overshoot penetration
• Prefer quadratic optimized smin/smax (fastest); use exponential version when extreme smoothness is needed
• k must not be zero (division by zero error); fall back to hard boolean when near zero
• For symmetric shapes, use abs() to fold coordinates and define only one side

Composition techniques:

• + Domain Repetition: mod()/fract() for infinite repetition of CSG shapes (mechanical arrays, railings)
• + Procedural Displacement: overlay noise displacement on SDF for surface detail
• + Procedural Texturing: use smin blend factor to simultaneously blend material ID / color
• + 2D SDF: equally applicable to 2D scenes (clouds, UI shape compositing)
• + Animation: bind k values, positions, and radii to iTime for dynamic deformation

Further Reading

Full step-by-step tutorials, mathematical derivations, and advanced usage in reference