// How nonlinear is a run? Compare the listener signal with the same strike
// scaled down by 1000 (then scaled back up). Usage: tsx tests/nonlin.ts script.m [duration] [v1 v2 ...]
import { readFileSync, readdirSync } from "fs";
import { runScript } from "../src/sim/runScript";
import { simulate } from "../src/sim/solver";
import { relDiff, spectrum } from "../src/sim/analysis";
const helpers = readdirSync("matlab/helpers").map(name => ({ name, source: readFileSync("matlab/helpers/" + name, "utf8") }));
const { sys } = runScript(readFileSync(process.argv[2], "utf8"), helpers, () => {});
sys.duration = Number(process.argv[3] ?? 0.3);
const scales = process.argv.slice(4).map(Number);
for (const sc of scales.length ? scales : [1]) {
  const hits = sys.hits.map(h => ({ ...h, velocity: h.velocity.map(v => v * sc) as [number, number, number], force: h.force.map(v => v * sc) as [number, number, number] }));
  const big = simulate({ ...sys, hits });
  const small = simulate({ ...sys, hits: hits.map(h => ({ ...h, velocity: h.velocity.map(v => v * 1e-3) as [number, number, number], force: h.force.map(v => v * 1e-3) as [number, number, number] })) });
  const lin = small.audio.map(v => v * 1e3);
  const fs = big.sampleRate;
  const win = [0.01, 0.03, 0.1, 0.3].filter(t => t <= sys.duration).map(t => `${t}s:${relDiff(lin, big.audio, Math.round(t * fs)).toFixed(3)}`);
  // spectral centroid, first 50 ms
  const cen = (x: Float32Array) => { const s = spectrum(x, 0, Math.round(0.05 * fs)); let a = 0, b = 0; s.db.forEach((d, i) => { const p = 10 ** (d / 10); a += p * i * s.df * fs; b += p; }); return a / b; };
  console.log(`scale ${sc}: dt ${big.dt.toExponential(2)} drift ${big.energyDrift.toExponential(1)} ${big.aborted ?? ""} | nonlin diff ${win.join(" ")} | centroid lin ${cen(lin).toFixed(0)} Hz, actual ${cen(big.audio).toFixed(0)} Hz | ${(big.elapsedMs/1000).toFixed(1)}s`);
}