ReMAS  1.5
Real-time Musical Accompaniment System
TempoFunctions.h
Go to the documentation of this file.
1 /**************************************************************************
2  * Copyright (C) 2017 by "Information Retrieval and Parallel Computing" *
3  * group (University of Oviedo, Spain), "Interdisciplinary Computation *
4  * and Communication" group (Polytechnic University of Valencia, Spain) *
5  * and "Signal Processing and Telecommunication Systems Research" group *
6  * (University of Jaen, Spain) *
7  * Contact: remaspack@gmail.com *
8  * *
9  * This program is free software; you can redistribute it and/or modify *
10  * it under the terms of the GNU General Public License as published by *
11  * the Free Software Foundation; either version 2 of the License, or *
12  * (at your option) any later version. *
13  * *
14  * This program is distributed in the hope that it will be useful, *
15  * but WITHOUT ANY WARRANTY; without even the implied warranty of *
16  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
17  * GNU General Public License for more details. *
18  * *
19  * You should have received a copy of the GNU General Public License *
20  * along with this program; if not, write to the *
21  * Free Software Foundation, Inc., *
22  * 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. *
23  **************************************************************************
24 */
35 #pragma once
36 
37 #ifndef ATEMPO_H
38 #define ATEMPO_H
39 
40 #include "defines.h"
41 
42 
43 
44 #ifdef OSC
45  #include <NetFunctions.h>
46  int ComputeTempoOSC (STempo *, int, int, int, int *, lo_address *, int, int);
47  int ComputeTempoOSCRL(STempoRL *, int, int, int, int *, lo_address *, int);
48 #endif
49 void ComputeTempo (STempo *, int, int, int, int *, int);
50 void ComputeTempoRL(STempoRL *, int, int, int, int *);
51 
52 //Size, weights and speeds (speed=1.0 at start time) for gaussian filter
53 //#define SizeGauss 3
54 //static const MyType WeighGauss[SizeGauss]={0.10, 0.36, 0.54};
55 //static MyType SpeedGauss[SizeGauss]={1.00, 1.00, 1.00};
56 #define SizeGauss 5
57 static const MyType WeighGauss[SizeGauss]={0.10, 0.15, 0.20, 0.25, 0.30};
58 static MyType SpeedGauss[SizeGauss]={1.00, 1.00, 1.00, 1.00, 1.00};
59 
60 
61 #ifdef OMP
62  #include <omp.h>
63  void pqsort(MyType *, int);
71  void pqsort(MyType *v, int size)
72  {
73  #pragma omp parallel
74  {
75  #pragma omp single nowait
76  { pqsort_inner(v, 0, size-1, size/omp_get_num_threads()); }
77  }
78  }
79 
80 
81  void pqsort_inner(MyType *, int, int, int);
91  void pqsort_inner(MyType *v, int left, int right, int cutoff)
92  {
93  int i=left, j=right;
94  MyType tmp, pivot=v[(left + right) / 2];
95 
96  while (i <= j)
97  {
98  while (v[i] < pivot) i++;
99  while (v[j] > pivot) j--;
100 
101  if (i <= j) { tmp = v[i]; v[i] = v[j]; v[j] = tmp; i++; j--; }
102  }
103 
104  if (((right-left)<cutoff)) {
105  if (left < j) { pqsort_inner(v, left, j, cutoff); }
106  if (i < right) { pqsort_inner(v, i, right, cutoff); }
107  } else {
108  #ifdef OMP4
109  #pragma omp task untied mergeable
110  { pqsort_inner(v, left, j, cutoff); }
111  #pragma omp task untied mergeable
112  { pqsort_inner(v, i, right, cutoff); }
113  #else
114  #pragma omp task untied
115  { pqsort_inner(v, left, j, cutoff); }
116  #pragma omp task untied
117  { pqsort_inner(v, i, right, cutoff); }
118  #endif
119  }
120  }
121 #else
122  int cmpfunc(const void *, const void *);
130  int cmpfunc(const void *a, const void *b)
131  {
132  if (*(MyType *)a < *(MyType *)b)
133  { return -1; }
134  else if (*(MyType *)a == *(MyType *)b)
135  { return 0; }
136  else
137  { return 1; }
138  }
139 #endif
140 
141 MyType TheilSenRegression(int *, int *, int, bool *);
151 MyType TheilSenRegression(int *x, int *y, int L, bool *outlier)
152 {
153  MyType a=.0, b=.0, lerr=.0, slopes[NUMAPPAIRS], err[NUMAPMAX];
154  int i, j, N, M, idx=0;
155 
156  if (L > NUMAPMAX) { N=NUMAPMAX; } else { N=L; }
157 
158  M=(N*(N-1))/2;
159 
160  /* slope between pairs */
161  #ifdef OMP
162  #pragma omp parallel for private(j, idx)
163  for (i=0; i<N; i++)
164  {
165  idx=i*N-((i+2)*(i+1)/2);
166  for (j=i+1; j<N; j++)
167  slopes[idx+j] = (MyType)(y[j]-y[i]) / (MyType)(x[j]-x[i]);
168  }
169  #else
170  for (i = 0; i < N; i++)
171  for (j = i+1; j < N; j++)
172  {
173  slopes[idx] = (MyType)(y[j]-y[i]) / (MyType)(x[j]-x[i]);
174  idx++;
175  }
176  #endif
177 
178  /* median (slow sorting method) */
179  #ifdef OMP
180  pqsort(slopes, M);
181  #else
182  qsort (slopes, M, sizeof(MyType), cmpfunc);
183  #endif
184 
185  if (M%2 == 0)
186  { a = (slopes[M/2-1] + slopes[M/2]) / 2; }
187  else
188  { a = slopes[(M-1)/2]; }
189  if ((a < 0.8) || (a > 1.2)) { *outlier = 1; return a; }
190 
191  #ifdef OMP
192  #pragma omp parallel for reduction(+: b)
193  #endif
194  for (i=0; i<N; i++)
195  b = b + (y[i] - a*x[i]);
196  b = b/N;
197 
198  /* fitting error for each point */
199  #ifdef OMP4
200  #pragma omp simd
201  #else
202  #ifdef OMP
203  #pragma omp parallel for reduction(+: b)
204  #endif
205  #endif
206  for (i = 0; i < N; i++)
207  #ifdef SIMPLE
208  err[i] = fabsf(y[i] - (a*x[i] + b));
209  #else
210  err[i] = fabs (y[i] - (a*x[i] + b));
211  #endif
212 
213  /* check if last point is an outlier */
214  lerr = err[(L-1)%NUMAPMAX];
215  #ifdef OMP
216  pqsort(err, N);
217  #else
218  qsort (err, N, sizeof(MyType), cmpfunc);
219  #endif
220  if (N%2 == 0)
221  { *outlier = (lerr > err[N/2-1]); }
222  else
223  { *outlier = (lerr > err[(N-1)/2]); }
224 
225  return a;
226 }
227 
228 
229 #ifdef OSC
230 
244  int ComputeTempoOSC(STempo *TEMPO, int current, int pos_min, int count_min, int *states_corr, lo_address *DirOSC, int NCli, int NStates)
245  {
246  int i;
247  MyType nextScoreTime, Speed;
248 
249  if ((count_min > TEMPO->PrevState) && (states_corr[count_min]==1) && (count_min < (NStates-1)))
250  {
251  /* Vmatch */
252  Speed=((pos_min - TEMPO->MidiFrame) * TrainJumpTime) / ((current - TEMPO->RealFrame) * RunJumpTime);
253 
254  TEMPO->SynthTime = ((current - TEMPO->RealFrame) * TrainJumpTime) * TEMPO->SynthSpeed + TEMPO->SynthTime;
255 
256  nextScoreTime = DelayTimeMidi * Speed + pos_min * TrainJumpTime;
257 
258  TEMPO->SynthSpeed = (nextScoreTime - TEMPO->SynthTime) / DelayTimeMidi;
259 
260  //TEMPO->SoloSpeed = Speed;
261  Speed=0.0;
262  for (i=0; i<SizeGauss-1; i++)
263  {
264  Speed += SpeedGauss[i]*WeighGauss[i];
265  SpeedGauss[i]= SpeedGauss[i+1];
266  }
267  Speed += SpeedGauss[SizeGauss-1]*WeighGauss[SizeGauss-1];
268 
269  if (TEMPO->SynthSpeed < 0)
270  { TEMPO->SynthSpeed = 0.25; }
271  else if (TEMPO->SynthSpeed > (1.2*Speed))
272  { TEMPO->SynthSpeed = (1.2*Speed); }
273  else if (TEMPO->SynthSpeed < (0.8*Speed))
274  { TEMPO->SynthSpeed = (0.8*Speed); }
275  SpeedGauss[SizeGauss-1]=TEMPO->SynthSpeed;
276 
277  TEMPO->MidiFrame=pos_min;
278  TEMPO->RealFrame=current;
279  TEMPO->PrevState=count_min;
280 
281  TEMPO->SoloSpeed = Speed;
282  TEMPO->NextFrame = current + round((TEMPO->SynthTime - pos_min*TrainJumpTime)/(TEMPO->SoloSpeed - TEMPO->SynthSpeed) / TrainJumpTime);
283  #ifdef TALK
284  printf("Vmatch: frame %d pos_min %d SynthSpeed %1.5f\n", current, pos_min, TEMPO->SynthSpeed);
285  #endif
286  for (i=0; i<NCli; i++) CHECKERR(SendTempo(DirOSC[i], TEMPO->SynthSpeed*100));
287  }
288  else if (TEMPO->NextFrame == current) {
289  /* Vtempo */
290  TEMPO->SynthSpeed = TEMPO->SoloSpeed;
291  #ifdef TALK
292  printf("Vtempo: matching with the soloist. Frame %d, pos_min %d SynthSpeed %1.5f\n", current, pos_min, TEMPO->SynthSpeed);
293  #endif
294  for (i=0; i<NCli; i++) CHECKERR(SendTempo(DirOSC[i], TEMPO->SynthSpeed*100));
295  }
296  return OK;
297  }
298 
299 
313  int ComputeTempoOSCRL(STempoRL *TEMPO, int current, int pos_min, int count_min, int *states_corr, lo_address *DirOSC, int NCli)
314  {
315  int i;
316  MyType SoloTime, TimeDiff, SpeedDiff, a;
317  bool outlier;
318 
319  /* Update synthesizer position */
320  TEMPO->SynthTime = TEMPO->SynthTime + TEMPO->SynthSpeed*(MyType)RunJumpTime;
321 
322  /* Check if anchor point */
323  if ((count_min > TEMPO->PrevState) && (states_corr[count_min]==1))
324  {
325  /* Store anchor point */
326  TEMPO->AudioTimeAP[TEMPO->numap % NUMAPMAX] = current;
327  TEMPO->ScoreTimeAP[TEMPO->numap % NUMAPMAX] = pos_min;
328  TEMPO->numap++;
329  TEMPO->PrevState = count_min;
330 
331  /* Performer position */
332  SoloTime = (MyType)pos_min * (MyType)TrainJumpTime;
333  #ifdef SIMPLE
334  TimeDiff = fabsf(TEMPO->SynthTime - SoloTime);
335  #else
336  TimeDiff = fabs (TEMPO->SynthTime - SoloTime);
337  #endif
338 
339  /* Performer speed */
340  if ((TEMPO->matched) && (TimeDiff>0.050))
341  {
342  a = TheilSenRegression(TEMPO->AudioTimeAP, TEMPO->ScoreTimeAP, TEMPO->numap, &outlier);
343  if (!outlier)
344  {
345  /* Vtempo */
346  TEMPO->SoloSpeed = a;
347 
348  /* Vmatch: synth speed to match performer in DelayTimeMidi score secs */
349  TEMPO->SynthSpeed = TEMPO->SoloSpeed * (SoloTime + DelayTimeMidi - TEMPO->SynthTime) / DelayTimeMidi;
350  if (TEMPO->SynthSpeed < (MyType)0.1)
351  TEMPO->SynthSpeed = (MyType)0.1;
352  else if (TEMPO->SynthSpeed > (MyType)3.0)
353  TEMPO->SynthSpeed = (MyType)3.0;
354 
355  /* Send Vmatch to the synthesizer (conviene bajar volumen) */
356  for (i=0; i<NCli; i++) CHECKERR(SendTempo(DirOSC[i], TEMPO->SynthSpeed*100));
357 
358  /* Schedule next speed command */
359  #ifdef SIMPLE
360  SpeedDiff = fabsf(TEMPO->SoloSpeed - TEMPO->SynthSpeed);
361  TEMPO->NextFrame = current + (int)roundf(TimeDiff / SpeedDiff / (MyType)RunJumpTime);
362  #else
363  SpeedDiff = fabs(TEMPO->SoloSpeed - TEMPO->SynthSpeed);
364  TEMPO->NextFrame = current + (int)round (TimeDiff / SpeedDiff / (MyType)RunJumpTime);
365  #endif
366  TEMPO->matched = 0;
367  #ifdef TALK
368  printf("Vmatch: frame %d pos_min %d SynthSpeed %1.5f\n", current, pos_min, TEMPO->SynthSpeed);
369  #endif
370  }
371  }
372  }
373  /* Check if NextFrame */
374  if (TEMPO->NextFrame == current)
375  {
376  /* Send Vtempo to the synthesizer */
377  TEMPO->SynthSpeed = TEMPO->SoloSpeed;
378  TEMPO->matched = 1;
379  for (i=0; i<NCli; i++) CHECKERR(SendTempo(DirOSC[i], TEMPO->SynthSpeed*100));
380  #ifdef TALK
381  printf("Vtempo: matching with the soloist. Frame %d, pos_min %d SynthSpeed %1.5f\n", current, pos_min, TEMPO->SynthSpeed);
382  #endif
383  }
384  return OK;
385  }
386 #endif
387 
388 
401 void ComputeTempo(STempo *TEMPO, int current, int pos_min, int count_min, int *states_corr, int NStates)
402 {
403  int i;
404  MyType nextScoreTime, Speed;
405 
406  if ((count_min > TEMPO->PrevState) && (states_corr[count_min]==1) && (count_min < (NStates-1)))
407  {
408  /* Vmatch */
409  Speed=((pos_min - TEMPO->MidiFrame) * TrainJumpTime) / ((current - TEMPO->RealFrame) * RunJumpTime);
410 
411  TEMPO->SynthTime = ((current - TEMPO->RealFrame) * TrainJumpTime) * TEMPO->SynthSpeed + TEMPO->SynthTime;
412 
413  nextScoreTime = DelayTimeMidi * Speed + pos_min * TrainJumpTime;
414 
415  TEMPO->SynthSpeed = (nextScoreTime - TEMPO->SynthTime) / DelayTimeMidi;
416 
417  //TEMPO->SoloSpeed = Speed;
418  Speed=0.0;
419  for (i=0; i<SizeGauss-1; i++)
420  {
421  Speed += SpeedGauss[i]*WeighGauss[i];
422  SpeedGauss[i]= SpeedGauss[i+1];
423  }
424  Speed += SpeedGauss[SizeGauss-1]*WeighGauss[SizeGauss-1];
425 
426  if (TEMPO->SynthSpeed < 0)
427  { TEMPO->SynthSpeed = 0.25; }
428  else if (TEMPO->SynthSpeed > (1.2*Speed))
429  { TEMPO->SynthSpeed = (1.2*Speed); }
430  else if (TEMPO->SynthSpeed < (0.8*Speed))
431  { TEMPO->SynthSpeed = (0.8*Speed); }
432  SpeedGauss[SizeGauss-1]=TEMPO->SynthSpeed;
433 
434  TEMPO->MidiFrame=pos_min;
435  TEMPO->RealFrame=current;
436  TEMPO->PrevState=count_min;
437 
438  TEMPO->SoloSpeed = Speed;
439  TEMPO->NextFrame = current + round((TEMPO->SynthTime - pos_min*TrainJumpTime)/(TEMPO->SoloSpeed - TEMPO->SynthSpeed) / TrainJumpTime);
440  #ifdef TALK
441  printf("Vmatch: frame %d pos_min %d SynthSpeed %1.5f\n", current, pos_min, TEMPO->SynthSpeed);
442  #endif
443  }
444  else if (TEMPO->NextFrame == current)
445  {
446  /* Vtempo */
447  TEMPO->SynthSpeed = TEMPO->SoloSpeed;
448  #ifdef TALK
449  printf("Vtempo: matching with the soloist. Frame %d, pos_min %d SynthSpeed %1.5f\n", current, pos_min, TEMPO->SynthSpeed);
450  #endif
451  }
452 }
453 
454 
466 void ComputeTempoRL(STempoRL *TEMPO, int current, int pos_min, int count_min, int *states_corr)
467 {
468  MyType SoloTime, TimeDiff, SpeedDiff, a;
469  bool outlier;
470 
471  /* Update synthesizer position */
472  TEMPO->SynthTime = TEMPO->SynthTime + TEMPO->SynthSpeed*(MyType)RunJumpTime;
473 
474  /* Check if anchor point */
475  if ((count_min > TEMPO->PrevState) && (states_corr[count_min]==1))
476  {
477  /* Store anchor point */
478  TEMPO->AudioTimeAP[TEMPO->numap % NUMAPMAX] = current;
479  TEMPO->ScoreTimeAP[TEMPO->numap % NUMAPMAX] = pos_min;
480  TEMPO->numap++;
481  TEMPO->PrevState = count_min;
482 
483  /* Performer position */
484  SoloTime = (MyType)pos_min * (MyType)TrainJumpTime;
485  #ifdef SIMPLE
486  TimeDiff = fabsf(TEMPO->SynthTime - SoloTime);
487  #else
488  TimeDiff = fabs (TEMPO->SynthTime - SoloTime);
489  #endif
490 
491  /* Performer speed */
492  if ((TEMPO->matched) && (TimeDiff>0.050))
493  {
494  a = TheilSenRegression(TEMPO->AudioTimeAP, TEMPO->ScoreTimeAP, TEMPO->numap, &outlier);
495  if (!outlier)
496  {
497  /* Vtempo */
498  TEMPO->SoloSpeed = a;
499 
500  /* Vmatch: synth speed to match performer in DelayTimeMidi score secs */
501  TEMPO->SynthSpeed = TEMPO->SoloSpeed * (SoloTime + DelayTimeMidi - TEMPO->SynthTime) / DelayTimeMidi;
502  if (TEMPO->SynthSpeed < (MyType)0.1)
503  TEMPO->SynthSpeed = (MyType)0.1;
504  else if (TEMPO->SynthSpeed > (MyType)3.0)
505  TEMPO->SynthSpeed = (MyType)3.0;
506 
507  /* Schedule next speed command */
508  #ifdef SIMPLE
509  SpeedDiff = fabsf(TEMPO->SoloSpeed - TEMPO->SynthSpeed);
510  TEMPO->NextFrame = current + (int)roundf(TimeDiff / SpeedDiff / (MyType)RunJumpTime);
511  #else
512  SpeedDiff = fabs(TEMPO->SoloSpeed - TEMPO->SynthSpeed);
513  TEMPO->NextFrame = current + (int)round (TimeDiff / SpeedDiff / (MyType)RunJumpTime);
514  #endif
515  TEMPO->matched = 0;
516  #ifdef TALK
517  printf("Vmatch: frame %d pos_min %d SynthSpeed %1.5f\n", current, pos_min, TEMPO->SynthSpeed);
518  #endif
519  }
520  }
521  }
522  /* Check if NextFrame */
523  if (TEMPO->NextFrame == current)
524  {
525  /* Send Vtempo to the synthesizer */
526  TEMPO->SynthSpeed = TEMPO->SoloSpeed;
527  TEMPO->matched = 1;
528  #ifdef TALK
529  printf("Vtempo: matching with the soloist. Frame %d, pos_min %d SynthSpeed %1.5f\n", current, pos_min, TEMPO->SynthSpeed);
530  #endif
531  }
532 }
533 
534 #endif
Struct for Compute tempos.
Definition: defines.h:292
File with functions used by ReMAS, both CPU and GPU, for management communications.
int cmpfunc(const void *, const void *)
cmpfunc compares two values.
void ComputeTempo(STempo *, int, int, int, int *, int)
ComputeTempo calculates the tempo for current frame using gaussian filter.
General header file with constants, defines, structs, etc. using by ReMAS, both CPU and GPU...
void ComputeTempoRL(STempoRL *, int, int, int, int *)
ComputeTempoRL calculates tempo and controls synthesizer speed using linear regression.
Struct for Compute tempos.
Definition: defines.h:276
MyType TheilSenRegression(int *, int *, int, bool *)
TheilSenRegression computes Theil-Sen linear regression.