annotate dsp/tempotracking/TempoTrack.cpp @ 321:f1e6be2de9a5

A threshold (delta) is added in the peak picking parameters structure (PPickParams). It is used as an offset when computing the smoothed detection function. A constructor for the structure PPickParams is also added to set the parameters to 0 when a structure instance is created. Hence programmes using the peak picking parameter structure and which do not set the delta parameter (e.g. QM Vamp note onset detector) won't be affected by the modifications. Functions modified: - dsp/onsets/PeakPicking.cpp - dsp/onsets/PeakPicking.h - dsp/signalconditioning/DFProcess.cpp - dsp/signalconditioning/DFProcess.h
author mathieub <mathieu.barthet@eecs.qmul.ac.uk>
date Mon, 20 Jun 2011 19:01:48 +0100
parents d5014ab8b0e5
children 2ae4ceb76ac3
rev   line source
c@264 1 /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*- vi:set ts=8 sts=4 sw=4: */
c@264 2
c@264 3 /*
c@264 4 QM DSP Library
c@264 5
c@264 6 Centre for Digital Music, Queen Mary, University of London.
c@309 7 This file copyright 2005-2006 Christian Landone.and Matthew Davies.
c@309 8
c@309 9 This program is free software; you can redistribute it and/or
c@309 10 modify it under the terms of the GNU General Public License as
c@309 11 published by the Free Software Foundation; either version 2 of the
c@309 12 License, or (at your option) any later version. See the file
c@309 13 COPYING included with this distribution for more information.
c@264 14 */
c@264 15
c@264 16 #include "TempoTrack.h"
c@264 17
c@264 18 #include "maths/MathAliases.h"
c@264 19 #include "maths/MathUtilities.h"
c@264 20
c@264 21 #include <iostream>
c@264 22
c@271 23 #include <cassert>
c@271 24
c@272 25 //#define DEBUG_TEMPO_TRACK 1
c@272 26
c@264 27
c@264 28 #define RAY43VAL
c@264 29
c@264 30 //////////////////////////////////////////////////////////////////////
c@264 31 // Construction/Destruction
c@264 32 //////////////////////////////////////////////////////////////////////
c@264 33
c@264 34 TempoTrack::TempoTrack( TTParams Params )
c@264 35 {
c@264 36 m_tempoScratch = NULL;
c@264 37 m_rawDFFrame = NULL;
c@264 38 m_smoothDFFrame = NULL;
c@264 39 m_frameACF = NULL;
c@264 40 m_smoothRCF = NULL;
c@264 41
c@264 42 m_dataLength = 0;
c@264 43 m_winLength = 0;
c@264 44 m_lagLength = 0;
c@264 45
c@264 46 m_rayparam = 0;
c@264 47 m_sigma = 0;
c@264 48 m_DFWVNnorm = 0;
c@264 49
c@264 50 initialise( Params );
c@264 51 }
c@264 52
c@264 53 TempoTrack::~TempoTrack()
c@264 54 {
c@264 55 deInitialise();
c@264 56 }
c@264 57
c@264 58 void TempoTrack::initialise( TTParams Params )
c@264 59 {
c@264 60 m_winLength = Params.winLength;
c@264 61 m_lagLength = Params.lagLength;
c@264 62
c@264 63 m_rayparam = 43.0;
c@264 64 m_sigma = sqrt(3.9017);
c@264 65 m_DFWVNnorm = exp( ( log( 2.0 ) / m_rayparam ) * ( m_winLength + 2 ) );
c@264 66
c@264 67 m_rawDFFrame = new double[ m_winLength ];
c@264 68 m_smoothDFFrame = new double[ m_winLength ];
c@264 69 m_frameACF = new double[ m_winLength ];
c@264 70 m_tempoScratch = new double[ m_lagLength ];
c@264 71 m_smoothRCF = new double[ m_lagLength ];
c@264 72
c@264 73
c@264 74 unsigned int winPre = Params.WinT.pre;
c@264 75 unsigned int winPost = Params.WinT.post;
c@264 76
c@264 77 m_DFFramer.configure( m_winLength, m_lagLength );
c@264 78
c@264 79 m_DFPParams.length = m_winLength;
c@264 80 m_DFPParams.AlphaNormParam = Params.alpha;
c@264 81 m_DFPParams.LPOrd = Params.LPOrd;
c@264 82 m_DFPParams.LPACoeffs = Params.LPACoeffs;
c@264 83 m_DFPParams.LPBCoeffs = Params.LPBCoeffs;
c@264 84 m_DFPParams.winPre = Params.WinT.pre;
c@264 85 m_DFPParams.winPost = Params.WinT.post;
c@264 86 m_DFPParams.isMedianPositive = true;
c@264 87
c@264 88 m_DFConditioning = new DFProcess( m_DFPParams );
c@264 89
c@264 90
c@264 91 // these are parameters for smoothing m_tempoScratch
c@264 92 m_RCFPParams.length = m_lagLength;
c@264 93 m_RCFPParams.AlphaNormParam = Params.alpha;
c@264 94 m_RCFPParams.LPOrd = Params.LPOrd;
c@264 95 m_RCFPParams.LPACoeffs = Params.LPACoeffs;
c@264 96 m_RCFPParams.LPBCoeffs = Params.LPBCoeffs;
c@264 97 m_RCFPParams.winPre = Params.WinT.pre;
c@264 98 m_RCFPParams.winPost = Params.WinT.post;
c@264 99 m_RCFPParams.isMedianPositive = true;
c@264 100
c@264 101 m_RCFConditioning = new DFProcess( m_RCFPParams );
c@264 102
c@264 103 }
c@264 104
c@264 105 void TempoTrack::deInitialise()
c@264 106 {
c@264 107 delete [] m_rawDFFrame;
c@264 108
c@264 109 delete [] m_smoothDFFrame;
c@264 110
c@264 111 delete [] m_smoothRCF;
c@264 112
c@264 113 delete [] m_frameACF;
c@264 114
c@264 115 delete [] m_tempoScratch;
c@264 116
c@264 117 delete m_DFConditioning;
c@264 118
c@264 119 delete m_RCFConditioning;
c@264 120
c@264 121 }
c@264 122
c@264 123 void TempoTrack::createCombFilter(double* Filter, unsigned int winLength, unsigned int TSig, double beatLag)
c@264 124 {
c@264 125 unsigned int i;
c@264 126
c@264 127 if( beatLag == 0 )
c@264 128 {
c@264 129 for( i = 0; i < winLength; i++ )
c@264 130 {
c@264 131 Filter[ i ] = ( ( i + 1 ) / pow( m_rayparam, 2.0) ) * exp( ( -pow(( i + 1 ),2.0 ) / ( 2.0 * pow( m_rayparam, 2.0))));
c@264 132 }
c@264 133 }
c@264 134 else
c@264 135 {
c@264 136 m_sigma = beatLag/4;
c@264 137 for( i = 0; i < winLength; i++ )
c@264 138 {
c@264 139 double dlag = (double)(i+1) - beatLag;
c@264 140 Filter[ i ] = exp(-0.5 * pow(( dlag / m_sigma), 2.0) ) / (sqrt( 2 * PI) * m_sigma);
c@264 141 }
c@264 142 }
c@264 143 }
c@264 144
c@264 145 double TempoTrack::tempoMM(double* ACF, double* weight, int tsig)
c@264 146 {
c@264 147
c@264 148 double period = 0;
c@264 149 double maxValRCF = 0.0;
c@264 150 unsigned int maxIndexRCF = 0;
c@264 151
c@264 152 double* pdPeaks;
c@264 153
c@264 154 unsigned int maxIndexTemp;
c@264 155 double maxValTemp;
c@264 156 unsigned int count;
c@264 157
c@264 158 unsigned int numelem,i,j;
c@264 159 int a, b;
c@264 160
c@264 161 for( i = 0; i < m_lagLength; i++ )
c@264 162 m_tempoScratch[ i ] = 0.0;
c@264 163
c@264 164 if( tsig == 0 )
c@264 165 {
c@264 166 //if time sig is unknown, use metrically unbiased version of Filterbank
c@264 167 numelem = 4;
c@264 168 }
c@264 169 else
c@264 170 {
c@264 171 numelem = tsig;
c@264 172 }
c@264 173
c@272 174 #ifdef DEBUG_TEMPO_TRACK
c@271 175 std::cerr << "tempoMM: m_winLength = " << m_winLength << ", m_lagLength = " << m_lagLength << ", numelem = " << numelem << std::endl;
c@272 176 #endif
c@271 177
c@264 178 for(i=1;i<m_lagLength-1;i++)
c@264 179 {
c@264 180 //first and last output values are left intentionally as zero
c@264 181 for (a=1;a<=numelem;a++)
c@264 182 {
c@264 183 for(b=(1-a);b<a;b++)
c@264 184 {
c@264 185 if( tsig == 0 )
c@264 186 {
c@264 187 m_tempoScratch[i] += ACF[a*(i+1)+b-1] * (1.0 / (2.0 * (double)a-1)) * weight[i];
c@264 188 }
c@264 189 else
c@264 190 {
c@264 191 m_tempoScratch[i] += ACF[a*(i+1)+b-1] * 1 * weight[i];
c@264 192 }
c@264 193 }
c@264 194 }
c@264 195 }
c@264 196
c@264 197
c@264 198 //////////////////////////////////////////////////
c@264 199 // MODIFIED BEAT PERIOD EXTRACTION //////////////
c@264 200 /////////////////////////////////////////////////
c@264 201
c@264 202 // find smoothed version of RCF ( as applied to Detection Function)
c@264 203 m_RCFConditioning->process( m_tempoScratch, m_smoothRCF);
c@264 204
c@264 205 if (tsig != 0) // i.e. in context dependent state
c@264 206 {
c@264 207 // NOW FIND MAX INDEX OF ACFOUT
c@271 208 for( i = 0; i < m_lagLength; i++)
c@271 209 {
c@271 210 if( m_tempoScratch[ i ] > maxValRCF)
c@271 211 {
c@271 212 maxValRCF = m_tempoScratch[ i ];
c@271 213 maxIndexRCF = i;
c@271 214 }
c@271 215 }
c@264 216 }
c@264 217 else // using rayleigh weighting
c@264 218 {
c@264 219 vector <vector<double> > rcfMat;
c@264 220
c@264 221 double sumRcf = 0.;
c@264 222
c@264 223 double maxVal = 0.;
c@264 224 // now find the two values which minimise rcfMat
c@264 225 double minVal = 0.;
c@264 226 int p_i = 1; // periodicity for row i;
c@264 227 int p_j = 1; //periodicity for column j;
c@264 228
c@264 229
c@264 230 for ( i=0; i<m_lagLength; i++)
c@264 231 {
c@264 232 m_tempoScratch[i] =m_smoothRCF[i];
c@264 233 }
c@264 234
c@264 235 // normalise m_tempoScratch so that it sums to zero.
c@264 236 for ( i=0; i<m_lagLength; i++)
c@264 237 {
c@264 238 sumRcf += m_tempoScratch[i];
c@264 239 }
c@264 240
c@264 241 for( i=0; i<m_lagLength; i++)
c@264 242 {
c@264 243 m_tempoScratch[i] /= sumRcf;
c@264 244 }
c@264 245
c@264 246 // create a matrix to store m_tempoScratchValues modified by log2 ratio
c@264 247 for ( i=0; i<m_lagLength; i++)
c@264 248 {
c@264 249 rcfMat.push_back ( vector<double>() ); // adds a new row...
c@264 250 }
c@264 251
c@264 252 for (i=0; i<m_lagLength; i++)
c@264 253 {
c@264 254 for (j=0; j<m_lagLength; j++)
c@264 255 {
c@264 256 rcfMat[i].push_back (0.);
c@264 257 }
c@264 258 }
c@264 259
c@264 260 // the 'i' and 'j' indices deliberately start from '1' and not '0'
c@264 261 for ( i=1; i<m_lagLength; i++)
c@264 262 {
c@264 263 for (j=1; j<m_lagLength; j++)
c@264 264 {
c@264 265 double log2PeriodRatio = log( static_cast<double>(i)/static_cast<double>(j) ) / log(2.0);
c@264 266 rcfMat[i][j] = ( abs(1.0-abs(log2PeriodRatio)) );
c@264 267 rcfMat[i][j] += ( 0.01*( 1./(m_tempoScratch[i]+m_tempoScratch[j]) ) );
c@264 268 }
c@264 269 }
c@264 270
c@264 271 // set diagonal equal to maximum value in rcfMat
c@264 272 // we don't want to pick one strong middle peak - we need a combination of two peaks.
c@264 273
c@264 274 for ( i=1; i<m_lagLength; i++)
c@264 275 {
c@264 276 for (j=1; j<m_lagLength; j++)
c@264 277 {
c@264 278 if (rcfMat[i][j] > maxVal)
c@264 279 {
c@264 280 maxVal = rcfMat[i][j];
c@264 281 }
c@264 282 }
c@264 283 }
c@264 284
c@264 285 for ( i=1; i<m_lagLength; i++)
c@264 286 {
c@264 287 rcfMat[i][i] = maxVal;
c@264 288 }
c@264 289
c@264 290 // now find the row and column number which minimise rcfMat
c@264 291 minVal = maxVal;
c@264 292
c@264 293 for ( i=1; i<m_lagLength; i++)
c@264 294 {
c@264 295 for ( j=1; j<m_lagLength; j++)
c@264 296 {
c@264 297 if (rcfMat[i][j] < minVal)
c@264 298 {
c@264 299 minVal = rcfMat[i][j];
c@264 300 p_i = i;
c@264 301 p_j = j;
c@264 302 }
c@264 303 }
c@264 304 }
c@264 305
c@264 306
c@264 307 // initially choose p_j (arbitrary) - saves on an else statement
c@264 308 int beatPeriod = p_j;
c@264 309 if (m_tempoScratch[p_i] > m_tempoScratch[p_j])
c@264 310 {
c@264 311 beatPeriod = p_i;
c@264 312 }
c@264 313
c@264 314 // now write the output
c@264 315 maxIndexRCF = static_cast<int>(beatPeriod);
c@264 316 }
c@264 317
c@264 318
c@264 319 double locked = 5168.f / maxIndexRCF;
c@264 320 if (locked >= 30 && locked <= 180) {
c@264 321 m_lockedTempo = locked;
c@264 322 }
c@264 323
c@272 324 #ifdef DEBUG_TEMPO_TRACK
c@272 325 std::cerr << "tempoMM: locked tempo = " << m_lockedTempo << std::endl;
c@272 326 #endif
c@272 327
c@264 328 if( tsig == 0 )
c@264 329 tsig = 4;
c@264 330
c@271 331
c@272 332 #ifdef DEBUG_TEMPO_TRACK
c@271 333 std::cerr << "tempoMM: maxIndexRCF = " << maxIndexRCF << std::endl;
c@272 334 #endif
c@264 335
c@264 336 if( tsig == 4 )
c@264 337 {
c@272 338 #ifdef DEBUG_TEMPO_TRACK
c@272 339 std::cerr << "tsig == 4" << std::endl;
c@272 340 #endif
c@272 341
c@264 342 pdPeaks = new double[ 4 ];
c@264 343 for( i = 0; i < 4; i++ ){ pdPeaks[ i ] = 0.0;}
c@264 344
c@264 345 pdPeaks[ 0 ] = ( double )maxIndexRCF + 1;
c@264 346
c@264 347 maxIndexTemp = 0;
c@264 348 maxValTemp = 0.0;
c@264 349 count = 0;
c@264 350
c@264 351 for( i = (2 * maxIndexRCF + 1) - 1; i < (2 * maxIndexRCF + 1) + 2; i++ )
c@264 352 {
c@264 353 if( ACF[ i ] > maxValTemp )
c@264 354 {
c@264 355 maxValTemp = ACF[ i ];
c@264 356 maxIndexTemp = count;
c@264 357 }
c@264 358 count++;
c@264 359 }
c@264 360 pdPeaks[ 1 ] = (double)( maxIndexTemp + 1 + ( (2 * maxIndexRCF + 1 ) - 2 ) + 1 )/2;
c@264 361
c@264 362 maxIndexTemp = 0;
c@264 363 maxValTemp = 0.0;
c@264 364 count = 0;
c@264 365
c@264 366 for( i = (3 * maxIndexRCF + 2 ) - 2; i < (3 * maxIndexRCF + 2 ) + 3; i++ )
c@264 367 {
c@264 368 if( ACF[ i ] > maxValTemp )
c@264 369 {
c@264 370 maxValTemp = ACF[ i ];
c@264 371 maxIndexTemp = count;
c@264 372 }
c@264 373 count++;
c@264 374 }
c@264 375 pdPeaks[ 2 ] = (double)( maxIndexTemp + 1 + ( (3 * maxIndexRCF + 2) - 4 ) + 1 )/3;
c@264 376
c@264 377 maxIndexTemp = 0;
c@264 378 maxValTemp = 0.0;
c@264 379 count = 0;
c@264 380
c@264 381 for( i = ( 4 * maxIndexRCF + 3) - 3; i < ( 4 * maxIndexRCF + 3) + 4; i++ )
c@264 382 {
c@264 383 if( ACF[ i ] > maxValTemp )
c@264 384 {
c@264 385 maxValTemp = ACF[ i ];
c@264 386 maxIndexTemp = count;
c@264 387 }
c@264 388 count++;
c@264 389 }
c@264 390 pdPeaks[ 3 ] = (double)( maxIndexTemp + 1 + ( (4 * maxIndexRCF + 3) - 9 ) + 1 )/4 ;
c@264 391
c@264 392
c@264 393 period = MathUtilities::mean( pdPeaks, 4 );
c@264 394 }
c@264 395 else
c@272 396 {
c@272 397 #ifdef DEBUG_TEMPO_TRACK
c@272 398 std::cerr << "tsig != 4" << std::endl;
c@272 399 #endif
c@272 400
c@264 401 pdPeaks = new double[ 3 ];
c@264 402 for( i = 0; i < 3; i++ ){ pdPeaks[ i ] = 0.0;}
c@264 403
c@264 404 pdPeaks[ 0 ] = ( double )maxIndexRCF + 1;
c@264 405
c@264 406 maxIndexTemp = 0;
c@264 407 maxValTemp = 0.0;
c@264 408 count = 0;
c@264 409
c@264 410 for( i = (2 * maxIndexRCF + 1) - 1; i < (2 * maxIndexRCF + 1) + 2; i++ )
c@264 411 {
c@264 412 if( ACF[ i ] > maxValTemp )
c@264 413 {
c@264 414 maxValTemp = ACF[ i ];
c@264 415 maxIndexTemp = count;
c@264 416 }
c@264 417 count++;
c@264 418 }
c@264 419 pdPeaks[ 1 ] = (double)( maxIndexTemp + 1 + ( (2 * maxIndexRCF + 1 ) - 2 ) + 1 )/2;
c@264 420
c@264 421 maxIndexTemp = 0;
c@264 422 maxValTemp = 0.0;
c@264 423 count = 0;
c@264 424
c@264 425 for( i = (3 * maxIndexRCF + 2 ) - 2; i < (3 * maxIndexRCF + 2 ) + 3; i++ )
c@264 426 {
c@264 427 if( ACF[ i ] > maxValTemp )
c@264 428 {
c@264 429 maxValTemp = ACF[ i ];
c@264 430 maxIndexTemp = count;
c@264 431 }
c@264 432 count++;
c@264 433 }
c@264 434 pdPeaks[ 2 ] = (double)( maxIndexTemp + 1 + ( (3 * maxIndexRCF + 2) - 4 ) + 1 )/3;
c@264 435
c@264 436
c@264 437 period = MathUtilities::mean( pdPeaks, 3 );
c@264 438 }
c@264 439
c@264 440 delete [] pdPeaks;
c@264 441
c@264 442 return period;
c@264 443 }
c@264 444
c@264 445 void TempoTrack::stepDetect( double* periodP, double* periodG, int currentIdx, int* flag )
c@264 446 {
c@264 447 double stepthresh = 1 * 3.9017;
c@264 448
c@264 449 if( *flag )
c@264 450 {
c@264 451 if(abs(periodG[ currentIdx ] - periodP[ currentIdx ]) > stepthresh)
c@264 452 {
c@264 453 // do nuffin'
c@264 454 }
c@264 455 }
c@264 456 else
c@264 457 {
c@264 458 if(fabs(periodG[ currentIdx ]-periodP[ currentIdx ]) > stepthresh)
c@264 459 {
c@264 460 *flag = 3;
c@264 461 }
c@264 462 }
c@264 463 }
c@264 464
c@264 465 void TempoTrack::constDetect( double* periodP, int currentIdx, int* flag )
c@264 466 {
c@264 467 double constthresh = 2 * 3.9017;
c@264 468
c@264 469 if( fabs( 2 * periodP[ currentIdx ] - periodP[ currentIdx - 1] - periodP[ currentIdx - 2] ) < constthresh)
c@264 470 {
c@264 471 *flag = 1;
c@264 472 }
c@264 473 else
c@264 474 {
c@264 475 *flag = 0;
c@264 476 }
c@264 477 }
c@264 478
c@264 479 int TempoTrack::findMeter(double *ACF, unsigned int len, double period)
c@264 480 {
c@264 481 int i;
c@264 482 int p = (int)MathUtilities::round( period );
c@264 483 int tsig;
c@264 484
c@264 485 double Energy_3 = 0.0;
c@264 486 double Energy_4 = 0.0;
c@264 487
c@264 488 double temp3A = 0.0;
c@264 489 double temp3B = 0.0;
c@264 490 double temp4A = 0.0;
c@264 491 double temp4B = 0.0;
c@264 492
c@264 493 double* dbf = new double[ len ]; int t = 0;
c@264 494 for( unsigned int u = 0; u < len; u++ ){ dbf[ u ] = 0.0; }
c@264 495
c@264 496 if( (double)len < 6 * p + 2 )
c@264 497 {
c@264 498 for( i = ( 3 * p - 2 ); i < ( 3 * p + 2 ) + 1; i++ )
c@264 499 {
c@264 500 temp3A += ACF[ i ];
c@264 501 dbf[ t++ ] = ACF[ i ];
c@264 502 }
c@264 503
c@264 504 for( i = ( 4 * p - 2 ); i < ( 4 * p + 2 ) + 1; i++ )
c@264 505 {
c@264 506 temp4A += ACF[ i ];
c@264 507 }
c@264 508
c@264 509 Energy_3 = temp3A;
c@264 510 Energy_4 = temp4A;
c@264 511 }
c@264 512 else
c@264 513 {
c@264 514 for( i = ( 3 * p - 2 ); i < ( 3 * p + 2 ) + 1; i++ )
c@264 515 {
c@264 516 temp3A += ACF[ i ];
c@264 517 }
c@264 518
c@264 519 for( i = ( 4 * p - 2 ); i < ( 4 * p + 2 ) + 1; i++ )
c@264 520 {
c@264 521 temp4A += ACF[ i ];
c@264 522 }
c@264 523
c@264 524 for( i = ( 6 * p - 2 ); i < ( 6 * p + 2 ) + 1; i++ )
c@264 525 {
c@264 526 temp3B += ACF[ i ];
c@264 527 }
c@264 528
c@264 529 for( i = ( 2 * p - 2 ); i < ( 2 * p + 2 ) + 1; i++ )
c@264 530 {
c@264 531 temp4B += ACF[ i ];
c@264 532 }
c@264 533
c@264 534 Energy_3 = temp3A + temp3B;
c@264 535 Energy_4 = temp4A + temp4B;
c@264 536 }
c@264 537
c@264 538 if (Energy_3 > Energy_4)
c@264 539 {
c@264 540 tsig = 3;
c@264 541 }
c@264 542 else
c@264 543 {
c@264 544 tsig = 4;
c@264 545 }
c@264 546
c@264 547
c@264 548 return tsig;
c@264 549 }
c@264 550
c@264 551 void TempoTrack::createPhaseExtractor(double *Filter, unsigned int winLength, double period, unsigned int fsp, unsigned int lastBeat)
c@264 552 {
c@264 553 int p = (int)MathUtilities::round( period );
c@264 554 int predictedOffset = 0;
c@264 555
c@272 556 #ifdef DEBUG_TEMPO_TRACK
c@271 557 std::cerr << "TempoTrack::createPhaseExtractor: period = " << period << ", p = " << p << std::endl;
c@272 558 #endif
c@271 559
c@272 560 if (p > 10000) {
c@272 561 std::cerr << "TempoTrack::createPhaseExtractor: WARNING! Highly implausible period value " << p << "!" << std::endl;
c@272 562 period = 5168 / 120;
c@272 563 }
c@271 564
c@272 565 double* phaseScratch = new double[ p*2 + 2 ];
c@272 566 for (int i = 0; i < p*2 + 2; ++i) phaseScratch[i] = 0.0;
c@264 567
c@264 568
c@264 569 if( lastBeat != 0 )
c@264 570 {
c@264 571 lastBeat = (int)MathUtilities::round((double)lastBeat );///(double)winLength);
c@264 572
c@272 573 predictedOffset = lastBeat + p - fsp;
c@264 574
c@272 575 if (predictedOffset < 0)
c@272 576 {
c@272 577 lastBeat = 0;
c@272 578 }
c@264 579 }
c@264 580
c@264 581 if( lastBeat != 0 )
c@264 582 {
c@264 583 int mu = p;
c@264 584 double sigma = (double)p/8;
c@264 585 double PhaseMin = 0.0;
c@264 586 double PhaseMax = 0.0;
c@264 587 unsigned int scratchLength = p*2;
c@264 588 double temp = 0.0;
c@264 589
c@264 590 for( int i = 0; i < scratchLength; i++ )
c@264 591 {
c@264 592 phaseScratch[ i ] = exp( -0.5 * pow( ( i - mu ) / sigma, 2 ) ) / ( sqrt( 2*PI ) *sigma );
c@264 593 }
c@264 594
c@264 595 MathUtilities::getFrameMinMax( phaseScratch, scratchLength, &PhaseMin, &PhaseMax );
c@264 596
c@264 597 for(int i = 0; i < scratchLength; i ++)
c@264 598 {
c@264 599 temp = phaseScratch[ i ];
c@264 600 phaseScratch[ i ] = (temp - PhaseMin)/PhaseMax;
c@264 601 }
c@264 602
c@272 603 #ifdef DEBUG_TEMPO_TRACK
c@271 604 std::cerr << "predictedOffset = " << predictedOffset << std::endl;
c@272 605 #endif
c@271 606
c@264 607 unsigned int index = 0;
c@272 608 for (int i = p - ( predictedOffset - 1); i < p + ( p - predictedOffset) + 1; i++)
c@264 609 {
c@272 610 #ifdef DEBUG_TEMPO_TRACK
c@272 611 std::cerr << "assigning to filter index " << index << " (size = " << p*2 << ")" << " value " << phaseScratch[i] << " from scratch index " << i << std::endl;
c@272 612 #endif
c@264 613 Filter[ index++ ] = phaseScratch[ i ];
c@264 614 }
c@264 615 }
c@264 616 else
c@264 617 {
c@264 618 for( int i = 0; i < p; i ++)
c@264 619 {
c@264 620 Filter[ i ] = 1;
c@264 621 }
c@264 622 }
c@264 623
c@264 624 delete [] phaseScratch;
c@264 625 }
c@264 626
c@264 627 int TempoTrack::phaseMM(double *DF, double *weighting, unsigned int winLength, double period)
c@264 628 {
c@264 629 int alignment = 0;
c@264 630 int p = (int)MathUtilities::round( period );
c@264 631
c@264 632 double temp = 0.0;
c@264 633
c@264 634 double* y = new double[ winLength ];
c@264 635 double* align = new double[ p ];
c@264 636
c@264 637 for( int i = 0; i < winLength; i++ )
c@264 638 {
c@264 639 y[ i ] = (double)( -i + winLength )/(double)winLength;
c@264 640 y[ i ] = pow(y [i ],2.0); // raise to power 2.
c@264 641 }
c@264 642
c@264 643 for( int o = 0; o < p; o++ )
c@264 644 {
c@264 645 temp = 0.0;
c@264 646 for(int i = 1 + (o - 1); i< winLength; i += (p + 1))
c@264 647 {
c@264 648 temp = temp + DF[ i ] * y[ i ];
c@264 649 }
c@264 650 align[ o ] = temp * weighting[ o ];
c@264 651 }
c@264 652
c@264 653
c@264 654 double valTemp = 0.0;
c@264 655 for(int i = 0; i < p; i++)
c@264 656 {
c@264 657 if( align[ i ] > valTemp )
c@264 658 {
c@264 659 valTemp = align[ i ];
c@264 660 alignment = i;
c@264 661 }
c@264 662 }
c@264 663
c@264 664 delete [] y;
c@264 665 delete [] align;
c@264 666
c@264 667 return alignment;
c@264 668 }
c@264 669
c@264 670 int TempoTrack::beatPredict(unsigned int FSP0, double alignment, double period, unsigned int step )
c@264 671 {
c@264 672 int beat = 0;
c@264 673
c@264 674 int p = (int)MathUtilities::round( period );
c@264 675 int align = (int)MathUtilities::round( alignment );
c@264 676 int FSP = (int)MathUtilities::round( FSP0 );
c@264 677
c@264 678 int FEP = FSP + ( step );
c@264 679
c@264 680 beat = FSP + align;
c@264 681
c@264 682 m_beats.push_back( beat );
c@264 683
c@264 684 while( beat + p < FEP )
c@264 685 {
c@264 686 beat += p;
c@264 687
c@264 688 m_beats.push_back( beat );
c@264 689 }
c@264 690
c@264 691 return beat;
c@264 692 }
c@264 693
c@264 694
c@264 695
c@264 696 vector<int> TempoTrack::process( vector <double> DF,
c@264 697 vector <double> *tempoReturn )
c@264 698 {
c@264 699 m_dataLength = DF.size();
c@264 700
c@264 701 m_lockedTempo = 0.0;
c@264 702
c@264 703 double period = 0.0;
c@264 704 int stepFlag = 0;
c@264 705 int constFlag = 0;
c@264 706 int FSP = 0;
c@264 707 int tsig = 0;
c@264 708 int lastBeat = 0;
c@264 709
c@264 710 vector <double> causalDF;
c@264 711
c@264 712 causalDF = DF;
c@264 713
c@264 714 //Prepare Causal Extension DFData
c@264 715 unsigned int DFCLength = m_dataLength + m_winLength;
c@264 716
c@264 717 for( unsigned int j = 0; j < m_winLength; j++ )
c@264 718 {
c@264 719 causalDF.push_back( 0 );
c@264 720 }
c@264 721
c@264 722
c@264 723 double* RW = new double[ m_lagLength ];
c@264 724 for( unsigned int clear = 0; clear < m_lagLength; clear++){ RW[ clear ] = 0.0;}
c@264 725
c@264 726 double* GW = new double[ m_lagLength ];
c@264 727 for(unsigned int clear = 0; clear < m_lagLength; clear++){ GW[ clear ] = 0.0;}
c@264 728
c@264 729 double* PW = new double[ m_lagLength ];
c@264 730 for(unsigned clear = 0; clear < m_lagLength; clear++){ PW[ clear ] = 0.0;}
c@264 731
c@264 732 m_DFFramer.setSource( &causalDF[0], m_dataLength );
c@264 733
c@264 734 unsigned int TTFrames = m_DFFramer.getMaxNoFrames();
c@272 735
c@272 736 #ifdef DEBUG_TEMPO_TRACK
c@272 737 std::cerr << "TTFrames = " << TTFrames << std::endl;
c@272 738 #endif
c@264 739
c@264 740 double* periodP = new double[ TTFrames ];
c@264 741 for(unsigned clear = 0; clear < TTFrames; clear++){ periodP[ clear ] = 0.0;}
c@264 742
c@264 743 double* periodG = new double[ TTFrames ];
c@264 744 for(unsigned clear = 0; clear < TTFrames; clear++){ periodG[ clear ] = 0.0;}
c@264 745
c@264 746 double* alignment = new double[ TTFrames ];
c@264 747 for(unsigned clear = 0; clear < TTFrames; clear++){ alignment[ clear ] = 0.0;}
c@264 748
c@264 749 m_beats.clear();
c@264 750
c@264 751 createCombFilter( RW, m_lagLength, 0, 0 );
c@264 752
c@264 753 int TTLoopIndex = 0;
c@264 754
c@264 755 for( unsigned int i = 0; i < TTFrames; i++ )
c@264 756 {
c@264 757 m_DFFramer.getFrame( m_rawDFFrame );
c@264 758
c@264 759 m_DFConditioning->process( m_rawDFFrame, m_smoothDFFrame );
c@264 760
c@264 761 m_correlator.doAutoUnBiased( m_smoothDFFrame, m_frameACF, m_winLength );
c@264 762
c@264 763 periodP[ TTLoopIndex ] = tempoMM( m_frameACF, RW, 0 );
c@264 764
c@264 765 if( GW[ 0 ] != 0 )
c@264 766 {
c@264 767 periodG[ TTLoopIndex ] = tempoMM( m_frameACF, GW, tsig );
c@264 768 }
c@264 769 else
c@264 770 {
c@264 771 periodG[ TTLoopIndex ] = 0.0;
c@264 772 }
c@264 773
c@264 774 stepDetect( periodP, periodG, TTLoopIndex, &stepFlag );
c@264 775
c@264 776 if( stepFlag == 1)
c@264 777 {
c@264 778 constDetect( periodP, TTLoopIndex, &constFlag );
c@264 779 stepFlag = 0;
c@264 780 }
c@264 781 else
c@264 782 {
c@264 783 stepFlag -= 1;
c@264 784 }
c@264 785
c@264 786 if( stepFlag < 0 )
c@264 787 {
c@264 788 stepFlag = 0;
c@264 789 }
c@264 790
c@264 791 if( constFlag != 0)
c@264 792 {
c@264 793 tsig = findMeter( m_frameACF, m_winLength, periodP[ TTLoopIndex ] );
c@264 794
c@264 795 createCombFilter( GW, m_lagLength, tsig, periodP[ TTLoopIndex ] );
c@264 796
c@264 797 periodG[ TTLoopIndex ] = tempoMM( m_frameACF, GW, tsig );
c@264 798
c@264 799 period = periodG[ TTLoopIndex ];
c@264 800
c@272 801 #ifdef DEBUG_TEMPO_TRACK
c@272 802 std::cerr << "TempoTrack::process: constFlag == " << constFlag << ", TTLoopIndex = " << TTLoopIndex << ", period from periodG = " << period << std::endl;
c@272 803 #endif
c@271 804
c@264 805 createPhaseExtractor( PW, m_winLength, period, FSP, 0 );
c@264 806
c@264 807 constFlag = 0;
c@264 808
c@264 809 }
c@264 810 else
c@264 811 {
c@264 812 if( GW[ 0 ] != 0 )
c@264 813 {
c@264 814 period = periodG[ TTLoopIndex ];
c@271 815
c@272 816 #ifdef DEBUG_TEMPO_TRACK
c@272 817 std::cerr << "TempoTrack::process: GW[0] == " << GW[0] << ", TTLoopIndex = " << TTLoopIndex << ", period from periodG = " << period << std::endl;
c@272 818 #endif
c@271 819
c@271 820 if (period > 10000) {
c@272 821 std::cerr << "TempoTrack::process: WARNING! Highly implausible period value " << period << "!" << std::endl;
c@271 822 std::cerr << "periodG contains (of " << TTFrames << " frames): " << std::endl;
c@271 823 for (int i = 0; i < TTLoopIndex + 3 && i < TTFrames; ++i) {
c@271 824 std::cerr << i << " -> " << periodG[i] << std::endl;
c@271 825 }
c@271 826 std::cerr << "periodP contains (of " << TTFrames << " frames): " << std::endl;
c@271 827 for (int i = 0; i < TTLoopIndex + 3 && i < TTFrames; ++i) {
c@271 828 std::cerr << i << " -> " << periodP[i] << std::endl;
c@271 829 }
c@272 830 period = 5168 / 120;
c@271 831 }
c@271 832
c@264 833 createPhaseExtractor( PW, m_winLength, period, FSP, lastBeat );
c@264 834
c@264 835 }
c@264 836 else
c@264 837 {
c@264 838 period = periodP[ TTLoopIndex ];
c@271 839
c@272 840 #ifdef DEBUG_TEMPO_TRACK
c@272 841 std::cerr << "TempoTrack::process: GW[0] == " << GW[0] << ", TTLoopIndex = " << TTLoopIndex << ", period from periodP = " << period << std::endl;
c@272 842 #endif
c@271 843
c@264 844 createPhaseExtractor( PW, m_winLength, period, FSP, 0 );
c@264 845 }
c@264 846 }
c@264 847
c@264 848 alignment[ TTLoopIndex ] = phaseMM( m_rawDFFrame, PW, m_winLength, period );
c@264 849
c@264 850 lastBeat = beatPredict(FSP, alignment[ TTLoopIndex ], period, m_lagLength );
c@264 851
c@264 852 FSP += (m_lagLength);
c@264 853
c@264 854 if (tempoReturn) tempoReturn->push_back(m_lockedTempo);
c@264 855
c@264 856 TTLoopIndex++;
c@264 857 }
c@264 858
c@264 859
c@264 860 delete [] periodP;
c@264 861 delete [] periodG;
c@264 862 delete [] alignment;
c@264 863
c@264 864 delete [] RW;
c@264 865 delete [] GW;
c@264 866 delete [] PW;
c@264 867
c@264 868 return m_beats;
c@264 869 }