% Render a cymbal crash with default parameters and print diagnostics.
% Set P.strike_velocity = 0.5 for a tap on the same instrument, or
% P.support_clamped = 1 to hear the rigid center clamp of the earlier
% literature instead of the washer. P.strike_count = 4 plays a repeated
% figure; the hand retracts between strokes, which the paper's
% constant-acceleration hand cannot do.
% Run with: numbl run cymbal2_demo.m

P = cymbal2_defaults();
tic;
[y, dinfo] = cymbal2_render(P);
elapsed = toc;

fprintf('grid       : %d radial x %d azimuthal, %.1f kHz simulation rate\n', ...
  P.n_radial, P.n_theta, dinfo.simRate / 1000);
fprintf('shell      : a = %.0f cm, rise %.0f mm -> Rs = %.2f m, washer at %.1f cm\n', ...
  P.radius_m * 100, P.rise_m * 1000, dinfo.RsM, P.washer_radius_m * 100);
f0 = cymbal2_modes(P, 0); f0 = f0(f0 > 1);
f1 = cymbal2_modes(P, 1); f1 = f1(f1 > 1);
f2 = cymbal2_modes(P, 2); f2 = f2(f2 > 1);
fprintf('low modes  : %.1f Hz (m=1 tilt), %.1f Hz (m=2), %.1f Hz (m=0 bounce)\n', ...
  f1(1), f2(1), f0(1));
fprintf('stroke     : hand %.1f m/s -> tip %.2f m/s at contact\n', ...
  dinfo.meanV, abs(dinfo.contactVel(1)));
fprintf('contact    : %.0f N peak, %.2f mm compression, %.1f ms over %d event(s)\n', ...
  dinfo.peakForceN, dinfo.peakXiMm, dinfo.contactMs, dinfo.contactEvents);
fprintf('shell peak : %.2f mm displacement = %.1f thicknesses\n', ...
  dinfo.peakDispM * 1000, dinfo.peakDispThickness);
fprintf('pickup     : %.2e m/s^2 peak acceleration at r = %.1f cm\n', ...
  dinfo.peakAccel, dinfo.pickupRadiusM * 100);
fprintf('length     : %.2f s @ %d Hz\n', dinfo.lengthS, P.sample_rate);
fprintf('render time: %.1f s\n', elapsed);