c@342: /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
c@337: 
c@337: #include "dsp/phasevocoder/PhaseVocoder.h"
c@337: 
c@337: #include "base/Window.h"
c@337: 
c@342: #include <iostream>
c@342: 
c@342: using std::cerr;
c@342: using std::endl;
c@342: 
c@337: #define BOOST_TEST_DYN_LINK
c@337: #define BOOST_TEST_MAIN
c@337: 
c@337: #include <boost/test/unit_test.hpp>
c@337: 
c@337: BOOST_AUTO_TEST_SUITE(TestFFT)
c@337: 
c@337: #define COMPARE_CONST(a, n) \
c@337:     for (int cmp_i = 0; cmp_i < (int)(sizeof(a)/sizeof(a[0])); ++cmp_i) { \
c@342:         BOOST_CHECK_SMALL(a[cmp_i] - n, 1e-7);				\
c@337:     }
c@337: 
c@337: #define COMPARE_ARRAY(a, b)						\
c@337:     for (int cmp_i = 0; cmp_i < (int)(sizeof(a)/sizeof(a[0])); ++cmp_i) { \
c@342:         BOOST_CHECK_SMALL(a[cmp_i] - b[cmp_i], 1e-7);			\
c@337:     }
c@337: 
c@337: #define COMPARE_ARRAY_EXACT(a, b)						\
c@337:     for (int cmp_i = 0; cmp_i < (int)(sizeof(a)/sizeof(a[0])); ++cmp_i) { \
c@337:         BOOST_CHECK_EQUAL(a[cmp_i], b[cmp_i]);			\
c@337:     }
c@337: 
c@337: BOOST_AUTO_TEST_CASE(fullcycle)
c@337: {
c@340:     // Cosine with one cycle exactly equal to pvoc hopsize. This is
c@340:     // pretty much the most trivial case -- in fact it's
c@340:     // indistinguishable from totally silent input (in the phase
c@340:     // values) because the measured phases are zero throughout.
c@340: 
c@340:     // We aren't windowing the input frame because (for once) it
c@340:     // actually *is* just a short part of a continuous infinite
c@340:     // sinusoid.
c@337: 
c@337:     double frame[] = { 1, 0, -1, 0, 1, 0, -1, 0 };
c@337: 
c@340:     PhaseVocoder pvoc(8, 4);
c@337: 
c@337:     // Make these arrays one element too long at each end, so as to
c@337:     // test for overruns. For frame size 8, we expect 8/2+1 = 5
c@337:     // mag/phase pairs.
c@337:     double mag[]   = { 999, 999, 999, 999, 999, 999, 999 };
c@337:     double phase[] = { 999, 999, 999, 999, 999, 999, 999 };
c@340:     double unw[]   = { 999, 999, 999, 999, 999, 999, 999 };
c@337: 
c@340:     pvoc.process(frame, mag + 1, phase + 1, unw + 1);
c@337: 
c@337:     double magExpected0[] = { 999, 0, 0, 4, 0, 0, 999 };
c@337:     COMPARE_ARRAY_EXACT(mag, magExpected0);
c@337: 
c@337:     double phaseExpected0[] = { 999, 0, 0, 0, 0, 0, 999 };
c@337:     COMPARE_ARRAY_EXACT(phase, phaseExpected0);
c@337: 
c@340:     double unwExpected0[] = { 999, 0, 0, 0, 0, 0, 999 };
c@340:     COMPARE_ARRAY(unw, unwExpected0);
c@340: 
c@340:     pvoc.process(frame, mag + 1, phase + 1, unw + 1);
c@337: 
c@337:     double magExpected1[] = { 999, 0, 0, 4, 0, 0, 999 };
c@337:     COMPARE_ARRAY_EXACT(mag, magExpected1);
c@337: 
c@340:     double phaseExpected1[] = { 999, 0, 0, 0, 0, 0, 999 };
c@337:     COMPARE_ARRAY(phase, phaseExpected1);
c@338: 
c@342:     // Derivation of unwrapped values:
c@340:     // 
c@340:     // * Bin 0 (DC) always has phase 0 and expected phase 0
c@340:     //
c@340:     // * Bin 1 has expected phase pi (the hop size is half a cycle at
c@340:     //   its frequency), but measured phase 0 (because there is no
c@340:     //   signal in that bin). So it has phase error -pi, which is
c@340:     //   mapped into (-pi,pi] range as pi, giving an unwrapped phase
c@340:     //   of 2*pi.
c@340:     //  
c@340:     // * Bin 2 has expected unwrapped phase 2*pi, measured phase 0,
c@340:     //   hence error 0 and unwrapped phase 2*pi.
c@340:     //
c@340:     // * Bin 3 is like bin 1: it has expected phase 3*pi, measured
c@340:     //   phase 0, so phase error -pi and unwrapped phase 4*pi.
c@340:     //
c@340:     // * Bin 4 (Nyquist) is like bin 2: expected phase 4*pi, measured
c@340:     //   phase 0, hence error 0 and unwrapped phase 4*pi.
c@340: 
c@340:     double unwExpected1[] = { 999, 0, 2*M_PI, 2*M_PI, 4*M_PI, 4*M_PI, 999 };
c@340:     COMPARE_ARRAY(unw, unwExpected1);
c@340: 
c@340:     pvoc.process(frame, mag + 1, phase + 1, unw + 1);
c@338: 
c@338:     double magExpected2[] = { 999, 0, 0, 4, 0, 0, 999 };
c@338:     COMPARE_ARRAY_EXACT(mag, magExpected2);
c@338: 
c@340:     double phaseExpected2[] = { 999, 0, 0, 0, 0, 0, 999 };
c@338:     COMPARE_ARRAY(phase, phaseExpected2);
c@340: 
c@340:     double unwExpected2[] = { 999, 0, 4*M_PI, 4*M_PI, 8*M_PI, 8*M_PI, 999 };
c@340:     COMPARE_ARRAY(unw, unwExpected2);
c@337: }
c@337: 
c@342: BOOST_AUTO_TEST_CASE(overlapping)
c@342: {
c@342:     // Sine (i.e. cosine starting at phase -pi/2) starting with the
c@342:     // first sample, introducing a cosine of half the frequency
c@342:     // starting at the fourth sample, i.e. the second hop. The cosine
c@342:     // is introduced "by magic", i.e. it doesn't appear in the second
c@342:     // half of the first frame (it would have quite strange effects on
c@342:     // the first frame if it did).
c@340: 
c@342:     double data[32] = { // 3 x 8-sample frames which we pretend are overlapping
c@342: 	0, 1, 0, -1, 0, 1, 0, -1,
c@342: 	1, 1.70710678, 0, -1.70710678, -1, 0.29289322, 0, -0.29289322,
c@342: 	-1, 0.29289322, 0, -0.29289322, 1, 1.70710678, 0, -1.70710678,
c@342:     };
c@342: 
c@342:     PhaseVocoder pvoc(8, 4);
c@342: 
c@342:     // Make these arrays one element too long at each end, so as to
c@342:     // test for overruns. For frame size 8, we expect 8/2+1 = 5
c@342:     // mag/phase pairs.
c@342:     double mag[]   = { 999, 999, 999, 999, 999, 999, 999 };
c@342:     double phase[] = { 999, 999, 999, 999, 999, 999, 999 };
c@342:     double unw[]   = { 999, 999, 999, 999, 999, 999, 999 };
c@342: 
c@342: cerr << "process 0" << endl;
c@342: 
c@342:     pvoc.process(data, mag + 1, phase + 1, unw + 1);
c@342: 
c@342:     double magExpected0[] = { 999, 0, 0, 4, 0, 0, 999 };
c@342:     COMPARE_ARRAY(mag, magExpected0);
c@342: 
c@342:     double phaseExpected0[] = { 999, 0, 0, -M_PI/2 , 0, 0, 999 };
c@342:     COMPARE_ARRAY(phase, phaseExpected0);
c@342: 
c@342:     double unwExpected0[] = { 999, 0, 0, -M_PI/2, 0, 0, 999 };
c@342:     COMPARE_ARRAY(unw, unwExpected0);
c@342: 
c@342: cerr << "process 1" << endl;
c@342: 
c@342:     pvoc.process(data + 8, mag + 1, phase + 1, unw + 1);
c@342: 
c@342:     double magExpected1[] = { 999, 0, 4, 4, 0, 0, 999 };
c@342:     COMPARE_ARRAY(mag, magExpected1);
c@342: 
c@342:     //!!! I don't know why [2] here is -M_PI and not M_PI; and I definitely don't know why [4] here is M_PI. Check these with care
c@342:     double phaseExpected1[] = { 999, 0, -M_PI, -M_PI/2, M_PI, 0, 999 };
c@342:     COMPARE_ARRAY(phase, phaseExpected1);
c@342: 
c@342:     double unwExpected1[] = { 999, 0, M_PI, 3*M_PI/2, 3*M_PI, 4*M_PI, 999 };
c@342:     COMPARE_ARRAY(unw, unwExpected1);
c@342: 
c@342: cerr << "process 2" << endl;
c@342: 
c@342:     pvoc.process(data + 16, mag + 1, phase + 1, unw + 1);
c@342: 
c@342:     double magExpected2[] = { 999, 0, 4, 4, 0, 0, 999 };
c@342:     COMPARE_ARRAY(mag, magExpected2);
c@342: 
c@342:     double phaseExpected2[] = { 999, 0, 0, -M_PI/2, 0, 0, 999 };
c@342:     COMPARE_ARRAY(phase, phaseExpected2);
c@342: 
c@342:     double unwExpected2[] = { 999, 0, 2*M_PI, 7*M_PI/2, 6*M_PI, 8*M_PI, 999 };
c@342:     COMPARE_ARRAY(unw, unwExpected2);
c@342: }
c@340: 
c@337: BOOST_AUTO_TEST_SUITE_END()
c@337: