Parametric score generation with IFS using WebGPU/WGSL

[gogins]

WebGPU makes it easier to write compute shaders that are far faster than CPU code. Back-of-the-envelope estimates indicate that parametric composition using IFS and implemented using WebGPU may, just possibly, become marginally tractable.

This project would consist of:

• IFS measures in dimensions {time, pitch, instrument} by random iteration, possibly using higher precision.
• Measures on IFS to CsoudAC scores, optionally using higher precision.
• Graphical displays:
  • The {time, pitch} dimensions of the IFS colored by {instrument -> hue, loudness -> value}.
  • Piano roll of generated score.
  • Two-dimensional parametric map of IFS using Hilbert indices on the unit square normalized to the graphics canvas and computational shortcuts for coloring points based their associated piano rolls:
    • Temporal regularity -> hue.
    • Tonal clarity -> saturation.
    • Voice-leading smoothness -> value.

[gogins]

Yes. Treat the piano roll as data (a binary/valued matrix time×pitch) and probe it for predictability, repetition, and structure. Below is a compact toolkit of measures you can compute and then combine into a single “non-randomness” score.

1) Set up

Let (P \in {0,1}^{T\times M}) (or velocities in ([0,1])). Also keep:

• Onset train (o_t=\mathbb{1}[\sum_p P_{t,p}>0])
• Active-pitch vector (a_t\in{0,1}^M) (row (t) of (P))

Use density (r=\frac{1}{TM}\sum_{t,p}P_{t,p}) to build “random” baselines (i.i.d. Bernoulli with rate (r), or Poisson onsets with rate (\lambda) matched to data).

2) Rhythm regularity & periodicity

• Autocorrelation peakiness: ( \mathrm{AC}(k)=\sum_t o_t o_{t-k} ). Take max peak (k≠0) / mean. More peak → more periodic.
• Spectral flatness (onsets): flatness of |FFT(o)|. Lower flatness → more metrical structure.
• Tempogram peak strength: peak of tempo-domain representation (e.g., Fourier tempogram). Strong peaks → regular meter.
• Runs test (Wald–Wolfowitz) on (o_t): rejects randomness if long runs or alternations occur.

3) Predictability / entropy

• Permutation entropy of (o_t) (and of row-sums (\sum_p P_{t,p})): lower = more order.
• n-gram/Markov cross-entropy: train an n-gram over bar-quantized onset patterns; evaluate per-token cross-entropy (H_n). Lower → more predictable (less random).
• Entropy rate via Lempel–Ziv complexity on serialized rows (e.g., interleave barwise). Smaller LZ size → more structure.
• Normalized Compression Distance (NCD):
  [
  \mathrm{NCD}(X,\text{noise})=\frac{C(X!\parallel!\text{noise})-\min{C(X),C(\text{noise})}}{\max{C(X),C(\text{noise})}}
  ]
  where (C) is compressed size (LZMA). Low distance from pure noise ⇒ random; high ⇒ structured.

4) Pitch organization & harmony

• Pitch-class histogram uniformity: KL divergence (D_{KL}(p(\text{pc})\parallel \text{uniform})). Higher → tonal bias.
• Krumhansl–Schmuckler key correlation: max correlation over 24 keys. Higher → clearer key center.
• Chord consonance / interval profile: proportion of intervals in {m3, M3, P5, m6, M6} vs dissonant set; or psychoacoustic roughness model—higher consonance share ⇒ less random.
• Voice-leading smoothness: mean absolute pitch step between successive active states of each voice (track voices greedily or by Hungarian matching). Smaller steps ⇒ structured lines.

5) Spatial/time–pitch dependencies

• Lagged mutual information:

  • Time: (I(a_t; a_{t-k})) for small k (e.g., 1–8 ticks). Larger MI ⇒ recurring patterns.
  • Pitch neighbors: (I(P_{t,p}; P_{t,p\pm1})) (texture smoothness).

• Block repetition: count unique bar-patterns and their reuse; Zipf-like skew indicates motif reuse.

6) Low-rank structure

