changeset 117:aef4b62069ba adaptive_diagonals

Merge from refactors branch
author Chris Cannam
date Fri, 05 Dec 2014 11:02:15 +0000
parents 1a422fa20b51 (diff) eed5f9594268 (current diff)
children 69b16c0bc903
files src/MatchVampPlugin.cpp
diffstat 6 files changed, 380 insertions(+), 59 deletions(-) [+]
line wrap: on
line diff
--- a/Makefile.linux	Fri Dec 05 10:05:31 2014 +0000
+++ b/Makefile.linux	Fri Dec 05 11:02:15 2014 +0000
@@ -1,6 +1,6 @@
 
-#CXXFLAGS	+= -fPIC -ffast-math -O3 -Wall -Werror -Wfloat-conversion
-CXXFLAGS	+= -fPIC -g -Wall -Werror -DPERFORM_ERROR_CHECKS=1
+CXXFLAGS	+= -fPIC -ffast-math -O3 -Wall -Werror
+#CXXFLAGS	+= -fPIC -g -Wall -Werror -Wfloat-conversion -DPERFORM_ERROR_CHECKS=1
 
 LDFLAGS		+= -shared -Wl,-Bstatic -lvamp-sdk -Wl,-Bdynamic -Wl,-Bsymbolic -Wl,-z,defs -lpthread -Wl,--version-script=vamp-plugin.map
 
--- a/src/Finder.cpp	Fri Dec 05 10:05:31 2014 +0000
+++ b/src/Finder.cpp	Fri Dec 05 11:02:15 2014 +0000
@@ -88,11 +88,16 @@
 }
 
 void
