From 49cf25d72483b7294ff531d9ed5b8ea241c0d0f9 Mon Sep 17 00:00:00 2001 From: Chris Cannam Date: Tue, 25 May 2010 21:49:05 +0100 Subject: [PATCH] * Further work on handling distinct analysis and synthesis window and FFT sizes --- Makefile.in | 13 ++--- src/StretcherImpl.cpp | 51 ++++++++++++++++--- src/StretcherImpl.h | 24 +++++++++ src/StretcherProcess.cpp | 67 +++++++++++++++++------- src/dsp/SincWindow.h | 107 +++++++++++++++++++++++++++++++++++++++ src/dsp/Window.h | 7 +-- src/system/VectorOps.h | 10 ++++ 7 files changed, 246 insertions(+), 33 deletions(-) create mode 100644 src/dsp/SincWindow.h diff --git a/Makefile.in b/Makefile.in index 7ab25ff..c8268b1 100644 --- a/Makefile.in +++ b/Makefile.in @@ -65,6 +65,7 @@ LIBRARY_INCLUDES := \ src/dsp/Resampler.h \ src/dsp/FFT.h \ src/dsp/MovingMedian.h \ + src/dsp/SincWindow.h \ src/dsp/Window.h \ src/system/Allocators.h \ src/system/Thread.h \ @@ -175,12 +176,12 @@ depend: src/rubberband-c.o: rubberband/rubberband-c.h src/rubberband-c.o: rubberband/RubberBandStretcher.h src/RubberBandStretcher.o: src/StretcherImpl.h -src/RubberBandStretcher.o: rubberband/RubberBandStretcher.h src/dsp/Window.h +src/RubberBandStretcher.o: rubberband/RubberBandStretcher.h src/dsp/Window.h src/dsp/SincWindow.h src/RubberBandStretcher.o: src/dsp/FFT.h src/base/RingBuffer.h src/RubberBandStretcher.o: src/base/Scavenger.h src/system/Thread.h src/RubberBandStretcher.o: src/system/Thread.h src/system/sysutils.h src/StretcherProcess.o: src/StretcherImpl.h rubberband/RubberBandStretcher.h -src/StretcherProcess.o: src/dsp/Window.h src/dsp/FFT.h src/base/RingBuffer.h +src/StretcherProcess.o: src/dsp/Window.h src/dsp/SincWindow.h src/dsp/FFT.h src/base/RingBuffer.h src/StretcherProcess.o: src/base/Scavenger.h src/system/Thread.h src/StretcherProcess.o: src/system/Thread.h src/system/sysutils.h src/StretcherProcess.o: src/dsp/PercussiveAudioCurve.h @@ -198,7 +199,7 @@ src/dsp/AudioCurveCalculator.o: src/system/sysutils.h src/dsp/SpectralDifferenceAudioCurve.o: src/dsp/SpectralDifferenceAudioCurve.h src/dsp/SpectralDifferenceAudioCurve.o: src/dsp/AudioCurveCalculator.h src/dsp/SpectralDifferenceAudioCurve.o: src/system/sysutils.h -src/dsp/SpectralDifferenceAudioCurve.o: src/dsp/Window.h +src/dsp/SpectralDifferenceAudioCurve.o: src/dsp/Window.h src/dsp/SincWindow.h src/dsp/SpectralDifferenceAudioCurve.o: src/system/VectorOps.h src/dsp/SpectralDifferenceAudioCurve.o: src/system/sysutils.h src/dsp/HighFrequencyAudioCurve.o: src/dsp/HighFrequencyAudioCurve.h @@ -223,19 +224,19 @@ src/system/Allocators.o: src/system/Allocators.h src/system/VectorOps.h src/system/Allocators.o: src/system/sysutils.h src/system/sysutils.o: src/system/sysutils.h src/StretcherChannelData.o: src/StretcherChannelData.h src/StretcherImpl.h -src/StretcherChannelData.o: rubberband/RubberBandStretcher.h src/dsp/Window.h +src/StretcherChannelData.o: rubberband/RubberBandStretcher.h src/dsp/Window.h src/dsp/SincWindow.h src/StretcherChannelData.o: src/dsp/FFT.h src/base/RingBuffer.h src/StretcherChannelData.o: src/base/Scavenger.h src/system/Thread.h src/StretcherChannelData.o: src/system/Thread.h src/system/sysutils.h src/StretcherChannelData.o: src/dsp/Resampler.h src/system/Allocators.h src/StretcherChannelData.o: src/system/VectorOps.h src/system/sysutils.h src/StretcherImpl.o: src/StretcherImpl.h rubberband/RubberBandStretcher.h -src/StretcherImpl.o: src/dsp/Window.h src/dsp/FFT.h src/base/RingBuffer.h +src/StretcherImpl.o: src/dsp/Window.h src/dsp/SincWindow.h src/dsp/FFT.h src/base/RingBuffer.h src/StretcherImpl.o: src/base/Scavenger.h src/system/Thread.h src/system/Thread.h src/StretcherImpl.o: src/system/sysutils.h src/dsp/PercussiveAudioCurve.h src/StretcherImpl.o: src/dsp/AudioCurveCalculator.h src/StretcherImpl.o: src/dsp/HighFrequencyAudioCurve.h -src/StretcherImpl.o: src/dsp/SpectralDifferenceAudioCurve.h src/dsp/Window.h +src/StretcherImpl.o: src/dsp/SpectralDifferenceAudioCurve.h src/dsp/Window.h src/dsp/SincWindow.h src/StretcherImpl.o: src/system/VectorOps.h src/system/sysutils.h src/StretcherImpl.o: src/dsp/SilentAudioCurve.h src/dsp/ConstantAudioCurve.h src/StretcherImpl.o: src/dsp/Resampler.h src/StretchCalculator.h diff --git a/src/StretcherImpl.cpp b/src/StretcherImpl.cpp index 26c7b63..39816f0 100644 --- a/src/StretcherImpl.cpp +++ b/src/StretcherImpl.cpp @@ -81,6 +81,7 @@ RubberBandStretcher::Impl::Impl(size_t sampleRate, m_debugLevel(m_defaultDebugLevel), m_mode(JustCreated), m_awindow(0), + m_asinc(0), m_swindow(0), m_studyFFT(0), m_spaceAvailable("space"), @@ -194,6 +195,10 @@ RubberBandStretcher::Impl::~Impl() i != m_windows.end(); ++i) { delete i->second; } + for (map *>::iterator i = m_sincs.begin(); + i != m_sincs.end(); ++i) { + delete i->second; + } } void @@ -487,7 +492,7 @@ RubberBandStretcher::Impl::calculateSizes() } } - // windowSize can be almost anything, but it can't be greater than + // m_fftSize can be almost anything, but it can't be greater than // 4 * m_baseFftSize unless ratio is less than 1/1024. m_fftSize = windowSize; @@ -515,7 +520,7 @@ RubberBandStretcher::Impl::calculateSizes() size_t (ceil(max (m_maxProcessSize / m_pitchScale, - m_sWindowSize * 2 * (m_timeRatio > 1.f ? m_timeRatio : 1.f)))); + m_maxProcessSize * 2 * (m_timeRatio > 1.f ? m_timeRatio : 1.f)))); if (m_realtime) { // This headroom is so as to try to avoid reallocation when @@ -586,8 +591,12 @@ RubberBandStretcher::Impl::configure() if (m_windows.find(*i) == m_windows.end()) { m_windows[*i] = new Window(HanningWindow, *i); } + if (m_sincs.find(*i) == m_sincs.end()) { + m_sincs[*i] = new SincWindow(*i, *i); + } } m_awindow = m_windows[m_aWindowSize]; + m_asinc = m_sincs[m_aWindowSize]; m_swindow = m_windows[m_sWindowSize]; if (m_debugLevel > 0) { @@ -728,15 +737,20 @@ RubberBandStretcher::Impl::reconfigure() std::cerr << "WARNING: reconfigure(): window allocation (size " << m_aWindowSize << ") required in RT mode" << std::endl; m_windows[m_aWindowSize] = new Window (HanningWindow, m_aWindowSize); + m_sincs[m_aWindowSize] = new SincWindow + (m_aWindowSize, m_aWindowSize); } if (m_windows.find(m_sWindowSize) == m_windows.end()) { std::cerr << "WARNING: reconfigure(): window allocation (size " << m_sWindowSize << ") required in RT mode" << std::endl; m_windows[m_sWindowSize] = new Window (HanningWindow, m_sWindowSize); + m_sincs[m_sWindowSize] = new SincWindow + (m_sWindowSize, m_sWindowSize); } m_awindow = m_windows[m_aWindowSize]; + m_asinc = m_sincs[m_aWindowSize]; m_swindow = m_windows[m_sWindowSize]; for (size_t c = 0; c < m_channels; ++c) { @@ -933,13 +947,36 @@ RubberBandStretcher::Impl::study(const float *const *input, size_t samples, bool size_t got = inbuf.peek(cd.accumulator, m_aWindowSize); assert(final || got == m_aWindowSize); - m_awindow->cut(cd.accumulator); + if (m_aWindowSize == m_fftSize) { - //!!! insert appropriate cunning here if m_aWindowSize != - //!!! m_fftSize + // We don't need the fftshift for studying, as we're + // only interested in magnitude. - // We don't need the fftshift for studying, as we're only - // interested in magnitude + m_awindow->cut(cd.accumulator); + + } else { + + // If we need to fold (i.e. if the window size is + // greater than the fft size so we are doing a + // time-aliased presum fft) or zero-pad, then we might + // as well use our standard function for it. This + // means we retain the m_asinc cut if folding as well, + // which is good for consistency with real-time mode. + // We get fftshift as well, which we don't want, but + // the penalty is nominal. + + // Note that we can't do this in-place. Pity + + float *tmp = (float *)alloca + (std::max(m_fftSize, m_aWindowSize) * sizeof(float)); + + if (m_aWindowSize > m_fftSize) { + m_asinc->cut(cd.accumulator); + } + + cutShiftAndFold(tmp, m_fftSize, cd.accumulator, m_awindow); + v_copy(cd.accumulator, tmp, m_fftSize); + } m_studyFFT->forwardMagnitude(cd.accumulator, cd.fltbuf); diff --git a/src/StretcherImpl.h b/src/StretcherImpl.h index e7c5201..eb5a545 100644 --- a/src/StretcherImpl.h +++ b/src/StretcherImpl.h @@ -18,6 +18,7 @@ #include "rubberband/RubberBandStretcher.h" #include "dsp/Window.h" +#include "dsp/SincWindow.h" #include "dsp/FFT.h" #include "dsp/CompoundAudioCurve.h" @@ -118,6 +119,27 @@ protected: size_t roundUp(size_t value); // to next power of two + template + void cutShiftAndFold(T *target, int targetSize, + float *src, // destructive to src + Window *window) { + window->cut(src); + const int windowSize = window->getSize(); + const int hs = targetSize / 2; + if (windowSize == targetSize) { + v_convert(target, src + hs, hs); + v_convert(target + hs, src, hs); + } else { + v_zero(target, targetSize); + int j = targetSize - windowSize/2; + while (j < 0) j += targetSize; + for (int i = 0; i < windowSize; ++i) { + target[j] += src[i]; + if (++j == targetSize) j = 0; + } + } + } + bool resampleBeforeStretching() const; double m_timeRatio; @@ -149,7 +171,9 @@ protected: ProcessMode m_mode; std::map *> m_windows; + std::map *> m_sincs; Window *m_awindow; + SincWindow *m_asinc; Window *m_swindow; FFT *m_studyFFT; diff --git a/src/StretcherProcess.cpp b/src/StretcherProcess.cpp index 8918812..46736fc 100644 --- a/src/StretcherProcess.cpp +++ b/src/StretcherProcess.cpp @@ -641,26 +641,17 @@ RubberBandStretcher::Impl::analyseChunk(size_t channel) double *const R__ dblbuf = cd.dblbuf; float *const R__ fltbuf = cd.fltbuf; - int sz = m_fftSize; - int hs = sz / 2; +//!!! int sz = m_fftSize; +// int hs = sz / 2; // cd.fltbuf is known to contain m_aWindowSize samples - m_awindow->cut(fltbuf); - - if (m_aWindowSize == m_fftSize) { - v_convert(dblbuf, fltbuf + hs, hs); - v_convert(dblbuf + hs, fltbuf, hs); - } else { - v_zero(dblbuf, sz); - int j = sz - m_aWindowSize/2; - while (j < 0) j += sz; - for (int i = 0; i < m_aWindowSize; ++i) { - dblbuf[j] += fltbuf[i]; - if (++j == sz) j = 0; - } + if (m_aWindowSize > m_fftSize) { + m_asinc->cut(fltbuf); } + cutShiftAndFold(dblbuf, m_fftSize, fltbuf, m_awindow); + cd.fft->forwardPolar(dblbuf, cd.mag, cd.phase); } @@ -902,12 +893,25 @@ RubberBandStretcher::Impl::synthesiseChunk(size_t channel) v_convert(fltbuf + hs, dblbuf, hs); } else { v_zero(fltbuf, m_sWindowSize); +#ifdef PARP +//!!! centre fft frame (effectively only if swindowsize == fftsize, else discontinuities) + const int count = std::min(m_fftSize, m_sWindowSize); + int i = m_sWindowSize/2 - count/2; + int j = m_fftSize - count/2; + for (int k = 0; k < count; ++k) { + fltbuf[i] += dblbuf[j]; + if (i++ == m_sWindowSize) i = 0; + if (j++ == m_fftSize) j = 0; + } +#else +//!!! unwrap int j = sz - m_sWindowSize/2; while (j < 0) j += sz; for (int i = 0; i < m_sWindowSize; ++i) { fltbuf[i] += dblbuf[j]; if (++j == sz) j = 0; } +#endif } } @@ -917,8 +921,37 @@ RubberBandStretcher::Impl::synthesiseChunk(size_t channel) cd.accumulatorFill = m_sWindowSize; - float fixed = m_swindow->getArea() * 1.5f; - m_swindow->add(windowAccumulator, fixed); +/* + std::cerr << "Note: synthesis window area = " << m_swindow->getArea() + << ", analysis window area = " << m_awindow->getArea() + << ", analysis sinc area = " << m_asinc->getArea() + << std::endl; +*/ + + //!!! + + float *tmp = (float *)alloca(m_sWindowSize * sizeof(float)); + v_set(tmp, 1.f, m_sWindowSize); + +//!!! m_awindow->cut(tmp); //!!! if (m_aWindowSize == m_sWindowSize) + m_swindow->cut(tmp); +/* + if (m_aWindowSize != m_fftSize) { + if (m_aWindowSize == m_sWindowSize) { + m_asinc->cut(tmp); + } else { + int sourceIndex = m_aWindowSize/2 - m_sWindowSize/2; + for (int i = 0; i < m_sWindowSize; ++i) { + if (sourceIndex >= 0 && sourceIndex < m_aWindowSize) { + tmp[i] *= m_asinc->getValue(sourceIndex); + } + ++sourceIndex; + } + } + } +*/ + v_add(windowAccumulator, tmp, m_sWindowSize); + } void diff --git a/src/dsp/SincWindow.h b/src/dsp/SincWindow.h new file mode 100644 index 0000000..1073872 --- /dev/null +++ b/src/dsp/SincWindow.h @@ -0,0 +1,107 @@ +/* -*- c-basic-offset: 4 indent-tabs-mode: nil -*- vi:set ts=8 sts=4 sw=4: */ + +/* + Rubber Band + An audio time-stretching and pitch-shifting library. + Copyright 2007-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 _RUBBERBAND_SINC_WINDOW_H_ +#define _RUBBERBAND_SINC_WINDOW_H_ + +#include +#include +#include +#include + +#include "system/sysutils.h" +#include "system/VectorOps.h" +#include "system/Allocators.h" + +namespace RubberBand { + +template +class SincWindow +{ +public: + /** + * Construct a sinc windower which produces a window of the given + * size containing the values of sinc(x) with x=0 at the centre, + * such that the distance from -pi to pi (the point at which the + * sinc function first crosses zero, for negative and positive + * arguments respectively) is p samples. + */ + SincWindow(int size, int p) : m_size(size), m_p(p) { encache(); } + SincWindow(const SincWindow &w) : m_size(w.m_size), m_p(w.m_p) { encache(); } + SincWindow &operator=(const SincWindow &w) { + if (&w == this) return *this; + m_size = w.m_size; + m_p = w.m_p; + encache(); + return *this; + } + virtual ~SincWindow() { delete[] m_cache; } + + inline void cut(T *const R__ block) const { + v_multiply(block, m_cache, m_size); + } + + inline void cut(const T *const R__ src, T *const R__ dst) const { + v_multiply(dst, src, m_cache, m_size); + } + + inline void add(T *const R__ dst, T scale) const { + v_add_with_gain(dst, m_cache, m_size, scale); + } + + inline T getArea() const { return m_area; } + inline T getValue(int i) const { return m_cache[i]; } + + inline int getSize() const { return m_size; } + inline int getP() const { return m_p; } + +protected: + int m_size; + int m_p; + T *R__ m_cache; + T m_area; + + void encache(); +}; + +template +void SincWindow::encache() +{ + const int n = m_size; + T *mult = allocate(n); + v_set(mult, T(1.0), n); + int i; + + for (i = 0; i < n; ++i) { + T extent = T(n)/2.; + T arg = (T(i) - extent) * (2. * M_PI) / m_p; + if (arg != 0.) { + mult[i] *= sin(arg) / arg; + } + } + + m_cache = mult; + + m_area = 0; + for (i = 0; i < n; ++i) { + std::cout << i << ":" << m_cache[i] << " "; + m_area += m_cache[i]; + } + std::cout << std::endl; + m_area /= n; +} + +} + +#endif diff --git a/src/dsp/Window.h b/src/dsp/Window.h index eb8313f..7cb8ae2 100644 --- a/src/dsp/Window.h +++ b/src/dsp/Window.h @@ -22,6 +22,7 @@ #include "system/sysutils.h" #include "system/VectorOps.h" +#include "system/Allocators.h" namespace RubberBand { @@ -86,10 +87,10 @@ protected: template void Window::encache() { - int n = int(m_size); - T *mult = new T[n]; + const int n = m_size; + T *mult = allocate(n); + v_set(mult, T(1.0), n); int i; - for (i = 0; i < n; ++i) mult[i] = 1.0; switch (m_type) { diff --git a/src/system/VectorOps.h b/src/system/VectorOps.h index 3a1d10c..d0cff53 100644 --- a/src/system/VectorOps.h +++ b/src/system/VectorOps.h @@ -51,6 +51,16 @@ inline void v_zero_channels(T *const R__ *const R__ ptr, } } +template +inline void v_set(T *const R__ ptr, + const T value, + const int count) +{ + for (int i = 0; i < count; ++i) { + ptr[i] = value; + } +} + template inline void v_copy(T *const R__ dst, const T *const R__ src,