Wykrywanie uderzeń i FFT

13

Pracuję nad platformówką, która zawiera muzykę z wykrywaniem beatów. Obecnie wykrywam uderzenia, sprawdzając, kiedy amplituda prądu przekracza próbkę historyczną. Nie działa to dobrze z gatunkami muzyki, takimi jak rock, które mają dość stałą amplitudę.

Spojrzałem więc dalej i znalazłem algorytmy dzielące dźwięk na wiele pasm za pomocą FFT ... potem znalazłem algorytm Cooley-Tukey FFt

Jedynym problemem, jaki mam, jest to, że jestem całkiem nowy w dziedzinie audio i nie mam pojęcia, jak go użyć, aby podzielić sygnał na wiele sygnałów.

Więc moje pytanie brzmi:

Jak używasz FFT do dzielenia sygnału na wiele pasm?

Dla zainteresowanych, to jest mój algorytm w języku c #:

// C = threshold, N = size of history buffer / 1024
    public void PlaceBeatMarkers(float C, int N)
    {
        List<float> instantEnergyList = new List<float>();
        short[] samples = soundData.Samples;

        float timePerSample = 1 / (float)soundData.SampleRate;
        int sampleIndex = 0;
        int nextSamples = 1024;

        // Calculate instant energy for every 1024 samples.
        while (sampleIndex + nextSamples < samples.Length)
        {

            float instantEnergy = 0;

            for (int i = 0; i < nextSamples; i++)
            {
                instantEnergy += Math.Abs((float)samples[sampleIndex + i]);
            }

            instantEnergy /= nextSamples;
            instantEnergyList.Add(instantEnergy);

            if(sampleIndex + nextSamples >= samples.Length)
                nextSamples = samples.Length - sampleIndex - 1;

            sampleIndex += nextSamples;
        }


        int index = N;
        int numInBuffer = index;
        float historyBuffer = 0;

        //Fill the history buffer with n * instant energy
        for (int i = 0; i < index; i++)
        {
            historyBuffer += instantEnergyList[i];
        }

        // If instantEnergy / samples in buffer < instantEnergy for the next sample then add beatmarker.
        while (index + 1 < instantEnergyList.Count)
        {
            if(instantEnergyList[index + 1] > (historyBuffer / numInBuffer) * C)
                beatMarkers.Add((index + 1) * 1024 * timePerSample); 
            historyBuffer -= instantEnergyList[index - numInBuffer];
            historyBuffer += instantEnergyList[index + 1];
            index++;
        }
    }
Quincy
źródło
Myślę, że dobrym punktem wyjścia są wpisy FFT i DSP wikipedii . Wpis dotyczący wykrywania uderzeń jest rzadki, ale zawiera linki do artykułu na gamedev.net
Tobias Kienzler

Odpowiedzi:

14

Cóż, jeśli twój sygnał wejściowy jest prawdziwy (jak w, każda próbka jest liczbą rzeczywistą), widmo będzie symetryczne i złożone. Wykorzystując symetrię, zwykle algorytmy FFT pakują wynik, zwracając tylko dodatnią połowę widma. Prawdziwa część każdego pasma znajduje się w próbkach parzystych, a część wyobrażona w próbkach nieparzystych. A czasem prawdziwe części są pakowane razem w pierwszej połowie odpowiedzi, a części urojone w drugiej połowie.

We wzorach, jeśli X [k] = FFT (x [n]), dajesz mu wektor i [n] = x [n] i otrzymujesz wynik o [m], a następnie

X[k] = o[2k] + j·o[2k+1]

(chociaż czasami otrzymujesz X [k] = o [k] + j · o [k + K / 2], gdzie K jest długością twojego okna, 1024 w twoim przykładzie). Nawiasem mówiąc, j jest urojoną jednostką sqrt (-1).

Wielkość prążka oblicza się jako pierwiastek iloczynu tego pasma z jego złożonym koniugatem:

|X[k]| = sqrt( X[k] · X[k]* )

A energia jest zdefiniowana jako kwadrat wielkości.

Jeśli nazwiemy a = o [2k] ib = o [2k + 1], otrzymamy

X[k] = a + j·b

w związku z tym

E[k] = |X[k]|^2 = (a+j·b)·(a-j·b) = a·a + b·b

Rozwijając całość, jeśli otrzymałeś o [m] jako wynik z algorytmu FFT, energia w paśmie k wynosi:

E[k] = o[2k] · o[2k] + o[2k+1] · o[2k+1]

(Uwaga: użyłem symbolu ·, aby wskazać mnożenie zamiast zwykłego *, aby uniknąć pomyłki z operatorem koniugacji)

