function [tstr, vstr, jstr, taustr, fstr, sstr, sinfo] = triangle_strikes(P)
% TRIANGLE_STRIKES Beater strike events: a single strike, or a seeded roll.
%
% strike_count = 1 gives one deterministic strike of strike_velocity at
% strike_pos at t = 0; no random numbers are drawn. strike_count > 1 gives a
% roll: nominally periodic hits at strike_rate_hz, alternating between
% strike_pos and roll_pos2 (the beater rattling in a corner strikes the two
% legs that form it in turn), with Gaussian timing jitter (timing_jitter_ms)
% and lognormal velocity spread (velocity_spread). The first hit stays at
% t = 0 so the onset is fixed. All randomness comes from triangle_rng seeded
% with P.seed, so a given seed reproduces the identical roll in the CLI,
% Node, and the browser.
%
% Each strike is a Hertz contact of the beater tip (a steel sphere of radius
% beater_radius_m) on the rod:
%   k_h  = (4/3) E* sqrt(R_eff),  R_eff = sqrt(r_b R''),
%          1/R'' = 1/r_b + 2/d      (sphere on cylinder: R_eff is the
%                                    geometric mean of the two relative
%                                    principal curvatures)
%   tau  = 3.218 (m*^2 / (k_h^2 v))^{1/5}    (Hertz contact time)
%   J    = (1 + restitution) m* v            (momentum transferred)
%   Fpk  = pi J / (2 tau)                    (half-sine pulse peak)
% The reduced mass m* pairs the beater with the rod mass engaged during the
% contact. The engaged length is estimated from the flexural group velocity
% at the pulse's own bandwidth, l = 2 sqrt(pi kappa cL tau), which depends
% on tau, so a short fixed-point iteration closes the loop. This is an
% estimate, not a contact simulation: what it captures is that a heavy
% beater lengthens the contact and darkens the strike, while a light beater
% approaches the free-sphere limit and tau ~ v^(-1/5).
%
% Returns row vectors sorted by time: tstr (s), vstr (m/s), jstr (N s),
% taustr (s), fstr (peak N), sstr (arclength fraction of the contact), and
% sinfo with summary statistics.

K = max(1, round(P.strike_count));
rate = max(1e-3, P.strike_rate_hz);
duration = (K - 1) / rate + P.tail_s;

E = P.youngs_modulus;
rho = P.density;
d = P.rod_diameter_m;
nub = 0.30;                          % beater Poisson ratio, fixed
Estar = 1 / ((1 - nub^2) / P.beater_youngs + (1 - P.poisson_ratio^2) / E);
rb = P.beater_radius_m;
Rpp = 1 / (1 / rb + 2 / d);
Reff = sqrt(rb * Rpp);
kh = (4 / 3) * Estar * sqrt(Reff);
kapcL = (d / 4) * sqrt(E / rho);     % kappa * cL: flexural dispersion constant
Aarea = pi * d^2 / 4;

state = mod(floor(abs(P.seed)), 4294967296);
z1 = zeros(K, 1);
z2 = zeros(K, 1);
if K > 1
  [u, state] = triangle_rng(state, 4 * K);
  U = reshape(u, K, 4);
  z1 = sqrt(-2 * log(U(:, 1))) .* cos(2 * pi * U(:, 2));
  z2 = sqrt(-2 * log(U(:, 3))) .* cos(2 * pi * U(:, 4));
end

tstr = zeros(1, K);
vstr = zeros(1, K);
jstr = zeros(1, K);
taustr = zeros(1, K);
fstr = zeros(1, K);
sstr = zeros(1, K);
for i = 1:K
  t = (i - 1) / rate;
  if i > 1
    t = max(0, t + P.timing_jitter_ms * 1e-3 * z1(i));
  end
  v = P.strike_velocity;
  if K > 1
    v = v * exp(P.velocity_spread * z2(i));
  end
  v = max(0.01, v);

  % engaged-mass fixed point for the Hertz contact time
  tau = 30e-6;
  mstar = P.beater_mass_kg;
  for it = 1:4
    leng = 2 * sqrt(pi * kapcL * tau);
    mloc = rho * Aarea * leng;
    mstar = P.beater_mass_kg * mloc / (P.beater_mass_kg + mloc);
    tau = 3.218 * (mstar^2 / (kh^2 * v))^(1 / 5);
  end

  tstr(i) = t;
  vstr(i) = v;
  jstr(i) = (1 + P.restitution) * mstar * v;
  taustr(i) = tau;
  fstr(i) = pi * jstr(i) / (2 * tau);
  if mod(i, 2) == 1
    sstr(i) = P.strike_pos;
  else
    sstr(i) = P.roll_pos2;
  end
end

[tstr, order] = sort(tstr);
tstr = reshape(tstr, 1, K);
vstr = reshape(vstr(order), 1, K);
jstr = reshape(jstr(order), 1, K);
taustr = reshape(taustr(order), 1, K);
fstr = reshape(fstr(order), 1, K);
sstr = reshape(sstr(order), 1, K);

sinfo.count = K;
sinfo.meanV = sum(vstr) / K;
sinfo.maxV = max(vstr);
sinfo.meanTauUs = sum(taustr) / K * 1e6;
sinfo.peakForceN = max(fstr);
sinfo.durationS = duration;
end