Treat (P) as a matrix:

• SVD/NMF reconstruction error: fit rank-k (e.g., k=4–8) and report ( |P-\hat P_k|_F / |P|_F ). Lower error ⇒ few latent patterns (motivic/meter grids).
• Matrix rank over (\mathbb{R}) (or approximate) relative to min(T,M): lower effective rank ⇒ regular.

7) Distributional tests (vs matched random)

• IOI distribution: compare inter-onset intervals to exponential (Poisson) or geometric via KS test; significant deviation ⇒ structure.
• Bar-position occupancy: χ² test of counts per 16th over a bar against uniform; strong peaks ⇒ metrical regularity.

---

A single “Non-Randomness Index” (NRI)

Compute features, z-score each against a matched random surrogate (same (T,M,r)) by Monte-Carlo shuffles. Then:
[
\mathrm{NRI}=\sum_i w_i, z_i
]
Suggested positive-polarity features (higher ⇒ more structure):
AC peakiness, 1–spectral-flatness, tempogram peak, −LZ size, −cross-entropy, KL(pc‖uniform), key correlation, consonance share, −VL step size, MI(lags), 1–NMF error, χ²(bar positions).

Weights (w_i): equal to start; tune by your notion of “musically listenable” (you’ve mentioned metric regularity, consonance, and degree of randomness—give those bigger weights).

---

Practical recipe (fast and robust)

1. Quantize to 16th or 24th notes; split into bars.

2. Build 100 surrogate rolls by time-scrambling within each bar and pitch-scrambling to preserve density/range.

3. Compute the 10–15 features above on the real roll and surrogates → z-scores.

4. Report NRI and a small dashboard:

   • Rhythm: AC peak, tempogram peak, χ²(bar), permutation entropy.
   • Pitch: KL(pc), key correlation, consonance share, VL smoothness.
   • Global: LZ size, NMF error, MI@lags.

---

Tiny pseudo-code sketch

def non_randomness_index(P, q=16, k_rank=6, S=100):
    feats = measure_features(P, q, k_rank)
    sur = [measure_features(scramble(P), q, k_rank) for _ in range(S)]
    mu, sigma = np.mean(sur, axis=0), np.std(sur, axis=0) + 1e-9
    z = (feats - mu) / sigma
    w = np.array([ # example weights (sum to 1)
        .12,.10,.10, .08, .10, .10, .10, .10, .10, .10
    ])
    return float(np.dot(w, z)), dict(zip(FEAT_NAMES, z))

---

Notes for edge cases

• Super-sparse rolls: rely more on rhythm (IOI/AC) and less on harmonic measures.
• Very dense textures: normalize by active notes per frame before consonance/voice-leading.
• Multiscale: compute features at two quantizations (e.g., 12th & 16th) and average.

If you want, I can turn this into a ready-to-run Python (or JS) script that ingests a MIDI, quantizes to a piano roll, and prints the NRI with a per-feature breakdown.

[gogins]

Here’s a self-contained ES module that computes the three indicators—rhythmic regularity, tonal clarity, and voice-leading smoothness—plus a combined musicality index.

• Input roll is a T × M matrix (array of arrays) with 0/1 or velocity values (time along rows, MIDI pitch 0..M-1 along columns).
• No dependencies. Fast enough for typical piano-roll sizes; you can swap the autocorr routine for an FFT later if needed.

// musicality.js  —  ES module
// Exported API:
//   rhythmicRegularity(roll, {assumeDurations=true})
//   tonalClarity(roll)
//   voiceLeadingSmoothness(roll, {bothDirections=true})
//   musicalityIndex(roll, {w1=1, w2=1, w3=1})

/** ------- utilities ------- **/

function shape(roll) {
  const T = roll.length;
  const M = T ? roll[0].length : 0;
  return { T, M };
}

function rowSum(arr) { return arr.reduce((a, b) => a + b, 0); }

function clamp01(x) { return Math.max(0, Math.min(1, x)); }

/** Make onset vector o_t (1 if any new note starts at t).
 * If assumeDurations=false, treats any activity at t as an onset. */
