summaryrefslogtreecommitdiff
path: root/libs/qm-dsp/dsp/tempotracking
diff options
context:
space:
mode:
Diffstat (limited to 'libs/qm-dsp/dsp/tempotracking')
-rw-r--r--libs/qm-dsp/dsp/tempotracking/DownBeat.cpp33
-rw-r--r--libs/qm-dsp/dsp/tempotracking/DownBeat.h11
-rw-r--r--libs/qm-dsp/dsp/tempotracking/TempoTrack.cpp1729
-rw-r--r--libs/qm-dsp/dsp/tempotracking/TempoTrack.h12
-rw-r--r--libs/qm-dsp/dsp/tempotracking/TempoTrackV2.cpp60
-rw-r--r--libs/qm-dsp/dsp/tempotracking/TempoTrackV2.h29
6 files changed, 971 insertions, 903 deletions
diff --git a/libs/qm-dsp/dsp/tempotracking/DownBeat.cpp b/libs/qm-dsp/dsp/tempotracking/DownBeat.cpp
index 1167bc6ea9..8570da6a17 100644
--- a/libs/qm-dsp/dsp/tempotracking/DownBeat.cpp
+++ b/libs/qm-dsp/dsp/tempotracking/DownBeat.cpp
@@ -44,7 +44,10 @@ DownBeat::DownBeat(float originalSampleRate,
// 16x decimation, which is our expected normal situation)
m_beatframesize = MathUtilities::nextPowerOfTwo
(int((m_rate / decimationFactor) * 1.3));
-// std::cerr << "rate = " << m_rate << ", bfs = " << m_beatframesize << std::endl;
+ if (m_beatframesize < 2) {
+ m_beatframesize = 2;
+ }
+// std::cerr << "rate = " << m_rate << ", dec = " << decimationFactor << ", bfs = " << m_beatframesize << std::endl;
m_beatframe = new double[m_beatframesize];
m_fftRealOut = new double[m_beatframesize];
m_fftImagOut = new double[m_beatframesize];
@@ -122,7 +125,7 @@ DownBeat::pushAudioBlock(const float *audio)
// std::cerr << "pushAudioBlock: rms in " << sqrt(rmsin) << ", out " << sqrt(rmsout) << std::endl;
m_buffill += m_increment / m_factor;
}
-
+
const float *
DownBeat::getBufferedAudio(size_t &length) const
{
@@ -192,9 +195,9 @@ DownBeat::findDownBeats(const float *audio,
}
// Now FFT beat frame
-
- m_fft->process(false, m_beatframe, m_fftRealOut, m_fftImagOut);
-
+
+ m_fft->forward(m_beatframe, m_fftRealOut, m_fftImagOut);
+
// Calculate magnitudes
for (size_t j = 0; j < m_beatframesize/2; ++j) {
@@ -257,7 +260,7 @@ DownBeat::measureSpecDiff(d_vec_t oldspec, d_vec_t newspec)
{
// JENSEN-SHANNON DIVERGENCE BETWEEN SPECTRAL FRAMES
- unsigned int SPECSIZE = 512; // ONLY LOOK AT FIRST 512 SAMPLES OF SPECTRUM.
+ unsigned int SPECSIZE = 512; // ONLY LOOK AT FIRST 512 SAMPLES OF SPECTRUM.
if (SPECSIZE > oldspec.size()/4) {
SPECSIZE = oldspec.size()/4;
}
@@ -266,37 +269,37 @@ DownBeat::measureSpecDiff(d_vec_t oldspec, d_vec_t newspec)
double sumnew = 0.;
double sumold = 0.;
-
+
for (unsigned int i = 0;i < SPECSIZE;i++)
{
newspec[i] +=EPS;
oldspec[i] +=EPS;
-
+
sumnew+=newspec[i];
sumold+=oldspec[i];
- }
-
+ }
+
for (unsigned int i = 0;i < SPECSIZE;i++)
{
newspec[i] /= (sumnew);
oldspec[i] /= (sumold);
-
+
// IF ANY SPECTRAL VALUES ARE 0 (SHOULDN'T BE ANY!) SET THEM TO 1
if (newspec[i] == 0)
{
newspec[i] = 1.;
}
-
+
if (oldspec[i] == 0)
{
oldspec[i] = 1.;
}
-
+
// JENSEN-SHANNON CALCULATION
- sd1 = 0.5*oldspec[i] + 0.5*newspec[i];
+ sd1 = 0.5*oldspec[i] + 0.5*newspec[i];
SD = SD + (-sd1*log(sd1)) + (0.5*(oldspec[i]*log(oldspec[i]))) + (0.5*(newspec[i]*log(newspec[i])));
}
-
+
return SD;
}
diff --git a/libs/qm-dsp/dsp/tempotracking/DownBeat.h b/libs/qm-dsp/dsp/tempotracking/DownBeat.h
index 3ef0d18127..fc1d7b2ec9 100644
--- a/libs/qm-dsp/dsp/tempotracking/DownBeat.h
+++ b/libs/qm-dsp/dsp/tempotracking/DownBeat.h
@@ -17,6 +17,7 @@
#define DOWNBEAT_H
#include <vector>
+#include <cstddef>
#include "dsp/rateconversion/Decimator.h"
@@ -28,7 +29,7 @@ class FFTReal;
* This class takes an input audio signal and a sequence of beat
* locations (calculated e.g. by TempoTrackV2) and estimates which of
* the beat locations are downbeats (first beat of the bar).
- *
+ *
* The input audio signal is expected to have been downsampled to a
* very low sampling rate (e.g. 2700Hz). A utility function for
* downsampling and buffering incoming block-by-block audio is
@@ -56,7 +57,7 @@ public:
/**
* Estimate which beats are down-beats.
- *
+ *
* audio contains the input audio stream after downsampling, and
* audioLength contains the number of samples in this downsampled
* stream.
@@ -83,18 +84,18 @@ public:
* and the region following it.
*/
void getBeatSD(vector<double> &beatsd) const;
-
+
/**
* For your downsampling convenience: call this function
* repeatedly with input audio blocks containing dfIncrement
* samples at the original sample rate, to decimate them to the
* downsampled rate and buffer them within the DownBeat class.
- *
+ *
* Call getBufferedAudio() to retrieve the results after all
* blocks have been processed.
*/
void pushAudioBlock(const float *audio);
-
+
/**
* Retrieve the accumulated audio produced by pushAudioBlock calls.
*/
diff --git a/libs/qm-dsp/dsp/tempotracking/TempoTrack.cpp b/libs/qm-dsp/dsp/tempotracking/TempoTrack.cpp
index c36385f250..389403edaa 100644
--- a/libs/qm-dsp/dsp/tempotracking/TempoTrack.cpp
+++ b/libs/qm-dsp/dsp/tempotracking/TempoTrack.cpp
@@ -1,869 +1,870 @@
-/* -*- c-basic-offset: 4 indent-tabs-mode: nil -*- vi:set ts=8 sts=4 sw=4: */
-
-/*
- QM DSP Library
-
- Centre for Digital Music, Queen Mary, University of London.
- This file copyright 2005-2006 Christian Landone.and Matthew Davies.
+/* -*- c-basic-offset: 4 indent-tabs-mode: nil -*- vi:set ts=8 sts=4 sw=4: */
+
+/*
+ QM DSP Library
+
+ Centre for Digital Music, Queen Mary, University of London.
+ This file copyright 2005-2006 Christian Landone.and Matthew Davies.
This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License as
published by the Free Software Foundation; either version 2 of the
License, or (at your option) any later version. See the file
- COPYING included with this distribution for more information.
-*/
-
-#include "TempoTrack.h"
-
-#include "maths/MathAliases.h"
-#include "maths/MathUtilities.h"
-
-#include <iostream>
-
-#include <cassert>
-
-//#define DEBUG_TEMPO_TRACK 1
-
-
-#define RAY43VAL
-
-//////////////////////////////////////////////////////////////////////
-// Construction/Destruction
-//////////////////////////////////////////////////////////////////////
-
-TempoTrack::TempoTrack( TTParams Params )
-{
- m_tempoScratch = NULL;
- m_rawDFFrame = NULL;
- m_smoothDFFrame = NULL;
- m_frameACF = NULL;
- m_smoothRCF = NULL;
-
- m_dataLength = 0;
- m_winLength = 0;
- m_lagLength = 0;
-
- m_rayparam = 0;
- m_sigma = 0;
- m_DFWVNnorm = 0;
-
- initialise( Params );
-}
-
-TempoTrack::~TempoTrack()
-{
- deInitialise();
-}
-
-void TempoTrack::initialise( TTParams Params )
-{
- m_winLength = Params.winLength;
- m_lagLength = Params.lagLength;
-
- m_rayparam = 43.0;
- m_sigma = sqrt(3.9017);
- m_DFWVNnorm = exp( ( log( 2.0 ) / m_rayparam ) * ( m_winLength + 2 ) );
-
- m_rawDFFrame = new double[ m_winLength ];
- m_smoothDFFrame = new double[ m_winLength ];
- m_frameACF = new double[ m_winLength ];
- m_tempoScratch = new double[ m_lagLength ];
- m_smoothRCF = new double[ m_lagLength ];
-
-
- unsigned int winPre = Params.WinT.pre;
- unsigned int winPost = Params.WinT.post;
-
- m_DFFramer.configure( m_winLength, m_lagLength );
-
- m_DFPParams.length = m_winLength;
- m_DFPParams.AlphaNormParam = Params.alpha;
- m_DFPParams.LPOrd = Params.LPOrd;
- m_DFPParams.LPACoeffs = Params.LPACoeffs;
- m_DFPParams.LPBCoeffs = Params.LPBCoeffs;
- m_DFPParams.winPre = Params.WinT.pre;
- m_DFPParams.winPost = Params.WinT.post;
- m_DFPParams.isMedianPositive = true;
-
- m_DFConditioning = new DFProcess( m_DFPParams );
-
-
- // these are parameters for smoothing m_tempoScratch
- m_RCFPParams.length = m_lagLength;
- m_RCFPParams.AlphaNormParam = Params.alpha;
- m_RCFPParams.LPOrd = Params.LPOrd;
- m_RCFPParams.LPACoeffs = Params.LPACoeffs;
- m_RCFPParams.LPBCoeffs = Params.LPBCoeffs;
- m_RCFPParams.winPre = Params.WinT.pre;
- m_RCFPParams.winPost = Params.WinT.post;
- m_RCFPParams.isMedianPositive = true;
-
- m_RCFConditioning = new DFProcess( m_RCFPParams );
-
-}
-
-void TempoTrack::deInitialise()
-{
- delete [] m_rawDFFrame;
-
- delete [] m_smoothDFFrame;
-
- delete [] m_smoothRCF;
-
- delete [] m_frameACF;
-
- delete [] m_tempoScratch;
-
- delete m_DFConditioning;
-
- delete m_RCFConditioning;
-
-}
-
-void TempoTrack::createCombFilter(double* Filter, unsigned int winLength, unsigned int TSig, double beatLag)
-{
- unsigned int i;
-
- if( beatLag == 0 )
- {
- for( i = 0; i < winLength; i++ )
- {
- Filter[ i ] = ( ( i + 1 ) / pow( m_rayparam, 2.0) ) * exp( ( -pow(( i + 1 ),2.0 ) / ( 2.0 * pow( m_rayparam, 2.0))));
- }
- }
- else
- {
- m_sigma = beatLag/4;
- for( i = 0; i < winLength; i++ )
- {
- double dlag = (double)(i+1) - beatLag;
- Filter[ i ] = exp(-0.5 * pow(( dlag / m_sigma), 2.0) ) / (sqrt( 2 * PI) * m_sigma);
- }
- }
-}
-
-double TempoTrack::tempoMM(double* ACF, double* weight, int tsig)
-{
-
- double period = 0;
- double maxValRCF = 0.0;
- unsigned int maxIndexRCF = 0;
-
- double* pdPeaks;
-
- unsigned int maxIndexTemp;
- double maxValTemp;
- unsigned int count;
-
- unsigned int numelem,i,j;
- int a, b;
-
- for( i = 0; i < m_lagLength; i++ )
- m_tempoScratch[ i ] = 0.0;
-
- if( tsig == 0 )
- {
- //if time sig is unknown, use metrically unbiased version of Filterbank
- numelem = 4;
- }
- else
- {
- numelem = tsig;
- }
-
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "tempoMM: m_winLength = " << m_winLength << ", m_lagLength = " << m_lagLength << ", numelem = " << numelem << std::endl;
-#endif
-
- for(i=1;i<m_lagLength-1;i++)
- {
- //first and last output values are left intentionally as zero
- for (a=1;a<=numelem;a++)
- {
- for(b=(1-a);b<a;b++)
- {
- if( tsig == 0 )
- {
- m_tempoScratch[i] += ACF[a*(i+1)+b-1] * (1.0 / (2.0 * (double)a-1)) * weight[i];
- }
- else
- {
- m_tempoScratch[i] += ACF[a*(i+1)+b-1] * 1 * weight[i];
- }
- }
- }
- }
-
-
- //////////////////////////////////////////////////
- // MODIFIED BEAT PERIOD EXTRACTION //////////////
- /////////////////////////////////////////////////
-
- // find smoothed version of RCF ( as applied to Detection Function)
- m_RCFConditioning->process( m_tempoScratch, m_smoothRCF);
-
- if (tsig != 0) // i.e. in context dependent state
- {
-// NOW FIND MAX INDEX OF ACFOUT
- for( i = 0; i < m_lagLength; i++)
- {
- if( m_tempoScratch[ i ] > maxValRCF)
- {
- maxValRCF = m_tempoScratch[ i ];
- maxIndexRCF = i;
- }
- }
- }
- else // using rayleigh weighting
- {
- vector <vector<double> > rcfMat;
-
- double sumRcf = 0.;
-
- double maxVal = 0.;
- // now find the two values which minimise rcfMat
- double minVal = 0.;
- int p_i = 1; // periodicity for row i;
- int p_j = 1; //periodicity for column j;
-
-
- for ( i=0; i<m_lagLength; i++)
- {
- m_tempoScratch[i] =m_smoothRCF[i];
- }
-
- // normalise m_tempoScratch so that it sums to zero.
- for ( i=0; i<m_lagLength; i++)
- {
- sumRcf += m_tempoScratch[i];
- }
-
- for( i=0; i<m_lagLength; i++)
- {
- m_tempoScratch[i] /= sumRcf;
- }
-
- // create a matrix to store m_tempoScratchValues modified by log2 ratio
- for ( i=0; i<m_lagLength; i++)
- {
- rcfMat.push_back ( vector<double>() ); // adds a new row...
- }
-
- for (i=0; i<m_lagLength; i++)
- {
- for (j=0; j<m_lagLength; j++)
- {
- rcfMat[i].push_back (0.);
- }
- }
-
- // the 'i' and 'j' indices deliberately start from '1' and not '0'
- for ( i=1; i<m_lagLength; i++)
- {
- for (j=1; j<m_lagLength; j++)
- {
- double log2PeriodRatio = log( static_cast<double>(i)/static_cast<double>(j) ) / log(2.0);
- rcfMat[i][j] = ( abs(1.0-abs(log2PeriodRatio)) );
- rcfMat[i][j] += ( 0.01*( 1./(m_tempoScratch[i]+m_tempoScratch[j]) ) );
- }
- }
-
- // set diagonal equal to maximum value in rcfMat
- // we don't want to pick one strong middle peak - we need a combination of two peaks.
-
- for ( i=1; i<m_lagLength; i++)
- {
- for (j=1; j<m_lagLength; j++)
- {
- if (rcfMat[i][j] > maxVal)
- {
- maxVal = rcfMat[i][j];
- }
- }
- }
-
- for ( i=1; i<m_lagLength; i++)
- {
- rcfMat[i][i] = maxVal;
- }
-
- // now find the row and column number which minimise rcfMat
- minVal = maxVal;
-
- for ( i=1; i<m_lagLength; i++)
- {
- for ( j=1; j<m_lagLength; j++)
- {
- if (rcfMat[i][j] < minVal)
- {
- minVal = rcfMat[i][j];
- p_i = i;
- p_j = j;
- }
- }
- }
-
-
- // initially choose p_j (arbitrary) - saves on an else statement
- int beatPeriod = p_j;
- if (m_tempoScratch[p_i] > m_tempoScratch[p_j])
- {
- beatPeriod = p_i;
- }
-
- // now write the output
- maxIndexRCF = static_cast<int>(beatPeriod);
- }
-
-
- double locked = 5168.f / maxIndexRCF;
- if (locked >= 30 && locked <= 180) {
- m_lockedTempo = locked;
- }
-
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "tempoMM: locked tempo = " << m_lockedTempo << std::endl;
-#endif
-
- if( tsig == 0 )
- tsig = 4;
-
-
-#ifdef DEBUG_TEMPO_TRACK
-std::cerr << "tempoMM: maxIndexRCF = " << maxIndexRCF << std::endl;
-#endif
-
- if( tsig == 4 )
- {
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "tsig == 4" << std::endl;
-#endif
-
- pdPeaks = new double[ 4 ];
- for( i = 0; i < 4; i++ ){ pdPeaks[ i ] = 0.0;}
-
- pdPeaks[ 0 ] = ( double )maxIndexRCF + 1;
-
- maxIndexTemp = 0;
- maxValTemp = 0.0;
- count = 0;
-
- for( i = (2 * maxIndexRCF + 1) - 1; i < (2 * maxIndexRCF + 1) + 2; i++ )
- {
- if( ACF[ i ] > maxValTemp )
- {
- maxValTemp = ACF[ i ];
- maxIndexTemp = count;
- }
- count++;
- }
- pdPeaks[ 1 ] = (double)( maxIndexTemp + 1 + ( (2 * maxIndexRCF + 1 ) - 2 ) + 1 )/2;
-
- maxIndexTemp = 0;
- maxValTemp = 0.0;
- count = 0;
-
- for( i = (3 * maxIndexRCF + 2 ) - 2; i < (3 * maxIndexRCF + 2 ) + 3; i++ )
- {
- if( ACF[ i ] > maxValTemp )
- {
- maxValTemp = ACF[ i ];
- maxIndexTemp = count;
- }
- count++;
- }
- pdPeaks[ 2 ] = (double)( maxIndexTemp + 1 + ( (3 * maxIndexRCF + 2) - 4 ) + 1 )/3;
-
- maxIndexTemp = 0;
- maxValTemp = 0.0;
- count = 0;
-
- for( i = ( 4 * maxIndexRCF + 3) - 3; i < ( 4 * maxIndexRCF + 3) + 4; i++ )
- {
- if( ACF[ i ] > maxValTemp )
- {
- maxValTemp = ACF[ i ];
- maxIndexTemp = count;
- }
- count++;
- }
- pdPeaks[ 3 ] = (double)( maxIndexTemp + 1 + ( (4 * maxIndexRCF + 3) - 9 ) + 1 )/4 ;
-
-
- period = MathUtilities::mean( pdPeaks, 4 );
- }
- else
- {
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "tsig != 4" << std::endl;
-#endif
-
- pdPeaks = new double[ 3 ];
- for( i = 0; i < 3; i++ ){ pdPeaks[ i ] = 0.0;}
-
- pdPeaks[ 0 ] = ( double )maxIndexRCF + 1;
-
- maxIndexTemp = 0;
- maxValTemp = 0.0;
- count = 0;
-
- for( i = (2 * maxIndexRCF + 1) - 1; i < (2 * maxIndexRCF + 1) + 2; i++ )
- {
- if( ACF[ i ] > maxValTemp )
- {
- maxValTemp = ACF[ i ];
- maxIndexTemp = count;
- }
- count++;
- }
- pdPeaks[ 1 ] = (double)( maxIndexTemp + 1 + ( (2 * maxIndexRCF + 1 ) - 2 ) + 1 )/2;
-
- maxIndexTemp = 0;
- maxValTemp = 0.0;
- count = 0;
-
- for( i = (3 * maxIndexRCF + 2 ) - 2; i < (3 * maxIndexRCF + 2 ) + 3; i++ )
- {
- if( ACF[ i ] > maxValTemp )
- {
- maxValTemp = ACF[ i ];
- maxIndexTemp = count;
- }
- count++;
- }
- pdPeaks[ 2 ] = (double)( maxIndexTemp + 1 + ( (3 * maxIndexRCF + 2) - 4 ) + 1 )/3;
-
-
- period = MathUtilities::mean( pdPeaks, 3 );
- }
-
- delete [] pdPeaks;
-
- return period;
-}
-
-void TempoTrack::stepDetect( double* periodP, double* periodG, int currentIdx, int* flag )
-{
- double stepthresh = 1 * 3.9017;
-
- if( *flag )
- {
- if(abs(periodG[ currentIdx ] - periodP[ currentIdx ]) > stepthresh)
- {
- // do nuffin'
- }
- }
- else
- {
- if(fabs(periodG[ currentIdx ]-periodP[ currentIdx ]) > stepthresh)
- {
- *flag = 3;
- }
- }
-}
-
-void TempoTrack::constDetect( double* periodP, int currentIdx, int* flag )
-{
- double constthresh = 2 * 3.9017;
-
- if( fabs( 2 * periodP[ currentIdx ] - periodP[ currentIdx - 1] - periodP[ currentIdx - 2] ) < constthresh)
- {
- *flag = 1;
- }
- else
- {
- *flag = 0;
- }
-}
-
-int TempoTrack::findMeter(double *ACF, unsigned int len, double period)
-{
- int i;
- int p = (int)MathUtilities::round( period );
- int tsig;
-
- double Energy_3 = 0.0;
- double Energy_4 = 0.0;
-
- double temp3A = 0.0;
- double temp3B = 0.0;
- double temp4A = 0.0;
- double temp4B = 0.0;
-
- double* dbf = new double[ len ]; int t = 0;
- for( unsigned int u = 0; u < len; u++ ){ dbf[ u ] = 0.0; }
-
- if( (double)len < 6 * p + 2 )
- {
- for( i = ( 3 * p - 2 ); i < ( 3 * p + 2 ) + 1; i++ )
- {
- temp3A += ACF[ i ];
- dbf[ t++ ] = ACF[ i ];
- }
-
- for( i = ( 4 * p - 2 ); i < ( 4 * p + 2 ) + 1; i++ )
- {
- temp4A += ACF[ i ];
- }
-
- Energy_3 = temp3A;
- Energy_4 = temp4A;
- }
- else
- {
- for( i = ( 3 * p - 2 ); i < ( 3 * p + 2 ) + 1; i++ )
- {
- temp3A += ACF[ i ];
- }
-
- for( i = ( 4 * p - 2 ); i < ( 4 * p + 2 ) + 1; i++ )
- {
- temp4A += ACF[ i ];
- }
-
- for( i = ( 6 * p - 2 ); i < ( 6 * p + 2 ) + 1; i++ )
- {
- temp3B += ACF[ i ];
- }
-
- for( i = ( 2 * p - 2 ); i < ( 2 * p + 2 ) + 1; i++ )
- {
- temp4B += ACF[ i ];
- }
-
- Energy_3 = temp3A + temp3B;
- Energy_4 = temp4A + temp4B;
- }
-
- if (Energy_3 > Energy_4)
- {
- tsig = 3;
- }
- else
- {
- tsig = 4;
- }
-
-
- return tsig;
-}
-
-void TempoTrack::createPhaseExtractor(double *Filter, unsigned int winLength, double period, unsigned int fsp, unsigned int lastBeat)
-{
- int p = (int)MathUtilities::round( period );
- int predictedOffset = 0;
-
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "TempoTrack::createPhaseExtractor: period = " << period << ", p = " << p << std::endl;
-#endif
-
- if (p > 10000) {
- std::cerr << "TempoTrack::createPhaseExtractor: WARNING! Highly implausible period value " << p << "!" << std::endl;
- period = 5168 / 120;
- }
-
- double* phaseScratch = new double[ p*2 + 2 ];
- for (int i = 0; i < p*2 + 2; ++i) phaseScratch[i] = 0.0;
-
-
- if( lastBeat != 0 )
- {
- lastBeat = (int)MathUtilities::round((double)lastBeat );///(double)winLength);
-
- predictedOffset = lastBeat + p - fsp;
-
- if (predictedOffset < 0)
- {
- lastBeat = 0;
- }
- }
-
- if( lastBeat != 0 )
- {
- int mu = p;
- double sigma = (double)p/8;
- double PhaseMin = 0.0;
- double PhaseMax = 0.0;
- unsigned int scratchLength = p*2;
- double temp = 0.0;
-
- for( int i = 0; i < scratchLength; i++ )
- {
- phaseScratch[ i ] = exp( -0.5 * pow( ( i - mu ) / sigma, 2 ) ) / ( sqrt( 2*PI ) *sigma );
- }
-
- MathUtilities::getFrameMinMax( phaseScratch, scratchLength, &PhaseMin, &PhaseMax );
-
- for(int i = 0; i < scratchLength; i ++)
- {
- temp = phaseScratch[ i ];
- phaseScratch[ i ] = (temp - PhaseMin)/PhaseMax;
- }
-
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "predictedOffset = " << predictedOffset << std::endl;
-#endif
-
- unsigned int index = 0;
- for (int i = p - ( predictedOffset - 1); i < p + ( p - predictedOffset) + 1; i++)
- {
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "assigning to filter index " << index << " (size = " << p*2 << ")" << " value " << phaseScratch[i] << " from scratch index " << i << std::endl;
-#endif
- Filter[ index++ ] = phaseScratch[ i ];
- }
- }
- else
- {
- for( int i = 0; i < p; i ++)
- {
- Filter[ i ] = 1;
- }
- }
-
- delete [] phaseScratch;
-}
-
-int TempoTrack::phaseMM(double *DF, double *weighting, unsigned int winLength, double period)
-{
- int alignment = 0;
- int p = (int)MathUtilities::round( period );
-
- double temp = 0.0;
-
- double* y = new double[ winLength ];
- double* align = new double[ p ];
-
- for( int i = 0; i < winLength; i++ )
- {
- y[ i ] = (double)( -i + winLength )/(double)winLength;
- y[ i ] = pow(y [i ],2.0); // raise to power 2.
- }
-
- for( int o = 0; o < p; o++ )
- {
- temp = 0.0;
- for(int i = 1 + (o - 1); i< winLength; i += (p + 1))
- {
- temp = temp + DF[ i ] * y[ i ];
- }
- align[ o ] = temp * weighting[ o ];
- }
-
-
- double valTemp = 0.0;
- for(int i = 0; i < p; i++)
- {
- if( align[ i ] > valTemp )
- {
- valTemp = align[ i ];
- alignment = i;
- }
- }
-
- delete [] y;
- delete [] align;
-
- return alignment;
-}
-
-int TempoTrack::beatPredict(unsigned int FSP0, double alignment, double period, unsigned int step )
-{
- int beat = 0;
-
- int p = (int)MathUtilities::round( period );
- int align = (int)MathUtilities::round( alignment );
- int FSP = (int)MathUtilities::round( FSP0 );
-
- int FEP = FSP + ( step );
-
- beat = FSP + align;
-
- m_beats.push_back( beat );
-
- while( beat + p < FEP )
- {
- beat += p;
-
- m_beats.push_back( beat );
- }
-
- return beat;
-}
-
-
-
-vector<int> TempoTrack::process( vector <double> DF,
- vector <double> *tempoReturn )
-{
- m_dataLength = DF.size();
-
- m_lockedTempo = 0.0;
-
- double period = 0.0;
- int stepFlag = 0;
- int constFlag = 0;
- int FSP = 0;
- int tsig = 0;
- int lastBeat = 0;
-
- vector <double> causalDF;
-
- causalDF = DF;
-
- //Prepare Causal Extension DFData
- unsigned int DFCLength = m_dataLength + m_winLength;
-
- for( unsigned int j = 0; j < m_winLength; j++ )
- {
- causalDF.push_back( 0 );
- }
-
-
- double* RW = new double[ m_lagLength ];
- for( unsigned int clear = 0; clear < m_lagLength; clear++){ RW[ clear ] = 0.0;}
-
- double* GW = new double[ m_lagLength ];
- for(unsigned int clear = 0; clear < m_lagLength; clear++){ GW[ clear ] = 0.0;}
-
- double* PW = new double[ m_lagLength ];
- for(unsigned clear = 0; clear < m_lagLength; clear++){ PW[ clear ] = 0.0;}
-
- m_DFFramer.setSource( &causalDF[0], m_dataLength );
-
- unsigned int TTFrames = m_DFFramer.getMaxNoFrames();
-
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "TTFrames = " << TTFrames << std::endl;
-#endif
-
- double* periodP = new double[ TTFrames ];
- for(unsigned clear = 0; clear < TTFrames; clear++){ periodP[ clear ] = 0.0;}
-
- double* periodG = new double[ TTFrames ];
- for(unsigned clear = 0; clear < TTFrames; clear++){ periodG[ clear ] = 0.0;}
-
- double* alignment = new double[ TTFrames ];
- for(unsigned clear = 0; clear < TTFrames; clear++){ alignment[ clear ] = 0.0;}
-
- m_beats.clear();
-
- createCombFilter( RW, m_lagLength, 0, 0 );
-
- int TTLoopIndex = 0;
-
- for( unsigned int i = 0; i < TTFrames; i++ )
- {
- m_DFFramer.getFrame( m_rawDFFrame );
-
- m_DFConditioning->process( m_rawDFFrame, m_smoothDFFrame );
-
- m_correlator.doAutoUnBiased( m_smoothDFFrame, m_frameACF, m_winLength );
-
- periodP[ TTLoopIndex ] = tempoMM( m_frameACF, RW, 0 );
-
- if( GW[ 0 ] != 0 )
- {
- periodG[ TTLoopIndex ] = tempoMM( m_frameACF, GW, tsig );
- }
- else
- {
- periodG[ TTLoopIndex ] = 0.0;
- }
-
- stepDetect( periodP, periodG, TTLoopIndex, &stepFlag );
-
- if( stepFlag == 1)
- {
- constDetect( periodP, TTLoopIndex, &constFlag );
- stepFlag = 0;
- }
- else
- {
- stepFlag -= 1;
- }
-
- if( stepFlag < 0 )
- {
- stepFlag = 0;
- }
-
- if( constFlag != 0)
- {
- tsig = findMeter( m_frameACF, m_winLength, periodP[ TTLoopIndex ] );
-
- createCombFilter( GW, m_lagLength, tsig, periodP[ TTLoopIndex ] );
-
- periodG[ TTLoopIndex ] = tempoMM( m_frameACF, GW, tsig );
-
- period = periodG[ TTLoopIndex ];
-
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "TempoTrack::process: constFlag == " << constFlag << ", TTLoopIndex = " << TTLoopIndex << ", period from periodG = " << period << std::endl;
-#endif
-
- createPhaseExtractor( PW, m_winLength, period, FSP, 0 );
-
- constFlag = 0;
-
- }
- else
- {
- if( GW[ 0 ] != 0 )
- {
- period = periodG[ TTLoopIndex ];
-
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "TempoTrack::process: GW[0] == " << GW[0] << ", TTLoopIndex = " << TTLoopIndex << ", period from periodG = " << period << std::endl;
-#endif
-
- if (period > 10000) {
- std::cerr << "TempoTrack::process: WARNING! Highly implausible period value " << period << "!" << std::endl;
- std::cerr << "periodG contains (of " << TTFrames << " frames): " << std::endl;
- for (int i = 0; i < TTLoopIndex + 3 && i < TTFrames; ++i) {
- std::cerr << i << " -> " << periodG[i] << std::endl;
- }
- std::cerr << "periodP contains (of " << TTFrames << " frames): " << std::endl;
- for (int i = 0; i < TTLoopIndex + 3 && i < TTFrames; ++i) {
- std::cerr << i << " -> " << periodP[i] << std::endl;
- }
- period = 5168 / 120;
- }
-
- createPhaseExtractor( PW, m_winLength, period, FSP, lastBeat );
-
- }
- else
- {
- period = periodP[ TTLoopIndex ];
-
-#ifdef DEBUG_TEMPO_TRACK
- std::cerr << "TempoTrack::process: GW[0] == " << GW[0] << ", TTLoopIndex = " << TTLoopIndex << ", period from periodP = " << period << std::endl;
-#endif
-
- createPhaseExtractor( PW, m_winLength, period, FSP, 0 );
- }
- }
-
- alignment[ TTLoopIndex ] = phaseMM( m_rawDFFrame, PW, m_winLength, period );
-
- lastBeat = beatPredict(FSP, alignment[ TTLoopIndex ], period, m_lagLength );
-
- FSP += (m_lagLength);
-
- if (tempoReturn) tempoReturn->push_back(m_lockedTempo);
-
- TTLoopIndex++;
- }
-
-
- delete [] periodP;
- delete [] periodG;
- delete [] alignment;
-
- delete [] RW;
- delete [] GW;
- delete [] PW;
-
- return m_beats;
-}
+ COPYING included with this distribution for more information.
+*/
+
+#include "TempoTrack.h"
+
+#include "maths/MathAliases.h"
+#include "maths/MathUtilities.h"
+
+#include <iostream>
+
+#include <cassert>
+
+//#define DEBUG_TEMPO_TRACK 1
+
+
+#define RAY43VAL
+
+//////////////////////////////////////////////////////////////////////
+// Construction/Destruction
+//////////////////////////////////////////////////////////////////////
+
+TempoTrack::TempoTrack( TTParams Params )
+{
+ m_tempoScratch = NULL;
+ m_rawDFFrame = NULL;
+ m_smoothDFFrame = NULL;
+ m_frameACF = NULL;
+ m_smoothRCF = NULL;
+
+ m_dataLength = 0;
+ m_winLength = 0;
+ m_lagLength = 0;
+
+ m_rayparam = 0;
+ m_sigma = 0;
+ m_DFWVNnorm = 0;
+
+ initialise( Params );
+}
+
+TempoTrack::~TempoTrack()
+{
+ deInitialise();
+}
+
+void TempoTrack::initialise( TTParams Params )
+{
+ m_winLength = Params.winLength;
+ m_lagLength = Params.lagLength;
+
+ m_rayparam = 43.0;
+ m_sigma = sqrt(3.9017);
+ m_DFWVNnorm = exp( ( log( 2.0 ) / m_rayparam ) * ( m_winLength + 2 ) );
+
+ m_rawDFFrame = new double[ m_winLength ];
+ m_smoothDFFrame = new double[ m_winLength ];
+ m_frameACF = new double[ m_winLength ];
+ m_tempoScratch = new double[ m_lagLength ];
+ m_smoothRCF = new double[ m_lagLength ];
+
+
+ unsigned int winPre = Params.WinT.pre;
+ unsigned int winPost = Params.WinT.post;
+
+ m_DFFramer.configure( m_winLength, m_lagLength );
+
+ m_DFPParams.length = m_winLength;
+ m_DFPParams.AlphaNormParam = Params.alpha;
+ m_DFPParams.LPOrd = Params.LPOrd;
+ m_DFPParams.LPACoeffs = Params.LPACoeffs;
+ m_DFPParams.LPBCoeffs = Params.LPBCoeffs;
+ m_DFPParams.winPre = Params.WinT.pre;
+ m_DFPParams.winPost = Params.WinT.post;
+ m_DFPParams.isMedianPositive = true;
+
+ m_DFConditioning = new DFProcess( m_DFPParams );
+
+
+ // these are parameters for smoothing m_tempoScratch
+ m_RCFPParams.length = m_lagLength;
+ m_RCFPParams.AlphaNormParam = Params.alpha;
+ m_RCFPParams.LPOrd = Params.LPOrd;
+ m_RCFPParams.LPACoeffs = Params.LPACoeffs;
+ m_RCFPParams.LPBCoeffs = Params.LPBCoeffs;
+ m_RCFPParams.winPre = Params.WinT.pre;
+ m_RCFPParams.winPost = Params.WinT.post;
+ m_RCFPParams.isMedianPositive = true;
+
+ m_RCFConditioning = new DFProcess( m_RCFPParams );
+
+}
+
+void TempoTrack::deInitialise()
+{
+ delete [] m_rawDFFrame;
+
+ delete [] m_smoothDFFrame;
+
+ delete [] m_smoothRCF;
+
+ delete [] m_frameACF;
+
+ delete [] m_tempoScratch;
+
+ delete m_DFConditioning;
+
+ delete m_RCFConditioning;
+
+}
+
+void TempoTrack::createCombFilter(double* Filter, unsigned int winLength, unsigned int TSig, double beatLag)
+{
+ unsigned int i;
+
+ if( beatLag == 0 )
+ {
+ for( i = 0; i < winLength; i++ )
+ {
+ Filter[ i ] = ( ( i + 1 ) / pow( m_rayparam, 2.0) ) * exp( ( -pow(( i + 1 ),2.0 ) / ( 2.0 * pow( m_rayparam, 2.0))));
+ }
+ }
+ else
+ {
+ m_sigma = beatLag/4;
+ for( i = 0; i < winLength; i++ )
+ {
+ double dlag = (double)(i+1) - beatLag;
+ Filter[ i ] = exp(-0.5 * pow(( dlag / m_sigma), 2.0) ) / (sqrt( 2 * PI) * m_sigma);
+ }
+ }
+}
+
+double TempoTrack::tempoMM(double* ACF, double* weight, int tsig)
+{
+
+ double period = 0;
+ double maxValRCF = 0.0;
+ unsigned int maxIndexRCF = 0;
+
+ double* pdPeaks;
+
+ unsigned int maxIndexTemp;
+ double maxValTemp;
+ unsigned int count;
+
+ unsigned int numelem,i,j;
+ int a, b;
+
+ for( i = 0; i < m_lagLength; i++ )
+ m_tempoScratch[ i ] = 0.0;
+
+ if( tsig == 0 )
+ {
+ //if time sig is unknown, use metrically unbiased version of Filterbank
+ numelem = 4;
+ }
+ else
+ {
+ numelem = tsig;
+ }
+
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "tempoMM: m_winLength = " << m_winLength << ", m_lagLength = " << m_lagLength << ", numelem = " << numelem << std::endl;
+#endif
+
+ for(i=1;i<m_lagLength-1;i++)
+ {
+ //first and last output values are left intentionally as zero
+ for (a=1;a<=numelem;a++)
+ {
+ for(b=(1-a);b<a;b++)
+ {
+ if( tsig == 0 )
+ {
+ m_tempoScratch[i] += ACF[a*(i+1)+b-1] * (1.0 / (2.0 * (double)a-1)) * weight[i];
+ }
+ else
+ {
+ m_tempoScratch[i] += ACF[a*(i+1)+b-1] * 1 * weight[i];
+ }
+ }
+ }
+ }
+
+
+ //////////////////////////////////////////////////
+ // MODIFIED BEAT PERIOD EXTRACTION //////////////
+ /////////////////////////////////////////////////
+
+ // find smoothed version of RCF ( as applied to Detection Function)
+ m_RCFConditioning->process( m_tempoScratch, m_smoothRCF);
+
+ if (tsig != 0) // i.e. in context dependent state
+ {
+// NOW FIND MAX INDEX OF ACFOUT
+ for( i = 0; i < m_lagLength; i++)
+ {
+ if( m_tempoScratch[ i ] > maxValRCF)
+ {
+ maxValRCF = m_tempoScratch[ i ];
+ maxIndexRCF = i;
+ }
+ }
+ }
+ else // using rayleigh weighting
+ {
+ vector <vector<double> > rcfMat;
+
+ double sumRcf = 0.;
+
+ double maxVal = 0.;
+ // now find the two values which minimise rcfMat
+ double minVal = 0.;
+ int p_i = 1; // periodicity for row i;
+ int p_j = 1; //periodicity for column j;
+
+
+ for ( i=0; i<m_lagLength; i++)
+ {
+ m_tempoScratch[i] =m_smoothRCF[i];
+ }
+
+ // normalise m_tempoScratch so that it sums to zero.
+ for ( i=0; i<m_lagLength; i++)
+ {
+ sumRcf += m_tempoScratch[i];
+ }
+
+ for( i=0; i<m_lagLength; i++)
+ {
+ m_tempoScratch[i] /= sumRcf;
+ }
+
+ // create a matrix to store m_tempoScratchValues modified by log2 ratio
+ for ( i=0; i<m_lagLength; i++)
+ {
+ rcfMat.push_back ( vector<double>() ); // adds a new row...
+ }
+
+ for (i=0; i<m_lagLength; i++)
+ {
+ for (j=0; j<m_lagLength; j++)
+ {
+ rcfMat[i].push_back (0.);
+ }
+ }
+
+ // the 'i' and 'j' indices deliberately start from '1' and not '0'
+ for ( i=1; i<m_lagLength; i++)
+ {
+ for (j=1; j<m_lagLength; j++)
+ {
+ double log2PeriodRatio = log( static_cast<double>(i)/static_cast<double>(j) ) / log(2.0);
+ rcfMat[i][j] = ( abs(1.0-abs(log2PeriodRatio)) );
+ rcfMat[i][j] += ( 0.01*( 1./(m_tempoScratch[i]+m_tempoScratch[j]) ) );
+ }
+ }
+
+ // set diagonal equal to maximum value in rcfMat
+ // we don't want to pick one strong middle peak - we need a combination of two peaks.
+
+ for ( i=1; i<m_lagLength; i++)
+ {
+ for (j=1; j<m_lagLength; j++)
+ {
+ if (rcfMat[i][j] > maxVal)
+ {
+ maxVal = rcfMat[i][j];
+ }
+ }
+ }
+
+ for ( i=1; i<m_lagLength; i++)
+ {
+ rcfMat[i][i] = maxVal;
+ }
+
+ // now find the row and column number which minimise rcfMat
+ minVal = maxVal;
+
+ for ( i=1; i<m_lagLength; i++)
+ {
+ for ( j=1; j<m_lagLength; j++)
+ {
+ if (rcfMat[i][j] < minVal)
+ {
+ minVal = rcfMat[i][j];
+ p_i = i;
+ p_j = j;
+ }
+ }
+ }
+
+
+ // initially choose p_j (arbitrary) - saves on an else statement
+ int beatPeriod = p_j;
+ if (m_tempoScratch[p_i] > m_tempoScratch[p_j])
+ {
+ beatPeriod = p_i;
+ }
+
+ // now write the output
+ maxIndexRCF = static_cast<int>(beatPeriod);
+ }
+
+
+ double locked = 5168.f / maxIndexRCF;
+ if (locked >= 30 && locked <= 180) {
+ m_lockedTempo = locked;
+ }
+
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "tempoMM: locked tempo = " << m_lockedTempo << std::endl;
+#endif
+
+ if( tsig == 0 )
+ tsig = 4;
+
+
+#ifdef DEBUG_TEMPO_TRACK
+std::cerr << "tempoMM: maxIndexRCF = " << maxIndexRCF << std::endl;
+#endif
+
+ if( tsig == 4 )
+ {
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "tsig == 4" << std::endl;
+#endif
+
+ pdPeaks = new double[ 4 ];
+ for( i = 0; i < 4; i++ ){ pdPeaks[ i ] = 0.0;}
+
+ pdPeaks[ 0 ] = ( double )maxIndexRCF + 1;
+
+ maxIndexTemp = 0;
+ maxValTemp = 0.0;
+ count = 0;
+
+ for( i = (2 * maxIndexRCF + 1) - 1; i < (2 * maxIndexRCF + 1) + 2; i++ )
+ {
+ if( ACF[ i ] > maxValTemp )
+ {
+ maxValTemp = ACF[ i ];
+ maxIndexTemp = count;
+ }
+ count++;
+ }
+ pdPeaks[ 1 ] = (double)( maxIndexTemp + 1 + ( (2 * maxIndexRCF + 1 ) - 2 ) + 1 )/2;
+
+ maxIndexTemp = 0;
+ maxValTemp = 0.0;
+ count = 0;
+
+ for( i = (3 * maxIndexRCF + 2 ) - 2; i < (3 * maxIndexRCF + 2 ) + 3; i++ )
+ {
+ if( ACF[ i ] > maxValTemp )
+ {
+ maxValTemp = ACF[ i ];
+ maxIndexTemp = count;
+ }
+ count++;
+ }
+ pdPeaks[ 2 ] = (double)( maxIndexTemp + 1 + ( (3 * maxIndexRCF + 2) - 4 ) + 1 )/3;
+
+ maxIndexTemp = 0;
+ maxValTemp = 0.0;
+ count = 0;
+
+ for( i = ( 4 * maxIndexRCF + 3) - 3; i < ( 4 * maxIndexRCF + 3) + 4; i++ )
+ {
+ if( ACF[ i ] > maxValTemp )
+ {
+ maxValTemp = ACF[ i ];
+ maxIndexTemp = count;
+ }
+ count++;
+ }
+ pdPeaks[ 3 ] = (double)( maxIndexTemp + 1 + ( (4 * maxIndexRCF + 3) - 9 ) + 1 )/4 ;
+
+
+ period = MathUtilities::mean( pdPeaks, 4 );
+ }
+ else
+ {
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "tsig != 4" << std::endl;
+#endif
+
+ pdPeaks = new double[ 3 ];
+ for( i = 0; i < 3; i++ ){ pdPeaks[ i ] = 0.0;}
+
+ pdPeaks[ 0 ] = ( double )maxIndexRCF + 1;
+
+ maxIndexTemp = 0;
+ maxValTemp = 0.0;
+ count = 0;
+
+ for( i = (2 * maxIndexRCF + 1) - 1; i < (2 * maxIndexRCF + 1) + 2; i++ )
+ {
+ if( ACF[ i ] > maxValTemp )
+ {
+ maxValTemp = ACF[ i ];
+ maxIndexTemp = count;
+ }
+ count++;
+ }
+ pdPeaks[ 1 ] = (double)( maxIndexTemp + 1 + ( (2 * maxIndexRCF + 1 ) - 2 ) + 1 )/2;
+
+ maxIndexTemp = 0;
+ maxValTemp = 0.0;
+ count = 0;
+
+ for( i = (3 * maxIndexRCF + 2 ) - 2; i < (3 * maxIndexRCF + 2 ) + 3; i++ )
+ {
+ if( ACF[ i ] > maxValTemp )
+ {
+ maxValTemp = ACF[ i ];
+ maxIndexTemp = count;
+ }
+ count++;
+ }
+ pdPeaks[ 2 ] = (double)( maxIndexTemp + 1 + ( (3 * maxIndexRCF + 2) - 4 ) + 1 )/3;
+
+
+ period = MathUtilities::mean( pdPeaks, 3 );
+ }
+
+ delete [] pdPeaks;
+
+ return period;
+}
+
+void TempoTrack::stepDetect( double* periodP, double* periodG, int currentIdx, int* flag )
+{
+ double stepthresh = 1 * 3.9017;
+
+ if( *flag )
+ {
+ if(abs(periodG[ currentIdx ] - periodP[ currentIdx ]) > stepthresh)
+ {
+ // do nuffin'
+ }
+ }
+ else
+ {
+ if(fabs(periodG[ currentIdx ]-periodP[ currentIdx ]) > stepthresh)
+ {
+ *flag = 3;
+ }
+ }
+}
+
+void TempoTrack::constDetect( double* periodP, int currentIdx, int* flag )
+{
+ double constthresh = 2 * 3.9017;
+
+ if( fabs( 2 * periodP[ currentIdx ] - periodP[ currentIdx - 1] - periodP[ currentIdx - 2] ) < constthresh)
+ {
+ *flag = 1;
+ }
+ else
+ {
+ *flag = 0;
+ }
+}
+
+int TempoTrack::findMeter(double *ACF, unsigned int len, double period)
+{
+ int i;
+ int p = (int)MathUtilities::round( period );
+ int tsig;
+
+ double Energy_3 = 0.0;
+ double Energy_4 = 0.0;
+
+ double temp3A = 0.0;
+ double temp3B = 0.0;
+ double temp4A = 0.0;
+ double temp4B = 0.0;
+
+ double* dbf = new double[ len ]; int t = 0;
+ for( unsigned int u = 0; u < len; u++ ){ dbf[ u ] = 0.0; }
+
+ if( (double)len < 6 * p + 2 )
+ {
+ for( i = ( 3 * p - 2 ); i < ( 3 * p + 2 ) + 1; i++ )
+ {
+ temp3A += ACF[ i ];
+ dbf[ t++ ] = ACF[ i ];
+ }
+
+ for( i = ( 4 * p - 2 ); i < ( 4 * p + 2 ) + 1; i++ )
+ {
+ temp4A += ACF[ i ];
+ }
+
+ Energy_3 = temp3A;
+ Energy_4 = temp4A;
+ }
+ else
+ {
+ for( i = ( 3 * p - 2 ); i < ( 3 * p + 2 ) + 1; i++ )
+ {
+ temp3A += ACF[ i ];
+ }
+
+ for( i = ( 4 * p - 2 ); i < ( 4 * p + 2 ) + 1; i++ )
+ {
+ temp4A += ACF[ i ];
+ }
+
+ for( i = ( 6 * p - 2 ); i < ( 6 * p + 2 ) + 1; i++ )
+ {
+ temp3B += ACF[ i ];
+ }
+
+ for( i = ( 2 * p - 2 ); i < ( 2 * p + 2 ) + 1; i++ )
+ {
+ temp4B += ACF[ i ];
+ }
+
+ Energy_3 = temp3A + temp3B;
+ Energy_4 = temp4A + temp4B;
+ }
+
+ if (Energy_3 > Energy_4)
+ {
+ tsig = 3;
+ }
+ else
+ {
+ tsig = 4;
+ }
+
+
+ return tsig;
+}
+
+void TempoTrack::createPhaseExtractor(double *Filter, unsigned int winLength, double period, unsigned int fsp, unsigned int lastBeat)
+{
+ int p = (int)MathUtilities::round( period );
+ int predictedOffset = 0;
+
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "TempoTrack::createPhaseExtractor: period = " << period << ", p = " << p << std::endl;
+#endif
+
+ if (p > 10000) {
+ std::cerr << "TempoTrack::createPhaseExtractor: WARNING! Highly implausible period value " << p << "!" << std::endl;
+ period = 5168 / 120;
+ }
+
+ double* phaseScratch = new double[ p*2 + 2 ];
+ for (int i = 0; i < p*2 + 2; ++i) phaseScratch[i] = 0.0;
+
+
+ if( lastBeat != 0 )
+ {
+ lastBeat = (int)MathUtilities::round((double)lastBeat );///(double)winLength);
+
+ predictedOffset = lastBeat + p - fsp;
+
+ if (predictedOffset < 0)
+ {
+ lastBeat = 0;
+ }
+ }
+
+ if( lastBeat != 0 )
+ {
+ int mu = p;
+ double sigma = (double)p/8;
+ double PhaseMin = 0.0;
+ double PhaseMax = 0.0;
+ unsigned int scratchLength = p*2;
+ double temp = 0.0;
+
+ for( int i = 0; i < scratchLength; i++ )
+ {
+ phaseScratch[ i ] = exp( -0.5 * pow( ( i - mu ) / sigma, 2 ) ) / ( sqrt( 2*PI ) *sigma );
+ }
+
+ MathUtilities::getFrameMinMax( phaseScratch, scratchLength, &PhaseMin, &PhaseMax );
+
+ for(int i = 0; i < scratchLength; i ++)
+ {
+ temp = phaseScratch[ i ];
+ phaseScratch[ i ] = (temp - PhaseMin)/PhaseMax;
+ }
+
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "predictedOffset = " << predictedOffset << std::endl;
+#endif
+
+ unsigned int index = 0;
+ for (int i = p - ( predictedOffset - 1); i < p + ( p - predictedOffset) + 1; i++)
+ {
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "assigning to filter index " << index << " (size = " << p*2 << ")" << " value " << phaseScratch[i] << " from scratch index " << i << std::endl;
+#endif
+ Filter[ index++ ] = phaseScratch[ i ];
+ }
+ }
+ else
+ {
+ for( int i = 0; i < p; i ++)
+ {
+ Filter[ i ] = 1;
+ }
+ }
+
+ delete [] phaseScratch;
+}
+
+int TempoTrack::phaseMM(double *DF, double *weighting, unsigned int winLength, double period)
+{
+ int alignment = 0;
+ int p = (int)MathUtilities::round( period );
+
+ double temp = 0.0;
+
+ double* y = new double[ winLength ];
+ double* align = new double[ p ];
+
+ for( int i = 0; i < winLength; i++ )
+ {
+ y[ i ] = (double)( -i + winLength )/(double)winLength;
+ y[ i ] = pow(y [i ],2.0); // raise to power 2.
+ }
+
+ for( int o = 0; o < p; o++ )
+ {
+ temp = 0.0;
+ for(int i = 1 + (o - 1); i< winLength; i += (p + 1))
+ {
+ temp = temp + DF[ i ] * y[ i ];
+ }
+ align[ o ] = temp * weighting[ o ];
+ }
+
+
+ double valTemp = 0.0;
+ for(int i = 0; i < p; i++)
+ {
+ if( align[ i ] > valTemp )
+ {
+ valTemp = align[ i ];
+ alignment = i;
+ }
+ }
+
+ delete [] y;
+ delete [] align;
+
+ return alignment;
+}
+
+int TempoTrack::beatPredict(unsigned int FSP0, double alignment, double period, unsigned int step )
+{
+ int beat = 0;
+
+ int p = (int)MathUtilities::round( period );
+ int align = (int)MathUtilities::round( alignment );
+ int FSP = (int)MathUtilities::round( FSP0 );
+
+ int FEP = FSP + ( step );
+
+ beat = FSP + align;
+
+ m_beats.push_back( beat );
+
+ while( beat + p < FEP )
+ {
+ beat += p;
+
+ m_beats.push_back( beat );
+ }
+
+ return beat;
+}
+
+
+
+vector<int> TempoTrack::process( vector <double> DF,
+ vector <double> *tempoReturn )
+{
+ m_dataLength = DF.size();
+
+ m_lockedTempo = 0.0;
+
+ double period = 0.0;
+ int stepFlag = 0;
+ int constFlag = 0;
+ int FSP = 0;
+ int tsig = 0;
+ int lastBeat = 0;
+
+ vector <double> causalDF;
+
+ causalDF = DF;
+
+ //Prepare Causal Extension DFData
+ unsigned int DFCLength = m_dataLength + m_winLength;
+
+ for( unsigned int j = 0; j < m_winLength; j++ )
+ {
+ causalDF.push_back( 0 );
+ }
+
+
+ double* RW = new double[ m_lagLength ];
+ for( unsigned int clear = 0; clear < m_lagLength; clear++){ RW[ clear ] = 0.0;}
+
+ double* GW = new double[ m_lagLength ];
+ for(unsigned int clear = 0; clear < m_lagLength; clear++){ GW[ clear ] = 0.0;}
+
+ double* PW = new double[ m_lagLength ];
+ for(unsigned clear = 0; clear < m_lagLength; clear++){ PW[ clear ] = 0.0;}
+
+ m_DFFramer.setSource( &causalDF[0], m_dataLength );
+
+ unsigned int TTFrames = m_DFFramer.getMaxNoFrames();
+
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "TTFrames = " << TTFrames << std::endl;
+#endif
+
+ double* periodP = new double[ TTFrames ];
+ for(unsigned clear = 0; clear < TTFrames; clear++){ periodP[ clear ] = 0.0;}
+
+ double* periodG = new double[ TTFrames ];
+ for(unsigned clear = 0; clear < TTFrames; clear++){ periodG[ clear ] = 0.0;}
+
+ double* alignment = new double[ TTFrames ];
+ for(unsigned clear = 0; clear < TTFrames; clear++){ alignment[ clear ] = 0.0;}
+
+ m_beats.clear();
+
+ createCombFilter( RW, m_lagLength, 0, 0 );
+
+ int TTLoopIndex = 0;
+
+ for( unsigned int i = 0; i < TTFrames; i++ )
+ {
+ m_DFFramer.getFrame( m_rawDFFrame );
+
+ m_DFConditioning->process( m_rawDFFrame, m_smoothDFFrame );
+
+ m_correlator.doAutoUnBiased( m_smoothDFFrame, m_frameACF, m_winLength );
+
+ periodP[ TTLoopIndex ] = tempoMM( m_frameACF, RW, 0 );
+
+ if( GW[ 0 ] != 0 )
+ {
+ periodG[ TTLoopIndex ] = tempoMM( m_frameACF, GW, tsig );
+ }
+ else
+ {
+ periodG[ TTLoopIndex ] = 0.0;
+ }
+
+ stepDetect( periodP, periodG, TTLoopIndex, &stepFlag );
+
+ if( stepFlag == 1)
+ {
+ constDetect( periodP, TTLoopIndex, &constFlag );
+ stepFlag = 0;
+ }
+ else
+ {
+ stepFlag -= 1;
+ }
+
+ if( stepFlag < 0 )
+ {
+ stepFlag = 0;
+ }
+
+ if( constFlag != 0)
+ {
+ tsig = findMeter( m_frameACF, m_winLength, periodP[ TTLoopIndex ] );
+
+ createCombFilter( GW, m_lagLength, tsig, periodP[ TTLoopIndex ] );
+
+ periodG[ TTLoopIndex ] = tempoMM( m_frameACF, GW, tsig );
+
+ period = periodG[ TTLoopIndex ];
+
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "TempoTrack::process: constFlag == " << constFlag << ", TTLoopIndex = " << TTLoopIndex << ", period from periodG = " << period << std::endl;
+#endif
+
+ createPhaseExtractor( PW, m_winLength, period, FSP, 0 );
+
+ constFlag = 0;
+
+ }
+ else
+ {
+ if( GW[ 0 ] != 0 )
+ {
+ period = periodG[ TTLoopIndex ];
+
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "TempoTrack::process: GW[0] == " << GW[0] << ", TTLoopIndex = " << TTLoopIndex << ", period from periodG = " << period << std::endl;
+#endif
+
+ if (period > 10000) {
+ std::cerr << "TempoTrack::process: WARNING! Highly implausible period value " << period << "!" << std::endl;
+ std::cerr << "periodG contains (of " << TTFrames << " frames): " << std::endl;
+ for (int i = 0; i < TTLoopIndex + 3 && i < TTFrames; ++i) {
+ std::cerr << i << " -> " << periodG[i] << std::endl;
+ }
+ std::cerr << "periodP contains (of " << TTFrames << " frames): " << std::endl;
+ for (int i = 0; i < TTLoopIndex + 3 && i < TTFrames; ++i) {
+ std::cerr << i << " -> " << periodP[i] << std::endl;
+ }
+ period = 5168 / 120;
+ }
+
+ createPhaseExtractor( PW, m_winLength, period, FSP, lastBeat );
+
+ }
+ else
+ {
+ period = periodP[ TTLoopIndex ];
+
+#ifdef DEBUG_TEMPO_TRACK
+ std::cerr << "TempoTrack::process: GW[0] == " << GW[0] << ", TTLoopIndex = " << TTLoopIndex << ", period from periodP = " << period << std::endl;
+#endif
+
+ createPhaseExtractor( PW, m_winLength, period, FSP, 0 );
+ }
+ }
+
+ alignment[ TTLoopIndex ] = phaseMM( m_rawDFFrame, PW, m_winLength, period );
+
+ lastBeat = beatPredict(FSP, alignment[ TTLoopIndex ], period, m_lagLength );
+
+ FSP += (m_lagLength);
+
+ if (tempoReturn) tempoReturn->push_back(m_lockedTempo);
+
+ TTLoopIndex++;
+ }
+
+
+ delete [] periodP;
+ delete [] periodG;
+ delete [] alignment;
+
+ delete [] RW;
+ delete [] GW;
+ delete [] PW;
+
+ return m_beats;
+}
+
diff --git a/libs/qm-dsp/dsp/tempotracking/TempoTrack.h b/libs/qm-dsp/dsp/tempotracking/TempoTrack.h
index 12eb977cd5..0c315717ba 100644
--- a/libs/qm-dsp/dsp/tempotracking/TempoTrack.h
+++ b/libs/qm-dsp/dsp/tempotracking/TempoTrack.h
@@ -5,11 +5,11 @@
Centre for Digital Music, Queen Mary, University of London.
This file 2005-2006 Christian Landone.
-
- This program is free software; you can redistribute it and/or
- modify it under the terms of the GNU General Public License as
- published by the Free Software Foundation; either version 2 of the
- License, or (at your option) any later version. See the file
+
+ This program is free software; you can redistribute it and/or
+ modify it under the terms of the GNU General Public License as
+ published by the Free Software Foundation; either version 2 of the
+ License, or (at your option) any later version. See the file
COPYING included with this distribution for more information.
*/
@@ -31,7 +31,7 @@ using std::vector;
struct WinThresh
{
unsigned int pre;
- unsigned int post;
+ unsigned int post;
};
struct TTParams
diff --git a/libs/qm-dsp/dsp/tempotracking/TempoTrackV2.cpp b/libs/qm-dsp/dsp/tempotracking/TempoTrackV2.cpp
index a2c2a24b0a..546693091e 100644
--- a/libs/qm-dsp/dsp/tempotracking/TempoTrackV2.cpp
+++ b/libs/qm-dsp/dsp/tempotracking/TempoTrackV2.cpp
@@ -91,10 +91,17 @@ TempoTrackV2::filter_df(d_vec_t &df)
}
+// MEPD 28/11/12
+// This function now allows for a user to specify an inputtempo (in BPM)
+// and a flag "constraintempo" which replaces the general rayleigh weighting for periodicities
+// with a gaussian which is centered around the input tempo
+// Note, if inputtempo = 120 and constraintempo = false, then functionality is
+// as it was before
void
TempoTrackV2::calculateBeatPeriod(const vector<double> &df,
vector<double> &beat_period,
- vector<double> &tempi)
+ vector<double> &tempi,
+ double inputtempo, bool constraintempo)
{
// to follow matlab.. split into 512 sample frames with a 128 hop size
// calculate the acf,
@@ -103,13 +110,42 @@ TempoTrackV2::calculateBeatPeriod(const vector<double> &df,
// and get best path
unsigned int wv_len = 128;
- double rayparam = 43.;
+
+ // MEPD 28/11/12
+ // the default value of inputtempo in the beat tracking plugin is 120
+ // so if the user specifies a different inputtempo, the rayparam will be updated
+ // accordingly.
+ // note: 60*44100/512 is a magic number
+ // this might (will?) break if a user specifies a different frame rate for the onset detection function
+ double rayparam = (60*44100/512)/inputtempo;
+
+ // these debug statements can be removed.
+// std::cerr << "inputtempo" << inputtempo << std::endl;
+// std::cerr << "rayparam" << rayparam << std::endl;
+// std::cerr << "constraintempo" << constraintempo << std::endl;
// make rayleigh weighting curve
d_vec_t wv(wv_len);
- for (unsigned int i=0; i<wv.size(); i++)
+
+ // check whether or not to use rayleigh weighting (if constraintempo is false)
+ // or use gaussian weighting it (constraintempo is true)
+ if (constraintempo)
{
- wv[i] = (static_cast<double> (i) / pow(rayparam,2.)) * exp((-1.*pow(-static_cast<double> (i),2.)) / (2.*pow(rayparam,2.)));
+ for (unsigned int i=0; i<wv.size(); i++)
+ {
+ // MEPD 28/11/12
+ // do a gaussian weighting instead of rayleigh
+ wv[i] = exp( (-1.*pow((static_cast<double> (i)-rayparam),2.)) / (2.*pow(rayparam/4.,2.)) );
+ }
+ }
+ else
+ {
+ for (unsigned int i=0; i<wv.size(); i++)
+ {
+ // MEPD 28/11/12
+ // standard rayleigh weighting over periodicities
+ wv[i] = (static_cast<double> (i) / pow(rayparam,2.)) * exp((-1.*pow(-static_cast<double> (i),2.)) / (2.*pow(rayparam,2.)));
+ }
}
// beat tracking frame size (roughly 6 seconds) and hop (1.5 seconds)
@@ -397,10 +433,14 @@ TempoTrackV2::normalise_vec(d_vec_t &df)
}
}
+// MEPD 28/11/12
+// this function has been updated to allow the "alpha" and "tightness" parameters
+// of the dynamic program to be set by the user
+// the default value of alpha = 0.9 and tightness = 4
void
TempoTrackV2::calculateBeats(const vector<double> &df,
const vector<double> &beat_period,
- vector<double> &beats)
+ vector<double> &beats, double alpha, double tightness)
{
if (df.empty() || beat_period.empty()) return;
@@ -414,8 +454,12 @@ TempoTrackV2::calculateBeats(const vector<double> &df,
backlink[i] = -1;
}
- double tightness = 4.;
- double alpha = 0.9;
+ //double tightness = 4.;
+ //double alpha = 0.9;
+ // MEPD 28/11/12
+ // debug statements that can be removed.
+// std::cerr << "alpha" << alpha << std::endl;
+// std::cerr << "tightness" << tightness << std::endl;
// main loop
for (unsigned int i=0; i<localscore.size(); i++)
@@ -462,7 +506,7 @@ TempoTrackV2::calculateBeats(const vector<double> &df,
int startpoint = get_max_ind(tmp_vec) + cumscore.size() - beat_period[beat_period.size()-1] ;
// can happen if no results obtained earlier (e.g. input too short)
- if (startpoint >= backlink.size()) startpoint = backlink.size()-1;
+ if (startpoint >= (int)backlink.size()) startpoint = backlink.size()-1;
// USE BACKLINK TO GET EACH NEW BEAT (TOWARDS THE BEGINNING OF THE FILE)
// BACKTRACKING FROM THE END TO THE BEGINNING.. MAKING SURE NOT TO GO BEFORE SAMPLE 0
diff --git a/libs/qm-dsp/dsp/tempotracking/TempoTrackV2.h b/libs/qm-dsp/dsp/tempotracking/TempoTrackV2.h
index 489c500bb2..5f30ba73b1 100644
--- a/libs/qm-dsp/dsp/tempotracking/TempoTrackV2.h
+++ b/libs/qm-dsp/dsp/tempotracking/TempoTrackV2.h
@@ -18,8 +18,7 @@
#define TEMPOTRACKV2_H
#include <vector>
-
-using std::vector;
+using namespace std;
//!!! Question: how far is this actually sample rate dependent? I
// think it does produce plausible results for e.g. 48000 as well as
@@ -40,15 +39,35 @@ public:
TempoTrackV2(float sampleRate, size_t dfIncrement);
~TempoTrackV2();
- // Returned beat periods are given in df increment units; tempi in bpm
+ // Returned beat periods are given in df increment units; inputtempo and tempi in bpm
+ void calculateBeatPeriod(const vector<double> &df,
+ vector<double> &beatPeriod,
+ vector<double> &tempi) {
+ calculateBeatPeriod(df, beatPeriod, tempi, 120.0, false);
+ }
+
+ // Returned beat periods are given in df increment units; inputtempo and tempi in bpm
+ // MEPD 28/11/12 Expose inputtempo and constraintempo parameters
+ // Note, if inputtempo = 120 and constraintempo = false, then functionality is as it was before
void calculateBeatPeriod(const vector<double> &df,
vector<double> &beatPeriod,
- vector<double> &tempi);
+ vector<double> &tempi,
+ double inputtempo, bool constraintempo);
+
+ // Returned beat positions are given in df increment units
+ void calculateBeats(const vector<double> &df,
+ const vector<double> &beatPeriod,
+ vector<double> &beats) {
+ calculateBeats(df, beatPeriod, beats, 0.9, 4.0);
+ }
// Returned beat positions are given in df increment units
+ // MEPD 28/11/12 Expose alpha and tightness parameters
+ // Note, if alpha = 0.9 and tightness = 4, then functionality is as it was before
void calculateBeats(const vector<double> &df,
const vector<double> &beatPeriod,
- vector<double> &beats);
+ vector<double> &beats,
+ double alpha, double tightness);
private:
typedef vector<int> i_vec_t;
generated by cgit v1.2.3 (git 2.39.1) at 2024-06-02 20:51:34 +0000