Substation/DrumDetector.cs

namespace Substation;

public class DrumDetector
{
    const int HistorySize = 50;

    readonly double[] _lowFluxHistory = new double[HistorySize];
    readonly double[] _midFluxHistory = new double[HistorySize];
    readonly double[] _highFluxHistory = new double[HistorySize];
    int _historyIndex;

    double[]? _prevLowMag;
    double[]? _prevMidMag;
    double[]? _prevHighMag;

    readonly bool[] _subKick = new bool[4];
    readonly bool[] _subSnare = new bool[4];
    readonly bool[] _subBrass = new bool[4];

    static readonly byte[] FreqKick = Freq.Compress4(new[] { 150, 150, 150, 150 });
    static readonly byte[] FreqSnare = Freq.Compress4(new[] { 50, 50, 50, 50 });
    static readonly byte[] FreqBrass = Freq.Compress4(new[] { 10, 10, 10, 10 });
    static readonly byte[] FreqSilent = Freq.Compress4(new[] { 10, 10, 10, 10 });

    public void ProcessWindow(double[] magnitude, int sampleRate, int windowSize, int subTickIndex)
    {
        if (subTickIndex < 0 || subTickIndex > 3) return;

        double binWidth = (double)sampleRate / windowSize;
        int lowLo = (int)(20 / binWidth), lowHi = (int)(150 / binWidth);
        int midLo = (int)(200 / binWidth), midHi = (int)(500 / binWidth);
        int highLo = (int)(5000 / binWidth), highHi = (int)(15000 / binWidth);

        var lowMag = ExtractBand(magnitude, lowLo, lowHi);
        var midMag = ExtractBand(magnitude, midLo, midHi);
        var highMag = ExtractBand(magnitude, highLo, highHi);

        double lowFlux = ComputeFlux(lowMag, _prevLowMag);
        double midFlux = ComputeFlux(midMag, _prevMidMag);
        double highFlux = ComputeFlux(highMag, _prevHighMag);

        _prevLowMag = lowMag;
        _prevMidMag = midMag;
        _prevHighMag = highMag;

        // Update history and compute adaptive thresholds
        _lowFluxHistory[_historyIndex] = lowFlux;
        _midFluxHistory[_historyIndex] = midFlux;
        _highFluxHistory[_historyIndex] = highFlux;
        _historyIndex = (_historyIndex + 1) % HistorySize;

        double lowThresh = AdaptiveThreshold(_lowFluxHistory);
        double midThresh = AdaptiveThreshold(_midFluxHistory);
        double highThresh = AdaptiveThreshold(_highFluxHistory);

        // Detect onsets
        bool kickOnset = lowFlux > lowThresh && lowFlux > 0.001;
        bool snareOnset = midFlux > midThresh && highFlux > highThresh && midFlux > 0.001;
        bool brassOnset = highFlux > highThresh && !snareOnset && highFlux > 0.001;

        if (kickOnset) _subKick[subTickIndex] = true;
        if (snareOnset) _subSnare[subTickIndex] = true;
        if (brassOnset) _subBrass[subTickIndex] = true;
    }

    public WaveFrame BuildFrame()
    {
        var intensity = new byte[4];
        var freqBytes = new byte[16]; // 4 sub-ticks × 4 bytes (but we use per-sub-tick freq)

        for (int i = 0; i < 4; i++)
        {
            byte freq;
            if (_subKick[i])
            {
                intensity[i] = 100;
                freq = Freq.Compress(150);
            }
            else if (_subSnare[i])
            {
                intensity[i] = 100;
                freq = Freq.Compress(50);
            }
            else if (_subBrass[i])
            {
                intensity[i] = 100;
                freq = Freq.Compress(10);
            }
            else
            {
                intensity[i] = 0;
                freq = 10;
            }

            freqBytes[i] = freq;
        }

        // Reset for next tick
        Array.Clear(_subKick, 0, 4);
        Array.Clear(_subSnare, 0, 4);
        Array.Clear(_subBrass, 0, 4);

        return new WaveFrame(freqBytes, intensity);
    }

    static double[] ExtractBand(double[] magnitude, int loBin, int hiBin)
    {
        if (hiBin <= loBin || hiBin >= magnitude.Length)
            return Array.Empty<double>();

        var band = new double[hiBin - loBin + 1];
        Array.Copy(magnitude, loBin, band, 0, band.Length);
        return band;
    }

    static double ComputeFlux(double[] current, double[]? previous)
    {
        if (previous == null || current.Length != previous.Length) return 0;

        double flux = 0;
        for (int i = 0; i < current.Length; i++)
        {
            double diff = current[i] - previous[i];
            if (diff > 0) flux += diff;
        }
        return flux;
    }

    static double AdaptiveThreshold(double[] history)
    {
        double sum = 0, sumSq = 0;
        int count = 0;
        foreach (var v in history)
        {
            if (v <= 0) continue;
            sum += v;
            sumSq += v * v;
            count++;
        }
        if (count == 0) return double.MaxValue;
        double mean = sum / count;
        double variance = sumSq / count - mean * mean;
        double stddev = variance > 0 ? Math.Sqrt(variance) : 0;
        return mean + 2 * stddev;
    }
}