cannam@0: /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
cannam@0: 
cannam@0: /*
cannam@0:     QM DSP Library
cannam@0: 
cannam@0:     Centre for Digital Music, Queen Mary, University of London.
Chris@84:     This file 2005-2006 Christian Landone.
Chris@84: 
Chris@84:     This program is free software; you can redistribute it and/or
Chris@84:     modify it under the terms of the GNU General Public License as
Chris@84:     published by the Free Software Foundation; either version 2 of the
Chris@84:     License, or (at your option) any later version.  See the file
Chris@84:     COPYING included with this distribution for more information.
cannam@0: */
cannam@0: 
cannam@0: #include "FiltFilt.h"
cannam@0: 
cannam@0: //////////////////////////////////////////////////////////////////////
cannam@0: // Construction/Destruction
cannam@0: //////////////////////////////////////////////////////////////////////
cannam@0: 
Chris@193: FiltFilt::FiltFilt(Filter::Parameters parameters) :
Chris@193:     m_filter(parameters)
cannam@0: {
Chris@193:     m_ord = m_filter.getOrder();
cannam@0: }
cannam@0: 
cannam@0: FiltFilt::~FiltFilt()
cannam@0: {
cannam@0: }
cannam@0: 
cannam@0: void FiltFilt::process(double *src, double *dst, unsigned int length)
cannam@0: {	
cannam@0:     unsigned int i;
cannam@0: 
cannam@58:     if (length == 0) return;
cannam@58: 
cannam@0:     unsigned int nFilt = m_ord + 1;
cannam@0:     unsigned int nFact = 3 * ( nFilt - 1);
cannam@0:     unsigned int nExt	= length + 2 * nFact;
cannam@0: 
Chris@193:     double *filtScratchIn = new double[ nExt ];
Chris@193:     double *filtScratchOut = new double[ nExt ];
cannam@0: 	
cannam@0:     for( i = 0; i< nExt; i++ ) 
cannam@0:     {
Chris@193: 	filtScratchIn[ i ] = 0.0;
Chris@193: 	filtScratchOut[ i ] = 0.0;
cannam@0:     }
cannam@0: 
cannam@0:     // Edge transients reflection
cannam@0:     double sample0 = 2 * src[ 0 ];
cannam@0:     double sampleN = 2 * src[ length - 1 ];
cannam@0: 
cannam@0:     unsigned int index = 0;
cannam@0:     for( i = nFact; i > 0; i-- )
cannam@0:     {
Chris@193: 	filtScratchIn[ index++ ] = sample0 - src[ i ];
cannam@0:     }
cannam@0:     index = 0;
cannam@0:     for( i = 0; i < nFact; i++ )
cannam@0:     {
Chris@193: 	filtScratchIn[ (nExt - nFact) + index++ ] = sampleN - src[ (length - 2) - i ];
cannam@0:     }
cannam@0: 
cannam@0:     index = 0;
cannam@0:     for( i = 0; i < length; i++ )
cannam@0:     {
Chris@193: 	filtScratchIn[ i + nFact ] = src[ i ];
cannam@0:     }
cannam@0: 	
cannam@0:     ////////////////////////////////
cannam@0:     // Do  0Ph filtering
Chris@193:     m_filter.process( filtScratchIn, filtScratchOut, nExt);
cannam@0: 	
cannam@0:     // reverse the series for FILTFILT 
cannam@0:     for ( i = 0; i < nExt; i++)
cannam@0:     { 
Chris@193: 	filtScratchIn[ i ] = filtScratchOut[ nExt - i - 1];
cannam@0:     }
cannam@0: 
cannam@0:     // do FILTER again 
Chris@193:     m_filter.process( filtScratchIn, filtScratchOut, nExt);
cannam@0: 	
cannam@0:     // reverse the series back 
cannam@0:     for ( i = 0; i < nExt; i++)
cannam@0:     {
Chris@193: 	filtScratchIn[ i ] = filtScratchOut[ nExt - i - 1 ];
cannam@0:     }
cannam@0:     for ( i = 0;i < nExt; i++)
cannam@0:     {
Chris@193: 	filtScratchOut[ i ] = filtScratchIn[ i ];
cannam@0:     }
cannam@0: 
cannam@0:     index = 0;
cannam@0:     for( i = 0; i < length; i++ )
cannam@0:     {
Chris@193: 	dst[ index++ ] = filtScratchOut[ i + nFact ];
cannam@0:     }	
cannam@0: 
Chris@193:     delete [] filtScratchIn;
Chris@193:     delete [] filtScratchOut;
cannam@0: 
cannam@0: }
cannam@0: 
cannam@0: void FiltFilt::reset()
cannam@0: {
cannam@0: 
cannam@0: }