git.sesse.net Git - kdenlive/blob - src/audioscopes/ffttools.cpp

   1 /***************************************************************************
   2  *   Copyright (C) 2010 by Simon Andreas Eugster (simon.eu@gmail.com)      *
   3  *   This file is part of kdenlive. See www.kdenlive.org.                  *
   4  *                                                                         *
   5  *   This program is free software; you can redistribute it and/or modify  *
   6  *   it under the terms of the GNU General Public License as published by  *
   7  *   the Free Software Foundation; either version 2 of the License, or     *
   8  *   (at your option) any later version.                                   *
   9  ***************************************************************************/
  10
  11 #include <math.h>
  12 #include <iostream>
  13
  14 #include <QString>
  15
  16 #include "ffttools.h"
  17
  18 //#define DEBUG_FFTTOOLS
  19 #ifdef DEBUG_FFTTOOLS
  20 #include <QDebug>
  21 #include <QTime>
  22 #endif
  23
  24 FFTTools::FFTTools() :
  25         m_fftCfgs(),
  26         m_windowFunctions()
  27 {
  28 }
  29 FFTTools::~FFTTools()
  30 {
  31     QHash<QString, kiss_fftr_cfg>::iterator i;
  32     for (i = m_fftCfgs.begin(); i != m_fftCfgs.end(); i++) {
  33         free(*i);
  34     }
  35 }
  36
  37 const QString FFTTools::windowSignature(const WindowType windowType, const int size, const float param)
  38 {
  39     return QString("s%1_t%2_p%3").arg(size).arg(windowType).arg(param, 0, 'f', 3);
  40 }
  41 const QString FFTTools::cfgSignature(const int size)
  42 {
  43     return QString("s%1").arg(size);
  44 }
  45
  46 // http://cplusplus.syntaxerrors.info/index.php?title=Cannot_declare_member_function_%E2%80%98static_int_Foo::bar%28%29%E2%80%99_to_have_static_linkage
  47 const QVector<float> FFTTools::window(const WindowType windowType, const int size, const float param)
  48 {
  49     Q_ASSERT(size > 0);
  50
  51     // Deliberately avoid converting size to a float
  52     // to keep mid an integer.
  53     float mid = (size-1)/2;
  54     float max = size-1;
  55     QVector<float> window;
  56
  57     switch (windowType) {
  58     case Window_Rect:
  59         return QVector<float>(size+1, 1);
  60         break;
  61     case Window_Triangle:
  62         window = QVector<float>(size+1);
  63
  64         for (int x = 0; x < mid; x++) {
  65             window[x] = x/mid + (mid-x)/mid*param;
  66         }
  67         for (int x = mid; x < size; x++) {
  68             window[x] = (x-mid)/(max-mid) * param + (max-x)/(max-mid);
  69         }
  70         window[size] = .5 + param/2;
  71
  72 #ifdef DEBUG_FFTTOOLS
  73         qDebug() << "Triangle window (factor " << window[size] << "):";
  74         for (int i = 0; i < size; i++) {
  75             qDebug() << window[i];
  76         }
  77         qDebug() << "Triangle window end.";
  78 #endif
  79
  80         return window;
  81         break;
  82     case Window_Hamming:
  83         // Use a quick version of the Hamming window here: Instead of
  84         // interpolating values between (-max/2) and (max/2)
  85         // we use integer values instead, ranging from -mid to (max-mid).
  86         window = QVector<float>(size+1);
  87
  88         for (int x = 0; x < size; x++) {
  89             window[x] = .54 + .46 * cos( 2*M_PI*(x-mid) / size );
  90         }
  91
  92         // Integrating the cosine over the window function results in
  93         // an area of 0; So only the constant factor 0.54 counts.
  94         window[size] = .54;
  95
  96 #ifdef DEBUG_FFTTOOLS
  97         qDebug() << "Hanning window (factor " << window[size] << "):";
  98         for (int i = 0; i < size; i++) {
  99             qDebug() << window[i];
 100         }
 101         qDebug() << "Hanning window end.";
 102 #endif
 103
 104         return window;
 105         break;
 106     }
 107     Q_ASSERT(false);
 108     return QVector<float>();
 109 }
 110
 111 void FFTTools::fftNormalized(const QVector<int16_t> audioFrame, const uint channel, const uint numChannels, float *freqSpectrum,
 112                              const WindowType windowType, const uint windowSize, const float param)
 113 {
 114 #ifdef DEBUG_FFTTOOLS
 115     QTime start = QTime::currentTime();
 116 #endif
 117
 118     const uint numSamples = audioFrame.size()/numChannels;
 119
 120     Q_ASSERT((windowSize & 1) == 0);
 121     Q_ASSERT(windowSize > 0);
 122
 123     const QString cfgSig = cfgSignature(windowSize);
 124     const QString winSig = windowSignature(windowType, windowSize, param);
 125
 126
 127     // Get the kiss_fft configuration from the config cache
 128     // or build a new configuration if the requested one is not available.
 129     kiss_fftr_cfg myCfg;
 130     if (m_fftCfgs.contains(cfgSig)) {
 131 #ifdef DEBUG_FFTTOOLS
 132         qDebug() << "Re-using FFT configuration with size " << windowSize;
 133 #endif
 134         myCfg = m_fftCfgs.value(cfgSig);
 135     } else {
 136 #ifdef DEBUG_FFTTOOLS
 137         qDebug() << "Creating FFT configuration with size " << windowSize;
 138 #endif
 139         myCfg = kiss_fftr_alloc(windowSize, 0,0,0);
 140         m_fftCfgs.insert(cfgSig, myCfg);
 141     }
 142
 143     // Get the window function from the cache
 144     // (except for a rectangular window; nothing to to there.
 145     QVector<float> window;
 146     float windowScaleFactor = 1;
 147     if (windowType != FFTTools::Window_Rect) {
 148
 149         if (m_windowFunctions.contains(winSig)) {
 150 #ifdef DEBUG_FFTTOOLS
 151             qDebug() << "Re-using window function with signature " << winSig;
 152 #endif
 153             window = m_windowFunctions.value(winSig);
 154         } else {
 155 #ifdef DEBUG_FFTTOOLS
 156             qDebug() << "Building new window function with signature " << winSig;
 157 #endif
 158             window = FFTTools::window(windowType, windowSize, 0);
 159             m_windowFunctions.insert(winSig, window);
 160         }
 161         windowScaleFactor = 1.0/window[windowSize];
 162     }
 163
 164
 165     // Prepare frequency space vector. The resulting FFT vector is only half as long.
 166     kiss_fft_cpx freqData[windowSize/2];
 167     float data[windowSize];
 168
 169     // Copy the first channel's audio into a vector for the FFT display;
 170     // Fill the data vector indices that cannot be covered with sample data with 0
 171     if (numSamples < windowSize) {
 172         std::fill(&data[numSamples], &data[windowSize-1], 0);
 173     }
 174     // Normalize signals to [0,1] to get correct dB values later on
 175     for (int i = 0; i < numSamples && i < windowSize; i++) {
 176         // Performance note: Benchmarking has shown that using the if/else inside the loop
 177         // does not do noticeable worse than keeping it outside (perhaps the branch predictor
 178         // is good enough), so it remains in there for better readability.
 179         if (windowType != FFTTools::Window_Rect) {
 180             data[i] = (float) audioFrame.data()[i*numChannels + channel] / 32767.0f * window[i];
 181         } else {
 182             data[i] = (float) audioFrame.data()[i*numChannels + channel] / 32767.0f;
 183         }
 184     }
 185
 186     // Calculate the Fast Fourier Transform for the input data
 187     kiss_fftr(myCfg, data, freqData);
 188
 189
 190     // Logarithmic scale: 20 * log ( 2 * magnitude / N ) with magnitude = sqrt(r² + i²)
 191     // with N = FFT size (after FFT, 1/2 window size)
 192     for (int i = 0; i < windowSize/2; i++) {
 193         // Logarithmic scale: 20 * log ( 2 * magnitude / N ) with magnitude = sqrt(r² + i²)
 194         // with N = FFT size (after FFT, 1/2 window size)
 195         freqSpectrum[i] = 20*log(pow(pow(fabs(freqData[i].r * windowScaleFactor),2) + pow(fabs(freqData[i].i * windowScaleFactor),2), .5)/((float)windowSize/2.0f))/log(10);;
 196     }
 197
 198 #ifdef DEBUG_FFTTOOLS
 199     qDebug() << "Calculated FFT in " << start.elapsed() << " ms.";
 200 #endif
 201 }
 202
 203 #ifdef DEBUG_FFTTOOLS
 204 #undef DEBUG_FFTTOOLS
 205 #endif