// MQ Extraction Algorithm
// McAulay-Quatieri sinusoidal analysis

// Extract partials from audio buffer
function extractPartials(params, stftCache) {
  const {fftSize, threshold, sampleRate} = params;
  const numFrames = stftCache.getNumFrames();

  const frames = [];
  for (let i = 0; i < numFrames; ++i) {
    const cachedFrame = stftCache.getFrameAtIndex(i);
    const squaredAmp = stftCache.getSquaredAmplitude(cachedFrame.time);
    const peaks = detectPeaks(squaredAmp, fftSize, sampleRate, threshold);
    frames.push({time: cachedFrame.time, peaks});
  }

  const partials = trackPartials(frames);

  for (const partial of partials) {
    partial.freqCurve = fitBezier(partial.times, partial.freqs);
    partial.ampCurve = fitBezier(partial.times, partial.amps);
  }

  return {partials, frames};
}

// Detect spectral peaks via local maxima + parabolic interpolation
// squaredAmp: pre-computed re*re+im*im per bin
function detectPeaks(squaredAmp, fftSize, sampleRate, thresholdDB) {
  const mag = new Float32Array(fftSize / 2);
  for (let i = 0; i < fftSize / 2; ++i) {
    mag[i] = 10 * Math.log10(Math.max(squaredAmp[i], 1e-20));
  }

  const peaks = [];
  for (let i = 2; i < mag.length - 2; ++i) {
    if (mag[i] > thresholdDB &&
        mag[i] > mag[i-1] && mag[i] > mag[i-2] &&
        mag[i] > mag[i+1] && mag[i] > mag[i+2]) {

      // Parabolic interpolation for sub-bin accuracy
      const alpha = mag[i-1];
      const beta  = mag[i];
      const gamma = mag[i+1];
      const p = 0.5 * (alpha - gamma) / (alpha - 2*beta + gamma);

      const freq  = (i + p) * sampleRate / fftSize;
      const ampDB = beta - 0.25 * (alpha - gamma) * p;
      peaks.push({freq, amp: Math.pow(10, ampDB / 20)});
    }
  }

  return peaks;
}

// Track partials across frames (birth/death/continuation)
function trackPartials(frames) {
  const partials        = [];
  const activePartials  = [];
  const candidates      = []; // pre-birth

  const trackingRatio   = 0.05; // 5% frequency tolerance
  const minTrackingHz   = 20;
  const birthPersistence = 3;   // frames before partial is born
  const deathAge        = 5;    // frames without match before death
  const minLength       = 10;   // frames required to keep partial

  for (const frame of frames) {
    const matched = new Set();

    // Continue active partials
    for (const partial of activePartials) {
      const lastFreq = partial.freqs[partial.freqs.length - 1];
      const tol = Math.max(lastFreq * trackingRatio, minTrackingHz);
      let bestIdx = -1, bestDist = Infinity;

      for (let i = 0; i < frame.peaks.length; ++i) {
        if (matched.has(i)) continue;
        const dist = Math.abs(frame.peaks[i].freq - lastFreq);
        if (dist < tol && dist < bestDist) { bestDist = dist; bestIdx = i; }
      }

      if (bestIdx >= 0) {
        const pk = frame.peaks[bestIdx];
        partial.times.push(frame.time);
        partial.freqs.push(pk.freq);
        partial.amps.push(pk.amp);
        partial.age = 0;
        matched.add(bestIdx);
      } else {
        partial.age++;
      }
    }

    // Advance candidates
    for (let i = candidates.length - 1; i >= 0; --i) {
      const cand = candidates[i];
      const lastFreq = cand.freqs[cand.freqs.length - 1];
      const tol = Math.max(lastFreq * trackingRatio, minTrackingHz);
      let bestIdx = -1, bestDist = Infinity;

      for (let j = 0; j < frame.peaks.length; ++j) {
        if (matched.has(j)) continue;
        const dist = Math.abs(frame.peaks[j].freq - lastFreq);
        if (dist < tol && dist < bestDist) { bestDist = dist; bestIdx = j; }
      }

      if (bestIdx >= 0) {
        const pk = frame.peaks[bestIdx];
        cand.times.push(frame.time);
        cand.freqs.push(pk.freq);
        cand.amps.push(pk.amp);
        matched.add(bestIdx);
        if (cand.times.length >= birthPersistence) {
          activePartials.push(cand);
          candidates.splice(i, 1);
        }
      } else {
        candidates.splice(i, 1);
      }
    }

    // Spawn candidates from unmatched peaks
    for (let i = 0; i < frame.peaks.length; ++i) {
      if (matched.has(i)) continue;
      const pk = frame.peaks[i];
      candidates.push({times: [frame.time], freqs: [pk.freq], amps: [pk.amp], age: 0});
    }

    // Kill aged-out partials
    for (let i = activePartials.length - 1; i >= 0; --i) {
      if (activePartials[i].age > deathAge) {
        if (activePartials[i].times.length >= minLength) partials.push(activePartials[i]);
        activePartials.splice(i, 1);
      }
    }
  }

  // Collect remaining active partials
  for (const partial of activePartials) {
    if (partial.times.length >= minLength) partials.push(partial);
  }

  return partials;
}

// Fit cubic bezier to trajectory using samples at ~1/3 and ~2/3 as control points
function fitBezier(times, values) {
  const n = times.length - 1;
  const t0 = times[0], v0 = values[0];
  const t3 = times[n], v3 = values[n];
  const dt = t3 - t0;

  if (dt <= 0 || n === 0) {
    return {t0, v0, t1: t0, v1: v0, t2: t3, v2: v3, t3, v3};
  }

  const v1 = values[Math.round(n / 3)];
  const v2 = values[Math.round(2 * n / 3)];

  return {t0, v0, t1: t0 + dt / 3, v1, t2: t0 + 2 * dt / 3, v2, t3, v3};
}