function onsetVector(roll, assumeDurations = true) {
  const { T, M } = shape(roll);
  const o = new Float32Array(T);
  if (T === 0) return o;
  for (let t = 0; t < T; t++) {
    let onset = 0;
    for (let p = 0; p < M; p++) {
      const v = roll[t][p] || 0;
      const prev = t > 0 ? (roll[t - 1][p] || 0) : 0;
      const starts = assumeDurations ? (v > 0 && prev <= 0) : (v > 0);
      if (starts) { onset = 1; break; }
    }
    o[t] = onset;
  }
  return o;
}

/** Unbiased autocorrelation (naive O(T^2)); returns AC[1..maxLag] */
function autocorrPeakiness(o, maxLag = null) {
  const T = o.length;
  if (T < 3) return 0;
  const mean = o.reduce((a, b) => a + b, 0) / T;
  let varsum = 0;
  for (let t = 0; t < T; t++) {
    const d = o[t] - mean;
    varsum += d * d;
  }
  if (varsum === 0) return 0;

  const L = maxLag ? Math.min(maxLag, Math.floor(T / 2)) : Math.floor(T / 2);
  const ac = new Float64Array(L + 1);
  ac[0] = 1; // normalized
  for (let k = 1; k <= L; k++) {
    let s = 0;
    for (let t = k; t < T; t++) {
      s += (o[t] - mean) * (o[t - k] - mean);
    }
    ac[k] = s / varsum;
  }
  // Peakiness: strongest non-zero AC peak relative to average remainder
  const rest = ac.slice(1);
  const maxPeak = Math.max(...rest);
  const meanRest = rest.reduce((a, b) => a + b, 0) / rest.length;
  const denom = Math.max(1e-12, Math.abs(meanRest));
  // Normalize to ~[0,1]: rel contrast, clipped
  return clamp01((maxPeak - meanRest) / (1 - meanRest));
}

/** Pitch-class histogram (12 bins), weighted by velocities or counts. */
function pitchClassHistogram(roll) {
  const { T, M } = shape(roll);
  const pc = new Float64Array(12);
  for (let t = 0; t < T; t++) {
    for (let p = 0; p < M; p++) {
      const v = roll[t][p] || 0;
      if (v > 0) pc[p % 12] += v;
    }
  }
  const sum = pc.reduce((a, b) => a + b, 0);
  if (sum > 0) for (let i = 0; i < 12; i++) pc[i] /= sum;
  return pc;
}

/** Pearson correlation between two equal-length arrays. */
function pearson(x, y) {
  const n = x.length;
  let sx = 0, sy = 0, sxx = 0, syy = 0, sxy = 0;
  for (let i = 0; i < n; i++) {
    const xi = x[i], yi = y[i];
    sx += xi; sy += yi;
    sxx += xi * xi; syy += yi * yi; sxy += xi * yi;
  }
  const mx = sx / n, my = sy / n;
  const vx = Math.max(1e-15, sxx - n * mx * mx);
  const vy = Math.max(1e-15, syy - n * my * my);
  return (sxy - n * mx * my) / Math.sqrt(vx * vy);
}

/** Minimal symmetric nearest-neighbor distance between pitch sets */
function meanNearestPitchStep(roll, bothDirections = true) {
  const { T, M } = shape(roll);
  // Extract active pitch lists per frame
  const act = new Array(T);
  for (let t = 0; t < T; t++) {
    const pts = [];
    for (let p = 0; p < M; p++) if ((roll[t][p] || 0) > 0) pts.push(p);
    act[t] = pts;
  }
  let total = 0, count = 0;

  const stepOneDir = (A, B) => {
    if (A.length === 0 || B.length === 0) return { sum: 0, n: 0 };
    let s = 0;
    for (const a of A) {
      let best = Infinity;
      // linear scan — M is small in practice
      for (const b of B) {
        const d = Math.abs(a - b);
        if (d < best) best = d;
        if (best === 0) break;
      }
      s += best;
    }
    return { sum: s, n: A.length };
  };

  for (let t = 0; t < T - 1; t++) {
    const A = act[t], B = act[t + 1];
    const fwd = stepOneDir(A, B);
    if (bothDirections) {
      const bwd = stepOneDir(B, A);
      total += fwd.sum + bwd.sum;
      count += fwd.n + bwd.n;
    } else {
      total += fwd.sum;
      count += fwd.n;
    }
  }
  if (count === 0) return 1; // degenerate: treat as large
  return total / count; // in semitones
}

