OSC:: needs the .h file too...

[ardour.git] / gtk2_ardour / fft.cc

diff --git a/gtk2_ardour/fft.cc b/gtk2_ardour/fft.cc
index 92f52bb64e665d85a6455da208e14affb204fde7..a4e34bf2aafedbb56116c372772ef017c2320506 100644 (file)
--- a/gtk2_ardour/fft.cc
+++ b/gtk2_ardour/fft.cc
@@ -23,19 +23,22 @@
 #include <string.h>
 #include <math.h>
 
+using namespace GTKArdour;
+
 FFT::FFT(uint32_t windowSize)
-       : _windowSize(windowSize),
-         _dataSize(_windowSize/2),
-       _iterations(0)
+       : _window_size(windowSize),
+         _data_size(_window_size/2),
+       _iterations(0),
+       _hann_window(0)
 {
-       _fftInput  = (float *) fftwf_malloc(sizeof(float) * _windowSize);
+       _fftInput  = (float *) fftwf_malloc(sizeof(float) * _window_size);
 
-       _fftOutput = (float *) fftwf_malloc(sizeof(float) * _windowSize);
+       _fftOutput = (float *) fftwf_malloc(sizeof(float) * _window_size);
 
-       _powerAtBin  = (float *) malloc(sizeof(float) * _dataSize);
-       _phaseAtBin  = (float *) malloc(sizeof(float) * _dataSize);
+       _power_at_bin  = (float *) malloc(sizeof(float) * _data_size);
+       _phase_at_bin  = (float *) malloc(sizeof(float) * _data_size);
 
-       _plan = fftwf_plan_r2r_1d(_windowSize, _fftInput, _fftOutput, FFTW_R2HC, FFTW_ESTIMATE);
+       _plan = fftwf_plan_r2r_1d(_window_size, _fftInput, _fftOutput, FFTW_R2HC, FFTW_ESTIMATE);
 
        reset();
 }
@@ -43,42 +46,49 @@ FFT::FFT(uint32_t windowSize)
 void
 FFT::reset()
 {
-       memset(_powerAtBin, 0, sizeof(float) * _dataSize);
-       memset(_phaseAtBin, 0, sizeof(float) * _dataSize);
-       
+       memset(_power_at_bin, 0, sizeof(float) * _data_size);
+       memset(_phase_at_bin, 0, sizeof(float) * _data_size);
+
        _iterations = 0;
 }
 
 void
-FFT::analyze(ARDOUR::Sample *input)
+FFT::analyze(ARDOUR::Sample *input, WindowingType windowing_type)
 {
        _iterations++;
 
-       memcpy(_fftInput, input, sizeof(float) * _windowSize);
+       memcpy(_fftInput, input, sizeof(float) * _window_size);
+
+       if (windowing_type == HANN) {
+               float *window = get_hann_window();
+               for (uint32_t i = 0; i < _window_size; i++) {
+                       _fftInput[i] *= window[i];
+               }
+       }
 
        fftwf_execute(_plan);
 
-       _powerAtBin[0] += _fftOutput[0] * _fftOutput[0];
-       _phaseAtBin[0] += 0.0;
+       _power_at_bin[0] += _fftOutput[0] * _fftOutput[0];
+       _phase_at_bin[0] += 0.0;
 
        float power;
        float phase;
 
 #define Re (_fftOutput[i])
-#define Im (_fftOutput[_windowSize-i])
-               for (uint32_t i=1; i < _dataSize - 1; i++) { 
+#define Im (_fftOutput[_window_size-i])
+               for (uint32_t i=1; i < _data_size - 1; i++) {
 
                power = (Re * Re) + (Im * Im);
                phase = atanf(Im / Re);
-       
+
                if (Re < 0.0 && Im > 0.0) {
                        phase += M_PI;
                } else if (Re < 0.0 && Im < 0.0) {
                        phase -= M_PI;
                }
 
-               _powerAtBin[i] += power;
-               _phaseAtBin[i] += phase;
+               _power_at_bin[i] += power;
+               _phase_at_bin[i] += phase;
        }
 #undef Re
 #undef Im
@@ -88,20 +98,48 @@ void
 FFT::calculate()
 {
        if (_iterations > 1) {
-               for (uint32_t i=0; i < _dataSize - 1; i++) { 
-                       _powerAtBin[i] /= (float)_iterations;
-                       _phaseAtBin[i] /= (float)_iterations;
+               for (uint32_t i=0; i < _data_size - 1; i++) {
+                       _power_at_bin[i] /= (float)_iterations;
+                       _phase_at_bin[i] /= (float)_iterations;
                }
                _iterations = 1;
        }
 }
 
+float *
+FFT::get_hann_window()
+{
+       if (_hann_window)
+               return _hann_window;
+
+
+        _hann_window = (float *) malloc(sizeof(float) * _window_size);
+
+       double sum = 0.0;
+
+        for (uint32_t i=0; i < _window_size; i++) {
+                _hann_window[i]=0.81f * ( 0.5f - (0.5f * (float) cos(2.0f * M_PI * (float)i / (float)(_window_size))));
+                sum += _hann_window[i];
+        }
+
+        double isum = 1.0 / sum;
+
+        for (uint32_t i=0; i < _window_size; i++) {
+                _hann_window[i] *= isum;
+        }
+
+       return _hann_window;
+}
+
 
 FFT::~FFT()
 {
+       if (_hann_window) {
+               free(_hann_window);
+       }
        fftwf_destroy_plan(_plan);
-       free(_powerAtBin);
-       free(_phaseAtBin);
+       free(_power_at_bin);
+       free(_phase_at_bin);
        free(_fftOutput);
        free(_fftInput);
 }