1 files changed, 30 insertions, 36 deletions

diff --git a/tools/mq_editor/mq_extract.js b/tools/mq_editor/mq_extract.js
index ff38d63..0940cf1 100644
--- a/tools/mq_editor/mq_extract.js
+++ b/tools/mq_editor/mq_extract.js
@@ -270,13 +270,12 @@ function expandPartialsLeft(partials, frames) {
   }
 }
 
-// Autodetect spread_above / spread_below from the spectrogram.
-// For each (subsampled) STFT frame within the partial, measures the
-// half-power (-3dB) width of the spectral peak above and below the center.
-// spread = half_bandwidth / f0  (fractional).
+// Autodetect spread from the spectrogram.
+// Measures -3dB half-bandwidth above and below the peak in each STFT frame,
+// returns spread = max(above, below) / f0 as a fractional frequency offset.
 function autodetectSpread(partial, stftCache, fftSize, sampleRate) {
   const curve = partial.freqCurve;
-  if (!curve || !stftCache) return {spread_above: 0.02, spread_below: 0.02};
+  if (!curve || !stftCache) return {spread: 0.02};
 
   const numFrames = stftCache.getNumFrames();
   const binHz    = sampleRate / fftSize;
@@ -331,9 +330,9 @@ function autodetectSpread(partial, stftCache, fftSize, sampleRate) {
     ++count;
   }
 
-  const spread_above = count > 0 ? Math.sqrt(sumAbove / count) : 0.01;
-  const spread_below = count > 0 ? Math.sqrt(sumBelow / count) : 0.01;
-  return {spread_above, spread_below};
+  const sa = count > 0 ? Math.sqrt(sumAbove / count) : 0.01;
+  const sb = count > 0 ? Math.sqrt(sumBelow / count) : 0.01;
+  return {spread: Math.max(sa, sb)};
 }
 
 // Track a single partial starting from a (time, freq) seed position.
@@ -428,7 +427,7 @@ function trackFromSeed(frames, seedTime, seedFreq, params) {
   return {
     times: allTimes, freqs: allFreqs, amps: allAmps, phases: allPhases,
     muted: false, freqCurve,
-    replicas: { decay_alpha: 0.1, jitter: 0.05, spread_above: 0.02, spread_below: 0.02 },
+    harmonics: { decay: 0.0, freq_mult: 2.0, jitter: 0.05, spread: 0.02 },
   };
 }
 
@@ -451,7 +450,7 @@ function trackIsoContour(stftCache, seedTime, seedFreq, params) {
   }
 
   const seedFrame = stftCache.getFrameAtIndex(seedFrameIdx);
-  const seedBin   = Math.max(1, Math.min(halfBins - 2, Math.round(seedFreq / binHz)));
+  const seedBin   = clamp(Math.round(seedFreq / binHz), 1, halfBins - 2);
   const targetSq  = seedFrame.squaredAmplitude[seedBin];
   if (targetSq <= 0) return null;
   const targetDB  = 10 * Math.log10(targetSq);
@@ -522,11 +521,14 @@ function trackIsoContour(stftCache, seedTime, seedFreq, params) {
     times: allTimes, freqs: allFreqs, amps: allAmps,
     phases: new Array(allTimes.length).fill(0),
     muted: false, freqCurve,
-    replicas: { decay_alpha: 0.1, jitter: 0.05, spread_above: 0.15, spread_below: 0.15 },
+    harmonics: { decay: 0.0, freq_mult: 2.0, jitter: 0.05, spread: 0.15 },
   };
 }
 
-// Fit cubic bezier to trajectory using least-squares for inner control points
+// Fit interpolating curve to trajectory via least-squares for inner control point values.
+// Inner knots fixed at u=1/3 and u=2/3 (t = t0+dt/3, t0+2*dt/3).
+// The curve passes through all 4 control points (Lagrange interpolation).
+// TODO: support arbitrary number of inner control points
 function fitBezier(times, values) {
   const n = times.length - 1;
   const t0 = times[0], v0 = values[0];
@@ -534,43 +536,35 @@ function fitBezier(times, values) {
   const dt = t3 - t0;
 
   if (dt <= 1e-9 || n < 2) {
-    // Linear fallback for too few points or zero duration
     return {t0, v0, t1: t0 + dt / 3, v1: v0 + (v3 - v0) / 3, t2: t0 + 2 * dt / 3, v2: v0 + 2 * (v3 - v0) / 3, t3, v3};
   }
 
-  // Least squares solve for v1, v2
-  // Bezier: B(u) = (1-u)^3*v0 + 3(1-u)^2*u*v1 + 3(1-u)*u^2*v2 + u^3*v3
-  // Target_i = val_i - (1-u)^3*v0 - u^3*v3
-  // Model_i = A_i*v1 + B_i*v2
-  // A_i = 3(1-u)^2*u
-  // B_i = 3(1-u)*u^2
+  // Lagrange basis with inner knots at u1=1/3, u2=2/3
+  // l1(u) = u*(u-2/3)*(u-1) / ((1/3)*(1/3-2/3)*(1/3-1))  = 13.5*u*(u-2/3)*(u-1)
+  // l2(u) = u*(u-1/3)*(u-1) / ((2/3)*(2/3-1/3)*(2/3-1))  = -13.5*u*(u-1/3)*(u-1)
+  // l0(u) = (u-1/3)*(u-2/3)*(u-1) / ((-1/3)*(-2/3)*(-1)) = -4.5*(u-1/3)*(u-2/3)*(u-1)
+  // l3(u) = u*(u-1/3)*(u-2/3) / ((2/3)*(1/3))             = 4.5*u*(u-1/3)*(u-2/3)
+  // Least-squares: minimize Σ(l1*v1 + l2*v2 - target_i)^2
+  // target_i = values[i] - l0*v0 - l3*v3
 
   let sA2 = 0, sB2 = 0, sAB = 0, sAT = 0, sBT = 0;
 
   for (let i = 0; i <= n; ++i) {
-    const u = (times[i] - t0) / dt;
-    const u2 = u * u;
-    const u3 = u2 * u;
-    const invU = 1.0 - u;
-    const invU2 = invU * invU;
-    const invU3 = invU2 * invU;
-
-    const A = 3 * invU2 * u;
-    const B = 3 * invU * u2;
-    const target = values[i] - (invU3 * v0 + u3 * v3);
-
-    sA2 += A * A;
-    sB2 += B * B;
-    sAB += A * B;
-    sAT += A * target;
-    sBT += B * target;
+    const u  = (times[i] - t0) / dt;
+    const l0 = -4.5 * (u - 1/3) * (u - 2/3) * (u - 1);
+    const l1 =  13.5 * u * (u - 2/3) * (u - 1);
+    const l2 = -13.5 * u * (u - 1/3) * (u - 1);
+    const l3 =   4.5 * u * (u - 1/3) * (u - 2/3);
+    const A = l1, B = l2;
+    const target = values[i] - l0 * v0 - l3 * v3;
+    sA2 += A * A; sB2 += B * B; sAB += A * B;
+    sAT += A * target; sBT += B * target;
   }
 
   const det = sA2 * sB2 - sAB * sAB;
   let v1, v2;
 
   if (Math.abs(det) < 1e-9) {
-    // Fallback to simple 1/3, 2/3 heuristic if matrix is singular
     const idx1 = Math.round(n / 3);
     const idx2 = Math.round(2 * n / 3);
     v1 = values[idx1];