/** Pitch range in semitones across the whole roll (for normalization). */
function pitchRange(roll) {
  const { T, M } = shape(roll);
  let lo = Infinity, hi = -Infinity;
  for (let t = 0; t < T; t++) {
    for (let p = 0; p < M; p++) {
      if ((roll[t][p] || 0) > 0) {
        if (p < lo) lo = p;
        if (p > hi) hi = p;
      }
    }
  }
  if (!Number.isFinite(lo)) return 1;
  return Math.max(1, hi - lo + 1);
}

/** ------- indicators ------- **/

/** 1) Rhythmic Regularity in [0,1] (higher = clearer beat) */
export function rhythmicRegularity(roll, { assumeDurations = true } = {}) {
  const o = onsetVector(roll, assumeDurations);
  return autocorrPeakiness(o);
}

/** 2) Tonal Clarity via Krumhansl–Schmuckler correlation in [0,1] */
export function tonalClarity(roll) {
  // Major/minor profiles (normalized)
  // Values from Krumhansl–Kessler (commonly used variant)
  const maj = [6.35, 2.23, 3.48, 2.33, 4.38, 4.09, 2.52, 5.19, 2.39, 3.66, 2.29, 2.88];
  const min = [6.33, 2.68, 3.52, 5.38, 2.60, 3.53, 2.54, 4.75, 3.98, 2.69, 3.34, 3.17];
  // Normalize profiles
  const norm = a => {
    const s = a.reduce((x, y) => x + y, 0);
    return a.map(v => v / s);
  };
  const majN = norm(maj), minN = norm(min);

  const pc = pitchClassHistogram(roll);
  // Consider all 12 rotations (keys)
  const rotate = (arr, k) => arr.map((_, i) => arr[(i - k + 12) % 12]);
  let best = 0;
  for (let k = 0; k < 12; k++) {
    const rMaj = pearson(pc, rotate(majN, k));
    const rMin = pearson(pc, rotate(minN, k));
    best = Math.max(best, rMaj, rMin);
  }
  // Map Pearson r [-1,1] to [0,1]
  return clamp01((best + 1) / 2);
}

/** 3) Voice-Leading Smoothness in [0,1] (higher = smoother) */
export function voiceLeadingSmoothness(roll, { bothDirections = true } = {}) {
  const meanStep = meanNearestPitchStep(roll, bothDirections); // semitones
  const pr = pitchRange(roll);
  // Normalize: 0 step -> 1.0; step ~ pr/2 -> ~0
  const norm = 1 - clamp01(meanStep / (pr / 2));
  return norm;
}

/** Combined Musicality Index (weighted average, returns [0,1]) */
export function musicalityIndex(
  roll,
  {
    w1 = 1, // rhythm
    w2 = 1, // tonal
    w3 = 1, // voice-leading
    assumeDurations = true,
    bothDirections = true
  } = {}
) {
  const r = rhythmicRegularity(roll, { assumeDurations });
  const t = tonalClarity(roll);
  const v = voiceLeadingSmoothness(roll, { bothDirections });
  const W = Math.max(1e-9, w1 + w2 + w3);
  return {
    rhythmicRegularity: r,
    tonalClarity: t,
    voiceLeadingSmoothness: v,
    musicality: (w1 * r + w2 * t + w3 * v) / W
  };
}

Usage

import {
  rhythmicRegularity,
  tonalClarity,
  voiceLeadingSmoothness,
  musicalityIndex
} from './musicality.js';

