1function [times, labels] = onset_knn(x, fs, training)
2% Detection and classification kept separate: events are found on a
3% smoothed energy envelope (rise above a robust threshold), then each
4% event is classified by distance-weighted k-nearest-neighbours over
5% spectral-shape features of the training snippets. Features are level
6% ratios, so the random per-event gain carries no information.
7nTypes = numel(training.snips);
8L = size(training.snips{1}, 1);
9pre = training.pre;
11% ---- detection ----------------------------------------------------------
12w = round(0.004 * fs);
13env = conv(x .* x, ones(w, 1) / w, 'same');
14env = sqrt(env);
15n = numel(env);
16thr = max(8 * median(env), 0.05 * max(env));
17refr = round(L / 2);
19starts = zeros(0, 1);
20i = 2;
21while i < n
22 if env(i) > thr && env(i-1) <= thr
23 % crossing; find the local peak, then walk back to the foot of the rise
24 j2 = min(n, i + refr - 1);
25 [pk, rel] = max(env(i:j2));
26 p = i + rel - 1;
27 s = p;
28 while s > 1 && env(s-1) < env(s) && env(s) > 0.05 * pk
29 s = s - 1;
30 end
31 starts(end+1, 1) = s;
32 % hold off until the envelope falls well below threshold (hysteresis),
33 % so ripple in a long ringing tail is not reported again
34 q = p;
35 while q < n && env(q) > 0.5 * thr
36 q = q + 1;
37 end
38 i = max(p + refr, q);
39 else
40 i = i + 1;
41 end
42end
44% ---- training features --------------------------------------------------
45nTrain = 0;
46for k = 1:nTypes
47 nTrain = nTrain + size(training.snips{k}, 2);
48end
49F = zeros(nTrain, 6);
50Fl = zeros(nTrain, 1);
51r = 0;
52for k = 1:nTypes
53 S = training.snips{k};
54 for j = 1:size(S, 2)
55 r = r + 1;
56 F(r, :) = event_features(S(pre+1:end, j), fs);
57 Fl(r) = k;
58 end
59end
60mu = mean(F, 1);
61sd = std(F, 0, 1) + 1e-9;
62Fz = (F - repmat(mu, nTrain, 1)) ./ repmat(sd, nTrain, 1);
64% ---- classify each detected event --------------------------------------
65kNN = min(5, nTrain);
66times = zeros(0, 1);
67labels = zeros(0, 1);
68for e = 1:numel(starts)
69 s = starts(e);
70 seg = x(s:min(n, s + (L - pre) - 1));
71 f = (event_features(seg, fs) - mu) ./ sd;
72 d = sqrt(sum((Fz - repmat(f, nTrain, 1)).^2, 2));
73 [ds, di] = sort(d);
74 votes = zeros(nTypes, 1);
75 for v = 1:kNN
76 votes(Fl(di(v))) = votes(Fl(di(v))) + 1 / (ds(v) + 1e-6);
77 end
78 [~, lab] = max(votes);
79 times(end+1, 1) = (s - 1) / fs;
80 labels(end+1, 1) = lab;
81end
82end
84function f = event_features(seg, fs)
85% Spectral shape of the event: four log band-energy ratios (mean removed,
86% so absolute level cancels), log spectral centroid, and zero-crossing
87% rate over the first 25 ms.
88nfft = 1024;
89m = min(numel(seg), nfft);
90s2 = seg(1:m);
91win = 0.5 - 0.5 * cos(2 * pi * (0:m-1)' / max(m - 1, 1));
92S = fft(s2 .* win, nfft);
93p = abs(S(1:nfft/2)).^2;
94hz = fs / nfft;
95edges = [0 200 800 2500 8000];
96b = zeros(1, 4);
97for q = 1:4
98 k0 = max(2, floor(edges(q) / hz) + 1);
99 k1 = min(nfft/2, floor(edges(q+1) / hz));
100 b(q) = log10(sum(p(k0:k1)) + 1e-12);
101end
102b = b - mean(b);
103kv = (0:nfft/2-1)' * hz;
104cen = sum(p .* kv) / (sum(p) + 1e-12);
105zn = min(numel(seg), round(0.025 * fs));
106zc = sum(abs(diff(seg(1:zn) >= 0)));
107f = [b, log10(cen + 1), zc / max(zn, 1) * fs / 1000];
108end