Improved beat detection routine so that it detects true base beat, even if Nth harmonic multiple sub-beat would appear slightly higher than the true base.

source/example/bpm/BPMDetect.cpp

@@ -297,7 +297,7 @@ void BPMDetect::init(int numChannels, int sampleRate)
 
 float BPMDetect::getBpm()
 {
-    float peakPos;
+    double peakPos;
     PeakFinder peakFinder;
 
     // find peak position
@@ -307,5 +307,5 @@ float BPMDetect::getBpm()
     if (peakPos < 1e-6) return 0.0; // detection failed.
 
     // calculate BPM
-    return 60.0f * (((float)sampleRate / (float)decimateBy) / peakPos);
+    return (float)(60.0 * (((double)sampleRate / (double)decimateBy) / peakPos));
 }

source/example/bpm/PeakFinder.cpp

@@ -44,6 +44,8 @@
 
 #include "PeakFinder.h"
 
+#define max(x, y) (((x) > (y)) ? (x) : (y))
+
 
 PeakFinder::PeakFinder()
 {
@@ -123,7 +125,7 @@ int PeakFinder::findCrossingLevel(const float *data, float level, int peakpos, i
 
 
 // Calculates the center of mass location of 'data' array items between 'firstPos' and 'lastPos'
-float PeakFinder::calcMassCenter(const float *data, int firstPos, int lastPos) const
+double PeakFinder::calcMassCenter(const float *data, int firstPos, int lastPos) const
 {
     int i;
     float sum;
@@ -140,38 +142,22 @@ float PeakFinder::calcMassCenter(const float *data, int firstPos, int lastPos) c
 }
 
 
-float PeakFinder::detectPeak(const float *data, int minPos, int maxPos) 
-{
-    #define max(x, y) (((x) > (y)) ? (x) : (y))
 
-    int i;
+/// get exact center of peak near given position by calculating local mass of center
-    int peakpos;                // position of peak level
+double PeakFinder::getPeakCenter(const float *data, int peakpos)
+{
     float peakLevel;            // peak level
     int crosspos1, crosspos2;   // position where the peak 'hump' crosses cutting level
     float cutLevel;             // cutting value
     float groundLevel;          // ground level of the peak
     int gp1, gp2;               // bottom positions of the peak 'hump'
 
-    this->minPos = minPos;
-    this->maxPos = maxPos;
-
-    // find absolute peak
-    peakpos = minPos;
-    peakLevel = data[minPos];
-    for (i = minPos + 1; i < maxPos; i ++)
-    {
-        if (data[i] > peakLevel) 
-        {
-            peakLevel = data[i];
-            peakpos = i;
-        }
-    }
-    
     // find ground positions.
     gp1 = findGround(data, peakpos, -1);
     gp2 = findGround(data, peakpos, 1);
 
     groundLevel = max(data[gp1], data[gp2]);
+    peakLevel = data[peakpos];
 
     if (groundLevel < 1e-6) return 0;                // ground level too small => detection failed
     if ((peakLevel / groundLevel) < 1.3) return 0;   // peak less than 30% of the ground level => no good peak detected
@@ -189,3 +175,60 @@ float PeakFinder::detectPeak(const float *data, int minPos, int maxPos)
 }
 
 
+
+double PeakFinder::detectPeak(const float *data, int minPos, int maxPos) 
+{
+
+    int i;
+    int peakpos;                // position of peak level
+    double highPeak, peak;
+
+    this->minPos = minPos;
+    this->maxPos = maxPos;
+
+    // find absolute peak
+    peakpos = minPos;
+    peak = data[minPos];
+    for (i = minPos + 1; i < maxPos; i ++)
+    {
+        if (data[i] > peak) 
+        {
+            peak = data[i];
+            peakpos = i;
+        }
+    }
+    
+    // Calculate exact location of the highest peak mass center
+    highPeak = getPeakCenter(data, peakpos);
+    peak = highPeak;
+
+    // Now check if the highest peak were in fact harmonic of the true base beat peak 
+    // - sometimes the highest peak can be Nth harmonic of the true base peak yet 
+    // just a slightly higher than the true base
+    for (i = 2; i < 10; i ++)
+    {
+        double peaktmp, tmp;
+        int i1,i2;
+
+        peakpos = (int)(highPeak / (double)i + 0.5f);
+        if (peakpos < minPos) break;
+
+        // calculate mass-center of possible base peak
+        peaktmp = getPeakCenter(data, peakpos);
+
+        // now compare to highest detected peak
+        i1 = (int)(highPeak + 0.5);
+        i2 = (int)(peaktmp + 0.5);
+        tmp = 2 * (data[i2] - data[i1]) / (data[i2] + data[i1]);
+        if (fabs(tmp) < 0.1)
+        {
+            // The highest peak is harmonic of almost as high base peak,
+            // thus use the base peak instead
+            peak = peaktmp;
+        }
+    }
+
+    return peak;
+}
+
+

source/example/bpm/PeakFinder.h

@@ -47,7 +47,7 @@ protected:
     int minPos, maxPos;
 
     /// Calculates the mass center between given vector items.
-    float calcMassCenter(const float *data, ///< Data vector.
+    double calcMassCenter(const float *data, ///< Data vector.
                          int firstPos,      ///< Index of first vector item beloging to the peak.
                          int lastPos        ///< Index of last vector item beloging to the peak.
                          ) const;
@@ -67,15 +67,18 @@ protected:
                      int direction          /// Direction where to proceed from the peak: 1 = right, -1 = left.
                      ) const;
 
+    /// get exact center of peak near given position by calculating local mass of center
+    double getPeakCenter(const float *data, int peakpos);
+
 public:
     /// Constructor. 
     PeakFinder();
 
     /// Detect exact peak position of the data vector by finding the largest peak 'hump'
-    /// and calculating the mass-center location of the peak hump. 
+    /// and calculating the mass-center location of the peak hump.
     ///
-    /// \return The exact mass-center location of the largest peak hump.
+    /// \return The location of the largest base harmonic peak hump.
-    float detectPeak(const float *data, /// Data vector to be analyzed. The data vector has
+    double detectPeak(const float *data, /// Data vector to be analyzed. The data vector has
                                         /// to be at least 'maxPos' items long.
                      int minPos,        ///< Min allowed peak location within the vector data.
                      int maxPos         ///< Max allowed peak location within the vector data.