libs/soundtouch/FIRFilter.cpp

   1 ////////////////////////////////////////////////////////////////////////////////
   2 ///
   3 /// General FIR digital filter routines with MMX optimization.
   4 ///
   5 /// Note : MMX optimized functions reside in a separate, platform-specific file,
   6 /// e.g. 'mmx_win.cpp' or 'mmx_gcc.cpp'
   7 ///
   8 /// Author        : Copyright (c) Olli Parviainen
   9 /// Author e-mail : oparviai @ iki.fi
  10 /// SoundTouch WWW: http://www.iki.fi/oparviai/soundtouch
  11 ///
  12 ////////////////////////////////////////////////////////////////////////////////
  13 //
  14 // Last changed  : $Date$
  15 // File revision : $Revision$
  16 //
  17 // $Id$
  18 //
  19 ////////////////////////////////////////////////////////////////////////////////
  20 //
  21 // License :
  22 //
  23 //  SoundTouch audio processing library
  24 //  Copyright (c) Olli Parviainen
  25 //
  26 //  This library is free software; you can redistribute it and/or
  27 //  modify it under the terms of the GNU Lesser General Public
  28 //  License as published by the Free Software Foundation; either
  29 //  version 2.1 of the License, or (at your option) any later version.
  30 //
  31 //  This library is distributed in the hope that it will be useful,
  32 //  but WITHOUT ANY WARRANTY; without even the implied warranty of
  33 //  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  34 //  Lesser General Public License for more details.
  35 //
  36 //  You should have received a copy of the GNU Lesser General Public
  37 //  License along with this library; if not, write to the Free Software
  38 //  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307  USA
  39 //
  40 ////////////////////////////////////////////////////////////////////////////////
  41
  42 #include <memory.h>
  43 #include <cassert>
  44 #include <cmath>
  45 #include <cstdlib>
  46 #include <stdexcept>
  47 #include "FIRFilter.h"
  48 #include "cpu_detect.h"
  49
  50 using namespace soundtouch;
  51
  52 /*****************************************************************************
  53  *
  54  * Implementation of the class 'FIRFilter'
  55  *
  56  *****************************************************************************/
  57
  58 FIRFilter::FIRFilter()
  59 {
  60     resultDivFactor = 0;
  61     length = 0;
  62     lengthDiv8 = 0;
  63     filterCoeffs = NULL;
  64 }
  65
  66
  67 FIRFilter::~FIRFilter()
  68 {
  69     delete[] filterCoeffs;
  70 }
  71
  72 // Usual C-version of the filter routine for stereo sound
  73 uint FIRFilter::evaluateFilterStereo(SAMPLETYPE *dest, const SAMPLETYPE *src, uint numSamples) const
  74 {
  75     uint i, j, end;
  76     LONG_SAMPLETYPE suml, sumr;
  77 #ifdef FLOAT_SAMPLES
  78     // when using floating point samples, use a scaler instead of a divider
  79     // because division is much slower operation than multiplying.
  80     double dScaler = 1.0 / (double)resultDivider;
  81 #endif
  82
  83     assert(length != 0);
  84
  85     end = 2 * (numSamples - length);
  86
  87     for (j = 0; j < end; j += 2)
  88     {
  89         const SAMPLETYPE *ptr;
  90
  91         suml = sumr = 0;
  92         ptr = src + j;
  93
  94         for (i = 0; i < length; i += 4)
  95         {
  96             // loop is unrolled by factor of 4 here for efficiency
  97             suml += ptr[2 * i + 0] * filterCoeffs[i + 0] +
  98                     ptr[2 * i + 2] * filterCoeffs[i + 1] +
  99                     ptr[2 * i + 4] * filterCoeffs[i + 2] +
 100                     ptr[2 * i + 6] * filterCoeffs[i + 3];
 101             sumr += ptr[2 * i + 1] * filterCoeffs[i + 0] +
 102                     ptr[2 * i + 3] * filterCoeffs[i + 1] +
 103                     ptr[2 * i + 5] * filterCoeffs[i + 2] +
 104                     ptr[2 * i + 7] * filterCoeffs[i + 3];
 105         }
 106
 107 #ifdef INTEGER_SAMPLES
 108         suml >>= resultDivFactor;
 109         sumr >>= resultDivFactor;
 110         // saturate to 16 bit integer limits
 111         suml = (suml < -32768) ? -32768 : (suml > 32767) ? 32767 : suml;
 112         // saturate to 16 bit integer limits
 113         sumr = (sumr < -32768) ? -32768 : (sumr > 32767) ? 32767 : sumr;
 114 #else
 115         suml *= dScaler;
 116         sumr *= dScaler;
 117 #endif // INTEGER_SAMPLES
 118         dest[j] = (SAMPLETYPE)suml;
 119         dest[j + 1] = (SAMPLETYPE)sumr;
 120     }
 121     return numSamples - length;
 122 }
 123
 124
 125
 126
 127 // Usual C-version of the filter routine for mono sound
 128 uint FIRFilter::evaluateFilterMono(SAMPLETYPE *dest, const SAMPLETYPE *src, uint numSamples) const
 129 {
 130     uint i, j, end;
 131     LONG_SAMPLETYPE sum;
 132 #ifdef FLOAT_SAMPLES
 133     // when using floating point samples, use a scaler instead of a divider
 134     // because division is much slower operation than multiplying.
 135     double dScaler = 1.0 / (double)resultDivider;
 136 #endif
 137
 138
 139     assert(length != 0);
 140
 141     end = numSamples - length;
 142     for (j = 0; j < end; j ++)
 143     {
 144         sum = 0;
 145         for (i = 0; i < length; i += 4)
 146         {
 147             // loop is unrolled by factor of 4 here for efficiency
 148             sum += src[i + 0] * filterCoeffs[i + 0] +
 149                    src[i + 1] * filterCoeffs[i + 1] +
 150                    src[i + 2] * filterCoeffs[i + 2] +
 151                    src[i + 3] * filterCoeffs[i + 3];
 152         }
 153 #ifdef INTEGER_SAMPLES
 154         sum >>= resultDivFactor;
 155         // saturate to 16 bit integer limits
 156         sum = (sum < -32768) ? -32768 : (sum > 32767) ? 32767 : sum;
 157 #else
 158         sum *= dScaler;
 159 #endif // INTEGER_SAMPLES
 160         dest[j] = (SAMPLETYPE)sum;
 161         src ++;
 162     }
 163     return end;
 164 }
 165
 166
 167 // Set filter coeffiecients and length.
 168 //
 169 // Throws an exception if filter length isn't divisible by 8
 170 void FIRFilter::setCoefficients(const SAMPLETYPE *coeffs, uint newLength, uint uResultDivFactor)
 171 {
 172     assert(newLength > 0);
 173     if (newLength % 8) throw std::runtime_error("FIR filter length not divisible by 8");
 174
 175     lengthDiv8 = newLength / 8;
 176     length = lengthDiv8 * 8;
 177     assert(length == newLength);
 178
 179     resultDivFactor = uResultDivFactor;
 180     resultDivider = (uint)pow(2, resultDivFactor);
 181
 182     delete[] filterCoeffs;
 183     filterCoeffs = new SAMPLETYPE[length];
 184     memcpy(filterCoeffs, coeffs, length * sizeof(SAMPLETYPE));
 185 }
 186
 187
 188 uint FIRFilter::getLength() const
 189 {
 190     return length;
 191 }
 192
 193
 194
 195 // Applies the filter to the given sequence of samples.
 196 //
 197 // Note : The amount of outputted samples is by value of 'filter_length'
 198 // smaller than the amount of input samples.
 199 uint FIRFilter::evaluate(SAMPLETYPE *dest, const SAMPLETYPE *src, uint numSamples, uint numChannels) const
 200 {
 201     assert(numChannels == 1 || numChannels == 2);
 202
 203     assert(length > 0);
 204     assert(lengthDiv8 * 8 == length);
 205     if (numSamples < length) return 0;
 206     assert(resultDivFactor >= 0);
 207     if (numChannels == 2)
 208     {
 209         return evaluateFilterStereo(dest, src, numSamples);
 210     } else {
 211         return evaluateFilterMono(dest, src, numSamples);
 212     }
 213 }
 214
 215 FIRFilter * FIRFilter::newInstance()
 216 {
 217     uint uExtensions;
 218
 219     uExtensions = detectCPUextensions();
 220
 221     // Check if MMX/SSE/3DNow! instruction set extensions supported by CPU
 222
 223 #ifdef ALLOW_MMX
 224     // MMX routines available only with integer sample types
 225     if (uExtensions & SUPPORT_MMX)
 226     {
 227         return ::new FIRFilterMMX;
 228     }
 229     else
 230 #endif // ALLOW_MMX
 231
 232 #ifdef ALLOW_SSE
 233     if (uExtensions & SUPPORT_SSE)
 234     {
 235         // SSE support
 236         return ::new FIRFilterSSE;
 237     }
 238     else
 239 #endif // ALLOW_SSE
 240
 241 #ifdef ALLOW_3DNOW
 242     if (uExtensions & SUPPORT_3DNOW)
 243     {
 244         // 3DNow! support
 245         return ::new FIRFilter3DNow;
 246     }
 247     else
 248 #endif // ALLOW_3DNOW
 249
 250     {
 251         // ISA optimizations not supported, use plain C version
 252         return ::new FIRFilter;
 253     }
 254 }