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