PageRenderTime 345ms CodeModel.GetById 161ms app.highlight 12ms RepoModel.GetById 169ms app.codeStats 0ms

/external/pysoundtouch14/libsoundtouch/AAFilter.cpp

http://echo-nest-remix.googlecode.com/
C++ | 184 lines | 88 code | 34 blank | 62 comment | 7 complexity | b2fb8d31d0f473d89fc2ded83c3578f2 MD5 | raw file
  1////////////////////////////////////////////////////////////////////////////////
  2///
  3/// FIR low-pass (anti-alias) filter with filter coefficient design routine and
  4/// MMX optimization. 
  5/// 
  6/// Anti-alias filter is used to prevent folding of high frequencies when 
  7/// transposing the sample rate with interpolation.
  8///
  9/// Author        : Copyright (c) Olli Parviainen
 10/// Author e-mail : oparviai 'at' iki.fi
 11/// SoundTouch WWW: http://www.surina.net/soundtouch
 12///
 13////////////////////////////////////////////////////////////////////////////////
 14//
 15// Last changed  : $Date: 2009-01-11 13:34:24 +0200 (Sun, 11 Jan 2009) $
 16// File revision : $Revision: 4 $
 17//
 18// $Id: AAFilter.cpp 45 2009-01-11 11:34:24Z oparviai $
 19//
 20////////////////////////////////////////////////////////////////////////////////
 21//
 22// License :
 23//
 24//  SoundTouch audio processing library
 25//  Copyright (c) Olli Parviainen
 26//
 27//  This library is free software; you can redistribute it and/or
 28//  modify it under the terms of the GNU Lesser General Public
 29//  License as published by the Free Software Foundation; either
 30//  version 2.1 of the License, or (at your option) any later version.
 31//
 32//  This library is distributed in the hope that it will be useful,
 33//  but WITHOUT ANY WARRANTY; without even the implied warranty of
 34//  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
 35//  Lesser General Public License for more details.
 36//
 37//  You should have received a copy of the GNU Lesser General Public
 38//  License along with this library; if not, write to the Free Software
 39//  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
 40//
 41////////////////////////////////////////////////////////////////////////////////
 42
 43#include <memory.h>

 44#include <assert.h>

 45#include <math.h>

 46#include <stdlib.h>

 47#include "AAFilter.h"

 48#include "FIRFilter.h"

 49
 50using namespace soundtouch;
 51
 52#define PI        3.141592655357989

 53#define TWOPI    (2 * PI)

 54
 55/*****************************************************************************
 56 *
 57 * Implementation of the class 'AAFilter'
 58 *
 59 *****************************************************************************/
 60
 61AAFilter::AAFilter(uint len)
 62{
 63    pFIR = FIRFilter::newInstance();
 64    cutoffFreq = 0.5;
 65    setLength(len);
 66}
 67
 68
 69
 70AAFilter::~AAFilter()
 71{
 72    delete pFIR;
 73}
 74
 75
 76
 77// Sets new anti-alias filter cut-off edge frequency, scaled to
 78// sampling frequency (nyquist frequency = 0.5).
 79// The filter will cut frequencies higher than the given frequency.
 80void AAFilter::setCutoffFreq(double newCutoffFreq)
 81{
 82    cutoffFreq = newCutoffFreq;
 83    calculateCoeffs();
 84}
 85
 86
 87
 88// Sets number of FIR filter taps
 89void AAFilter::setLength(uint newLength)
 90{
 91    length = newLength;
 92    calculateCoeffs();
 93}
 94
 95
 96
 97// Calculates coefficients for a low-pass FIR filter using Hamming window
 98void AAFilter::calculateCoeffs()
 99{
100    uint i;
101    double cntTemp, temp, tempCoeff,h, w;
102    double fc2, wc;
103    double scaleCoeff, sum;
104    double *work;
105    SAMPLETYPE *coeffs;
106
107    assert(length >= 2);
108    assert(length % 4 == 0);
109    assert(cutoffFreq >= 0);
110    assert(cutoffFreq <= 0.5);
111
112    work = new double[length];
113    coeffs = new SAMPLETYPE[length];
114
115    fc2 = 2.0 * cutoffFreq; 
116    wc = PI * fc2;
117    tempCoeff = TWOPI / (double)length;
118
119    sum = 0;
120    for (i = 0; i < length; i ++) 
121    {
122        cntTemp = (double)i - (double)(length / 2);
123
124        temp = cntTemp * wc;
125        if (temp != 0) 
126        {
127            h = fc2 * sin(temp) / temp;                     // sinc function
128        } 
129        else 
130        {
131            h = 1.0;
132        }
133        w = 0.54 + 0.46 * cos(tempCoeff * cntTemp);       // hamming window
134
135        temp = w * h;
136        work[i] = temp;
137
138        // calc net sum of coefficients 
139        sum += temp;
140    }
141
142    // ensure the sum of coefficients is larger than zero
143    assert(sum > 0);
144
145    // ensure we've really designed a lowpass filter...
146    assert(work[length/2] > 0);
147    assert(work[length/2 + 1] > -1e-6);
148    assert(work[length/2 - 1] > -1e-6);
149
150    // Calculate a scaling coefficient in such a way that the result can be
151    // divided by 16384
152    scaleCoeff = 16384.0f / sum;
153
154    for (i = 0; i < length; i ++) 
155    {
156        // scale & round to nearest integer
157        temp = work[i] * scaleCoeff;
158        temp += (temp >= 0) ? 0.5 : -0.5;
159        // ensure no overfloods
160        assert(temp >= -32768 && temp <= 32767);
161        coeffs[i] = (SAMPLETYPE)temp;
162    }
163
164    // Set coefficients. Use divide factor 14 => divide result by 2^14 = 16384
165    pFIR->setCoefficients(coeffs, length, 14);
166
167    delete[] work;
168    delete[] coeffs;
169}
170
171
172// Applies the filter to the given sequence of samples. 
173// Note : The amount of outputted samples is by value of 'filter length' 
174// smaller than the amount of input samples.
175uint AAFilter::evaluate(SAMPLETYPE *dest, const SAMPLETYPE *src, uint numSamples, uint numChannels) const
176{
177    return pFIR->evaluate(dest, src, numSamples, numChannels);
178}
179
180
181uint AAFilter::getLength() const
182{
183    return pFIR->getLength();
184}