summaryrefslogtreecommitdiff
path: root/libs/qm-dsp/dsp/phasevocoder
diff options
context:
space:
mode:
Diffstat (limited to 'libs/qm-dsp/dsp/phasevocoder')
-rw-r--r--libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.cpp116
-rw-r--r--libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.h62
2 files changed, 136 insertions, 42 deletions
diff --git a/libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.cpp b/libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.cpp
index cd52f5240a..7ee0be2c96 100644
--- a/libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.cpp
+++ b/libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.cpp
@@ -4,7 +4,7 @@
QM DSP Library
Centre for Digital Music, Queen Mary, University of London.
- This file 2005-2006 Christian Landone.
+ This file 2005-2006 Christian Landone, copyright 2013 QMUL.
This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License as
@@ -15,30 +15,47 @@
#include "PhaseVocoder.h"
#include "dsp/transforms/FFT.h"
+#include "maths/MathUtilities.h"
#include <math.h>
-//////////////////////////////////////////////////////////////////////
-// Construction/Destruction
-//////////////////////////////////////////////////////////////////////
+#include <cassert>
-PhaseVocoder::PhaseVocoder(unsigned int n) :
- m_n(n)
+#include <iostream>
+using std::cerr;
+using std::endl;
+
+PhaseVocoder::PhaseVocoder(int n, int hop) :
+ m_n(n),
+ m_hop(hop)
{
m_fft = new FFTReal(m_n);
- m_realOut = new double[m_n];
- m_imagOut = new double[m_n];
+ m_time = new double[m_n];
+ m_real = new double[m_n];
+ m_imag = new double[m_n];
+ m_phase = new double[m_n/2 + 1];
+ m_unwrapped = new double[m_n/2 + 1];
+
+ for (int i = 0; i < m_n/2 + 1; ++i) {
+ m_phase[i] = 0.0;
+ m_unwrapped[i] = 0.0;
+ }
+
+ reset();
}
PhaseVocoder::~PhaseVocoder()
{
- delete [] m_realOut;
- delete [] m_imagOut;
+ delete[] m_unwrapped;
+ delete[] m_phase;
+ delete[] m_real;
+ delete[] m_imag;
+ delete[] m_time;
delete m_fft;
}
-void PhaseVocoder::FFTShift(unsigned int size, double *src)
+void PhaseVocoder::FFTShift(double *src)
{
- const int hs = size/2;
+ const int hs = m_n/2;
for (int i = 0; i < hs; ++i) {
double tmp = src[i];
src[i] = src[i + hs];
@@ -46,34 +63,73 @@ void PhaseVocoder::FFTShift(unsigned int size, double *src)
}
}
-void PhaseVocoder::process(double *src, double *mag, double *theta)
+void PhaseVocoder::processTimeDomain(const double *src,
+ double *mag, double *theta,
+ double *unwrapped)
{
- FFTShift( m_n, src);
-
- m_fft->process(0, src, m_realOut, m_imagOut);
+ for (int i = 0; i < m_n; ++i) {
+ m_time[i] = src[i];
+ }
+ FFTShift(m_time);
+ m_fft->forward(m_time, m_real, m_imag);
+ getMagnitudes(mag);
+ getPhases(theta);
+ unwrapPhases(theta, unwrapped);
+}
- getMagnitude( m_n/2, mag, m_realOut, m_imagOut);
- getPhase( m_n/2, theta, m_realOut, m_imagOut);
+void PhaseVocoder::processFrequencyDomain(const double *reals,
+ const double *imags,
+ double *mag, double *theta,
+ double *unwrapped)
+{
+ for (int i = 0; i < m_n/2 + 1; ++i) {
+ m_real[i] = reals[i];
+ m_imag[i] = imags[i];
+ }
+ getMagnitudes(mag);
+ getPhases(theta);
+ unwrapPhases(theta, unwrapped);
}
-void PhaseVocoder::getMagnitude(unsigned int size, double *mag, double *real, double *imag)
+void PhaseVocoder::reset()
{
- unsigned int j;
+ for (int i = 0; i < m_n/2 + 1; ++i) {
+ // m_phase stores the "previous" phase, so set to one step
+ // behind so that a signal with initial phase at zero matches
+ // the expected values. This is completely unnecessary for any
+ // analytical purpose, it's just tidier.
+ double omega = (2 * M_PI * m_hop * i) / m_n;
+ m_phase[i] = -omega;
+ m_unwrapped[i] = -omega;
+ }
+}
- for( j = 0; j < size; j++)
- {
- mag[ j ] = sqrt( real[ j ] * real[ j ] + imag[ j ] * imag[ j ]);
+void PhaseVocoder::getMagnitudes(double *mag)
+{
+ for (int i = 0; i < m_n/2 + 1; i++) {
+ mag[i] = sqrt(m_real[i] * m_real[i] + m_imag[i] * m_imag[i]);
}
}
-void PhaseVocoder::getPhase(unsigned int size, double *theta, double *real, double *imag)
+void PhaseVocoder::getPhases(double *theta)
+{
+ for (int i = 0; i < m_n/2 + 1; i++) {
+ theta[i] = atan2(m_imag[i], m_real[i]);
+ }
+}
+
+void PhaseVocoder::unwrapPhases(double *theta, double *unwrapped)
{
- unsigned int k;
+ for (int i = 0; i < m_n/2 + 1; ++i) {
+
+ double omega = (2 * M_PI * m_hop * i) / m_n;
+ double expected = m_phase[i] + omega;
+ double error = MathUtilities::princarg(theta[i] - expected);
- // Phase Angle "matlab" style
- //Watch out for quadrant mapping !!!
- for( k = 0; k < size; k++)
- {
- theta[ k ] = atan2( -imag[ k ], real[ k ]);
+ unwrapped[i] = m_unwrapped[i] + omega + error;
+
+ m_phase[i] = theta[i];
+ m_unwrapped[i] = unwrapped[i];
}
}
+
diff --git a/libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.h b/libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.h
index 8a76335d41..20242175f0 100644
--- a/libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.h
+++ b/libs/qm-dsp/dsp/phasevocoder/PhaseVocoder.h
@@ -4,7 +4,7 @@
QM DSP Library
Centre for Digital Music, Queen Mary, University of London.
- This file 2005-2006 Christian Landone.
+ This file 2005-2006 Christian Landone, copyright 2013 QMUL.
This program is free software; you can redistribute it and/or
modify it under the terms of the GNU General Public License as
@@ -18,25 +18,63 @@
class FFTReal;
-class PhaseVocoder
+class PhaseVocoder
{
public:
- PhaseVocoder( unsigned int size );
+ PhaseVocoder(int size, int hop);
virtual ~PhaseVocoder();
- void process( double* src, double* mag, double* theta);
+ /**
+ * Given one frame of time-domain samples, FFT and return the
+ * magnitudes, instantaneous phases, and unwrapped phases.
+ *
+ * src must have size values (where size is the frame size value
+ * as passed to the PhaseVocoder constructor), and should have
+ * been windowed as necessary by the caller (but not fft-shifted).
+ *
+ * mag, phase, and unwrapped must each be non-NULL and point to
+ * enough space for size/2 + 1 values. The redundant conjugate
+ * half of the output is not returned.
+ */
+ void processTimeDomain(const double *src,
+ double *mag, double *phase, double *unwrapped);
+
+ /**
+ * Given one frame of frequency-domain samples, return the
+ * magnitudes, instantaneous phases, and unwrapped phases.
+ *
+ * reals and imags must each contain size/2+1 values (where size
+ * is the frame size value as passed to the PhaseVocoder
+ * constructor).
+ *
+ * mag, phase, and unwrapped must each be non-NULL and point to
+ * enough space for size/2+1 values.
+ */
+ void processFrequencyDomain(const double *reals, const double *imags,
+ double *mag, double *phase, double *unwrapped);
+
+ /**
+ * Reset the stored phases to zero. Note that this may be
+ * necessary occasionally (depending on the application) to avoid
+ * loss of floating-point precision in the accumulated unwrapped
+ * phase values as they grow.
+ */
+ void reset();
protected:
- void getPhase(unsigned int size, double *theta, double *real, double *imag);
-// void coreFFT( unsigned int NumSamples, double *RealIn, double* ImagIn, double *RealOut, double *ImagOut);
- void getMagnitude( unsigned int size, double* mag, double* real, double* imag);
- void FFTShift( unsigned int size, double* src);
+ void FFTShift(double *src);
+ void getMagnitudes(double *mag);
+ void getPhases(double *theta);
+ void unwrapPhases(double *theta, double *unwrapped);
- unsigned int m_n;
+ int m_n;
+ int m_hop;
FFTReal *m_fft;
- double *m_imagOut;
- double *m_realOut;
-
+ double *m_time;
+ double *m_imag;
+ double *m_real;
+ double *m_phase;
+ double *m_unwrapped;
};
#endif