// Example: roll is T×M array, e.g., 128 time steps × 88 MIDI pitches
const roll = Array.from({ length: 128 }, () => Array(88).fill(0));
// ... populate roll with 0/1 or velocities ...

const r = rhythmicRegularity(roll);
const t = tonalClarity(roll);
const v = voiceLeadingSmoothness(roll);
const combo = musicalityIndex(roll, { w1: 1, w2: 1, w3: 1 });

console.log({ r, t, v, combo });

Notes / tweaks

• Rhythm: for very long rolls, you can cap maxLag inside autocorrPeakiness (e.g., to 1–2 bars in ticks) for speed and relevance.
• Tonal clarity: if you work with microtonal or pitch-rich material, replace the 12-bin histogram with your own scale template(s).
• Voice-leading: this uses nearest-neighbor across frames. If you want stricter polyphonic tracking, swap in a Hungarian assignment per frame pair.

If you want a variant that ingests MIDI files directly (Node or browser with Web MIDI/MIDIFile), I can add a small loader and quantizer around this module.

[gogins]

Evaluation of the project plan:

Here’s a straight, engineering-minded take on feasibility + a concrete plan to make this tractable.

Feasibility (short verdict)

• Yes, marginally tractable—if you: (1) keep GPU kernels branch-light, (2) tile histograms locally to avoid global atomic contention, (3) reserve “higher precision” for only the numerically touchy steps, and (4) render the parametric map progressively with adaptive sampling + caching.
• Expect 10–100× over a well-tuned CPU baseline for your IFS random-iteration and feature extraction, with discrete GPUs nearer the top end and integrated GPUs nearer the bottom. This is enough to make exploratory, progressive maps practical.

---

Architecture at a glance

1. GPU pass A — IFS random iteration → measures
   Input: an IFS definition (affine/similarity maps + probs) and the parameter point(s).
   Output: per-point streaming measures sufficient to (a) render the {time,pitch,instrument} scatter, (b) synthesize a score, (c) compute the three musicality indicators without pulling huge traces back to the CPU.

2. GPU pass B — feature extraction (on GPU)
   Compute approximations of the three indicators directly from streaming tallies:

   • Rhythmic regularity ≈ lag-autocorrelation on onset counts (few lags).
   • Tonal clarity ≈ pitch-class histogram → best key correlation (12 rotations).
   • Voice-leading smoothness ≈ mean nearest-neighbor step (1-frame lookahead).

3. CPU/JS (optional) post
   If you need exact versions (e.g., your JS module for musicality), transfer compact summaries (histograms, lag sums) rather than full rolls.

4. Visualization

   • {time,pitch} scatter plot colored by {instrument→hue, loudness→value}.
   • Piano-roll: render from compact event buffers.
   • Parametric Hilbert map: progressive quadtree tiles; each tile caches point color = (temporal regularity hue, tonal clarity sat, voice-leading value).

---

Bullet-by-bullet evaluation

1) IFS measures in {time, pitch, instrument} by random iteration, possibly higher precision

Tractable.

• Kernel shape: each thread advances its own IFS random walk (a few hundred–thousand steps), emitting events into thread-local accumulators. Periodically reduce into workgroup-local bins; at the end, one atomic add per bin to a global buffer. This avoids hammering global atomics each step.
• Coordinate model: treat {time, pitch} as continuous in [0,1]^2 and map to discrete bins for measures/visuals; instrument as a small int or categorical index; loudness as value in [0,1].
• RNG: use a per-thread counter-based generator (e.g., Philox-like) or a hashed LCG with Weyl sequences; entirely branchless and seedable by (pointId, threadId, step).
• Higher precision: only for (a) deep param zoom parameters, (b) accumulation of numerically sensitive stats. Start with f32; use compensated sums (Kahan/Neumaier) and fixed-point for histogram totals. If/when you hit precision limits (e.g., extremely deep zooms or small probability maps), switch specific inner products to the double-double routines (the 2×f32 WGSL you already have). Do not blanket everything in DD—keep it isolated (e.g., when updating z^2 + c-style maps or applying very ill-conditioned transforms).

