cannam@0: /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
cannam@0: 
cannam@0: /*
cannam@0:     Vamp feature extraction plugin using the MATCH audio alignment
cannam@0:     algorithm.
cannam@0: 
cannam@0:     Centre for Digital Music, Queen Mary, University of London.
cannam@0:     This file copyright 2007 Simon Dixon, Chris Cannam and QMUL.
cannam@0:     
cannam@0:     This program is free software; you can redistribute it and/or
cannam@0:     modify it under the terms of the GNU General Public License as
cannam@0:     published by the Free Software Foundation; either version 2 of the
cannam@0:     License, or (at your option) any later version.  See the file
cannam@0:     COPYING included with this distribution for more information.
cannam@0: */
cannam@0: 
cannam@0: #include "Finder.h"
cannam@0: 
Chris@30: #include "Path.h"
Chris@30: 
Chris@30: #include <algorithm>
Chris@92: #include <iomanip>
Chris@30: 
Chris@72: using namespace std;
cannam@0: 
Chris@140: //#define DEBUG_FINDER 1
Chris@140: //#define PERFORM_ERROR_CHECKS 1
Chris@140: 
Chris@72: Finder::Finder(Matcher *pm)
cannam@0: {
Chris@72:     m_m = pm;
Chris@72:     m_duration1 = -1;
Chris@72:     m_duration2 = -1;
cannam@0: } // constructor
cannam@0: 
cannam@0: Finder::~Finder()
cannam@0: {
cannam@0: }
cannam@0: 
Chris@60: void
Chris@60: Finder::setDurations(int d1, int d2)
Chris@60: {
Chris@140: #ifdef DEBUG_FINDER
Chris@140:     cerr << "*** setDurations: " << d1 << ", " << d2 << endl;
Chris@140: #endif
Chris@72:     m_duration1 = d1;
Chris@72:     m_duration2 = d2;
Chris@60: }
Chris@60: 
Chris@45: Matcher::Advance
cannam@0: Finder::getExpandDirection(int row, int col)
cannam@0: {
Chris@72:     double min = m_m->getPathCost(row, col);
Chris@72:     
Chris@72:     int bestRow = row;
Chris@72:     int bestCol = col;
Chris@72: 
Chris@72:     pair<int, int> rowRange = m_m->getRowRange(col);
Chris@72:     if (rowRange.second > row+1) {
Chris@72:         rowRange.second = row+1;	// don't cheat by looking at future :)
Chris@72:     }
Chris@72:     for (int index = rowRange.first; index < rowRange.second; index++) {
Chris@135:         double tmp = m_m->getNormalisedPathCost(index, col);
cannam@0:         if (tmp < min) {
cannam@0:             min = tmp;
cannam@0:             bestRow = index;
cannam@0:         }
cannam@0:     }
Chris@72: 
Chris@72:     pair<int, int> colRange = m_m->getColRange(row);
Chris@72:     if (colRange.second > col+1) {
Chris@72:         colRange.second = col+1;	// don't cheat by looking at future :)
Chris@72:     }
Chris@72:     for (int index = colRange.first; index < colRange.second; index++) {
Chris@135:         double tmp = m_m->getNormalisedPathCost(row, index);
cannam@0:         if (tmp < min) {
cannam@0:             min = tmp;
cannam@0:             bestCol = index;
cannam@0:             bestRow = row;
cannam@0:         }
cannam@0:     }
Chris@72: 
Chris@135: //    cerr << "at [" << row << "," << col << "] (cost " << m_m->getPathCost(row, col) << ") blocksize = " << m_m->getBlockSize() << " best is [" << bestRow << "," << bestCol << "] (cost " << min << ")" << endl;
Chris@135:     
Chris@45:     if (bestRow == row) {
Chris@45:         if (bestCol == col) {
Chris@45:             return Matcher::AdvanceBoth;
Chris@45:         } else {
Chris@45:             return Matcher::AdvanceThis;
Chris@45:         }
Chris@45:     } else if (bestCol == col) {
Chris@45:         return Matcher::AdvanceOther;
Chris@45:     } else {
Chris@46:         return Matcher::AdvanceNone;
Chris@45:     }
Chris@73: }
cannam@0: 
cannam@0: void
cannam@0: Finder::recalculatePathCostMatrix(int r1, int c1, int r2, int c2) 
cannam@0: {
Chris@72:     int prevRowStart = 0, prevRowStop = 0;
Chris@72: 
Chris@83:     float diagonalWeight = m_m->getDiagonalWeight();
Chris@83:     
Chris@72:     for (int r = r1; r <= r2; r++) {
Chris@72: 
Chris@72:         pair<int, int> colRange = m_m->getColRange(r);
Chris@72: 
Chris@72:         int rowStart = max(c1, colRange.first);
Chris@72:         int rowStop = min(c2 + 1, colRange.second);
Chris@72:         
Chris@72:         for (int c = rowStart; c < rowStop; c++) {
Chris@72: 
Chris@72:             float newCost = m_m->getDistance(r, c);
Chris@72:             Matcher::Advance dir = Matcher::AdvanceNone;
Chris@72: 
Chris@72:             if (r > r1) {	// not first row
Chris@72:                 double min = -1;
Chris@72:                 if ((c > prevRowStart) && (c <= prevRowStop)) {
Chris@72:                     // diagonal from (r-1,c-1)
Chris@83:                     min = m_m->getPathCost(r-1, c-1) + newCost * diagonalWeight;
Chris@72:                     dir = Matcher::AdvanceBoth;
Chris@72:                 }
Chris@72:                 if ((c >= prevRowStart) && (c < prevRowStop)) {
Chris@72:                     // vertical from (r-1,c)
Chris@72:                     double cost = m_m->getPathCost(r-1, c) + newCost;
Chris@72:                     if ((min < 0) || (cost < min)) {
Chris@72:                         min = cost;
Chris@72:                         dir = Matcher::AdvanceThis;
Chris@72:                     }
Chris@72:                 }
Chris@72:                 if (c > rowStart) {
Chris@72:                     // horizontal from (r,c-1)
Chris@72:                     double cost = m_m->getPathCost(r, c-1) + newCost;
Chris@72:                     if ((min < 0) || (cost < min)) {
Chris@72:                         min = cost;
Chris@72:                         dir = Matcher::AdvanceOther;
Chris@72:                     }
Chris@72:                 }
Chris@72:                 
Chris@72:                 m_m->setPathCost(r, c, dir, min);
Chris@72: 
Chris@72:             } else if (c > rowStart) {	// first row
Chris@72:                 // horizontal from (r,c-1)
Chris@72:                 m_m->setPathCost(r, c, Matcher::AdvanceOther,
Chris@72:                                    m_m->getPathCost(r, c-1) + newCost);
Chris@72:             }
Chris@72:         }
Chris@72: 
Chris@72:         prevRowStart = rowStart;
Chris@72:         prevRowStop = rowStop;
cannam@0:     }
Chris@72: } 
Chris@30: 
Chris@82: #ifdef PERFORM_ERROR_CHECKS
Chris@81: Finder::ErrorPosition
Chris@81: Finder::checkPathCostMatrix() 
Chris@81: {
Chris@81:     ErrorPosition err;
Chris@81: 
Chris@81:     int r1 = 0;
Chris@81:     int c1 = 0;
Chris@81:     int r2 = m_m->getFrameCount() - 1;
Chris@81:     int c2 = m_m->getOtherFrameCount() - 1;
Chris@81: 
Chris@81:     if (r2 < r1 || c2 < c1) {
Chris@81:         return err;
Chris@81:     }
Chris@81: 
Chris@81:     int prevRowStart = 0, prevRowStop = 0;
Chris@81: 
Chris@83:     float diagonalWeight = m_m->getDiagonalWeight();
Chris@83:     
Chris@81:     for (int r = r1; r <= r2; r++) {
Chris@81: 
Chris@81:         pair<int, int> colRange = m_m->getColRange(r);
Chris@81: 
Chris@81:         int rowStart = max(c1, colRange.first);
Chris@81:         int rowStop = min(c2 + 1, colRange.second);
Chris@81:         
Chris@81:         for (int c = rowStart; c < rowStop; c++) {
Chris@81: 
Chris@81:             float newCost = m_m->getDistance(r, c);
Chris@81:             double updateTo = -1.0;
Chris@81:             Matcher::Advance dir = Matcher::AdvanceNone;
Chris@81: 
Chris@95:             if (r > r1) { // not first row
Chris@81:                 double min = -1;
Chris@81:                 if ((c > prevRowStart) && (c <= prevRowStop)) {
Chris@81:                     // diagonal from (r-1,c-1)
Chris@83:                     min = m_m->getPathCost(r-1, c-1) + newCost * diagonalWeight;
Chris@81:                     err.prevCost = m_m->getPathCost(r-1, c-1);
Chris@83:                     err.distance = newCost * diagonalWeight;
Chris@81:                     dir = Matcher::AdvanceBoth;
Chris@81:                 }
Chris@81:                 if ((c >= prevRowStart) && (c < prevRowStop)) {
Chris@81:                     // vertical from (r-1,c)
Chris@81:                     double cost = m_m->getPathCost(r-1, c) + newCost;
Chris@81:                     if ((min < 0) || (cost < min)) {
Chris@81:                         min = cost;
Chris@81:                         err.prevCost = m_m->getPathCost(r-1, c);
Chris@81:                         err.distance = newCost;
Chris@81:                         dir = Matcher::AdvanceThis;
Chris@81:                     }
Chris@81:                 }
Chris@81:                 if (c > rowStart) {
Chris@81:                     // horizontal from (r,c-1)
Chris@81:                     double cost = m_m->getPathCost(r, c-1) + newCost;
Chris@81:                     if ((min < 0) || (cost < min)) {
Chris@81:                         min = cost;
Chris@81:                         err.prevCost = m_m->getPathCost(r, c-1);
Chris@81:                         err.distance = newCost;
Chris@81:                         dir = Matcher::AdvanceOther;
Chris@81:                     }
Chris@81:                 }
Chris@81: 
Chris@81:                 updateTo = min;
Chris@81: 
Chris@82:             } else { // first row
Chris@82: 
Chris@82:                 if (c > rowStart) {
Chris@82:                     // horizontal from (r,c-1)
Chris@83:                     updateTo = m_m->getPathCost(r, c-1) + newCost;
Chris@83:                     err.prevCost = m_m->getPathCost(r, c-1);
Chris@83:                     err.distance = newCost;
Chris@82:                     dir = Matcher::AdvanceOther;
Chris@82:                 }
Chris@81:             }
Chris@81: 
Chris@82:             if (dir != Matcher::AdvanceNone) {
Chris@86:                 if (m_m->getAdvance(r, c) != dir) {
Chris@86:                     err.type = ErrorPosition::WrongAdvance;
Chris@86:                     err.r = r;
Chris@86:                     err.c = c;
Chris@86:                     err.costWas = m_m->getPathCost(r, c);
Chris@86:                     err.costShouldBe = updateTo;
Chris@86:                     err.advanceWas = m_m->getAdvance(r, c);
Chris@86:                     err.advanceShouldBe = dir;
Chris@86:                     return err;
Chris@86:                 }
Chris@84:                 if (m_m->getPathCost(r, c) != updateTo) {
Chris@84:                     err.type = ErrorPosition::WrongCost;
Chris@84:                     err.r = r;
Chris@84:                     err.c = c;
Chris@84:                     err.costWas = m_m->getPathCost(r, c);
Chris@84:                     err.costShouldBe = updateTo;
Chris@84:                     err.advanceWas = m_m->getAdvance(r, c);
Chris@84:                     err.advanceShouldBe = dir;
Chris@82:                     return err;
Chris@82:                 }
Chris@82:             } else {
Chris@82:                 // AdvanceNone should occur only at r = r1, c = c1
Chris@82:                 if (r != r1 || c != c1) {
Chris@82:                     err.type = ErrorPosition::NoAdvance;
Chris@82:                     err.r = r;
Chris@82:                     err.c = c;
Chris@82:                     err.costWas = m_m->getPathCost(r, c);
Chris@82:                     err.costShouldBe = updateTo;
Chris@84:                     err.advanceWas = m_m->getAdvance(r, c);
Chris@84:                     err.advanceShouldBe = dir;
Chris@82:                     return err;
Chris@82:                 }
Chris@81:             }
Chris@81:         }
Chris@81: 
Chris@81:         prevRowStart = rowStart;
Chris@81:         prevRowStop = rowStop;
Chris@81:     }
Chris@81: 
Chris@81:     return err;
Chris@82: }
Chris@81: 
Chris@92: void
Chris@92: Finder::checkAndReport()
Chris@30: {
Chris@92:     cerr << "Finder: Checking path-cost matrix..." << endl;
Chris@82:     ErrorPosition err = checkPathCostMatrix();
Chris@92:     if (err.type == ErrorPosition::NoError) {
Chris@92:         cerr << "No errors found" << endl;
Chris@92:     } else {
Chris@82:         cerr << "\nWARNING: Checking path-cost matrix returned mismatch:" << endl;
Chris@92:         cerr << "Type: " << err.type << ": ";
Chris@92:         switch (err.type) {
Chris@92:         case ErrorPosition::NoError: break;
Chris@92:         case ErrorPosition::WrongCost: cerr << "WrongCost"; break;
Chris@92:         case ErrorPosition::WrongAdvance: cerr << "WrongAdvance"; break;
Chris@92:         case ErrorPosition::NoAdvance: cerr << "NoAdvance"; break;
Chris@92:         }
Chris@92:         cerr << endl;
Chris@84:         cerr << "At row " << err.r << ", column " << err.c
Chris@84:              << "\nShould be advancing "
Chris@84:              << Matcher::advanceToString(err.advanceShouldBe)
Chris@84:              << ", advance in matrix is "
Chris@84:              << Matcher::advanceToString(err.advanceWas)
Chris@83:              << "\nPrev cost " << err.prevCost
Chris@82:              << " plus distance " << err.distance << " gives "
Chris@84:              << err.costShouldBe << ", matrix contains " << err.costWas
Chris@83:              << endl;
Chris@83:         cerr << "Note: diagonal weight = " << m_m->getDiagonalWeight() << endl;
Chris@83:         cerr << endl;
Chris@92: 
Chris@95:         int w(4);
Chris@95:         int ww(15);
Chris@92: 
Chris@92:         cerr << "Distance matrix leading up to this point:" << endl;
Chris@95:         cerr << setprecision(12) << setw(w) << "";
Chris@92:         for (int i = -4; i <= 0; ++i) {
Chris@95:             cerr << setw(ww) << i;
Chris@92:         }
Chris@92:         cerr << endl;
Chris@92:         for (int j = -4; j <= 0; ++j) {
Chris@92:             cerr << setw(w) << j;
Chris@92:             for (int i = -4; i <= 0; ++i) {
Chris@95:                 cerr << setw(ww) << m_m->getDistance(err.r + j, err.c + i);
Chris@92:             }
Chris@92:             cerr << endl;
Chris@92:         }
Chris@92:         cerr << endl;
Chris@92: 
Chris@92:         cerr << "Cost matrix leading up to this point:" << endl;
Chris@92:         cerr << setw(w) << "";
Chris@92:         for (int i = -4; i <= 0; ++i) {
Chris@95:             cerr << setw(ww) << i;
Chris@92:         }
Chris@92:         cerr << endl;
Chris@92:         for (int j = -4; j <= 0; ++j) {
Chris@92:             cerr << setw(w) << j;
Chris@92:             for (int i = -4; i <= 0; ++i) {
Chris@95:                 cerr << setw(ww) << m_m->getPathCost(err.r + j, err.c + i);
Chris@92:             }
Chris@92:             cerr << endl;
Chris@92:         }
Chris@92:         cerr << endl;
Chris@82:     }
Chris@92: }
Chris@92: #endif
Chris@92: 
Chris@92: int
Chris@92: Finder::retrievePath(bool smooth, vector<int> &pathx, vector<int> &pathy)
Chris@92: {
Chris@92:     pathx.clear();
Chris@92:     pathy.clear();
Chris@92: 
Chris@92: #ifdef PERFORM_ERROR_CHECKS
Chris@92:     checkAndReport();
Chris@82: #endif
Chris@82: 
Chris@72:     int ex = m_m->getOtherFrameCount() - 1;
Chris@72:     int ey = m_m->getFrameCount() - 1;
Chris@69: 
Chris@69:     if (ex < 0 || ey < 0) {
Chris@69:         return 0;
Chris@69:     }
Chris@66:     
Chris@66:     int x = ex;
Chris@66:     int y = ey;
Chris@30: 
Chris@140: #ifdef DEBUG_FINDER
Chris@140:     cerr << "*** retrievePath: smooth = " << smooth << endl;
Chris@140:     cerr << "*** retrievePath: before: x = " << x << ", y = " << y << endl;
Chris@140: #endif
Chris@140: 
Chris@72:     if (m_duration2 > 0 && m_duration2 < m_m->getOtherFrameCount()) {
Chris@72:         x = m_duration2 - 1;
Chris@60:     }
Chris@72:     if (m_duration1 > 0 && m_duration1 < m_m->getFrameCount()) {
Chris@72:         y = m_duration1 - 1;
Chris@60:     }
Chris@60: 
Chris@72:     if (!m_m->isAvailable(y, x)) {
Chris@66:         // Path did not pass through the expected end point --
Chris@66:         // probably means the pieces are substantially different in
Chris@66:         // the later bits. Reset the expected end point to the end of
Chris@66:         // both files including any trailing silence.
Chris@66:         cerr << "NOTE: Path did not pass through expected end point, inputs are probably significantly different" << endl;
Chris@66:         x = ex;
Chris@66:         y = ey;
Chris@66:     }
Chris@66: 
Chris@55:     recalculatePathCostMatrix(0, 0, y, x);
Chris@55: 
Chris@140: #ifdef DEBUG_FINDER
Chris@140:     cerr << "*** retrievePath: start: x = " << x << ", y = " << y << endl;
Chris@140: #endif
Chris@66:     
Chris@72:     while (m_m->isAvailable(y, x) && (x > 0 || y > 0)) {
Chris@30: 
Chris@33: //        cerr << "x = " << x << ", y = " << y;
Chris@33:         
Chris@30:         pathx.push_back(x);
Chris@30:         pathy.push_back(y);
Chris@30: 
Chris@72:         switch (m_m->getAdvance(y, x)) {
Chris@45:         case Matcher::AdvanceThis:
Chris@70: //            cerr << ", going down (dist = " << getDistance() << ")" << endl;
Chris@33:             y--;
Chris@33:             break;
Chris@45:         case Matcher::AdvanceOther:
Chris@70: //            cerr << ", going left (dist = " << getDistance() << ")" << endl;
Chris@33:             x--;
Chris@33:             break;
Chris@45:         case Matcher::AdvanceBoth:
Chris@70: //            cerr << ", going diag (dist = " << getDistance() << ")" << endl;
Chris@33:             x--;
Chris@33:             y--;
Chris@33:             break;
Chris@45:         case Matcher::AdvanceNone: // this would indicate a bug, but we wouldn't want to hang
Chris@69:             cerr << "WARNING: Neither matcher advanced in path backtrack at (" << x << "," << y << ")" << endl;
Chris@33:             if (x > y) {
Chris@33:                 x--;
Chris@33:             } else {
Chris@33:                 y--;
Chris@33:             }
Chris@33:             break;
Chris@30:         }
Chris@30:     }
Chris@30: 
Chris@72:     if (x > 0 || y > 0) {
Chris@72:         cerr << "WARNING: Ran out of available path at (" << y << "," << x
Chris@72:              << ")!" << endl;
Chris@72:     }
Chris@72:     
Chris@72:     reverse(pathx.begin(), pathx.end());
Chris@72:     reverse(pathy.begin(), pathy.end());
Chris@30: 
Chris@31:     if (smooth) {
Chris@31:         int smoothedLen = Path().smooth(pathx, pathy, pathx.size());
Chris@31:         return smoothedLen;
Chris@31:     } else {
Chris@31:         return pathx.size();
Chris@31:     }
Chris@30: }
Chris@30: 
Chris@30: