Chris@26: /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
Chris@26: 
Chris@26: /*
Chris@26:     Vamp feature extraction plugin using the MATCH audio alignment
Chris@26:     algorithm.
Chris@26: 
Chris@26:     Centre for Digital Music, Queen Mary, University of London.
Chris@236:     Copyright (c) 2007-2020 Simon Dixon, Chris Cannam, and Queen Mary
Chris@230:     University of London, Copyright (c) 2014-2015 Tido GmbH.
Chris@26:     
Chris@26:     This program is free software; you can redistribute it and/or
Chris@26:     modify it under the terms of the GNU General Public License as
Chris@26:     published by the Free Software Foundation; either version 2 of the
Chris@26:     License, or (at your option) any later version.  See the file
Chris@26:     COPYING included with this distribution for more information.
Chris@26: */
Chris@26: 
Chris@26: #include "DistanceMetric.h"
Chris@26: 
Chris@26: #include <cassert>
Chris@26: #include <cmath>
Chris@133: #include <iostream>
Chris@26: 
Chris@133: using namespace std;
Chris@26: 
Chris@196: //#define DEBUG_DISTANCE_METRIC 1
Chris@140: 
Chris@185: template <> uint8_t
Chris@185: DistanceMetric::scaleIntoRange(double distance)
Chris@185: {
Chris@186:     double scaled = m_params.scale * distance;
Chris@190:     if (scaled < 0) {
Chris@190:         scaled = 0;
Chris@190:     }
Chris@212:     if (scaled > DISTANCE_MAX) {
Chris@212:         scaled = DISTANCE_MAX;
Chris@190:         ++m_overcount;
Chris@190:     }
Chris@186:     return uint8_t(scaled);
Chris@185: }
Chris@185: 
Chris@185: template <> float
Chris@185: DistanceMetric::scaleIntoRange(double distance)
Chris@185: {
Chris@185:     return float(distance);
Chris@185: }
Chris@185: 
Chris@185: template <> double
Chris@185: DistanceMetric::scaleIntoRange(double distance)
Chris@185: {
Chris@185:     return distance;
Chris@185: }
Chris@185: 
Chris@143: DistanceMetric::DistanceMetric(Parameters params) :
Chris@190:     m_params(params),
Chris@190:     m_max(0),
Chris@190:     m_overcount(0)
Chris@140: {
Chris@140: #ifdef DEBUG_DISTANCE_METRIC
Chris@190:     cerr << "*** DistanceMetric: metric = " << m_params.metric
Chris@190:          << ", norm = " << m_params.norm
Chris@190:          << ", noise = " << m_params.noise
Chris@190:          << ", scale = " << m_params.scale
Chris@143:          << endl;
Chris@140: #endif
Chris@140: }
Chris@140: 
Chris@190: DistanceMetric::~DistanceMetric()
Chris@190: {
Chris@190: #ifdef DEBUG_DISTANCE_METRIC
Chris@190:     cerr << "*** DistanceMetric::~DistanceMetric: metric = " << m_params.metric
Chris@190:          << ", norm = " << m_params.norm
Chris@190:          << ", noise = " << m_params.noise;
Chris@190: #ifdef USE_COMPACT_TYPES
Chris@190:     cerr << ", scale = " << m_params.scale;
Chris@190:     cerr << "\n*** DistanceMetric::~DistanceMetric: max scaled value = "
Chris@190:          << distance_print_t(m_max)
Chris@190:          << ", " << m_overcount << " clipped" << endl;
Chris@190: #else
Chris@190:     cerr << ", no scaling";
Chris@190:     cerr << "\n*** DistanceMetric::~DistanceMetric: max value = "
Chris@190:          << distance_print_t(m_max)
Chris@190:          << endl;
Chris@190: #endif
Chris@190: #endif
Chris@190: }
Chris@190: 
Chris@183: distance_t
Chris@186: DistanceMetric::scaleValueIntoDistanceRange(double value)
Chris@186: {
Chris@186:     return scaleIntoRange<distance_t>(value);
Chris@186: }
Chris@190: 
Chris@190: distance_t
Chris@190: DistanceMetric::scaleAndTally(double value)
Chris@190: {
Chris@190:     distance_t dist = scaleIntoRange<distance_t>(value);
Chris@190:     if (dist > m_max) m_max = dist;
Chris@190:     return dist;
Chris@190: }
Chris@190: 
Chris@186: distance_t
Chris@183: DistanceMetric::calcDistance(const feature_t &f1,
Chris@183: 			     const feature_t &f2)
Chris@26: {
Chris@26:     double d = 0;
Chris@26:     double sum = 0;
Chris@156:     double eps = 1e-16;
Chris@26: 
Chris@180:     assert(f2.size() == f1.size());
Chris@180:     int featureSize = static_cast<int>(f1.size());
Chris@145: 
Chris@197:     double minNoise = 0.0;
Chris@197: #ifdef USE_COMPACT_TYPES
Chris@197:     minNoise = 1.0 / m_params.scale;
Chris@197: #endif
Chris@197:     
Chris@156:     if (m_params.metric == Cosine) {
Chris@156: 
Chris@156:         double num = 0, denom1 = 0, denom2 = 0;
Chris@156:         
Chris@156:         for (int i = 0; i < featureSize; ++i) {
Chris@156:             num += f1[i] * f2[i];
Chris@156:             denom1 += f1[i] * f1[i];
Chris@156:             denom2 += f2[i] * f2[i];
Chris@156:         }
Chris@156: 
Chris@156:         d = 1.0 - (num / (eps + sqrt(denom1 * denom2)));
Chris@156: 
Chris@156:         if (m_params.noise == AddNoise) {
Chris@197:             double noise = 1e-2;
Chris@197:             if (noise < minNoise) noise = minNoise;
Chris@197:             d += noise;
Chris@156:         }
Chris@156:         if (d > 1.0) d = 1.0;
Chris@156:         
Chris@190:         return scaleAndTally(d); // normalisation param ignored
Chris@156:     }
Chris@156: 
Chris@157:     if (m_params.metric == Manhattan) {
Chris@157:         for (int i = 0; i < featureSize; i++) {
Chris@157:             d += fabs(f1[i] - f2[i]);
Chris@157:             sum += fabs(f1[i]) + fabs(f2[i]);
Chris@157:         }
Chris@157:     } else {
Chris@157:         // Euclidean
Chris@157:         for (int i = 0; i < featureSize; i++) {
Chris@157:             d += (f1[i] - f2[i]) * (f1[i] - f2[i]);
Chris@157:             sum += fabs(f1[i]) + fabs(f2[i]);
Chris@157:         }
Chris@157:         d = sqrt(d);
Chris@26:     }
Chris@26: 
Chris@150:     if (m_params.noise == AddNoise) {
Chris@197:         double noise = 1e-3 * featureSize;
Chris@197:         if (noise < minNoise) noise = minNoise;
Chris@150:         d += noise;
Chris@150:         sum += noise;
Chris@150:     }
Chris@145:     
Chris@143:     if (sum == 0) {
Chris@190:         return scaleAndTally(0);
Chris@143:     }
Chris@26: 
Chris@143:     double distance = 0;
Chris@26: 
Chris@143:     if (m_params.norm == NormaliseDistanceToSum) {
Chris@143: 
Chris@143:         distance = d / sum; // 0 <= d/sum <= 2
Chris@143: 
Chris@143:     } else if (m_params.norm == NormaliseDistanceToLogSum) {
Chris@143: 
Chris@143:         // note if this were to be restored, it would have to use
Chris@143:         // totalEnergies vector instead of f1[freqMapSize] which used to
Chris@143:         // store the total energy:
Chris@143:         //	double weight = (5 + Math.log(f1[freqMapSize] + f2[freqMapSize]))/10.0;
Chris@143: 
Chris@143:         double weight = (8 + log(sum)) / 10.0;
Chris@133:     
Chris@143:         if (weight < 0) weight = 0;
Chris@143:         else if (weight > 1) weight = 1;
Chris@26: 
Chris@143:         distance = d / sum * weight;
Chris@143: 
Chris@143:     } else {
Chris@143: 
Chris@143:         distance = d;
Chris@143:     }
Chris@143:     
Chris@190:     return scaleAndTally(distance);
Chris@26: }