comparison src/BTrack.cpp @ 96:c58f01834337

Merge branch 'release/1.0.4'
author Adam Stark <adamstark.uk@gmail.com>
date Sat, 18 Jun 2016 10:50:06 +0100
parents 4aa362058011
children 6a4dd7478954
comparison
equal deleted inserted replaced
87:496d12635af8 96:c58f01834337
21 21
22 #include <cmath> 22 #include <cmath>
23 #include <algorithm> 23 #include <algorithm>
24 #include "BTrack.h" 24 #include "BTrack.h"
25 #include "samplerate.h" 25 #include "samplerate.h"
26 26 #include <iostream>
27 //======================================================================= 27
28 BTrack::BTrack() : odf(512,1024,ComplexSpectralDifferenceHWR,HanningWindow) 28 //=======================================================================
29 { 29 BTrack::BTrack()
30 initialise(512, 1024); 30 : odf (512, 1024, ComplexSpectralDifferenceHWR, HanningWindow)
31 } 31 {
32 32 initialise (512, 1024);
33 //======================================================================= 33 }
34 BTrack::BTrack(int hopSize_) : odf(hopSize_,2*hopSize_,ComplexSpectralDifferenceHWR,HanningWindow) 34
35 //=======================================================================
36 BTrack::BTrack (int hopSize_)
37 : odf(hopSize_, 2*hopSize_, ComplexSpectralDifferenceHWR, HanningWindow)
35 { 38 {
36 initialise(hopSize_, 2*hopSize_); 39 initialise (hopSize_, 2*hopSize_);
37 } 40 }
38 41
39 //======================================================================= 42 //=======================================================================
40 BTrack::BTrack(int hopSize_,int frameSize_) : odf(hopSize_,frameSize_,ComplexSpectralDifferenceHWR,HanningWindow) 43 BTrack::BTrack (int hopSize_, int frameSize_)
41 { 44 : odf (hopSize_, frameSize_, ComplexSpectralDifferenceHWR, HanningWindow)
42 initialise(hopSize_, frameSize_); 45 {
43 } 46 initialise (hopSize_, frameSize_);
44 47 }
45 //======================================================================= 48
46 double BTrack::getBeatTimeInSeconds(long frameNumber,int hopSize,int fs) 49 //=======================================================================
50 BTrack::~BTrack()
51 {
52 #ifdef USE_FFTW
53 // destroy fft plan
54 fftw_destroy_plan (acfForwardFFT);
55 fftw_destroy_plan (acfBackwardFFT);
56 fftw_free (complexIn);
57 fftw_free (complexOut);
58 #endif
59
60 #ifdef USE_KISS_FFT
61 free (cfgForwards);
62 free (cfgBackwards);
63 delete [] fftIn;
64 delete [] fftOut;
65 #endif
66 }
67
68 //=======================================================================
69 double BTrack::getBeatTimeInSeconds (long frameNumber, int hopSize, int fs)
47 { 70 {
48 double hop = (double) hopSize; 71 double hop = (double) hopSize;
49 double samplingFrequency = (double) fs; 72 double samplingFrequency = (double) fs;
50 double frameNum = (double) frameNumber; 73 double frameNum = (double) frameNumber;
51 74
52 return ((hop / samplingFrequency) * frameNum); 75 return ((hop / samplingFrequency) * frameNum);
53 } 76 }
54 77
55 //======================================================================= 78 //=======================================================================
56 double BTrack::getBeatTimeInSeconds(int frameNumber,int hopSize,int fs) 79 double BTrack::getBeatTimeInSeconds (int frameNumber, int hopSize, int fs)
57 { 80 {
58 long frameNum = (long) frameNumber; 81 long frameNum = (long) frameNumber;
59 82
60 return getBeatTimeInSeconds(frameNum, hopSize, fs); 83 return getBeatTimeInSeconds (frameNum, hopSize, fs);
61 } 84 }
62 85
63 86
64 87
65 //======================================================================= 88 //=======================================================================
66 void BTrack::initialise(int hopSize_, int frameSize_) 89 void BTrack::initialise (int hopSize_, int frameSize_)
67 { 90 {
68 double rayparam = 43; 91 double rayparam = 43;
69 double pi = 3.14159265; 92 double pi = 3.14159265;
70 93
71 94
81 104
82 beatDueInFrame = false; 105 beatDueInFrame = false;
83 106
84 107
85 // create rayleigh weighting vector 108 // create rayleigh weighting vector
86 for (int n = 0;n < 128;n++) 109 for (int n = 0; n < 128; n++)
87 { 110 {
88 weightingVector[n] = ((double) n / pow(rayparam,2)) * exp((-1*pow((double)-n,2)) / (2*pow(rayparam,2))); 111 weightingVector[n] = ((double) n / pow(rayparam,2)) * exp((-1*pow((double)-n,2)) / (2*pow(rayparam,2)));
89 } 112 }
90 113
91 // initialise prev_delta 114 // initialise prev_delta
92 for (int i = 0;i < 41;i++) 115 for (int i = 0; i < 41; i++)
93 { 116 {
94 prevDelta[i] = 1; 117 prevDelta[i] = 1;
95 } 118 }
96 119
97 double t_mu = 41/2; 120 double t_mu = 41/2;
116 // in case it is requested before any processing takes place 139 // in case it is requested before any processing takes place
117 latestCumulativeScoreValue = 0; 140 latestCumulativeScoreValue = 0;
118 141
119 // initialise algorithm given the hopsize 142 // initialise algorithm given the hopsize
120 setHopSize(hopSize_); 143 setHopSize(hopSize_);
121 } 144
122 145
123 //======================================================================= 146 // Set up FFT for calculating the auto-correlation function
124 void BTrack::setHopSize(int hopSize_) 147 FFTLengthForACFCalculation = 1024;
148
149 #ifdef USE_FFTW
150 complexIn = (fftw_complex*) fftw_malloc (sizeof(fftw_complex) * FFTLengthForACFCalculation); // complex array to hold fft data
151 complexOut = (fftw_complex*) fftw_malloc (sizeof(fftw_complex) * FFTLengthForACFCalculation); // complex array to hold fft data
152
153 acfForwardFFT = fftw_plan_dft_1d (FFTLengthForACFCalculation, complexIn, complexOut, FFTW_FORWARD, FFTW_ESTIMATE); // FFT plan initialisation
154 acfBackwardFFT = fftw_plan_dft_1d (FFTLengthForACFCalculation, complexOut, complexIn, FFTW_BACKWARD, FFTW_ESTIMATE); // FFT plan initialisation
155 #endif
156
157 #ifdef USE_KISS_FFT
158 fftIn = new kiss_fft_cpx[FFTLengthForACFCalculation];
159 fftOut = new kiss_fft_cpx[FFTLengthForACFCalculation];
160 cfgForwards = kiss_fft_alloc (FFTLengthForACFCalculation, 0, 0, 0);
161 cfgBackwards = kiss_fft_alloc (FFTLengthForACFCalculation, 1, 0, 0);
162 #endif
163 }
164
165 //=======================================================================
166 void BTrack::setHopSize (int hopSize_)
125 { 167 {
126 hopSize = hopSize_; 168 hopSize = hopSize_;
127 onsetDFBufferSize = (512*512)/hopSize; // calculate df buffer size 169 onsetDFBufferSize = (512*512)/hopSize; // calculate df buffer size
128 170
129 beatPeriod = round(60/((((double) hopSize)/44100)*tempo)); 171 beatPeriod = round(60/((((double) hopSize)/44100)*tempo));
130 172
131 // set size of onset detection function buffer 173 // set size of onset detection function buffer
132 onsetDF.resize(onsetDFBufferSize); 174 onsetDF.resize (onsetDFBufferSize);
133 175
134 // set size of cumulative score buffer 176 // set size of cumulative score buffer
135 cumulativeScore.resize(onsetDFBufferSize); 177 cumulativeScore.resize (onsetDFBufferSize);
136 178
137 // initialise df_buffer to zeros 179 // initialise df_buffer to zeros
138 for (int i = 0;i < onsetDFBufferSize;i++) 180 for (int i = 0; i < onsetDFBufferSize; i++)
139 { 181 {
140 onsetDF[i] = 0; 182 onsetDF[i] = 0;
141 cumulativeScore[i] = 0; 183 cumulativeScore[i] = 0;
142 184
143
144 if ((i % ((int) round(beatPeriod))) == 0) 185 if ((i % ((int) round(beatPeriod))) == 0)
145 { 186 {
146 onsetDF[i] = 1; 187 onsetDF[i] = 1;
147 } 188 }
148 } 189 }
149 } 190 }
150 191
151 //======================================================================= 192 //=======================================================================
152 void BTrack::updateHopAndFrameSize(int hopSize_,int frameSize_) 193 void BTrack::updateHopAndFrameSize (int hopSize_, int frameSize_)
153 { 194 {
154 // update the onset detection function object 195 // update the onset detection function object
155 odf.initialise(hopSize_, frameSize_); 196 odf.initialise (hopSize_, frameSize_);
156 197
157 // update the hop size being used by the beat tracker 198 // update the hop size being used by the beat tracker
158 setHopSize(hopSize_); 199 setHopSize (hopSize_);
159 } 200 }
160 201
161 //======================================================================= 202 //=======================================================================
162 bool BTrack::beatDueInCurrentFrame() 203 bool BTrack::beatDueInCurrentFrame()
163 { 204 {
181 { 222 {
182 return latestCumulativeScoreValue; 223 return latestCumulativeScoreValue;
183 } 224 }
184 225
185 //======================================================================= 226 //=======================================================================
186 void BTrack::processAudioFrame(double *frame) 227 void BTrack::processAudioFrame (double* frame)
187 { 228 {
188 // calculate the onset detection function sample for the frame 229 // calculate the onset detection function sample for the frame
189 double sample = odf.calculateOnsetDetectionFunctionSample(frame); 230 double sample = odf.calculateOnsetDetectionFunctionSample (frame);
190
191
192 231
193 // process the new onset detection function sample in the beat tracking algorithm 232 // process the new onset detection function sample in the beat tracking algorithm
194 processOnsetDetectionFunctionSample(sample); 233 processOnsetDetectionFunctionSample (sample);
195 } 234 }
196 235
197 //======================================================================= 236 //=======================================================================
198 void BTrack::processOnsetDetectionFunctionSample(double newSample) 237 void BTrack::processOnsetDetectionFunctionSample (double newSample)
199 { 238 {
200 // we need to ensure that the onset 239 // we need to ensure that the onset
201 // detection function sample is positive 240 // detection function sample is positive
202 newSample = fabs(newSample); 241 newSample = fabs (newSample);
203 242
204 // add a tiny constant to the sample to stop it from ever going 243 // add a tiny constant to the sample to stop it from ever going
205 // to zero. this is to avoid problems further down the line 244 // to zero. this is to avoid problems further down the line
206 newSample = newSample + 0.0001; 245 newSample = newSample + 0.0001;
207 246
208 m0--; 247 m0--;
209 beatCounter--; 248 beatCounter--;
210 beatDueInFrame = false; 249 beatDueInFrame = false;
211 250
212 // move all samples back one step
213 for (int i=0;i < (onsetDFBufferSize-1);i++)
214 {
215 onsetDF[i] = onsetDF[i+1];
216 }
217
218 // add new sample at the end 251 // add new sample at the end
219 onsetDF[onsetDFBufferSize-1] = newSample; 252 onsetDF.addSampleToEnd (newSample);
220 253
221 // update cumulative score 254 // update cumulative score
222 updateCumulativeScore(newSample); 255 updateCumulativeScore (newSample);
223 256
224 // if we are halfway between beats 257 // if we are halfway between beats
225 if (m0 == 0) 258 if (m0 == 0)
226 { 259 {
227 predictBeat(); 260 predictBeat();
237 calculateTempo(); 270 calculateTempo();
238 } 271 }
239 } 272 }
240 273
241 //======================================================================= 274 //=======================================================================
242 void BTrack::setTempo(double tempo) 275 void BTrack::setTempo (double tempo)
243 { 276 {
244 277
245 /////////// TEMPO INDICATION RESET ////////////////// 278 /////////// TEMPO INDICATION RESET //////////////////
246 279
247 // firstly make sure tempo is between 80 and 160 bpm.. 280 // firstly make sure tempo is between 80 and 160 bpm..
304 // offbeat is half of new beat period away 337 // offbeat is half of new beat period away
305 m0 = (int) round(((double) new_bperiod)/2); 338 m0 = (int) round(((double) new_bperiod)/2);
306 } 339 }
307 340
308 //======================================================================= 341 //=======================================================================
309 void BTrack::fixTempo(double tempo) 342 void BTrack::fixTempo (double tempo)
310 { 343 {
311 // firstly make sure tempo is between 80 and 160 bpm.. 344 // firstly make sure tempo is between 80 and 160 bpm..
312 while (tempo > 160) 345 while (tempo > 160)
313 { 346 {
314 tempo = tempo/2; 347 tempo = tempo/2;
344 377
345 //======================================================================= 378 //=======================================================================
346 void BTrack::resampleOnsetDetectionFunction() 379 void BTrack::resampleOnsetDetectionFunction()
347 { 380 {
348 float output[512]; 381 float output[512];
382
349 float input[onsetDFBufferSize]; 383 float input[onsetDFBufferSize];
350 384
351 for (int i = 0;i < onsetDFBufferSize;i++) 385 for (int i = 0;i < onsetDFBufferSize;i++)
352 { 386 {
353 input[i] = (float) onsetDF[i]; 387 input[i] = (float) onsetDF[i];
354 } 388 }
355 389
356 double src_ratio = 512.0/((double) onsetDFBufferSize); 390 double src_ratio = 512.0/((double) onsetDFBufferSize);
357 int BUFFER_LEN = onsetDFBufferSize; 391 int BUFFER_LEN = onsetDFBufferSize;
358 int output_len; 392 int output_len;
359 SRC_DATA src_data ; 393 SRC_DATA src_data ;
360 394
361 //output_len = (int) floor (((double) BUFFER_LEN) * src_ratio) ; 395 //output_len = (int) floor (((double) BUFFER_LEN) * src_ratio) ;
362 output_len = 512; 396 output_len = 512;
363 397
364 src_data.data_in = input; 398 src_data.data_in = input;
365 src_data.input_frames = BUFFER_LEN; 399 src_data.input_frames = BUFFER_LEN;
366 400
367 src_data.src_ratio = src_ratio; 401 src_data.src_ratio = src_ratio;
368 402
369 src_data.data_out = output; 403 src_data.data_out = output;
370 src_data.output_frames = output_len; 404 src_data.output_frames = output_len;
371 405
372 src_simple (&src_data, SRC_SINC_BEST_QUALITY, 1); 406 src_simple (&src_data, SRC_SINC_BEST_QUALITY, 1);
373 407
374 for (int i = 0;i < output_len;i++) 408 for (int i = 0;i < output_len;i++)
375 { 409 {
376 resampledOnsetDF[i] = (double) src_data.data_out[i]; 410 resampledOnsetDF[i] = (double) src_data.data_out[i];
377 } 411 }
378 } 412 }
379 413
380 //======================================================================= 414 //=======================================================================
381 void BTrack::calculateTempo() 415 void BTrack::calculateTempo()
382 { 416 {
383 // adaptive threshold on input 417 // adaptive threshold on input
384 adaptiveThreshold(resampledOnsetDF,512); 418 adaptiveThreshold (resampledOnsetDF,512);
385 419
386 // calculate auto-correlation function of detection function 420 // calculate auto-correlation function of detection function
387 calculateBalancedACF(resampledOnsetDF); 421 calculateBalancedACF (resampledOnsetDF);
388 422
389 // calculate output of comb filterbank 423 // calculate output of comb filterbank
390 calculateOutputOfCombFilterBank(); 424 calculateOutputOfCombFilterBank();
391 425
392
393 // adaptive threshold on rcf 426 // adaptive threshold on rcf
394 adaptiveThreshold(combFilterBankOutput,128); 427 adaptiveThreshold (combFilterBankOutput,128);
395 428
396 429
397 int t_index; 430 int t_index;
398 int t_index2; 431 int t_index2;
399 // calculate tempo observation vector from beat period observation vector 432 // calculate tempo observation vector from beat period observation vector
400 for (int i = 0;i < 41;i++) 433 for (int i = 0;i < 41;i++)
401 { 434 {
402 t_index = (int) round(tempoToLagFactor / ((double) ((2*i)+80))); 435 t_index = (int) round (tempoToLagFactor / ((double) ((2*i)+80)));
403 t_index2 = (int) round(tempoToLagFactor / ((double) ((4*i)+160))); 436 t_index2 = (int) round (tempoToLagFactor / ((double) ((4*i)+160)));
404 437
405 438
406 tempoObservationVector[i] = combFilterBankOutput[t_index-1] + combFilterBankOutput[t_index2-1]; 439 tempoObservationVector[i] = combFilterBankOutput[t_index-1] + combFilterBankOutput[t_index2-1];
407 } 440 }
408 441
423 for (int j=0;j < 41;j++) 456 for (int j=0;j < 41;j++)
424 { 457 {
425 maxval = -1; 458 maxval = -1;
426 for (int i = 0;i < 41;i++) 459 for (int i = 0;i < 41;i++)
427 { 460 {
428 curval = prevDelta[i]*tempoTransitionMatrix[i][j]; 461 curval = prevDelta[i] * tempoTransitionMatrix[i][j];
429 462
430 if (curval > maxval) 463 if (curval > maxval)
431 { 464 {
432 maxval = curval; 465 maxval = curval;
433 } 466 }
434 } 467 }
435 468
436 delta[j] = maxval*tempoObservationVector[j]; 469 delta[j] = maxval * tempoObservationVector[j];
437 } 470 }
438 471
439 472
440 normaliseArray(delta,41); 473 normaliseArray(delta,41);
441 474
451 } 484 }
452 485
453 prevDelta[j] = delta[j]; 486 prevDelta[j] = delta[j];
454 } 487 }
455 488
456 beatPeriod = round((60.0*44100.0)/(((2*maxind)+80)*((double) hopSize))); 489 beatPeriod = round ((60.0*44100.0)/(((2*maxind)+80)*((double) hopSize)));
457 490
458 if (beatPeriod > 0) 491 if (beatPeriod > 0)
459 { 492 {
460 estimatedTempo = 60.0/((((double) hopSize) / 44100.0)*beatPeriod); 493 estimatedTempo = 60.0/((((double) hopSize) / 44100.0) * beatPeriod);
461 } 494 }
462 } 495 }
463 496
464 //======================================================================= 497 //=======================================================================
465 void BTrack::adaptiveThreshold(double *x,int N) 498 void BTrack::adaptiveThreshold (double* x, int N)
466 { 499 {
467 int i = 0; 500 int i = 0;
468 int k,t = 0; 501 int k,t = 0;
469 double x_thresh[N]; 502 double x_thresh[N];
470 503
474 t = std::min(N,p_post); // what is smaller, p_post of df size. This is to avoid accessing outside of arrays 507 t = std::min(N,p_post); // what is smaller, p_post of df size. This is to avoid accessing outside of arrays
475 508
476 // find threshold for first 't' samples, where a full average cannot be computed yet 509 // find threshold for first 't' samples, where a full average cannot be computed yet
477 for (i = 0;i <= t;i++) 510 for (i = 0;i <= t;i++)
478 { 511 {
479 k = std::min((i+p_pre),N); 512 k = std::min ((i+p_pre),N);
480 x_thresh[i] = calculateMeanOfArray(x,1,k); 513 x_thresh[i] = calculateMeanOfArray (x,1,k);
481 } 514 }
482 // find threshold for bulk of samples across a moving average from [i-p_pre,i+p_post] 515 // find threshold for bulk of samples across a moving average from [i-p_pre,i+p_post]
483 for (i = t+1;i < N-p_post;i++) 516 for (i = t+1;i < N-p_post;i++)
484 { 517 {
485 x_thresh[i] = calculateMeanOfArray(x,i-p_pre,i+p_post); 518 x_thresh[i] = calculateMeanOfArray (x,i-p_pre,i+p_post);
486 } 519 }
487 // for last few samples calculate threshold, again, not enough samples to do as above 520 // for last few samples calculate threshold, again, not enough samples to do as above
488 for (i = N-p_post;i < N;i++) 521 for (i = N-p_post;i < N;i++)
489 { 522 {
490 k = std::max((i-p_post),1); 523 k = std::max ((i-p_post),1);
491 x_thresh[i] = calculateMeanOfArray(x,k,N); 524 x_thresh[i] = calculateMeanOfArray (x,k,N);
492 } 525 }
493 526
494 // subtract the threshold from the detection function and check that it is not less than 0 527 // subtract the threshold from the detection function and check that it is not less than 0
495 for (i = 0;i < N;i++) 528 for (i = 0; i < N; i++)
496 { 529 {
497 x[i] = x[i] - x_thresh[i]; 530 x[i] = x[i] - x_thresh[i];
498 if (x[i] < 0) 531 if (x[i] < 0)
499 { 532 {
500 x[i] = 0; 533 x[i] = 0;
512 combFilterBankOutput[i] = 0; 545 combFilterBankOutput[i] = 0;
513 } 546 }
514 547
515 numelem = 4; 548 numelem = 4;
516 549
517 for (int i = 2;i <= 127;i++) // max beat period 550 for (int i = 2; i <= 127; i++) // max beat period
518 { 551 {
519 for (int a = 1;a <= numelem;a++) // number of comb elements 552 for (int a = 1; a <= numelem; a++) // number of comb elements
520 { 553 {
521 for (int b = 1-a;b <= a-1;b++) // general state using normalisation of comb elements 554 for (int b = 1-a; b <= a-1; b++) // general state using normalisation of comb elements
522 { 555 {
523 combFilterBankOutput[i-1] = combFilterBankOutput[i-1] + (acf[(a*i+b)-1]*weightingVector[i-1])/(2*a-1); // calculate value for comb filter row 556 combFilterBankOutput[i-1] = combFilterBankOutput[i-1] + (acf[(a*i+b)-1]*weightingVector[i-1])/(2*a-1); // calculate value for comb filter row
524 } 557 }
525 } 558 }
526 } 559 }
527 } 560 }
528 561
529 //======================================================================= 562 //=======================================================================
530 void BTrack::calculateBalancedACF(double *onsetDetectionFunction) 563 void BTrack::calculateBalancedACF (double* onsetDetectionFunction)
531 { 564 {
532 int l, n = 0; 565 int onsetDetectionFunctionLength = 512;
533 double sum, tmp; 566
534 567 #ifdef USE_FFTW
535 // for l lags from 0-511 568 // copy into complex array and zero pad
536 for (l = 0;l < 512;l++) 569 for (int i = 0;i < FFTLengthForACFCalculation;i++)
537 { 570 {
538 sum = 0; 571 if (i < onsetDetectionFunctionLength)
539 572 {
540 // for n samples from 0 - (512-lag) 573 complexIn[i][0] = onsetDetectionFunction[i];
541 for (n = 0;n < (512-l);n++) 574 complexIn[i][1] = 0.0;
542 { 575 }
543 tmp = onsetDetectionFunction[n] * onsetDetectionFunction[n+l]; // multiply current sample n by sample (n+l) 576 else
544 sum = sum + tmp; // add to sum 577 {
545 } 578 complexIn[i][0] = 0.0;
546 579 complexIn[i][1] = 0.0;
547 acf[l] = sum / (512-l); // weight by number of mults and add to acf buffer 580 }
548 } 581 }
549 } 582
550 583 // perform the fft
551 //======================================================================= 584 fftw_execute (acfForwardFFT);
552 double BTrack::calculateMeanOfArray(double *array,int startIndex,int endIndex) 585
586 // multiply by complex conjugate
587 for (int i = 0;i < FFTLengthForACFCalculation;i++)
588 {
589 complexOut[i][0] = complexOut[i][0]*complexOut[i][0] + complexOut[i][1]*complexOut[i][1];
590 complexOut[i][1] = 0.0;
591 }
592
593 // perform the ifft
594 fftw_execute (acfBackwardFFT);
595
596 #endif
597
598 #ifdef USE_KISS_FFT
599 // copy into complex array and zero pad
600 for (int i = 0;i < FFTLengthForACFCalculation;i++)
601 {
602 if (i < onsetDetectionFunctionLength)
603 {
604 fftIn[i].r = onsetDetectionFunction[i];
605 fftIn[i].i = 0.0;
606 }
607 else
608 {
609 fftIn[i].r = 0.0;
610 fftIn[i].i = 0.0;
611 }
612 }
613
614 // execute kiss fft
615 kiss_fft (cfgForwards, fftIn, fftOut);
616
617 // multiply by complex conjugate
618 for (int i = 0;i < FFTLengthForACFCalculation;i++)
619 {
620 fftOut[i].r = fftOut[i].r * fftOut[i].r + fftOut[i].i * fftOut[i].i;
621 fftOut[i].i = 0.0;
622 }
623
624 // perform the ifft
625 kiss_fft (cfgBackwards, fftOut, fftIn);
626
627 #endif
628
629 double lag = 512;
630
631 for (int i = 0; i < 512; i++)
632 {
633 #ifdef USE_FFTW
634 // calculate absolute value of result
635 double absValue = sqrt (complexIn[i][0]*complexIn[i][0] + complexIn[i][1]*complexIn[i][1]);
636 #endif
637
638 #ifdef USE_KISS_FFT
639 // calculate absolute value of result
640 double absValue = sqrt (fftIn[i].r * fftIn[i].r + fftIn[i].i * fftIn[i].i);
641 #endif
642 // divide by inverse lad to deal with scale bias towards small lags
643 acf[i] = absValue / lag;
644
645 // this division by 1024 is technically unnecessary but it ensures the algorithm produces
646 // exactly the same ACF output as the old time domain implementation. The time difference is
647 // minimal so I decided to keep it
648 acf[i] = acf[i] / 1024.;
649
650 lag = lag - 1.;
651 }
652 }
653
654 //=======================================================================
655 double BTrack::calculateMeanOfArray (double* array, int startIndex, int endIndex)
553 { 656 {
554 int i; 657 int i;
555 double sum = 0; 658 double sum = 0;
556 659
557 int length = endIndex - startIndex; 660 int length = endIndex - startIndex;
558 661
559 // find sum 662 // find sum
560 for (i = startIndex;i < endIndex;i++) 663 for (i = startIndex; i < endIndex; i++)
561 { 664 {
562 sum = sum + array[i]; 665 sum = sum + array[i];
563 } 666 }
564 667
565 if (length > 0) 668 if (length > 0)
571 return 0; 674 return 0;
572 } 675 }
573 } 676 }
574 677
575 //======================================================================= 678 //=======================================================================
576 void BTrack::normaliseArray(double *array,int N) 679 void BTrack::normaliseArray (double* array, int N)
577 { 680 {
578 double sum = 0; 681 double sum = 0;
579 682
580 for (int i = 0;i < N;i++) 683 for (int i = 0; i < N; i++)
581 { 684 {
582 if (array[i] > 0) 685 if (array[i] > 0)
583 { 686 {
584 sum = sum + array[i]; 687 sum = sum + array[i];
585 } 688 }
586 } 689 }
587 690
588 if (sum > 0) 691 if (sum > 0)
589 { 692 {
590 for (int i = 0;i < N;i++) 693 for (int i = 0; i < N; i++)
591 { 694 {
592 array[i] = array[i] / sum; 695 array[i] = array[i] / sum;
593 } 696 }
594 } 697 }
595 } 698 }
596 699
597 //======================================================================= 700 //=======================================================================
598 void BTrack::updateCumulativeScore(double odfSample) 701 void BTrack::updateCumulativeScore (double odfSample)
599 { 702 {
600 int start, end, winsize; 703 int start, end, winsize;
601 double max; 704 double max;
602 705
603 start = onsetDFBufferSize - round(2*beatPeriod); 706 start = onsetDFBufferSize - round (2 * beatPeriod);
604 end = onsetDFBufferSize - round(beatPeriod/2); 707 end = onsetDFBufferSize - round (beatPeriod / 2);
605 winsize = end-start+1; 708 winsize = end-start+1;
606 709
607 double w1[winsize]; 710 double w1[winsize];
608 double v = -2*beatPeriod; 711 double v = -2*beatPeriod;
609 double wcumscore; 712 double wcumscore;
610 713
611
612 // create window 714 // create window
613 for (int i = 0;i < winsize;i++) 715 for (int i = 0; i < winsize; i++)
614 { 716 {
615 w1[i] = exp((-1*pow(tightness*log(-v/beatPeriod),2))/2); 717 w1[i] = exp((-1 * pow (tightness * log (-v / beatPeriod), 2)) / 2);
616 v = v+1; 718 v = v+1;
617 } 719 }
618 720
619 // calculate new cumulative score value 721 // calculate new cumulative score value
620 max = 0; 722 max = 0;
621 int n = 0; 723 int n = 0;
622 for (int i=start;i <= end;i++) 724 for (int i=start; i <= end; i++)
623 { 725 {
624 wcumscore = cumulativeScore[i]*w1[n]; 726 wcumscore = cumulativeScore[i]*w1[n];
625 727
626 if (wcumscore > max) 728 if (wcumscore > max)
627 { 729 {
628 max = wcumscore; 730 max = wcumscore;
629 } 731 }
630 n++; 732 n++;
631 } 733 }
632 734
633 735 latestCumulativeScoreValue = ((1 - alpha) * odfSample) + (alpha * max);
634 // shift cumulative score back one 736
635 for (int i = 0;i < (onsetDFBufferSize-1);i++) 737 cumulativeScore.addSampleToEnd (latestCumulativeScoreValue);
636 {
637 cumulativeScore[i] = cumulativeScore[i+1];
638 }
639
640 // add new value to cumulative score
641 cumulativeScore[onsetDFBufferSize-1] = ((1-alpha)*odfSample) + (alpha*max);
642
643 latestCumulativeScoreValue = cumulativeScore[onsetDFBufferSize-1];
644
645 } 738 }
646 739
647 //======================================================================= 740 //=======================================================================
648 void BTrack::predictBeat() 741 void BTrack::predictBeat()
649 { 742 {
650 int windowSize = (int) beatPeriod; 743 int windowSize = (int) beatPeriod;
651 double futureCumulativeScore[onsetDFBufferSize + windowSize]; 744 double futureCumulativeScore[onsetDFBufferSize + windowSize];
652 double w2[windowSize]; 745 double w2[windowSize];
746
653 // copy cumscore to first part of fcumscore 747 // copy cumscore to first part of fcumscore
654 for (int i = 0;i < onsetDFBufferSize;i++) 748 for (int i = 0;i < onsetDFBufferSize;i++)
655 { 749 {
656 futureCumulativeScore[i] = cumulativeScore[i]; 750 futureCumulativeScore[i] = cumulativeScore[i];
657 } 751 }
658 752
659 // create future window 753 // create future window
660 double v = 1; 754 double v = 1;
661 for (int i = 0;i < windowSize;i++) 755 for (int i = 0; i < windowSize; i++)
662 { 756 {
663 w2[i] = exp((-1*pow((v - (beatPeriod/2)),2)) / (2*pow((beatPeriod/2) ,2))); 757 w2[i] = exp((-1*pow((v - (beatPeriod/2)),2)) / (2*pow((beatPeriod/2) ,2)));
664 v++; 758 v++;
665 } 759 }
666 760
675 { 769 {
676 w1[i] = exp((-1*pow(tightness*log(-v/beatPeriod),2))/2); 770 w1[i] = exp((-1*pow(tightness*log(-v/beatPeriod),2))/2);
677 v = v+1; 771 v = v+1;
678 } 772 }
679 773
680
681
682 // calculate future cumulative score 774 // calculate future cumulative score
683 double max; 775 double max;
684 int n; 776 int n;
685 double wcumscore; 777 double wcumscore;
686 for (int i = onsetDFBufferSize;i < (onsetDFBufferSize+windowSize);i++) 778 for (int i = onsetDFBufferSize; i < (onsetDFBufferSize + windowSize); i++)
687 { 779 {
688 start = i - round(2*beatPeriod); 780 start = i - round (2*beatPeriod);
689 end = i - round(beatPeriod/2); 781 end = i - round (beatPeriod/2);
690 782
691 max = 0; 783 max = 0;
692 n = 0; 784 n = 0;
693 for (int k=start;k <= end;k++) 785 for (int k=start;k <= end;k++)
694 { 786 {
702 } 794 }
703 795
704 futureCumulativeScore[i] = max; 796 futureCumulativeScore[i] = max;
705 } 797 }
706 798
707
708 // predict beat 799 // predict beat
709 max = 0; 800 max = 0;
710 n = 0; 801 n = 0;
711 802
712 for (int i = onsetDFBufferSize;i < (onsetDFBufferSize+windowSize);i++) 803 for (int i = onsetDFBufferSize; i < (onsetDFBufferSize + windowSize); i++)
713 { 804 {
714 wcumscore = futureCumulativeScore[i]*w2[n]; 805 wcumscore = futureCumulativeScore[i]*w2[n];
715 806
716 if (wcumscore > max) 807 if (wcumscore > max)
717 { 808 {
721 812
722 n++; 813 n++;
723 } 814 }
724 815
725 // set next prediction time 816 // set next prediction time
726 m0 = beatCounter+round(beatPeriod/2); 817 m0 = beatCounter + round (beatPeriod / 2);
727 818 }
728
729 }