1 /*
2     SuperCollider real time audio synthesis system
3  Copyright (c) 2002 James McCartney. All rights reserved.
4     http://www.audiosynth.com
5 
6  This program is free software; you can redistribute it and/or modify
7  it under the terms of the GNU General Public License as published by
8  the Free Software Foundation; either version 2 of the License, or
9  (at your option) any later version.
10 
11  This program is distributed in the hope that it will be useful,
12  but WITHOUT ANY WARRANTY; without even the implied warranty of
13  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
14  GNU General Public License for more details.
15 
16  You should have received a copy of the GNU General Public License
17  along with this program; if not, write to the Free Software
18  Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301  USA
19  */
20 
21 // BeatTrack2 UGen implemented by Nick Collins (http://www.informatics.sussex.ac.uk/users/nc81/)
22 // 6 Nov 2007
23 
24 #include "ML.h"
25 
26 // need to add bestgroove option to store groove, else remove output which is currently always straight 16ths
27 
28 
29 static const int g_numtempi = 120;
30 static float g_periods[g_numtempi] = { 1,
31                                        0.98360655737705,
32                                        0.96774193548387,
33                                        0.95238095238095,
34                                        0.9375,
35                                        0.92307692307692,
36                                        0.90909090909091,
37                                        0.8955223880597,
38                                        0.88235294117647,
39                                        0.8695652173913,
40                                        0.85714285714286,
41                                        0.84507042253521,
42                                        0.83333333333333,
43                                        0.82191780821918,
44                                        0.81081081081081,
45                                        0.8,
46                                        0.78947368421053,
47                                        0.77922077922078,
48                                        0.76923076923077,
49                                        0.75949367088608,
50                                        0.75,
51                                        0.74074074074074,
52                                        0.73170731707317,
53                                        0.72289156626506,
54                                        0.71428571428571,
55                                        0.70588235294118,
56                                        0.69767441860465,
57                                        0.68965517241379,
58                                        0.68181818181818,
59                                        0.67415730337079,
60                                        0.66666666666667,
61                                        0.65934065934066,
62                                        0.65217391304348,
63                                        0.64516129032258,
64                                        0.63829787234043,
65                                        0.63157894736842,
66                                        0.625,
67                                        0.61855670103093,
68                                        0.61224489795918,
69                                        0.60606060606061,
70                                        0.6,
71                                        0.59405940594059,
72                                        0.58823529411765,
73                                        0.58252427184466,
74                                        0.57692307692308,
75                                        0.57142857142857,
76                                        0.56603773584906,
77                                        0.5607476635514,
78                                        0.55555555555556,
79                                        0.55045871559633,
80                                        0.54545454545455,
81                                        0.54054054054054,
82                                        0.53571428571429,
83                                        0.53097345132743,
84                                        0.52631578947368,
85                                        0.52173913043478,
86                                        0.51724137931034,
87                                        0.51282051282051,
88                                        0.50847457627119,
89                                        0.50420168067227,
90                                        0.5,
91                                        0.49586776859504,
92                                        0.49180327868852,
93                                        0.48780487804878,
94                                        0.48387096774194,
95                                        0.48,
96                                        0.47619047619048,
97                                        0.47244094488189,
98                                        0.46875,
99                                        0.46511627906977,
100                                        0.46153846153846,
101                                        0.45801526717557,
102                                        0.45454545454545,
103                                        0.45112781954887,
104                                        0.44776119402985,
105                                        0.44444444444444,
106                                        0.44117647058824,
107                                        0.43795620437956,
108                                        0.43478260869565,
109                                        0.43165467625899,
110                                        0.42857142857143,
111                                        0.42553191489362,
112                                        0.42253521126761,
113                                        0.41958041958042,
114                                        0.41666666666667,
115                                        0.41379310344828,
116                                        0.41095890410959,
117                                        0.40816326530612,
118                                        0.40540540540541,
119                                        0.40268456375839,
120                                        0.4,
121                                        0.39735099337748,
122                                        0.39473684210526,
123                                        0.3921568627451,
124                                        0.38961038961039,
125                                        0.38709677419355,
126                                        0.38461538461538,
127                                        0.38216560509554,
128                                        0.37974683544304,
129                                        0.37735849056604,
130                                        0.375,
131                                        0.37267080745342,
132                                        0.37037037037037,
133                                        0.3680981595092,
134                                        0.36585365853659,
135                                        0.36363636363636,
136                                        0.36144578313253,
137                                        0.35928143712575,
138                                        0.35714285714286,
139                                        0.35502958579882,
140                                        0.35294117647059,
141                                        0.35087719298246,
142                                        0.34883720930233,
143                                        0.34682080924855,
144                                        0.3448275862069,
145                                        0.34285714285714,
146                                        0.34090909090909,
147                                        0.33898305084746,
148                                        0.33707865168539,
149                                        0.33519553072626 };
150 // float g_tempoweight[g_numtempi]= { 0.8, 0.82581988897472, 0.83651483716701, 0.84472135955, 0.85163977794943,
151 // 0.85773502691896, 0.86324555320337, 0.8683130051064, 0.87302967433402, 0.87745966692415, 0.88164965809277,
152 // 0.88563488385777, 0.88944271909999, 0.89309493362513, 0.89660917830793, 0.9, 0.90327955589886, 0.90645812948448,
153 // 0.90954451150103, 0.91254628677423, 0.91547005383793, 0.91832159566199, 0.9211060141639, 0.92382783747338,
154 // 0.92649110640674, 0.92909944487358, 0.93165611772088, 0.93416407864999, 0.93662601021279, 0.93904435743076,
155 // 0.94142135623731, 0.94375905768565, 0.94605934866804, 0.94832396974191, 0.95055453054182, 0.95275252316519,
156 // 0.9549193338483, 0.95705625319186, 0.95916448515084, 0.96124515496597, 0.96329931618555, 0.96532795690183,
157 // 0.96733200530682, 0.969312334656, 0.97126976771554, 0.97320508075689, 0.97511900715418, 0.97701224063136,
158 // 0.97888543819998, 0.98073922282301, 0.98257418583506, 0.98439088914586, 0.98618986725025, 0.98797162906496,
159 // 0.9897366596101, 0.99148542155127, 0.99321835661586, 0.99493588689618, 0.99663841605004, 0.99832633040858, 1,
160 // 0.99832633040858, 0.99663841605004, 0.99493588689618, 0.99321835661586, 0.99148542155127, 0.9897366596101,
161 // 0.98797162906496, 0.98618986725025, 0.98439088914586, 0.98257418583506, 0.98073922282301, 0.97888543819998,
162 // 0.97701224063136, 0.97511900715418, 0.97320508075689, 0.97126976771554, 0.969312334656, 0.96733200530682,
163 // 0.96532795690183, 0.96329931618555, 0.96124515496597, 0.95916448515084, 0.95705625319186, 0.9549193338483,
164 // 0.95275252316519, 0.95055453054182, 0.94832396974191, 0.94605934866804, 0.94375905768565, 0.94142135623731,
165 // 0.93904435743076, 0.93662601021279, 0.93416407864999, 0.93165611772088, 0.92909944487358, 0.92649110640674,
166 // 0.92382783747338, 0.9211060141639, 0.91832159566199, 0.91547005383793, 0.91254628677423, 0.90954451150103,
167 // 0.90645812948448, 0.90327955589886, 0.9, 0.89660917830793, 0.89309493362513, 0.88944271909999, 0.88563488385777,
168 // 0.88164965809277, 0.87745966692415, 0.87302967433402, 0.8683130051064, 0.86324555320337, 0.85773502691896,
169 // 0.85163977794943, 0.84472135955, 0.83651483716701, 0.82581988897472 }; const float g_groove = 0.32;
170 
171 
172 static float g_sep[8] = { 0.0, 0.25, 0.5, 0.75, 0.0, 0.32, 0.5, 0.82 };
173 // weight for particular step
174 static float g_weight[4] = { 1.0, 0.5, 0.9, 0.6 };
175 // weight for blurring feature envelope locally
176 static float g_weight2[9] = { 0.05, 0.1, 0.3, 0.7, 1.0, 0.7, 0.3, 0.1, 0.05 };
177 
178 // void BeatTrack2_dofft(BeatTrack2 *unit, uint32);
179 static void calculatetemplate(BeatTrack2* unit, int which, int j);
180 static void finaldecision(BeatTrack2* unit);
181 
BeatTrack2_Ctor(BeatTrack2 * unit)182 void BeatTrack2_Ctor(BeatTrack2* unit) {
183     // unit->m_srate = unit->mWorld->mFullRate.mSampleRate;
184     float kblocklength = unit->mWorld->mFullRate.mBufDuration; // seconds per control block
185     unit->m_krlength = kblocklength;
186     // N features per block over numphases*2 variants for one of 120 tempi, so need at least 120 blocks to complete
187 
188     unit->m_phaseaccuracy = ZIN0(3); // 0.02; //20 msec resolution; could be argument of UGen
189 
190     unit->m_numphases = (int*)RTAlloc(unit->mWorld, g_numtempi * sizeof(int));
191     // unit->m_phases = (float**)RTAlloc(unit->mWorld, g_numtempi * sizeof(float*));
192 
193     for (int j = 0; j < g_numtempi; ++j) {
194         float period = g_periods[j];
195 
196         int num = (int)(period / unit->m_phaseaccuracy); // maximum will be 1.0/0.02 = 50
197 
198         unit->m_numphases[j] = num;
199 
200         //
201         //	unit->m_phases[j]= (float*)RTAlloc(unit->mWorld, unit->m_numphases[j] * sizeof(float));
202         //
203         //	float phase=0.0;
204         //
205         //	for (i=0; i<num; ++i) {
206         //	unit->m_phases[j][i] = phase;
207         //	phase += unit->m_phaseaccuracy;
208         //	}
209     }
210 
211     unit->m_numfeatures = (int)(ZIN0(1) + 0.001);
212 
213 
214     // for efficiency
215     unit->m_scores = (float*)RTAlloc(unit->mWorld, (2 * unit->m_numfeatures) * sizeof(float));
216 
217     unit->m_temporalwindowsize =
218         ZIN0(2); // typically small, 2 seconds for fast reactions compared to 6 secs for BeatTrack
219 
220     unit->m_fullwindowsize = unit->m_temporalwindowsize + 1.0
221         + 0.1; // plus one to cover all phases of the 60bpm based period, and a further 0.1 for indexing safety; ie
222                // looking at areas around the point you're interested in
223 
224     unit->m_buffersize = (int)(unit->m_fullwindowsize / unit->m_krlength); // in control blocks
225 
226 
227     // printf("loading test blocklength %f numfeatures %d temporal %f full %f blocks %d \n",unit->m_krlength,
228     // unit->m_numfeatures, unit->m_temporalwindowsize, unit->m_fullwindowsize, unit->m_buffersize);
229 
230 
231     // float ** m_pastfeatures;  //for each feature, a trail of last m_workingmemorysize values
232     unit->m_pastfeatures = (float**)RTAlloc(unit->mWorld, unit->m_numfeatures * sizeof(float*));
233 
234     for (int j = 0; j < unit->m_numfeatures; ++j) {
235         unit->m_pastfeatures[j] = (float*)RTAlloc(unit->mWorld, unit->m_buffersize * sizeof(float));
236 
237         Clear(unit->m_buffersize, unit->m_pastfeatures[j]); // set all to zero at first
238 
239         // for (i=0; i<unit->m_buffersize; ++i) {
240         //	unit->m_pastfeatures[j][i] = 0.0;
241         //	}
242         //
243     }
244 
245     // main counter
246     unit->m_counter = 0;
247 
248     // could avoid allocation by having a hard limit on
249     unit->bestscore = (float*)RTAlloc(unit->mWorld, 4 * unit->m_numfeatures * sizeof(float));
250     unit->bestphase = (int*)RTAlloc(unit->mWorld, 4 * unit->m_numfeatures * sizeof(int));
251     unit->besttempo = (int*)RTAlloc(unit->mWorld, 4 * unit->m_numfeatures * sizeof(int));
252     unit->bestgroove = (int*)RTAlloc(unit->mWorld, 4 * unit->m_numfeatures * sizeof(int));
253 
254     for (int i = 0; i < 4; ++i) {
255         int basepos = i * unit->m_numfeatures;
256 
257         for (int j = 0; j < unit->m_numfeatures; ++j) {
258             unit->bestscore[basepos + j] = -9999.0;
259             unit->bestphase[basepos + j] = 0;
260             unit->besttempo[basepos + j] = 60;
261             unit->bestgroove[basepos + j] = 0;
262         }
263     }
264 
265     unit->m_phase = 0.0;
266     unit->m_period = 0.5;
267     unit->m_groove = 0;
268     unit->m_currtempo = 2;
269     unit->m_phaseperblock = unit->m_krlength / unit->m_period;
270 
271     unit->m_predictphase = 0.4f;
272     unit->m_predictperiod = 0.3f;
273 
274 
275     unit->m_outputphase = unit->m_phase;
276     unit->m_outputtempo = unit->m_currtempo;
277     unit->m_outputgroove = unit->m_groove;
278     unit->m_outputphaseperblock = unit->m_phaseperblock;
279 
280 
281     unit->m_calculationperiod = 0.5; // every half second; could also be additional argument to UGen
282     unit->m_calculationschedule = 0.0;
283 
284     // printf("srate %f conversion factor %f frame period %f \n", unit->m_srate, unit->m_srateconversion,
285     // unit->m_frameperiod);
286 
287 
288     int bufnum = (int)(ZIN0(5) + 0.001f);
289     if (bufnum >= unit->mWorld->mNumSndBufs)
290         bufnum = 0;
291 
292     if (bufnum < 0)
293         unit->m_weightingscheme = bufnum < 2 ? 0 : 1;
294     else {
295         SndBuf* buf = unit->mWorld->mSndBufs + bufnum;
296         unit->m_tempoweights = buf;
297         unit->m_weightingscheme = 2;
298     }
299 
300     // printf("bufnum %d weightingscheme %d check %f %f\n", bufnum, unit->m_weightingscheme, unit->m_tempoweights[0],
301     // unit->m_tempoweights[119]);
302 
303 
304     unit->halftrig = 0;
305     unit->q1trig = 0;
306     unit->q2trig = 0;
307 
308 
309     unit->mCalcFunc = (UnitCalcFunc)&BeatTrack2_next;
310 
311     // initialize outputs
312     ZOUT0(0) = 0.0;
313     ZOUT0(1) = 0.0;
314     ZOUT0(2) = 0.0;
315     ZOUT0(3) = unit->m_outputtempo;
316     ZOUT0(4) = unit->m_outputphase;
317     ZOUT0(5) = unit->m_outputgroove;
318 }
319 
320 
BeatTrack2_Dtor(BeatTrack2 * unit)321 void BeatTrack2_Dtor(BeatTrack2* unit) {
322     RTFree(unit->mWorld, unit->m_numphases);
323 
324     RTFree(unit->mWorld, unit->m_scores);
325 
326     RTFree(unit->mWorld, unit->bestscore);
327     RTFree(unit->mWorld, unit->bestphase);
328     RTFree(unit->mWorld, unit->besttempo);
329 
330     for (int j = 0; j < unit->m_numfeatures; ++j)
331         RTFree(unit->mWorld, unit->m_pastfeatures[j]);
332 
333     RTFree(unit->mWorld, unit->m_pastfeatures);
334 }
335 
336 
337 // over phases and for each groove
calculatetemplate(BeatTrack2 * unit,int which,int j)338 void calculatetemplate(BeatTrack2* unit, int which, int j) {
339     int tmpindex;
340 
341     int startcounter = unit->m_startcounter;
342 
343     int numphases = unit->m_numphases[which];
344 
345     float period = g_periods[which];
346 
347     float blockconvert = unit->m_krlength;
348 
349     float windowsize = unit->m_temporalwindowsize;
350 
351     int buffersize = unit->m_buffersize; // unit->m_fullwindowsize/unit->m_krlength; //in control blocks
352 
353     float** pastfeatures = unit->m_pastfeatures;
354     // unit->m_pastfeatures = (float**)RTAlloc(unit->mWorld, unit->m_numfeatures * sizeof(float*));
355 
356     int beatsfit = (int)(windowsize / period); // complete beats only, or also fit as many as possible?
357 
358     float weight; // compensation for number of events matched; may alter equation later
359 
360     switch (unit->m_weightingscheme) {
361     case 0:
362         weight = 1.0f; // flat
363         break;
364     case 1:
365         weight = 1.0f / (beatsfit * 4); // compensate for number of time points tested
366         break;
367     case 2:
368         SndBuf* buf = unit->m_tempoweights;
369         if (buf->data)
370             weight = buf->data[which]; // user defined temmpo biases (usually a mask on allowed tempi)
371         else
372             weight = 1.f;
373         break;
374     }
375 
376     int numfeatures = unit->m_numfeatures;
377 
378     float* scores = unit->m_scores; //[2*numfeatures];
379 
380     float* bestscore = unit->bestscore;
381     int* bestphase = unit->bestphase;
382     int* besttempo = unit->besttempo;
383     int* bestgroove = unit->bestgroove;
384 
385     for (int i = 0; i < numphases; ++i) {
386         // initialise scores
387         // for (j=0; j<2; ++j)
388         for (int k = 0; k < numfeatures; ++k)
389             scores[2 * k + j] = 0.0;
390 
391         float phaseadd = i * unit->m_phaseaccuracy;
392 
393         // calculation for a particular phase of template
394         // for (j=0; j<2; ++j) {
395 
396         for (int h = 0; h < beatsfit; ++h) {
397             for (int l = 0; l < 4; ++l) {
398                 float sep = phaseadd + (h * period) + ((g_sep[j * 4 + l]) * period);
399                 float weight = g_weight[l];
400 
401                 int blocks = (int)((sep / blockconvert) + 0.5); // round to nearest
402 
403                 // convert sep to control periods and find appropriate point in source data
404                 int index = (startcounter + buffersize - blocks) % (buffersize);
405 
406 
407                 // widen over four either side
408                 for (int m = (-4); m < 5; ++m) {
409                     int actualindex = (index + buffersize + m) % (buffersize);
410 
411                     for (int k = 0; k < numfeatures; ++k) {
412                         int scoreindexnow = 2 * k + j;
413 
414                         // could widen this value here, even based on cubic interpolation etc
415                         scores[scoreindexnow] += weight * (g_weight2[m + 4]) * (pastfeatures[k][actualindex]);
416                     }
417 
418                     // scores[2*k+j] += weight * (pastfeatures[k][index]);
419                 }
420             }
421         }
422         //}
423 
424         // update any winners from scores
425         // for (j=0; j<2; ++j) {
426 
427         for (int k = 0; k < numfeatures; ++k) {
428             float scorenow = (scores[2 * k + j]) * weight;
429 
430             // NEED TO STORE J IF PRESERVING SENSE OF GROOVE
431 
432             if (scorenow > bestscore[k]) {
433                 tmpindex = numfeatures + k;
434                 // shift up to make room
435                 bestscore[tmpindex] = bestscore[k];
436                 bestphase[tmpindex] = bestphase[k];
437                 besttempo[tmpindex] = besttempo[k];
438                 bestgroove[tmpindex] = bestgroove[k];
439 
440                 bestscore[k] = scorenow;
441                 bestphase[k] = i;
442                 besttempo[k] = which;
443                 bestgroove[k] = j;
444 
445                 // printf("bestscore %f bestphase %d besttempo %d bestgroove %d \n",
446                 // bestscore[k],bestphase[k],besttempo[k], bestgroove[k]);
447             } else if (scorenow > bestscore[numfeatures + k]) {
448                 tmpindex = numfeatures + k;
449                 bestscore[tmpindex] = scorenow;
450                 bestphase[tmpindex] = i;
451                 besttempo[tmpindex] = which;
452                 bestgroove[tmpindex] = j;
453             }
454         }
455         //}
456     }
457 }
458 
459 
460 // a winner must appear at least twice, across features, and be superior to the secondbest in those features too by some
461 // margins a consistency check could also run to look at change from last time to this
finaldecision(BeatTrack2 * unit)462 void finaldecision(BeatTrack2* unit) {
463     // int foundgood = 0;
464     int bestcandidate = 0;
465     int bestpreviousmatchsum = 0; //(-1);  //should be 0, but allowing different for now
466     float excess; //, consistency;
467                   // int exactmatches, closematches;  //can be out by a few indices on period; could match on tempo but
468                   // not phase etc combine these four factors in one master score?
469 
470     for (int i = 0; i < unit->m_numfeatures; ++i) {
471         int matchsum = 0;
472 
473         float secondbest = unit->bestscore[unit->m_numfeatures + i];
474         excess = (secondbest != 0) ? (unit->bestscore[i] / secondbest) : unit->bestscore[i];
475         int tempo = unit->besttempo[i];
476 
477         // could check consistency too by looking at phase update from last prediction in same feature
478 
479         for (int j = 0; j < unit->m_numfeatures; ++j) {
480             if (j != i) {
481                 if (abs(unit->besttempo[j] - tempo) < 5)
482                     matchsum++;
483             }
484 
485             // check over all previous features
486             if (abs(unit->besttempo[2 * unit->m_numfeatures + j] - tempo) < 5)
487                 matchsum++;
488         }
489 
490         // printf("i %d matchsum %d excess %f \n",i, matchsum, excess);
491 
492         if (secondbest != 0)
493             matchsum += (int)excess;
494 
495         // so must have at least one match //&& (excess>1.03)
496         if ((matchsum > bestpreviousmatchsum)) {
497             bestcandidate = i;
498             bestpreviousmatchsum = matchsum;
499             // foundgood = 1;
500         }
501     }
502 
503 
504     // consistency: could require it to win twice; have a candidatepending which makes a phase prediction; only let
505     // through if prediction fulfilled
506 
507     // unit->m_amortlength will be numtempi *2  = 240
508 
509     float bestphase =
510         fmod(((unit->bestphase[bestcandidate] * unit->m_phaseaccuracy) + (unit->m_krlength * (unit->m_amortlength)))
511                  / (unit->m_period),
512              (float)1.0);
513 
514     // if(unit->m_prediction) {
515 
516     if ((fabs(bestphase - unit->m_predictphase) < ((2 * (unit->m_phaseaccuracy)) / unit->m_predictperiod))
517         && (fabs((g_periods[unit->besttempo[bestcandidate]]) - unit->m_predictperiod) < 0.04)) {
518         unit->m_period = unit->m_predictperiod;
519         // time elapsed since a known beat is phase of winner in seconds, to calculation start point, plus time for
520         // calculation (120 control blocks) divided by period, modulo 1.0
521         unit->m_phase = bestphase;
522         unit->m_currtempo = 1.f / unit->m_period;
523         unit->m_phaseperblock = unit->m_krlength / unit->m_period;
524     }
525 
526     //}
527 
528     // unit->m_prediction=false;
529 
530 
531     // if(foundgood) {
532     // if clear winner
533 
534     unit->m_predictperiod = g_periods[unit->besttempo[bestcandidate]];
535 
536     // time elapsed since a known beat is phase of winner in seconds, to calculation start point, plus time for
537     // calculation (120 control blocks) divided by period, modulo 1.0
538     unit->m_predictphase = fmod(((unit->bestphase[bestcandidate] * unit->m_phaseaccuracy)
539                                  + (unit->m_krlength * (unit->m_amortlength)) + unit->m_calculationperiod)
540                                     / (unit->m_period),
541                                 (float)1.0);
542 
543 
544     // if(foundgood) {
545     ////if clear winner
546     //
547     // unit->m_period = g_periods[unit->besttempo[bestcandidate]];
548     ////time elapsed since a known beat is phase of winner in seconds, to calculation start point, plus time for
549     /// calculation (120 control blocks) divided by period, modulo 1.0
550     // unit->m_phase= fmod( ((unit->bestphase[bestcandidate] * unit->m_phaseaccuracy)  + (unit->m_krlength *
551     // 120))/(unit->m_period), 1.0);
552     //
553     // unit->m_currtempo = 1.0/unit->m_period;
554     // unit->m_phaseperblock = unit->m_krlength/unit->m_period;
555     //}
556 }
557 
558 
BeatTrack2_next(BeatTrack2 * unit,int wrongNumSamples)559 void BeatTrack2_next(BeatTrack2* unit, int wrongNumSamples) {
560     // keep updating feature memories
561     unit->m_counter = (unit->m_counter + 1) % (unit->m_buffersize);
562 
563     int busnum = (int)(ZIN0(0) + 0.001f);
564 
565     // unit->m_features = unit->mWorld->mControlBus + busnum;
566 
567     float* features = unit->mWorld->mControlBus + busnum;
568 
569     // hmm, is this pointer guaranteed to stay the same? may have to update each time...
570     for (int j = 0; j < unit->m_numfeatures; ++j) {
571         unit->m_pastfeatures[j][unit->m_counter] = features[j]; // unit->m_features[j];
572     }
573 
574     unit->m_calculationschedule += unit->m_krlength;
575 
576     // check for new calculation round
577     if (unit->m_calculationschedule > unit->m_calculationperiod) {
578         unit->m_calculationschedule -= unit->m_calculationperiod;
579 
580         // reset best scores and move old to previous slots
581         for (int i = 0; i < 2; ++i) {
582             int pos1 = (2 + i) * unit->m_numfeatures;
583             int pos2 = i * unit->m_numfeatures;
584 
585             for (int j = 0; j < unit->m_numfeatures; ++j) {
586                 unit->bestscore[pos1 + j] = unit->bestscore[pos2 + j];
587                 unit->bestscore[pos2 + j] = -9999.0;
588                 unit->bestphase[pos1 + j] = unit->bestphase[pos2 + j];
589                 unit->bestphase[pos2 + j] = 0;
590                 unit->besttempo[pos1 + j] = unit->besttempo[pos2 + j];
591                 unit->besttempo[pos2 + j] = 60;
592             }
593         }
594 
595         // state 0 is do nothing
596         unit->m_amortisationstate = 1;
597         unit->m_amortcount = 0;
598         unit->m_amortlength = g_numtempi * 2; //
599                                               // unit->m_amortisationsteps=0;
600 
601         // store essential data
602         unit->m_startcounter = unit->m_counter;
603 
604         unit->m_currphase = unit->m_phase;
605     }
606 
607 
608     // keeps incrementing but will be reset with each calculation run
609     // unit->m_amortisationsteps=unit->m_amortisationsteps+1;
610 
611     // if state nonzero do something
612     switch (unit->m_amortisationstate) {
613     case 0:
614         break; // do nothing case
615     case 1: // calculate acf
616         calculatetemplate(unit, unit->m_amortcount >> 1, unit->m_amortcount % 2);
617 
618         unit->m_amortcount = unit->m_amortcount + 1;
619 
620         if (unit->m_amortcount == unit->m_amortlength) {
621             unit->m_amortisationstate = 2;
622             // unit->m_amortlength=1;
623             // unit->m_amortcount=0;
624         }
625         break;
626     case 2: // done calculating template matches, now decide whether to follow through
627         finaldecision(unit);
628         unit->m_amortisationstate = 0;
629         break;
630 
631     default:
632         break;
633     }
634 
635 
636     // test if impulse to output
637     unit->m_phase += unit->m_phaseperblock;
638 
639     // if(unit->m_counter%400==0) printf("phase %f period %f\n", unit->m_phase, unit->m_period);
640 
641     // if not locked, update output phase from model phase, else keep a separate output phase
642 
643     float lock = ZIN0(4);
644     // printf("lock %f \n",lock);
645 
646     if (lock < 0.5f) {
647         unit->m_outputphase = unit->m_phase;
648         unit->m_outputtempo = unit->m_currtempo;
649         unit->m_outputgroove = unit->m_groove;
650         unit->m_outputphaseperblock = unit->m_phaseperblock;
651     } else {
652         unit->m_outputphase += unit->m_outputphaseperblock;
653     }
654 
655     if (unit->m_phase >= 1.f) {
656         unit->m_phase -= 1.f;
657     }
658 
659     // 0 is beat, 1 is quaver, 2 is semiquaver, 3 is actual current tempo in bps
660     // so no audio accuracy with beats, just asap, may as well be control rate
661     ZOUT0(0) = 0.0;
662     ZOUT0(1) = 0.0;
663     ZOUT0(2) = 0.0;
664     ZOUT0(3) = unit->m_outputtempo; //*0.016666667;
665     ZOUT0(4) = unit->m_outputphase;
666     ZOUT0(5) = unit->m_outputgroove;
667 
668     // output beat
669     if (unit->m_outputphase >= 1.f) {
670         // printf("beat \n");
671 
672         unit->m_outputphase -= 1.f;
673         ZOUT0(0) = 1.0;
674         ZOUT0(1) = 1.0;
675         ZOUT0(2) = 1.0;
676         unit->halftrig = 0;
677         unit->q1trig = 0;
678         unit->q2trig = 0;
679     }
680 
681     if (unit->m_outputphase >= 0.5 && unit->halftrig == 0) {
682         ZOUT0(1) = 1.0;
683         ZOUT0(2) = 1.0;
684         unit->halftrig = 1;
685     }
686 
687     float groove = unit->m_outputgroove * 0.07;
688 
689     if (unit->m_outputphase >= (0.25 + groove) && unit->q1trig == 0) {
690         ZOUT0(2) = 1.0;
691         unit->q1trig = 1;
692     }
693 
694     if (unit->m_outputphase >= (0.75 + groove) && unit->q2trig == 0) {
695         ZOUT0(2) = 1.0;
696         unit->q2trig = 1;
697     }
698 }
699