libs/soundtouch/AAFilter.cpp

   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 @ iki.fi
  11 /// SoundTouch WWW: http://www.iki.fi/oparviai/soundtouch
  12 ///
  13 ////////////////////////////////////////////////////////////////////////////////
  14 //
  15 // Last changed  : $Date$
  16 // File revision : $Revision$
  17 //
  18 // $Id$
  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
  50 using 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
  61 AAFilter::AAFilter(const uint length)
  62 {
  63     pFIR = FIRFilter::newInstance();
  64     cutoffFreq = 0.5;
  65     setLength(length);
  66 }
  67
  68
  69
  70 AAFilter::~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.
  80 void AAFilter::setCutoffFreq(const double newCutoffFreq)
  81 {
  82     cutoffFreq = newCutoffFreq;
  83     calculateCoeffs();
  84 }
  85
  86
  87
  88 // Sets number of FIR filter taps
  89 void AAFilter::setLength(const uint newLength)
  90 {
  91     length = newLength;
  92     calculateCoeffs();
  93 }
  94
  95
  96
  97 // Calculates coefficients for a low-pass FIR filter using Hamming window
  98 void 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 > 0);
 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.
 175 uint AAFilter::evaluate(SAMPLETYPE *dest, const SAMPLETYPE *src, uint numSamples, uint numChannels) const
 176 {
 177     return pFIR->evaluate(dest, src, numSamples, numChannels);
 178 }
 179
 180
 181 uint AAFilter::getLength() const
 182 {
 183     return pFIR->getLength();
 184 }