-Finder::recalculatePathCostMatrix(int r1, int c1, int r2, int c2) 
+Finder::recalculatePathCostMatrix(int r1, int c1, int r2, int c2)
+{
+    recalculatePathCostMatrix(r1, c1, r2, c2, m_m->getDiagonalWeight());
+}
+
+void
+Finder::recalculatePathCostMatrix(int r1, int c1, int r2, int c2,
+                                  float diagonalWeight) 
 {
     int prevRowStart = 0, prevRowStop = 0;
-
-    float diagonalWeight = m_m->getDiagonalWeight();
     
     for (int r = r1; r <= r2; r++) {
 
@@ -327,11 +332,14 @@
 }
 #endif
 
-int
-Finder::retrievePath(bool smooth, vector<int> &pathx, vector<int> &pathy)
+void
+Finder::retrievePath(vector<int> &pathx,
+                     vector<int> &pathy,
+                     vector<float> &distances)
 {
     pathx.clear();
     pathy.clear();
+    distances.clear();
 
 #ifdef PERFORM_ERROR_CHECKS
     checkAndReport();
@@ -341,7 +349,7 @@
     int ey = m_m->getFrameCount() - 1;
 
     if (ex < 0 || ey < 0) {
-        return 0;
+        return;
     }
     
     int x = ex;
@@ -366,7 +374,7 @@
         y = ey;
     }
 
-    recalculatePathCostMatrix(0, 0, y, x);
+//    recalculatePathCostMatrix(0, 0, y, x);
 
 //    cerr << "start: x = " << x << ", y = " << y << endl;
     
@@ -376,6 +384,7 @@
         
         pathx.push_back(x);
         pathy.push_back(y);
+        distances.push_back(m_m->getDistance(y, x));
 
         switch (m_m->getAdvance(y, x)) {
         case Matcher::AdvanceThis:
@@ -409,12 +418,35 @@
     
     reverse(pathx.begin(), pathx.end());
     reverse(pathy.begin(), pathy.end());
+    reverse(distances.begin(), distances.end());
+}
 
-    if (smooth) {
-        int smoothedLen = Path().smooth(pathx, pathy, pathx.size());
-        return smoothedLen;
-    } else {
-        return pathx.size();
+void
+Finder::smoothWithPinPoints(const map<int, int> &pinpoints)
+{
+    cerr << "Pin points are:" << endl;
+
+    typedef map<int, int> PPMap;
+
+    for (PPMap::const_iterator i = pinpoints.begin(); i != pinpoints.end(); ++i) {
+        cerr << "[" << i->first << "," << i->second << "] ";
+    }
+    
+    cerr << endl;
+
+    if (pinpoints.size() < 2) return;
+
+    pair<int, int> prev = *pinpoints.begin();
+
+    for (PPMap::const_iterator i = pinpoints.begin(); i != pinpoints.end(); ++i) {
+        if (i == pinpoints.begin()) continue;
+
+        pair<int, int> curr = *i;
+
+        recalculatePathCostMatrix(prev.second, prev.first,
+                                  curr.second, curr.first, 1.0);
+        
+        prev = curr;
     }
 }
 
--- a/src/Finder.h	Fri Dec 05 10:05:31 2014 +0000
+++ b/src/Finder.h	Fri Dec 05 11:02:15 2014 +0000
@@ -18,6 +18,7 @@
 #define _FINDER_H_
 
 #include <vector>
+#include <map>
 #include <iostream>
 
 #include "Matcher.h"
@@ -57,14 +58,16 @@
      * Track back after all of the matchers have been fed in order to
      * obtain the lowest cost path available. Path x and y coordinate
      * pairs are returned in corresponding elements of pathx and
-     * pathy. Return value is the length of the returned path: only
-     * this many elements from pathx and pathy are valid (any
-     * subsequent ones may be spurious).
-     *
-     * @param smooth whether to smooth the path before returning it
+     * pathy, and the corresponding elements of the returned distance
+     * vector contain the distance between the pair of features
+     * selected at each step in the path.
      */
-    int retrievePath(bool smooth, std::vector<int> &pathx, std::vector<int> &pathy);
+    void retrievePath(std::vector<int> &pathx,
+                      std::vector<int> &pathy,
+                      std::vector<float> &distances);
 
+    void smoothWithPinPoints(const std::map<int, int> &);
+    
 protected:
 #ifdef PERFORM_ERROR_CHECKS
     struct ErrorPosition {
@@ -83,6 +86,9 @@
     ErrorPosition checkPathCostMatrix();
     void checkAndReport();
 #endif
+
+    void recalculatePathCostMatrix(int r1, int c1, int r2, int c2,
+                                   float diagonalWeight);
     
     Matcher *m_m;
     int m_duration1;
--- a/src/MatchVampPlugin.cpp	Fri Dec 05 10:05:31 2014 +0000
+++ b/src/MatchVampPlugin.cpp	Fri Dec 05 11:02:15 2014 +0000
@@ -27,6 +27,9 @@
 
 #include <vector>
 #include <algorithm>
+#include <map>
+
+using namespace std;
 
 //static int extant = 0;
 
@@ -66,10 +69,10 @@
     m_defaultFcParams()
 {
     if (inputSampleRate < sampleRateMin) {
-        std::cerr << "MatchVampPlugin::MatchVampPlugin: input sample rate "
+        cerr << "MatchVampPlugin::MatchVampPlugin: input sample rate "
                   << inputSampleRate << " < min supported rate "
                   << sampleRateMin << ", plugin will refuse to initialise"
-                  << std::endl;
+                  << endl;
     }
 
     if (!m_serialisingMutexInitialised) {
@@ -82,12 +85,11 @@
     }
 
     m_pipeline = 0;
-//    std::cerr << "MatchVampPlugin::MatchVampPlugin(" << this << "): extant = " << ++extant << std::endl;
 }
 
 MatchVampPlugin::~MatchVampPlugin()
 {
-//    std::cerr << "MatchVampPlugin::~MatchVampPlugin(" << this << "): extant = " << --extant << std::endl;
+//    cerr << "MatchVampPlugin::~MatchVampPlugin(" << this << "): extant = " << --extant << endl;
 
     delete m_pipeline;
 
@@ -250,7 +252,7 @@
 }
 
 float
-MatchVampPlugin::getParameter(std::string name) const
+MatchVampPlugin::getParameter(string name) const
 {
     if (name == "serialise") {
         return m_serialise ? 1.0 : 0.0; 
@@ -276,7 +278,7 @@
 }
 
 void
-MatchVampPlugin::setParameter(std::string name, float value)
+MatchVampPlugin::setParameter(string name, float value)
 {
     if (name == "serialise") {
         m_serialise = (value > 0.5);
@@ -324,9 +326,9 @@
 MatchVampPlugin::initialise(size_t channels, size_t stepSize, size_t blockSize)
 {
     if (m_inputSampleRate < sampleRateMin) {
-        std::cerr << "MatchVampPlugin::MatchVampPlugin: input sample rate "
+        cerr << "MatchVampPlugin::MatchVampPlugin: input sample rate "
                   << m_inputSampleRate << " < min supported rate "
-                  << sampleRateMin << std::endl;
+                  << sampleRateMin << endl;
         return false;
     }
     if (channels < getMinChannelCount() ||
@@ -351,6 +353,10 @@
     delete m_pipeline;
     m_pipeline = 0;
     m_frameNo = 0;
+
+    m_mag1.clear();
+    m_mag2.clear();
+    
     createMatchers();
     m_begin = true;
     m_locked = false;
@@ -458,9 +464,72 @@
     m_bFeaturesOutNo = list.size();
     list.push_back(desc);
 
+    desc.identifier = "feature_distance";
+    desc.name = "Feature Distance";
+    desc.description = "Value of distance metric between features at each point-in-A along the chosen alignment path";
+    desc.unit = "";
+    desc.hasFixedBinCount = true;
+    desc.binCount = 1;
+    desc.hasKnownExtents = false;
+    desc.isQuantized = false;
+    desc.sampleType = OutputDescriptor::FixedSampleRate;
+    desc.sampleRate = outRate;
+    m_distOutNo = list.size();
+    list.push_back(desc);
+
+    desc.identifier = "feature_mag";
+    desc.name = "Feature Magnitudes";
+    desc.description = "Sum of magnitudes of feature pairs for each point-in-A along the chosen alignment path";
+    desc.unit = "";
+    desc.hasFixedBinCount = true;
+    desc.binCount = 1;
+    desc.hasKnownExtents = false;
+    desc.isQuantized = false;
+    desc.sampleType = OutputDescriptor::FixedSampleRate;
+    desc.sampleRate = outRate;
+    m_magOutNo = list.size();
+    list.push_back(desc);
+
+    desc.identifier = "confidence";
+    desc.name = "Confidence";
+    desc.description = "Confidence metric for the quality of match at each point-in-A along the chosen alignment path";
+    desc.unit = "";
+    desc.hasFixedBinCount = true;
+    desc.binCount = 1;
+    desc.hasKnownExtents = false;
+    desc.isQuantized = false;
+    desc.sampleType = OutputDescriptor::FixedSampleRate;
+    desc.sampleRate = outRate;
+    m_confidenceOutNo = list.size();
+    list.push_back(desc);
+
+    desc.identifier = "confidence_peaks";
+    desc.name = "Confidence Peaks";
+    desc.description = "Peak locations for the confidence metric";
+    desc.unit = "";
+    desc.hasFixedBinCount = true;
+    desc.binCount = 0;
+    desc.hasKnownExtents = false;
+    desc.isQuantized = false;
+    desc.sampleType = OutputDescriptor::FixedSampleRate;
+    desc.sampleRate = outRate;
+    m_confPeakOutNo = list.size();
+    list.push_back(desc);
+
     return list;
 }
 
+static double
+magOf(const vector<double> &f)
+{
+    double mag = 0.0;
+    for (int j = 0; j < (int)f.size(); ++j) {
+        mag += f[j] * f[j];
+    }
+    mag = sqrt(mag);
+    return mag;
+}
+
 MatchVampPlugin::FeatureSet
 MatchVampPlugin::process(const float *const *inputBuffers,
                          Vamp::RealTime timestamp)
@@ -478,12 +547,19 @@
         m_begin = false;
     }
     
-//    std::cerr << timestamp.toString();
+//    cerr << timestamp.toString();
 
     m_pipeline->feedFrequencyDomainAudio(inputBuffers[0], inputBuffers[1]);
 
-    vector<double> f1, f2;
-    m_pipeline->extractConditionedFeatures(f1, f2);
+    vector<double> f1, f2, c1, c2;
+
+    m_pipeline->extractFeatures(f1, f2);
+    m_pipeline->extractConditionedFeatures(c1, c2);    
+
+    m_mag1.push_back(magOf(f1));
+    m_mag2.push_back(magOf(f2));
+    m_cmag1.push_back(magOf(c1));
+    m_cmag2.push_back(magOf(c2));
 
     FeatureSet returnFeatures;
 
@@ -491,19 +567,19 @@
     f.hasTimestamp = false;
 
     f.values.clear();
-    for (int j = 0; j < (int)f1.size(); ++j) {
-        f.values.push_back(float(f1[j]));
+    for (int j = 0; j < (int)c1.size(); ++j) {
+        f.values.push_back(float(c1[j]));
     }
     returnFeatures[m_aFeaturesOutNo].push_back(f);
 
     f.values.clear();
-    for (int j = 0; j < (int)f2.size(); ++j) {
-        f.values.push_back(float(f2[j]));
+    for (int j = 0; j < (int)c2.size(); ++j) {
+        f.values.push_back(float(c2[j]));
     }
     returnFeatures[m_bFeaturesOutNo].push_back(f);
 
-//    std::cerr << ".";
-//    std::cerr << std::endl;
+//    cerr << ".";
+//    cerr << endl;
 
     ++m_frameNo;
     
@@ -518,18 +594,86 @@
     FeatureSet returnFeatures;
     
     Finder *finder = m_pipeline->getFinder();
-    std::vector<int> pathx;
-    std::vector<int> pathy;
-    int len = finder->retrievePath(m_smooth, pathx, pathy);
+    vector<int> pathx;
+    vector<int> pathy;
+    vector<float> distances;
+    finder->retrievePath(pathx, pathy, distances);
 
     int prevx = 0;
     int prevy = 0;
+    int len = pathx.size();
+
+//!!!
+//  m_smooth = true;
+    
+    if (m_smooth) {
+    
+        vector<float> confidence;
+
+        for (int i = 0; i < len; ++i) {
+            int x = pathx[i];
+            int y = pathy[i];
+
+            double magSum = (m_cmag1[y] * m_mag1[y]) + (m_cmag2[x] * m_mag2[x]);
+            double distance = distances[i];
+            float c = magSum - distance;
+            confidence.push_back(c);
+
+            if (x != prevx) {
+                Feature f;
+                f.values.push_back(c);
+                returnFeatures[m_confidenceOutNo].push_back(f);
+
+                f.values.clear();
+                f.values.push_back(magSum);
+                returnFeatures[m_magOutNo].push_back(f);
+
+                f.values.clear();
+                f.values.push_back(distance);
+                returnFeatures[m_distOutNo].push_back(f);
+            }
+        }
+        
+        map<int, int> pinpoints;
+            
+        vector<float> csorted = confidence;
+        sort(csorted.begin(), csorted.end());
+        float thresh = csorted[int(csorted.size() * 0.7)]; // 70th percentile
+
+        for (int i = 1; i + 1 < len; ++i) {
+
+            int x = pathx[i];
+            int y = pathy[i];
+
+            if (confidence[i] > thresh &&
+                confidence[i] > confidence[i-1] &&
+                confidence[i] >= confidence[i+1]) {
+
+                pinpoints[x] = y;
+                
+                Vamp::RealTime xt = Vamp::RealTime::frame2RealTime
+                    (x * m_stepSize, lrintf(m_inputSampleRate));
+                Feature feature;
+                feature.hasTimestamp = true;
+                feature.timestamp = m_startTime + xt;
+                returnFeatures[m_confPeakOutNo].push_back(feature);
+            }
+        }
+
+        finder->smoothWithPinPoints(pinpoints);
+
+        pathx.clear();
+        pathy.clear();
+        distances.clear();
+        finder->retrievePath(pathx, pathy, distances);
+        len = pathx.size();
+    }    
 
     for (int i = 0; i < len; ++i) {
 
         int x = pathx[i];
         int y = pathy[i];
-
+            
         Vamp::RealTime xt = Vamp::RealTime::frame2RealTime
             (x * m_stepSize, lrintf(m_inputSampleRate));
         Vamp::RealTime yt = Vamp::RealTime::frame2RealTime
@@ -600,31 +744,31 @@
 
 /*
     for (int i = 0; i < len; ++i) {
-        std::cerr << i << ": [" << pathx[i] << "," << pathy[i] << "]" << std::endl;
+        cerr << i << ": [" << pathx[i] << "," << pathy[i] << "]" << endl;
     }
 
-    std::cerr << std::endl;
-    std::cerr << "File: A" << std::endl;
-    std::cerr << "Marks: -1" << std::endl;
-    std::cerr << "FixedPoints: true 0" << std::endl;
-    std::cerr << "0" << std::endl;
-    std::cerr << "0" << std::endl;
-    std::cerr << "0" << std::endl;
-    std::cerr << "0" << std::endl;
-    std::cerr << "File: B" << std::endl;
-    std::cerr << "Marks: 0" << std::endl;
-    std::cerr << "FixedPoints: true 0" << std::endl;
-    std::cerr << "0.02" << std::endl;
-    std::cerr << "0.02" << std::endl;
+    cerr << endl;
+    cerr << "File: A" << endl;
+    cerr << "Marks: -1" << endl;
+    cerr << "FixedPoints: true 0" << endl;
+    cerr << "0" << endl;
+    cerr << "0" << endl;
+    cerr << "0" << endl;
+    cerr << "0" << endl;
+    cerr << "File: B" << endl;
+    cerr << "Marks: 0" << endl;
+    cerr << "FixedPoints: true 0" << endl;
+    cerr << "0.02" << endl;
+    cerr << "0.02" << endl;
 
-    std::cerr << len << std::endl;
+    cerr << len << endl;
     for (int i = 0; i < len; ++i) {
-        std::cerr << pathx[i] << std::endl;
+        cerr << pathx[i] << endl;
     }
 
-    std::cerr << len << std::endl;
+    cerr << len << endl;
     for (int i = 0; i < len; ++i) {
-        std::cerr << pathy[i] << std::endl;
+        cerr << pathy[i] << endl;
     }
 */
 }
--- a/src/MatchVampPlugin.h	Fri Dec 05 10:05:31 2014 +0000
+++ b/src/MatchVampPlugin.h	Fri Dec 05 11:02:15 2014 +0000
@@ -90,6 +90,11 @@
     FeatureConditioner::Parameters m_fcParams;
     FeatureConditioner::Parameters m_defaultFcParams;
 
+    std::vector<float> m_mag1;
+    std::vector<float> m_mag2;
+    std::vector<float> m_cmag1;
+    std::vector<float> m_cmag2;
+    
     mutable int m_pathOutNo;
     mutable int m_abOutNo;
     mutable int m_baOutNo;
@@ -97,6 +102,10 @@
     mutable int m_abRatioOutNo;
     mutable int m_aFeaturesOutNo;
     mutable int m_bFeaturesOutNo;
+    mutable int m_magOutNo;
+    mutable int m_distOutNo;
+    mutable int m_confidenceOutNo;
+    mutable int m_confPeakOutNo;
 
 #ifdef _WIN32
     static HANDLE m_serialisingMutex;
--- /dev/null	Thu Jan 01 00:00:00 1970 +0000
+++ b/src/MedianFilter.h	Fri Dec 05 11:02:15 2014 +0000
@@ -0,0 +1,130 @@
+/* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
+
+/*
+    QM DSP Library
+
+    Centre for Digital Music, Queen Mary, University of London.
+    This file Copyright 2010 Chris Cannam.
+
+    This program is free software; you can redistribute it and/or
+    modify it under the terms of the GNU General Public License as
+    published by the Free Software Foundation; either version 2 of the
+    License, or (at your option) any later version.  See the file
+    COPYING included with this distribution for more information.
+*/
+
+#ifndef MEDIAN_FILTER_H
+#define MEDIAN_FILTER_H
+
+#include <algorithm>
+#include <cassert>
+#include <cmath>
+#include <iostream>
+#include <vector>
+
+template <typename T>
+class MedianFilter
+{
+public:
+    MedianFilter(int size, float percentile = 50.f) :
+	m_size(size),
+	m_frame(new T[size]),
+	m_sorted(new T[size]),
+	m_sortend(m_sorted + size - 1) {
+        setPercentile(percentile);
+	reset();
+    }
+
+    ~MedianFilter() { 
+	delete[] m_frame;
+	delete[] m_sorted;
+    }
+
+    void setPercentile(float p) {
+        m_index = int((m_size * p) / 100.f);
+        if (m_index >= m_size) m_index = m_size-1;
+        if (m_index < 0) m_index = 0;
+    }
+
+    void push(T value) {
+        if (value != value) {
+            std::cerr << "WARNING: MedianFilter::push: attempt to push NaN" << std::endl;
+            return; // nan
+        }
+	drop(m_frame[0]);
+	const int sz1 = m_size-1;
+	for (int i = 0; i < sz1; ++i) m_frame[i] = m_frame[i+1];
+	m_frame[m_size-1] = value;
+	put(value);
+    }
+
+    T get() const {
+	return m_sorted[m_index];
+    }
+
+    int getSize() const {
+        return m_size; 
+    }
+
+    T getAt(float percentile) {
+	int ix = int((m_size * percentile) / 100.f);
+        if (ix >= m_size) ix = m_size-1;
+        if (ix < 0) ix = 0;
+	return m_sorted[ix];
+    }
+
+    void reset() {
+	for (int i = 0; i < m_size; ++i) m_frame[i] = 0;
+	for (int i = 0; i < m_size; ++i) m_sorted[i] = 0;
+    }
+
+    static std::vector<T> filter(int size, const std::vector<T> &in) {
+        std::vector<T> out;
+        MedianFilter<T> f(size);
+        for (int i = 0; i < int(in.size()); ++i) {
+            f.push(in[i]);
+            T median = f.get();
+            if (i >= size/2) out.push_back(median);
+        }
+        while (out.size() < in.size()) {
+            f.push(T());
+            out.push_back(f.get());
+        }
+        return out;
+    }
+
+private:
+    const int m_size;
+    T *const m_frame;
+    T *const m_sorted;
+    T *const m_sortend;
+    int m_index;
+
+    void put(T value) {
+	// precondition: m_sorted contains m_size-1 values, packed at start
+	// postcondition: m_sorted contains m_size values, one of which is value
+	T *point = std::lower_bound(m_sorted, m_sortend, value);
+	const int n = m_sortend - point;
+	for (int i = n; i > 0; --i) point[i] = point[i-1];
+	*point = value;
+    }
+
+    void drop(T value) {
+	// precondition: m_sorted contains m_size values, one of which is value
+	// postcondition: m_sorted contains m_size-1 values, packed at start
+	T *point = std::lower_bound(m_sorted, m_sortend + 1, value);
+	if (*point != value) {
+	    std::cerr << "WARNING: MedianFilter::drop: *point is " << *point
+		      << ", expected " << value << std::endl;
+	}
+	const int n = m_sortend - point;
+	for (int i = 0; i < n; ++i) point[i] = point[i+1];
+	*m_sortend = T(0);
+    }
+
+    MedianFilter(const MedianFilter &); // not provided
+    MedianFilter &operator=(const MedianFilter &); // not provided
+};
+
+#endif
+