Częstotliwość pasma k, przy założeniu częstotliwości próbkowania 44,1 Khz i okna 1024 próbek, wynosi

freq(k) = k / 1024 * 44100 [Hz]

Na przykład, twoje pierwsze pasmo k = 0 oznacza 0 Hz, k = 1 to 43 Hz, a ostatnie k = 511 wynosi 22 kHz (częstotliwość Nyquista).

Mam nadzieję, że to odpowiada na twoje pytanie o to, jak uzyskać energię sygnału na pasmo za pomocą FFT.

Dodatek : Odpowiadając na twoje pytanie w komentarzu i zakładając, że używasz kodu z linku zamieszczonego w pytaniu (algorytm Cooley-Tukey w C): Powiedzmy, że masz dane wejściowe jako wektor krótkich ints:

// len is 1024 in this example.  It MUST be a power of 2
// centerFreq is given in Hz, for example 43.0
double EnergyForBand( short *input, int len, double centerFreq)
{
  int i;
  int band;
  complex *xin;
  complex *xout;
  double magnitude;
  double samplingFreq = 44100.0; 

  // 1. Get the input as a vector of complex samples
  xin = (complex *)malloc(sizeof(struct complex_t) * len);

  for (i=0;i<len;i++) {
    xin[i].re = (double)input[i];
    xin[i].im = 0;
  }

  // 2. Transform the signal
  xout = FFT_simple(xin, len);

  // 3. Find the band ( Note: floor(x+0.5) = round(x) )
  band = (int) floor(centerFreq * len / samplingFreq + 0.5); 

  // 4. Get the magnitude
  magnitude = complex_magnitude( xout[band] );

  // 5. Don't leak memory
  free( xin );
  free( xout );

  // 6. Return energy
  return magnitude * magnitude;
}

Moje C jest trochę zardzewiałe (obecnie głównie koduję w C ++), ale mam nadzieję, że nie popełniłem żadnego dużego błędu w tym kodzie. Oczywiście, jeśli interesuje Cię energia innych pasm, nie ma sensu przekształcać całego okna dla każdego z nich, co byłoby stratą czasu procesora. W takim przypadku wykonaj raz transformację i uzyskaj wszystkie potrzebne wartości od xout.

CeeJay
źródło
Och, właśnie spojrzałem na kod, który podłączyłeś, już daje wyniki w postaci „złożonej”, a nawet zapewnia funkcję obliczania wielkości liczby zespolonej. Wtedy musiałbyś tylko obliczyć kwadrat tej wielkości dla każdego elementu wektora wyjściowego, nie musisz się martwić sortowaniem wyników.
CeeJay,
Na przykład, jeśli mam wszystkie 1024 próbki z okna 0-1024 i otrzymałem je jako wartości rzeczywiste, więc nie mam złożonej części. i chcę obliczyć tam energię w paśmie częstotliwości 43 Hz. Jak w takim razie miałbym to zintegrować? (Potrzebuję tylko prawdziwej części z powrotem, części dodatniej) Jeśli mógłbyś to zrobić w jakimś pseudokodzie, będę dogłębnie w tobie na zawsze, a wtedy rzeczywiście mogę trochę zrozumieć pojęcie :)
Quincy
Kod, który napisałem, korzysta z biblioteki C, którą podłączyłeś, która zawiera już „złożoną” strukturę. To sprawia, że ​​rozpakowywanie, które opisałem w moim pytaniu, jest niepotrzebne (a kod to odzwierciedla)
CeeJay
0

Sam tego nie zrobiłem ani nie czytałem dużo na ten temat, ale moje pierwsze ujęcie wygląda mniej więcej tak:

Przede wszystkim musisz zastosować funkcję okna, aby uzyskać widmo zależne od czasu za pomocą FFT. Bicie zwykle leży w niższych częstotliwościach, więc zastosuj kolejny FFT z większym oknem czasowym na intensywności niektórych z tych częstotliwości (dla uproszczenia zacznij od tylko 1 przy np. 100 Hz i sprawdź, czy to jest wystarczająco niezawodne). Znajdź szczyt w tym spektrum, a ta częstotliwość jest przypuszczeniem rytmu.

Tobias Kienzler
źródło
Nie mam właściwie problemu z wykrywaniem rytmu, ale rozumiem, jak działa FFT. Jestem naprawdę nowy w przetwarzaniu sygnałów i takie rzeczy jak: „zastosuj funkcję okna, aby uzyskać widmo zależne od czasu z FFT”, nie ma dla mnie żadnego sensu. W każdym razie dzięki :)
Quincy