Mercurial > hg > svapp
changeset 780:8fa98f89eda8
Add subsequence DTW (not yet in use)
| author | Chris Cannam |
|---|---|
| date | Wed, 01 Jul 2020 15:34:46 +0100 |
| parents | 5de2b710cfae |
| children | b651dc5ff555 |
| files | align/DTW.h align/TransformDTWAligner.cpp |
| diffstat | 2 files changed, 131 insertions(+), 77 deletions(-) [+] |
line wrap: on
line diff
--- a/align/DTW.h Wed Jul 01 11:41:07 2020 +0100 +++ b/align/DTW.h Wed Jul 01 15:34:46 2020 +0100 @@ -27,54 +27,22 @@ DTW(std::function<double(const Value &, const Value &)> distanceMetric) : m_metric(distanceMetric) { } - std::vector<size_t> alignSeries(std::vector<Value> s1, - std::vector<Value> s2) { + /** + * Align the sequence s2 against the whole of the sequence s1, + * returning the index into s1 for each element in s2. + */ + std::vector<size_t> alignSequences(std::vector<Value> s1, + std::vector<Value> s2) { + return align(s1, s2, false); + } - // Return the index into s2 for each element in s1 - - std::vector<size_t> alignment(s1.size(), 0); - - if (s1.empty() || s2.empty()) { - return alignment; - } - - auto costs = costSeries(s1, s2); - -#ifdef DEBUG_DTW - SVCERR << "Cost matrix:" << endl; - for (auto v: costs) { - for (auto x: v) { - SVCERR << x << " "; - } - SVCERR << "\n"; - } -#endif - - size_t j = s1.size() - 1; - size_t i = s2.size() - 1; - - while (j > 0 && i > 0) { - - alignment[j] = i; - - cost_t a = costs[j-1][i]; - cost_t b = costs[j][i-1]; - cost_t both = costs[j-1][i-1]; - - if (a < b) { - --j; - if (both <= a) { - --i; - } - } else { - --i; - if (both <= b) { - --j; - } - } - } - - return alignment; + /** + * Align the sequence sub against the best-matching subsequence of + * s, returning the index into s for each element in sub. + */ + std::vector<size_t> alignSubsequence(std::vector<Value> s, + std::vector<Value> sub) { + return align(s, sub, true); } private: @@ -102,8 +70,9 @@ } } - std::vector<std::vector<cost_t>> costSeries(std::vector<Value> s1, - std::vector<Value> s2) { + std::vector<std::vector<cost_t>> costSequences(std::vector<Value> s1, + std::vector<Value> s2, + bool subsequence) { std::vector<std::vector<cost_t>> costs (s1.size(), std::vector<cost_t>(s2.size(), 0.0)); @@ -111,22 +80,97 @@ for (size_t j = 0; j < s1.size(); ++j) { for (size_t i = 0; i < s2.size(); ++i) { cost_t c = m_metric(s1[j], s2[i]); - costs[j][i] = choose - ( - { j > 0, - j > 0 ? c + costs[j-1][i] : 0.0 - }, - { i > 0, - i > 0 ? c + costs[j][i-1] : 0.0 - }, - { j > 0 && i > 0, - j > 0 && i > 0 ? c + costs[j-1][i-1] : 0.0 - }); + if (i == 0 && subsequence) { + costs[j][i] = c; + } else { + costs[j][i] = choose + ( + { j > 0, + j > 0 ? c + costs[j-1][i] : 0.0 + }, + { i > 0, + i > 0 ? c + costs[j][i-1] : 0.0 + }, + { j > 0 && i > 0, + j > 0 && i > 0 ? c + costs[j-1][i-1] : 0.0 + }); + } } } return costs; } + + std::vector<size_t> align(const std::vector<Value> &s1, + const std::vector<Value> &s2, + bool subsequence) { + + // Return the index into s1 for each element in s2 + + std::vector<size_t> alignment(s2.size(), 0); + + if (s1.empty() || s2.empty()) { + return alignment; + } + + auto costs = costSequences(s1, s2, subsequence); + +#ifdef DEBUG_DTW + SVCERR << "Cost matrix:" << endl; + for (auto v: cost) { + for (auto x: v) { + SVCERR << x << " "; + } + SVCERR << "\n"; + } +#endif + + size_t j = s1.size() - 1; + size_t i = s2.size() - 1; + + if (subsequence) { + cost_t min = 0.0; + size_t minidx = 0; + for (size_t j = 0; j < s1.size(); ++j) { + if (j == 0 || costs[j][i] < min) { + min = costs[j][i]; + minidx = j; + } + } + j = minidx; +#ifdef DEBUG_DTW + SVCERR << "Lowest cost at end of subsequence = " << min + << " at index " << j << ", tracking back from there" << endl; +#endif + } + + while (i > 0 && (j > 0 || subsequence)) { + + alignment[i] = j; + + cost_t a = costs[j-1][i]; + cost_t b = costs[j][i-1]; + cost_t both = costs[j-1][i-1]; + + if (a < b) { + --j; + if (both <= a) { + --i; + } + } else { + --i; + if (both <= b) { + --j; + } + } + } + + if (subsequence) { + alignment[0] = j; + } + + return alignment; + } }; class MagnitudeDTW @@ -134,9 +178,14 @@ public: MagnitudeDTW() : m_dtw(metric) { } - std::vector<size_t> alignSeries(std::vector<double> s1, - std::vector<double> s2) { - return m_dtw.alignSeries(s1, s2); + std::vector<size_t> alignSequences(std::vector<double> s1, + std::vector<double> s2) { + return m_dtw.alignSequences(s1, s2); + } + + std::vector<size_t> alignSubsequence(std::vector<double> s, + std::vector<double> sub) { + return m_dtw.alignSubsequence(s, sub); } private: @@ -163,9 +212,14 @@ RiseFallDTW() : m_dtw(metric) { } - std::vector<size_t> alignSeries(std::vector<Value> s1, - std::vector<Value> s2) { - return m_dtw.alignSeries(s1, s2); + std::vector<size_t> alignSequences(std::vector<Value> s1, + std::vector<Value> s2) { + return m_dtw.alignSequences(s1, s2); + } + + std::vector<size_t> alignSubsequence(std::vector<Value> s, + std::vector<Value> sub) { + return m_dtw.alignSubsequence(s, sub); } private:
--- a/align/TransformDTWAligner.cpp Wed Jul 01 11:41:07 2020 +0100 +++ b/align/TransformDTWAligner.cpp Wed Jul 01 15:34:46 2020 +0100 @@ -297,19 +297,19 @@ for (int i = 0; in_range_for(alignment, i); ++i) { - // DTW returns "the index into s2 for each element in s1" - sv_frame_t refFrame = refFrames[i]; - - if (!in_range_for(otherFrames, alignment[i])) { + // DTW returns "the index into s1 for each element in s2" + sv_frame_t alignedFrame = otherFrames[i]; + + if (!in_range_for(refFrames, alignment[i])) { SVCERR << "TransformDTWAligner::makePath: Internal error: " - << "DTW maps index " << i << " in reference frame vector " - << "(size " << refFrames.size() << ") onto index " - << alignment[i] << " in other frame vector " - << "(only size " << otherFrames.size() << ")" << endl; + << "DTW maps index " << i << " in other frame vector " + << "(size " << otherFrames.size() << ") onto index " + << alignment[i] << " in ref frame vector " + << "(only size " << refFrames.size() << ")" << endl; continue; } - sv_frame_t alignedFrame = otherFrames[alignment[i]]; + sv_frame_t refFrame = refFrames[alignment[i]]; path.add(PathPoint(alignedFrame, refFrame)); } @@ -361,7 +361,7 @@ << "]: serialising DTW to avoid over-allocation" << endl; #endif QMutexLocker locker(&m_dtwMutex); - alignment = dtw.alignSeries(s1, s2); + alignment = dtw.alignSequences(s1, s2); } #ifdef DEBUG_TRANSFORM_DTW_ALIGNER @@ -450,7 +450,7 @@ << "]: serialising DTW to avoid over-allocation" << endl; #endif QMutexLocker locker(&m_dtwMutex); - alignment = dtw.alignSeries(s1, s2); + alignment = dtw.alignSequences(s1, s2); } #ifdef DEBUG_TRANSFORM_DTW_ALIGNER