c@225: /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
c@225: 
c@225: /*
c@225:     QM DSP Library
c@225: 
c@225:     Centre for Digital Music, Queen Mary, University of London.
c@309:     This file 2005-2006 Christian Landone.
c@309: 
mathieu@321:     Modifications:
mathieu@321: 
mathieu@321:     - delta threshold
mathieu@321:     Description: add delta threshold used as offset in the smoothed
mathieu@321:     detection function
mathieu@321:     Author: Mathieu Barthet
mathieu@321:     Date: June 2010
mathieu@321: 
c@309:     This program is free software; you can redistribute it and/or
c@309:     modify it under the terms of the GNU General Public License as
c@309:     published by the Free Software Foundation; either version 2 of the
c@309:     License, or (at your option) any later version.  See the file
mathieu@321:     COPYING included with this distribution for more information.    
c@225: */
c@225: 
c@225: #include "PeakPicking.h"
c@241: #include "maths/Polyfit.h"
c@225: 
c@237: #include <iostream>
c@272: #include <cstring>
c@272: 
c@225: 
c@225: //////////////////////////////////////////////////////////////////////
c@225: // Construction/Destruction
c@225: //////////////////////////////////////////////////////////////////////
c@225: 
c@225: PeakPicking::PeakPicking( PPickParams Config )
c@225: {
c@225:     m_workBuffer = NULL;
c@225:     initialise( Config );
c@225: }
c@225: 
c@225: PeakPicking::~PeakPicking()
c@225: {
c@225:     deInitialise();
c@225: }
c@225: 
c@225: void PeakPicking::initialise( PPickParams Config )
c@225: {
c@225:     m_DFLength = Config.length ;
c@225:     Qfilta = Config.QuadThresh.a ;
c@225:     Qfiltb = Config.QuadThresh.b ;
c@225:     Qfiltc = Config.QuadThresh.c ;
cannam@483:         
c@225:     m_DFProcessingParams.length = m_DFLength; 
c@225:     m_DFProcessingParams.LPOrd = Config.LPOrd; 
c@225:     m_DFProcessingParams.LPACoeffs = Config.LPACoeffs; 
c@225:     m_DFProcessingParams.LPBCoeffs = Config.LPBCoeffs; 
c@225:     m_DFProcessingParams.winPre  = Config.WinT.pre;
c@225:     m_DFProcessingParams.winPost = Config.WinT.post; 
c@225:     m_DFProcessingParams.AlphaNormParam = Config.alpha;
c@225:     m_DFProcessingParams.isMedianPositive = false;
c@410:     m_DFProcessingParams.delta = Config.delta; //add the delta threshold as an adjustable parameter
mathieu@321: 
c@225:     m_DFSmoothing = new DFProcess( m_DFProcessingParams );
c@225: 
c@225:     m_workBuffer = new double[ m_DFLength ];
c@225:     memset( m_workBuffer, 0, sizeof(double)*m_DFLength);
c@225: }
c@225: 
c@225: void PeakPicking::deInitialise()
c@225: {
c@225:     delete [] m_workBuffer;
c@225:     delete m_DFSmoothing;
c@225:     m_workBuffer = NULL;
c@225: }
c@225: 
cannam@499: void PeakPicking::process( double* src, int len, vector<int> &onsets )
c@225: {
c@283:     if (len < 4) return;
c@283: 
cannam@483:     vector <double> m_maxima;   
c@225: 
c@225:     // Signal conditioning 
c@225:     m_DFSmoothing->process( src, m_workBuffer );
cannam@483:         
cannam@499:     for (int i = 0; i < len; i++) {
cannam@499:         m_maxima.push_back( m_workBuffer[ i ] );                
c@225:     }
cannam@483:         
c@225:     quadEval( m_maxima, onsets );
c@225: 
cannam@499:     for( int b = 0; b < (int)m_maxima.size(); b++) {
cannam@483:         src[ b ] = m_maxima[ b ];
c@225:     }
c@225: }
c@225: 
c@225: int PeakPicking::quadEval( vector<double> &src, vector<int> &idx )
c@225: {
cannam@499:     int maxLength;
c@225: 
c@225:     vector <int> m_maxIndex;
c@225:     vector <int> m_onsetPosition;
cannam@483:         
c@225:     vector <double> m_maxFit;
c@225:     vector <double> m_poly;
c@225:     vector <double> m_err;
c@225: 
c@225:     m_poly.push_back(0);
c@225:     m_poly.push_back(0);
c@225:     m_poly.push_back(0);
c@225: 
cannam@499:     for (int t = -2; t < 3; t++) {
cannam@483:         m_err.push_back( (double)t );
c@225:     }
cannam@483: 
cannam@499:     for (int i = 2; i < int(src.size()) - 2; i++) {
cannam@499:         if ((src[i] > src[i-1]) && (src[i] > src[i+1]) && (src[i] > 0) ) {
c@237:             m_maxIndex.push_back(i);
cannam@483:         }
c@225:     }
c@225: 
cannam@499:     maxLength = int(m_maxIndex.size());
c@225: 
c@225:     double selMax = 0;
c@225: 
cannam@499:     for (int j = 0; j < maxLength ; j++) {
cannam@483:         for (int k = -2; k <= 2; ++k) {
cannam@483:             selMax = src[ m_maxIndex[j] + k ] ;
cannam@483:             m_maxFit.push_back(selMax);                 
cannam@483:         }
c@225: 
cannam@483:         TPolyFit::PolyFit2(m_err, m_maxFit, m_poly);
c@225: 
cannam@483:         double f = m_poly[0];
cannam@483:         double h = m_poly[2];
c@225: 
cannam@483:         if (h < -Qfilta || f > Qfiltc) {
cannam@483:             idx.push_back(m_maxIndex[j]);
cannam@483:         }
cannam@483:                 
cannam@483:         m_maxFit.clear();
c@225:     }
c@225: 
c@225:     return 1;
c@225: }