Key data structures (global):

• histTP[timeBins][pitchBins] (u32/atomic) — for quick scatter and onset counts.
• pcHist[12] (u32/atomic) — pitch-class usage.
• onsetLag[kMax] (u64/atomic) — ∑ o(t)*o(t−k).
• activePerFrame[timeBins] (u32/atomic) — #notes per frame (to derive onsets).
• stepSum, stepCount (u64/atomic) — voice-leading distance totals.
• Optional compact event ring buffer per workgroup for piano-roll rendering.

2) Measures → CsoundAC scores, optionally higher precision

Tractable with a compact event buffer.

• Emit only quantized onsets, pitches, velocities into a GPU->CPU staging buffer (struct of {t,p,vel,inst}). Use a fixed-size ring per workgroup; if overflow, (a) increase work quota next pass, or (b) probabilistically downsample (since you’re aiming for perceptual stats more than exact micro-detail).
• On CPU, transform events into CsoundAC score objects. If needed, refine timing/velocity with higher precision arithmetic on the CPU side—cheap, because the event count is already pruned.

3) Graphics displays

{time,pitch} colored by {instrument→hue, loudness→value}

• Render from histTP or from the event buffer. A simple compute → storage texture write per bin is enough. Hue map by instrument index, value by normalized loudness.

Piano roll

• GPU rasterize events into a 2D storage texture (time × pitch). Or draw instanced quads in a render pass. The same event buffer serves both the roll and the CsoundAC step.

Parametric map via Hilbert indices on the unit square

• Use a progressive tile renderer:

  • Queue parameter points by Hilbert order to keep spatial coherence (nice perceived convergence).
  • For each point, run Pass A+B with a budget (N steps). Store its triple {regularity, clarity, smoothness}.
  • Color mapping: hue=regularity, sat=clarity, value=smoothness (your chosen scheme).

• Shortcuts / early stopping: after every M steps, check convergence of the three features (delta below ε for a few consecutive checks) → stop early for “easy” points; reallocate steps to “hard” points.

• Caching: memoize (parameter → features) in an LRU keyed by (IFS params rounded to a sensible grid).

---

Concrete kernel sketches (WGSL-ish pseudocode)

Random iteration + local tallies

struct Params {
  timeBins: u32, pitchBins: u32, kMax: u32, stepsPerThread: u32;
  // IFS transforms, probs, etc…
};

@group(0) @binding(0) var<uniform> params: Params;
@group(0) @binding(1) var<storage, read_write> histTP: array<atomic<u32>>;
@group(0) @binding(2) var<storage, read_write> pcHist: array<atomic<u32>>;
@group(0) @binding(3) var<storage, read_write> onsetLag: array<atomic<u64>>;
@group(0) @binding(4) var<storage, read_write> activePerFrame: array<atomic<u32>>;
@group(0) @binding(5) var<storage, read_write> stepSum: atomic<u64>;
@group(0) @binding(6) var<storage, read_write> stepCount: atomic<u64>;
// … event buffer bindings if you want a piano roll

var<workgroup> tpTile: array<u32, TILE_T * TILE_P>;
var<workgroup> pcLocal: array<u32, 12>;
var<workgroup> lagLocal: array<u64, KMAX_LOCAL>;

@compute @workgroup_size(256)
fn iterate(@builtin(global_invocation_id) gid: vec3<u32>) {
  // init local tiles
  // walk a per-thread RNG; advance IFS state; compute (t,p,inst,vel)
  // accumulate:
  //   - tpTile[bin]++
  //   - pcLocal[p % 12] += 1
  //   - activePerFrame[t]++ (atomic global is okay: once per step)
  //   - for small k=1..K: track recent onsets in a shift register and update lagLocal[k]
  //   - voice-leading: nearest-neighbor to prev frame's active pitches in a tiny scratch buffer; update stepSum/stepCount (use local first, then global)

  // end: one wave of atomic adds from local → global per bin
}

Why this works: 99% of updates hit workgroup-local memory; you pay global atomics only once per bin per workgroup, which scales.

