// The contact law and what it does on impact. Left: f = K [xi]+^alpha, the
// Hertz-type compression law of the mallet head. Right: the head thrown at a
// small plate-like target, five modal oscillators sharing the contact point,
// integrated through and past the contact. The same simulation renders a
// second of the plate's response as audio, so mallet hardness can be heard,
// not just plotted: a soft head excites only the low modes, a hard head all
// of them.

import { useMemo, useState } from 'react';
import { useCanvas, drawFrame, PLOT_COLORS } from '../../components/plots/canvasUtils';
import { Figure, FigControls, FigSlider, usePlayer } from '../figparts';

const HEIGHT = 260;
const M_HEAD = 0.028;
// A plate-like ladder of point-coupled modes: frequency (Hz), modal mass
// (kg), and amplitude decay rate (1/s).
const MODES = [
  { f: 180, m: 0.35, g: 2.5 },
  { f: 476, m: 0.22, g: 3.5 },
  { f: 894, m: 0.16, g: 5.0 },
  { f: 1520, m: 0.12, g: 7.0 },
  { f: 2350, m: 0.10, g: 10.0 },
];
const PLOT_T = 0.012;
const AUDIO_SR = 25000;
const DT = 1e-6; // 40 substeps per audio sample

interface StrikeState {
  x1: number; v1: number;
  q: Float64Array; qd: Float64Array;
}

function makeState(v0: number): StrikeState {
  return { x1: -1e-5, v1: v0, q: new Float64Array(MODES.length), qd: new Float64Array(MODES.length) };
}

// One symplectic-Euler substep; returns the contact force and the summed
// modal acceleration (the collocated pickup).
function step(s: StrikeState, K: number, alpha: number): { fc: number; acc: number } {
  let ub = 0;
  for (let i = 0; i < MODES.length; i++) ub += s.q[i];
  const xi = s.x1 - ub;
  const fc = xi > 0 ? K * Math.pow(xi, alpha) : 0;
  s.v1 += (-fc / M_HEAD) * DT;
  s.x1 += s.v1 * DT;
  let acc = 0;
  for (let i = 0; i < MODES.length; i++) {
    const md = MODES[i];
    const w = 2 * Math.PI * md.f;
    const ai = (fc - md.m * w * w * s.q[i] - 2 * md.m * md.g * s.qd[i]) / md.m;
    s.qd[i] += ai * DT;
    s.q[i] += s.qd[i] * DT;
    acc += ai;
  }
  return { fc, acc };
}

function simulate(K: number, alpha: number, v0: number) {
  const s = makeState(v0);
  const n = Math.round(PLOT_T / DT);
  const keep = 600;
  const stride = Math.max(1, Math.floor(n / keep));
  const t: number[] = [];
  const f: number[] = [];
  let peakF = 0;
  let peakXi = 0;
  let contacts = 0;
  let inC = false;
  for (let k = 0; k < n; k++) {
    let ub = 0;
    for (let i = 0; i < MODES.length; i++) ub += s.q[i];
    const xi = s.x1 - ub;
    if (xi > peakXi) peakXi = xi;
    const { fc } = step(s, K, alpha);
    if (fc > 0 && !inC) contacts++;
    inC = fc > 0;
    if (fc > peakF) peakF = fc;
    if (k % stride === 0) {
      t.push(k * DT);
      f.push(fc);
    }
  }
  return { t, f, peakF, peakXi, contacts };
}

export function renderStrikeAudio(K: number, alpha: number, v0: number): Float32Array<ArrayBuffer> {
  const dur = 1.2;
  const nOut = Math.round(AUDIO_SR * dur);
  const sub = Math.round(1 / (AUDIO_SR * DT));
  const s = makeState(v0);
  const y = new Float32Array(nOut);
  for (let k = 0; k < nOut; k++) {
    let acc = 0;
    for (let j = 0; j < sub; j++) acc = step(s, K, alpha).acc;
    y[k] = acc;
  }
  let pk = 0;
  for (let k = 0; k < nOut; k++) pk = Math.max(pk, Math.abs(y[k]));
  const fadeOut = Math.round(0.1 * AUDIO_SR);
  for (let k = 0; k < nOut; k++) {
    let g = 0.5 / (pk || 1);
    if (k > nOut - fadeOut) g *= (nOut - k) / fadeOut;
    y[k] *= g;
  }
  return y;
}

