c@0: //
c@0: //  FIRFilter.cpp
c@0: //  Tempogram
c@0: //
c@0: //  Created by Carl Bussey on 25/06/2014.
c@0: //  Copyright (c) 2014 Carl Bussey. All rights reserved.
c@0: //
c@0: 
c@0: #include "FIRFilter.h"
c@7: 
c@0: using namespace std;
c@0: using Vamp::FFT;
c@0: 
c@13: FIRFilter::FIRFilter(const size_t &lengthInput, const size_t &numberOfCoefficients) :
c@8:     m_lengthInput(lengthInput),
c@8:     m_numberOfCoefficients(numberOfCoefficients),
c@13:     m_pFftInput(0),
c@13:     m_pFftCoefficients(0),
c@13:     m_pFftReal1(0),
c@13:     m_pFftImag1(0),
c@13:     m_pFftReal2(0),
c@13:     m_pFftImag2(0),
c@13:     m_pFftFilteredReal(0),
c@13:     m_pFftFilteredImag(0),
c@13:     m_pFftOutputReal(0),
c@13:     m_pFftOutputImag(0)
c@0: {
c@0:     initialise();
c@0: }
c@0: 
c@0: FIRFilter::~FIRFilter()
c@0: {
c@0:     cleanup();
c@0: }
c@0: 
c@9: //allocate memory
c@0: void
c@0: FIRFilter::initialise()
c@0: {
c@9:     //next power of 2
c@8:     m_lengthFIRFFT = pow(2,(ceil(log2(m_lengthInput+m_numberOfCoefficients-1))));
c@0:     
c@13:     m_pFftInput = new double[m_lengthFIRFFT];
c@13:     m_pFftCoefficients = new double[m_lengthFIRFFT];
c@13:     m_pFftReal1 = new double[m_lengthFIRFFT];
c@13:     m_pFftImag1 = new double[m_lengthFIRFFT];
c@13:     m_pFftReal2 = new double[m_lengthFIRFFT];
c@13:     m_pFftImag2 = new double[m_lengthFIRFFT];
c@13:     m_pFftFilteredReal = new double[m_lengthFIRFFT];
c@13:     m_pFftFilteredImag = new double[m_lengthFIRFFT];
c@13:     m_pFftOutputReal = new double[m_lengthFIRFFT];
c@13:     m_pFftOutputImag = new double[m_lengthFIRFFT];
c@0:     
c@20:     for(int i = 0; i < (int)m_lengthFIRFFT; i++){
c@13:         m_pFftInput[i] = m_pFftCoefficients[i] = m_pFftReal1[i] = m_pFftImag1[i] = m_pFftReal2[i] = m_pFftImag2[i] = m_pFftFilteredReal[i] = m_pFftFilteredImag[i] = m_pFftOutputReal[i] = m_pFftOutputImag[i] = 0.0;
c@0:     }
c@0: }
c@0: 
c@0: void
c@15: FIRFilter::process(const float* pInput, const float* pCoefficients, float* pOutput, OutputTypeArgument outputType)
c@0: {
c@15:     
c@9:     //Copy to same length FFT buffers
c@20:     for(int i = 0; i < (int)m_lengthFIRFFT; i++){
c@20:         m_pFftInput[i] = i < (int)m_lengthInput ? pInput[i] : 0.0;
c@20:         m_pFftCoefficients[i] = i < (int)m_numberOfCoefficients ? pCoefficients[i] : 0.0;
c@0:     }
c@0:     
c@13:     FFT::forward(m_lengthFIRFFT, m_pFftInput, 0, m_pFftReal1, m_pFftImag1);
c@13:     FFT::forward(m_lengthFIRFFT, m_pFftCoefficients, 0, m_pFftReal2, m_pFftImag2);
c@7:     
c@9:     //Multiply FFT coefficients. Multiplication in freq domain is convolution in time domain.
c@20:     for (int i = 0; i < (int)m_lengthFIRFFT; i++){
c@13:         m_pFftFilteredReal[i] = (m_pFftReal1[i] * m_pFftReal2[i]) - (m_pFftImag1[i] * m_pFftImag2[i]);
c@13:         m_pFftFilteredImag[i] = (m_pFftReal1[i] * m_pFftImag2[i]) + (m_pFftReal2[i] * m_pFftImag1[i]);
c@0:     }
c@13:     
c@13:     FFT::inverse(m_lengthFIRFFT, m_pFftFilteredReal, m_pFftFilteredImag, m_pFftOutputReal, m_pFftOutputImag);
c@0:     
c@9:     //copy to output
c@15:     int offset = 0;
c@20:     int outputLength = m_lengthInput;
c@20:     if (outputType == all) outputLength = (int)m_lengthFIRFFT;
c@20:     else if (outputType == middle) offset = floor(m_numberOfCoefficients/2.0f);
c@15:     else if (outputType != first) cerr << "FIRFilter::process(params) - " << outputType << " is not a valid argument. outputType is set to first." << endl;
c@15:     
c@20:     for (int i = 0; i < outputLength; i++){
c@14:         pOutput[i] = m_pFftOutputReal[i + offset];
c@0:     }
c@0: }
c@0: 
c@9: //remove memory allocations
c@0: void
c@0: FIRFilter::cleanup()
c@0: {
c@13:     delete []m_pFftInput;
c@13:     delete []m_pFftCoefficients;
c@13:     delete []m_pFftReal1;
c@13:     delete []m_pFftImag1;
c@13:     delete []m_pFftReal2;
c@13:     delete []m_pFftImag2;
c@13:     delete []m_pFftFilteredReal;
c@13:     delete []m_pFftFilteredImag;
c@13:     delete []m_pFftOutputReal;
c@13:     delete []m_pFftOutputImag;
c@13:     m_pFftInput = m_pFftCoefficients = m_pFftReal1 = m_pFftImag1 = m_pFftReal2 = m_pFftImag2 = m_pFftFilteredReal = m_pFftFilteredImag = m_pFftOutputReal = m_pFftOutputImag = 0;
Chris@10: }