Feature extraction on GPU (from tallies)

• Rhythmic regularity: from onsetLag[k] and activePerFrame, compute normalized autocorr coefficients and a simple peakiness metric.
• Tonal clarity: from pcHist[12], test 12 rotations × {maj,min} templates; take max correlation.
• Voice-leading: stepSum/stepCount → mean step; normalize by pitch range (keep track of minPitch/maxPitch during iteration).

Write the triple to a compact buffer that the CPU reads when it colors the parametric canvas. This is a few dozen numbers—cheap to copy.

---

Precision plan

• Default: f32 + careful accumulation.

• Compensation: use Neumaier/Kahan in local sums (WGSL friendly, ~2–3 ops more).

• Fixed-point: for very long counts, use u64 atomics where needed (onsetLag, stepSum) to avoid overflow.

• Double-double: keep your 2×f32 routines handy for:

  • numerically sensitive transforms in the IFS (rare),
  • deep zoom coordinate mapping when you quantize time,pitch (less rare if you zoom absurdly deep).
    Use them only at those touch points; the rest stays f32.

---

Performance knobs (so you can keep it “marginally tractable”)

• Steps per point: start with ~1e5–1e6 total steps shared across threads for a param point; adapt by convergence of the three features.
• Lags: keep kMax tiny (e.g., 16 or 32 ticks). You don’t need full FFTs to detect pulse.
• Bins: don’t exceed what you can render meaningfully. E.g., 1024×128 for {time,pitch} is already detailed.
• Workgroup tiling: choose tiles that fit in workgroup memory (e.g., 32×16 bins).
• Scheduling: traverse the param map by Hilbert order and subdivide tiles quadtree-style; this makes progressive previews pleasant.

---

Risk & mitigation

• Atomic hot-spots: solved by workgroup tiling and batched reduction.
• RNG quality vs cost: counter-based or 64-bit hashed LCG is fine; test a few seeds for banding.
• Precision drift: add unit tests comparing GPU tallies vs CPU reference on small cases; enable a “DD guardrail” path behind a flag.
• Data transfer: never bring full traces back—only compact summaries and optional event samples.

---

Minimal prototype (2–3 days of focused work)

1. Compute pass: IFS random iteration with tiled histTP, pcHist, onsetLag[8], stepSum/count.
2. Feature pass: compute (regularity, clarity, smoothness) per param point.
3. Canvas: draw param map with your hue/sat/value mapping; progressive tiles.
4. Roll: simple raster from an event sample buffer (downsampled if needed).
5. Toggle: CPU verification mode for a tiny instance; print deltas.

---

“Go / No-Go” for your statement

• “WebGPU makes it easier to write compute shaders that are far faster than CPU code.” — Go. The memory model and WGSL ergonomics are good enough; the key is to minimize atomics and branch divergence.
• “Parametric composition via IFS may become marginally tractable.” — Go (with conditions). With progressive rendering, local tiling, and selective precision, you can explore sizeable swaths of parameter space interactively.

If you want, I can drop in a skeleton WGSL + JS host that (a) runs the tiled random-iteration kernel, (b) computes the three features on-GPU, and (c) paints the param map progressively.

[gogins]

Dropped a working single-file prototype in the canvas. It:

• runs one WebGPU compute invocation per pixel (param point),
• does a placeholder IFS random walk and computes the three features on-GPU,
• colors the param map with hue = rhythmic regularity, saturation = tonal clarity, value = voice-leading smoothness, and
• blits to the canvas via a tiny render pass.

Where to plug in your domain logic:

• IFS mapping: replace ifs_step() with your {time, pitch, instrument, loudness} update.
• Onset detection: swap the placeholder (s & 7u)==0u with your event rule.
• Pitch mapping: replace tp_bin() and MIDI/PC handling with your quantizer.
• Bins & lags: tune pitchBins, lagCount, and steps/px in the UI.

If you want a second pass that emits a piano roll texture and an event buffer suitable for CsoundAC (plus a toggle for double-double precision in the hot spots), I can extend this skeleton accordingly.