export default function HertzFig() {
  const [alpha, setAlpha] = useState(2.54);
  const [K, setK] = useState(1.6e8);
  const [v0, setV0] = useState(5);
  const sim = useMemo(() => simulate(K, alpha, v0), [K, alpha, v0]);
  const play = usePlayer();

  const ref = useCanvas(
    (ctx, w, h) => {
      drawFrame(ctx, w, h);
      const pw = (w - 54) / 2;
      const ph = h - 62;
      const y0 = 14;

      // ---- left: force vs compression ----
      const xiMax = Math.max(1.5 * sim.peakXi, 1e-4);
      const fAtXiMax = K * Math.pow(xiMax, alpha);
      const x0 = 34;
      ctx.strokeStyle = PLOT_COLORS.grid;
      ctx.strokeRect(x0, y0, pw, ph);
      ctx.strokeStyle = PLOT_COLORS.trace;
      ctx.lineWidth = 2;
      ctx.beginPath();
      for (let i = 0; i <= 120; i++) {
        const xi = (i / 120) * xiMax;
        const fv = K * Math.pow(xi, alpha);
        const x = x0 + (xi / xiMax) * pw;
        const y = y0 + ph - (fv / fAtXiMax) * ph;
        if (i === 0) ctx.moveTo(x, y);
        else ctx.lineTo(x, y);
      }
      ctx.stroke();
      ctx.lineWidth = 1;
      // operating point
      if (sim.peakXi > 0) {
        ctx.fillStyle = PLOT_COLORS.accent;
        ctx.beginPath();
        ctx.arc(
          x0 + (sim.peakXi / xiMax) * pw,
          y0 + ph - (sim.peakF / fAtXiMax) * ph,
          4, 0, 2 * Math.PI,
        );
        ctx.fill();
      }
      ctx.fillStyle = PLOT_COLORS.axis;
      ctx.font = '11px system-ui';
      ctx.textAlign = 'center';
      ctx.fillText(
        `f = K·ξ^${alpha.toFixed(2)}   (compression 0…${(xiMax * 1000).toFixed(2)} mm)`,
        x0 + pw / 2,
        y0 + ph + 18,
      );
      ctx.textAlign = 'left';
      ctx.fillText('force vs compression', x0 + 6, y0 + 14);

      // ---- right: force vs time ----
      const x1p = x0 + pw + 20;
      const fMax = Math.max(1, sim.peakF) * 1.08;
      ctx.strokeStyle = PLOT_COLORS.grid;
      ctx.strokeRect(x1p, y0, pw, ph);
      ctx.strokeStyle = PLOT_COLORS.accent;
      ctx.lineWidth = 2;
      ctx.beginPath();
      for (let i = 0; i < sim.t.length; i++) {
        const x = x1p + (sim.t[i] / PLOT_T) * pw;
        const y = y0 + ph - (sim.f[i] / fMax) * ph;
        if (i === 0) ctx.moveTo(x, y);
        else ctx.lineTo(x, y);
      }
      ctx.stroke();
      ctx.lineWidth = 1;
      ctx.fillStyle = PLOT_COLORS.axis;
      ctx.textAlign = 'center';
      ctx.fillText(`force vs time (0…${(PLOT_T * 1000).toFixed(0)} ms)`, x1p + pw / 2, y0 + ph + 18);
      ctx.textAlign = 'left';
      ctx.fillText(
        `peak ${sim.peakF.toFixed(0)} N, ${sim.contacts} contact${sim.contacts === 1 ? '' : 's'}`,
        x1p + 6,
        y0 + 14,
      );
    },
    [sim, K, alpha],
    HEIGHT,
  );

  return (
    <Figure
      caption={
        <>
          Left: the contact law, a power law in the compression with exponent &alpha;. Hertz
          theory gives &alpha; = 1.5 for a hard sphere on a plate; the paper's measured felt
          mallet is much more strongly hardening, &alpha; = 2.54, so the felt is soft at
          grazing contact and stiffens as it packs. Right: the force the head actually
          delivers when thrown at a small plate-like target (five modes sharing the contact
          point). Soften the head and the pulse lengthens; stiffen it and the strike becomes
          a brief spike; at some settings the target rebounds into the head and the pulse
          splits, the double peak Samejima's Fig. 12 shows for the full stick model. The play
          button renders a second of the target's response to this exact strike: a long soft
          pulse cannot excite the upper modes (the thud of felt), a short hard one excites
          them all (the ping of a wood tip). The excitation spectrum is the pulse's Fourier
          transform, and you are hearing its width change.
        </>
      }
    >
      <canvas ref={ref} />
      <FigControls>
        <FigSlider
          label="exponent α"
          value={alpha}
          min={1.2}
          max={3.2}
          step={0.02}
          format={(v) => v.toFixed(2)}
          onChange={setAlpha}
        />
        <FigSlider
          label="stiffness K"
          value={K}
          min={1e6}
          max={1e10}
          log
          format={(v) => v.toExponential(1)}
          onChange={setK}
        />
        <FigSlider
          label="impact speed"
          value={v0}
          min={0.5}
          max={8}
          step={0.1}
          format={(v) => `${v.toFixed(1)} m/s`}
          onChange={setV0}
        />
        <button onClick={() => play(renderStrikeAudio(K, alpha, v0), AUDIO_SR)}>
          ▶ Hear this strike
        </button>
      </FigControls>
    </Figure>
  );
}