libs/qm-dsp/dsp/signalconditioning/FiltFilt.cpp

   1 /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
   2
   3 /*
   4     QM DSP Library
   5
   6     Centre for Digital Music, Queen Mary, University of London.
   7     This file 2005-2006 Christian Landone.
   8
   9     This program is free software; you can redistribute it and/or
  10     modify it under the terms of the GNU General Public License as
  11     published by the Free Software Foundation; either version 2 of the
  12     License, or (at your option) any later version.  See the file
  13     COPYING included with this distribution for more information.
  14 */
  15
  16 #include "FiltFilt.h"
  17
  18 //////////////////////////////////////////////////////////////////////
  19 // Construction/Destruction
  20 //////////////////////////////////////////////////////////////////////
  21
  22 FiltFilt::FiltFilt( FiltFiltConfig Config )
  23 {
  24     m_filtScratchIn = NULL;
  25     m_filtScratchOut = NULL;
  26     m_ord = 0;
  27
  28     initialise( Config );
  29 }
  30
  31 FiltFilt::~FiltFilt()
  32 {
  33     deInitialise();
  34 }
  35
  36 void FiltFilt::initialise( FiltFiltConfig Config )
  37 {
  38     m_ord = Config.ord;
  39     m_filterConfig.ord = Config.ord;
  40     m_filterConfig.ACoeffs = Config.ACoeffs;
  41     m_filterConfig.BCoeffs = Config.BCoeffs;
  42
  43     m_filter = new Filter( m_filterConfig );
  44 }
  45
  46 void FiltFilt::deInitialise()
  47 {
  48     delete m_filter;
  49 }
  50
  51
  52 void FiltFilt::process(double *src, double *dst, unsigned int length)
  53 {
  54     unsigned int i;
  55
  56     if (length == 0) return;
  57
  58     unsigned int nFilt = m_ord + 1;
  59     unsigned int nFact = 3 * ( nFilt - 1);
  60     unsigned int nExt   = length + 2 * nFact;
  61
  62     m_filtScratchIn = new double[ nExt ];
  63     m_filtScratchOut = new double[ nExt ];
  64
  65
  66     for( i = 0; i< nExt; i++ )
  67     {
  68         m_filtScratchIn[ i ] = 0.0;
  69         m_filtScratchOut[ i ] = 0.0;
  70     }
  71
  72     // Edge transients reflection
  73     double sample0 = 2 * src[ 0 ];
  74     double sampleN = 2 * src[ length - 1 ];
  75
  76     unsigned int index = 0;
  77     for( i = nFact; i > 0; i-- )
  78     {
  79         m_filtScratchIn[ index++ ] = sample0 - src[ i ];
  80     }
  81     index = 0;
  82     for( i = 0; i < nFact; i++ )
  83     {
  84         m_filtScratchIn[ (nExt - nFact) + index++ ] = sampleN - src[ (length - 2) - i ];
  85     }
  86
  87     index = 0;
  88     for( i = 0; i < length; i++ )
  89     {
  90         m_filtScratchIn[ i + nFact ] = src[ i ];
  91     }
  92
  93     ////////////////////////////////
  94     // Do  0Ph filtering
  95     m_filter->process( m_filtScratchIn, m_filtScratchOut, nExt);
  96
  97     // reverse the series for FILTFILT
  98     for ( i = 0; i < nExt; i++)
  99     {
 100         m_filtScratchIn[ i ] = m_filtScratchOut[ nExt - i - 1];
 101     }
 102
 103     // do FILTER again
 104     m_filter->process( m_filtScratchIn, m_filtScratchOut, nExt);
 105
 106     // reverse the series back
 107     for ( i = 0; i < nExt; i++)
 108     {
 109         m_filtScratchIn[ i ] = m_filtScratchOut[ nExt - i - 1 ];
 110     }
 111     for ( i = 0;i < nExt; i++)
 112     {
 113         m_filtScratchOut[ i ] = m_filtScratchIn[ i ];
 114     }
 115
 116     index = 0;
 117     for( i = 0; i < length; i++ )
 118     {
 119         dst[ index++ ] = m_filtScratchOut[ i + nFact ];
 120     }
 121
 122     delete [] m_filtScratchIn;
 123     delete [] m_filtScratchOut;
 124
 125 }
 126
 127 void FiltFilt::reset()
 128 {
 129
 130 }