Essentia  2.1-beta5-dev
bpmutil.h
Go to the documentation of this file.
1 
2 /*
3  * Copyright (C) 2006-2016 Music Technology Group - Universitat Pompeu Fabra
4  *
5  * This file is part of Essentia
6  *
7  * Essentia is free software: you can redistribute it and/or modify it under
8  * the terms of the GNU Affero General Public License as published by the Free
9  * Software Foundation (FSF), either version 3 of the License, or (at your
10  * option) any later version.
11  *
12  * This program is distributed in the hope that it will be useful, but WITHOUT
13  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14  * FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
15  * details.
16  *
17  * You should have received a copy of the Affero GNU General Public License
18  * version 3 along with this program. If not, see http://www.gnu.org/licenses/
19  */
20 
21 #ifndef ESSENTIA_BPMUTILS_H
22 #define ESSENTIA_BPMUTILS_H
23 
24 #include "../essentiamath.h"
25 #include <cassert>
26 
27 namespace essentia {
28 
29 inline
30 Real lagToBpm(Real lag, Real sampleRate, Real hopSize) {
31  return 60.0*sampleRate/lag/hopSize;
32 }
33 
34 inline
35 Real bpmToLag(Real bpm, Real sampleRate, Real hopSize) {
36  return lagToBpm(bpm, sampleRate, hopSize);
37 }
38 
39 inline
40 int longestChain(const std::vector<Real>& dticks, int startpos, Real period, Real tolerance) {
41  int pos = startpos;
42  Real ubound = period*(1+tolerance);
43  Real lbound = period*(1-tolerance);
44  while ((pos < (int)dticks.size()) &&
45  (lbound < dticks[pos] && dticks[pos] < ubound)) {
46  pos++;
47  }
48 
49  return pos - startpos;
50 }
51 
52 inline
53 void bpmDistance(Real x, Real y, Real& error, Real& ratio) {
54  ratio = x/y;
55  error = -1;
56  if (ratio < 1) {
57  ratio = round(1./ratio);
58  error=(x*ratio-y)/std::min(y,Real(x*ratio))*100;
59  }
60  else {
61  ratio = round(ratio);
62  error = (x-y*ratio)/std::min(x,Real(y*ratio))*100;
63  }
64 }
65 
66 inline
67 bool areEqual(Real a, Real b, Real tolerance) {
68  //return fabs(a-b) <= epsilon;
69  Real error=0;
70  Real ratio=0;
71  bpmDistance(a,b,error,ratio);
72  return (std::fabs(error)<tolerance) && (int(ratio)==1);
73 }
74 
75 inline
76 bool areHarmonics(Real x, Real y, Real epsilon, bool bPower2) {
77  // epsilon must be in %. a strict choice could be 3
78  Real ratio = 0;
79  Real error = 0;
80  bpmDistance(x, y, error, ratio);
81  error = std::fabs(error);
82  if (error <= epsilon) {
83  if (bPower2) return isPowerTwo(int(fabs(ratio)));
84  return true;
85  }
86  return false;
87 }
88 
89 inline
91  // epsilon must be in %. a strict choice could be 3
92  if (x<y) return greatestCommonDivisor(y,x,epsilon);
93  if (x==0) return 0;
94  Real error = std::numeric_limits<int>::max(),
95  ratio=std::numeric_limits<int>::max();
96  bpmDistance(x,y,error,ratio);
97  if (fabs(error)<epsilon) return y;
98  int a = int(x+0.5);
99  int b = int(y+0.5);
100  while (fabs(error) > epsilon) {
101  bpmDistance(a,b,error,ratio);
102  int remainder = a%b;
103  a=b;
104  b=remainder;
105  //if(x<1) return 1;
106  }
107  return a;
108 }
109 
110 
111 inline
112 std::vector<Real> roundBpms(const std::vector<Real>& bpms) {
113  Real epsilon = 3; // 3%
114  Real mainBpm=bpms[0];
115  std::vector<Real> harmonicBpms;
116  harmonicBpms.reserve(bpms.size());
117  for (int i=0; i<int(bpms.size()); i++) {
118  Real ratio=bpms[0]/mainBpm;
119  if (ratio < Real(1.0)) ratio = 1.0/ratio;
120  ratio = round(ratio*10.)/10.; // rounding to 1 float pos
121  int iRatio = int(ratio);
122  if (ratio-iRatio <= 0.100001) { // allow 2.9, 3.1 be considered as 3
123  harmonicBpms.push_back(bpms[i]);
124  }
125  else {
126  if ((ratio-iRatio) == 0.5) { // only interested in pure halfs
127  harmonicBpms.push_back(bpms[i]);
128  harmonicBpms.push_back(greatestCommonDivisor(bpms[i], mainBpm,epsilon));
129  }
130  }
131  }
132  return harmonicBpms;
133 }
134 
135 
136 // original postprocessticks from essentia 1.0
137 inline
138 std::vector<Real> postProcessTicks(const std::vector<Real>& origticks) {
139  if (origticks.size() < 3) return origticks;
140 
141  // find the most likely beat period
142  const int nticks = origticks.size();
143  std::vector<Real> dticks(nticks-1);
144 
145  for (int i=0; i<nticks-1; i++) dticks[i] = origticks[i+1] - origticks[i];
146 
147  // we might have had 6 secs frames during which we didn't find any beat, in which
148  // case we'll have one huge dtick value, which we actually want to prune
149  for (int i=0; i<(int)dticks.size(); i++) {
150  if (dticks[i] > 2.) {
151  dticks.erase(dticks.begin() + i);
152  i--;
153  }
154  }
155 
156  const int nbins = 100;
157  std::vector<int> dist(nbins);
158  std::vector<Real> distx(nbins);
159 
160  hist(&dticks[0], nticks-1, &dist[0], &distx[0], nbins);
161 
162  int maxidx = max_element(dist.begin(), dist.end()) - dist.begin();
163  Real maxbinCenter = distx[maxidx];
164 
165  // find the longest sequence of beats which has a fixed period of the previously
166  // found value; use a tolerance of about 10%
167  // Note: this favors high BPMs, because they will have more beats in the same amount of time
168  int maxl = 0;
169  int idx = 0;
170  Real period = maxbinCenter;
171 
172  for (int startpos = 0; startpos < nticks-1; startpos++) {
173  int l = longestChain(dticks, startpos, period, 0.1);
174  if (l > maxl) {
175  maxl = l;
176  idx = startpos;
177  }
178  }
179 
180  if (idx == 0 && maxl == 0) {
181  std::cerr << "Internal error while processing the beats, returning the original ones" << std::endl;
182  return origticks;
183  }
184 
185  // let's assume those beats are correct, and try to replace all the other ones
186  // with respect to the fixed period we have and the old positions of the beats
187  std::deque<Real> ticks(origticks.begin() + idx,
188  origticks.begin() + idx + maxl + 1);
189 
190  // take this value as the period for the whole track
191  Real targetPeriod = mean(dticks, idx, idx+maxl);
192  // 0.15, because 0.1 might be too strict, while 0.2 will allow false positives more easily
193  Real tolerance = 0.15 * targetPeriod;
194 
195 
196  // do the beats after the current beat base
197  Real cpos = ticks.back() + targetPeriod;
198  std::deque<Real> remaining(origticks.begin() + idx + maxl + 1,
199  origticks.end());
200 
201  while (!remaining.empty()) {
202  Real nbeat = remaining.front();
203 
204  if (nbeat < cpos - tolerance) {
205  // too far before, drop next beat
206  remaining.pop_front();
207  }
208  else {
209  // right in our expected zone, adjust the estimated beat to the one
210  // we actually found (NB: if we're too far away in front, just keep the
211  // beat as is)
212  if (nbeat < cpos + tolerance) {
213  cpos = nbeat;
214  remaining.pop_front();
215  }
216 
217  // in any case, mark the beat and jump on the next one
218  ticks.push_back(cpos);
219  cpos += targetPeriod;
220  }
221  }
222 
223  // do the beats before the current beat base
224  cpos = ticks.front() - targetPeriod;
225  remaining = std::deque<Real>(origticks.begin(),
226  origticks.begin() + idx);
227 
228  while (!remaining.empty()) {
229  Real nbeat = remaining.back();
230 
231  if (nbeat > cpos + tolerance) {
232  // too far after, drop beat
233  remaining.pop_back();
234  }
235  else {
236  // right in our expected zone, adjust the estimated beat to the one
237  // we actually found
238  if (nbeat > cpos - tolerance) {
239  cpos = nbeat;
240  remaining.pop_back();
241  }
242 
243  // in any case, mark the beat and jump on the next one
244  ticks.push_front(cpos);
245  cpos -= targetPeriod;
246  }
247  }
248 
249 
250  return std::vector<Real>(ticks.begin(), ticks.end());
251 }
252 
253  // modified version of the postprocessticks from tempotapticks, so it does not
254 // tend to favour fast bpms
255 inline
256 std::vector<Real> postProcessTicks(const std::vector<Real>& origticks,
257  const std::vector<Real>& ticksAmplitudes,
258  const Real& preferredPeriod) {
259  if (origticks.size() < 3) return origticks;
260 
261  // find the most likely beat period
262  const int nticks = origticks.size();
263  std::vector<Real> dticks(nticks-1);
264 
265  for (int i=0; i<nticks-1; i++) dticks[i] = origticks[i+1] - origticks[i];
266 
267  // we might have had 6 secs frames during which we didn't find any beat, in which
268  // case we'll have one huge dtick value, which we actually want to prune
269  for (int i=0; i<(int)dticks.size(); i++) {
270  if (dticks[i] > 2.) {
271  dticks.erase(dticks.begin() + i);
272  i--;
273  }
274  }
275 
276  const int nbins = 100;
277  std::vector<int> dist(nbins);
278  std::vector<Real> distx(nbins);
279 
280  hist(&dticks[0], dticks.size(), &dist[0], &distx[0], nbins);
281 
282  int maxidx = max_element(dist.begin(), dist.end()) - dist.begin();
283  Real maxbinCenter = distx[maxidx];
284 
285  // find the longest sequence of beats which has a fixed period of the previously
286  // found value; use a tolerance of about 10%
287  // Note: this favors high BPMs, because they will have more beats in the same amount of time
288  int maxl = 0;
289  int idx = 0;
290  Real period = maxbinCenter;
291 
292  //std::cout << "period: " << period << std::endl;
293  for (int startpos = 0; startpos < nticks-1; startpos++) {
294  int l = longestChain(dticks, startpos, period, 0.05);
295  if (l > maxl) {
296  maxl = l;
297  idx = startpos;
298  }
299  }
300 
301  Real targetPeriod = preferredPeriod;
302  if (idx ==0 && maxl==0) {
303  idx = argmax(ticksAmplitudes);
304  }
305  // take this value as the period for the whole track
306  else targetPeriod = mean(dticks, idx, idx+maxl);
307 
308  targetPeriod = (targetPeriod+preferredPeriod)/2.0;
309 
310  // if the targetPeriod is differs too much from the preferred period we
311  // search for the tick with max amplitude and take that as the reference tick
312  if (targetPeriod < 0.95*preferredPeriod || targetPeriod > 1.05*preferredPeriod) {
313  idx = idx + argmax(std::vector<Real>(ticksAmplitudes.begin()+idx, ticksAmplitudes.begin()+idx+maxl+1));
314  maxl = 0;
315  targetPeriod = preferredPeriod;
316  //std::cout << "Targets differ too much!. New target period will be the preferred one " << targetPeriod << std::endl;
317  }
318 
319  Real originalTargetPeriod = targetPeriod;
320 
321  // let's assume those beats are correct, and try to replace all the other ones
322  // with respect to the fixed period we have and the old positions of the beats
323  std::deque<Real> ticks(origticks.begin() + idx,
324  origticks.begin() + idx + maxl + 1);
325 
326  //fix tolerance at no more than 30ms:
327  Real tolerance =0.03;
328 
329  if (targetPeriod < 0.05) {
330  std::cout << "PostProcessTicks: target Period too short. Returning the original ticks" << std::endl;
331  return origticks;
332  }
333  Real cummulatedPeriod = targetPeriod;
334  int nAccumulations = 1;
335 
336 
337  // do the beats after the current beat base
338  Real cpos = ticks.back() + targetPeriod;
339  std::deque<Real> remaining(origticks.begin() + idx + maxl + 1,
340  origticks.end());
341 
342  while (!remaining.empty()) {
343  Real nbeat = remaining.front();
344 
345  if (nbeat < cpos - tolerance) {
346  // too far before, drop next beat
347  remaining.pop_front();
348  }
349  else {
350  // right in our expected zone, adjust the estimated beat to the one
351  // we actually found (NB: if we're too far away in front, just keep the
352  // beat as is)
353  if (nbeat < cpos + tolerance) {
354  cummulatedPeriod += (nbeat - (cpos - targetPeriod));
355  nAccumulations++;
356  targetPeriod = cummulatedPeriod/nAccumulations;
357  //std::cout << "new target Period: " << targetPeriod << " bpm: " << 60./targetPeriod << std::endl;
358  //std::cout << " \tbeat at: " << nbeat << " belongs to [" << cpos-tolerance << ", " << cpos+tolerance <<"], cpos: " << cpos << std::endl;
359  cpos = nbeat;
360  remaining.pop_front();
361  }
362 
363  // in any case, mark the beat and jump on the next one
364  ticks.push_back(cpos);
365  cpos += targetPeriod;
366  }
367  }
368 
369  // do the beats before the current beat base
370  cpos = ticks.front() - targetPeriod;
371  remaining = std::deque<Real>(origticks.begin(),
372  origticks.begin() + idx);
373 
374  targetPeriod = originalTargetPeriod;
375  cummulatedPeriod = targetPeriod;
376  nAccumulations = 1;
377 
378  while (!remaining.empty()) {
379  Real nbeat = remaining.back();
380 
381  if (nbeat > cpos + tolerance) {
382  // too far after, drop beat
383  remaining.pop_back();
384  }
385  else {
386  // right in our expected zone, adjust the estimated beat to the one
387  // we actually found
388  if (nbeat > cpos - tolerance) {
389  cummulatedPeriod += ((cpos + targetPeriod)-nbeat);
390  nAccumulations++;
391  targetPeriod = cummulatedPeriod/nAccumulations;
392  //std::cout << "new target Period: " << targetPeriod << " bpm: " << 60./targetPeriod << std::endl;
393  //std::cout << " \tbeat at: " << nbeat << " belongs to [" << cpos-tolerance << ", " << cpos+tolerance <<"], cpos: " << cpos << std::endl;
394  cpos = nbeat;
395  remaining.pop_back();
396  }
397 
398  // in any case, mark the beat and jump on the next one
399  ticks.push_front(cpos);
400  cpos -= targetPeriod;
401  }
402  }
403 
404  return std::vector<Real>(ticks.begin(), ticks.end());
405 }
406 } // namespace essentia
407 
408 #endif
bool areEqual(Real a, Real b, Real tolerance)
Definition: bpmutil.h:67
T round(const T value)
Definition: essentiamath.h:519
Real bpmToLag(Real bpm, Real sampleRate, Real hopSize)
Definition: bpmutil.h:35
std::vector< Real > roundBpms(const std::vector< Real > &bpms)
Definition: bpmutil.h:112
void bpmDistance(Real x, Real y, Real &error, Real &ratio)
Definition: bpmutil.h:53
float Real
Definition: types.h:68
Real lagToBpm(Real lag, Real sampleRate, Real hopSize)
Definition: bpmutil.h:30
bool isPowerTwo(T n)
Definition: essentiamath.h:45
Definition: algorithm.h:28
std::vector< Real > postProcessTicks(const std::vector< Real > &origticks)
Definition: bpmutil.h:138
T mean(const std::vector< T > &array, int start, int end)
Definition: essentiamath.h:142
bool areHarmonics(Real x, Real y, Real epsilon, bool bPower2)
Definition: bpmutil.h:76
Real greatestCommonDivisor(Real x, Real y, Real epsilon)
Definition: bpmutil.h:90
int longestChain(const std::vector< Real > &dticks, int startpos, Real period, Real tolerance)
Definition: bpmutil.h:40
void hist(const T *array, uint n, int *n_array, T *x_array, uint n_bins)
Definition: essentiamath.h:794
int argmax(const std::vector< Real > &input)
Definition: essentiamath.h:676