/ proof-of-concept / audio-spike-sorting
Sign in
proof-of-concept / audio-spike-sorting
audio-spike-sorting / src / algorithms / naive_threshold.m
47 lines · 1.3 KBBlameHistoryRaw
1function [times, labels] = naive_threshold(x, fs, training)
2% Deliberately weak floor: detect where the smoothed amplitude envelope
3% peaks above a threshold, and classify every event by its duration alone
4% (time the envelope stays above a fraction of the peak), nearest to the
5% mean duration of each training type. Scores from this baseline are the
6% floor real methods should clear.
7nTypes = numel(training.snips);
8L = size(training.snips{1}, 1);
10w = round(0.005 * fs);
11kern = ones(w, 1) / w;
12env = conv(abs(x), kern, 'same');
13n = numel(env);
14thr = max(8 * median(env), 0.05 * max(env));
15refr = round(L / 2);
17% mean above-threshold duration per type, from the training snippets
18durs = zeros(nTypes, 1);
19for k = 1:nTypes
20 S = training.snips{k};
21 dk = 0;
22 for j = 1:size(S, 2)
23 ev = conv(abs(S(:, j)), kern, 'same');
24 dk = dk + sum(ev > 0.2 * max(ev));
25 end
26 durs(k) = dk / max(size(S, 2), 1);
27end
29times = zeros(0, 1);
30labels = zeros(0, 1);
31i = 2;
32while i < n
33 if env(i) > thr && env(i-1) <= thr
34 j2 = min(n, i + refr - 1);
35 [pk, rel] = max(env(i:j2));
36 p = i + rel - 1;
37 j3 = min(n, i + L - 1);
38 dur = sum(env(i:j3) > 0.2 * pk);
39 [~, lab] = min(abs(durs - dur));
40 times(end+1, 1) = (i - 1) / fs;
41 labels(end+1, 1) = lab;
42 i = p + refr;
43 else
44 i = i + 1;
45 end
46end
47end
↑↓ move↵openescclose