PluginInputDomainAdapter.cpp

Go to the documentation of this file.
00001 /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
00002 
00003 /*
00004     Vamp
00005 
00006     An API for audio analysis and feature extraction plugins.
00007 
00008     Centre for Digital Music, Queen Mary, University of London.
00009     Copyright 2006-2007 Chris Cannam and QMUL.
00010   
00011     This file is based in part on Don Cross's public domain FFT
00012     implementation.
00013 
00014     Permission is hereby granted, free of charge, to any person
00015     obtaining a copy of this software and associated documentation
00016     files (the "Software"), to deal in the Software without
00017     restriction, including without limitation the rights to use, copy,
00018     modify, merge, publish, distribute, sublicense, and/or sell copies
00019     of the Software, and to permit persons to whom the Software is
00020     furnished to do so, subject to the following conditions:
00021 
00022     The above copyright notice and this permission notice shall be
00023     included in all copies or substantial portions of the Software.
00024 
00025     THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
00026     EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
00027     MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
00028     NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS BE LIABLE FOR
00029     ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF
00030     CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
00031     WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
00032 
00033     Except as contained in this notice, the names of the Centre for
00034     Digital Music; Queen Mary, University of London; and Chris Cannam
00035     shall not be used in advertising or otherwise to promote the sale,
00036     use or other dealings in this Software without prior written
00037     authorization.
00038 */
00039 
00040 #include "PluginInputDomainAdapter.h"
00041 
00042 #include <cmath>
00043 
00044 
00068 #ifdef HAVE_FFTW3
00069 #include <fftw3.h>
00070 #endif
00071 
00072 
00073 namespace Vamp {
00074 
00075 namespace HostExt {
00076 
00077 class PluginInputDomainAdapter::Impl
00078 {
00079 public:
00080     Impl(Plugin *plugin, float inputSampleRate);
00081     ~Impl();
00082     
00083     bool initialise(size_t channels, size_t stepSize, size_t blockSize);
00084 
00085     size_t getPreferredStepSize() const;
00086     size_t getPreferredBlockSize() const;
00087 
00088     FeatureSet process(const float *const *inputBuffers, RealTime timestamp);
00089 
00090 protected:
00091     Plugin *m_plugin;
00092     float m_inputSampleRate;
00093     int m_channels;
00094     int m_blockSize;
00095     float **m_freqbuf;
00096 
00097     double *m_ri;
00098     double *m_window;
00099 
00100 #ifdef HAVE_FFTW3
00101     fftw_plan m_plan;
00102     fftw_complex *m_cbuf;
00103 #else
00104     double *m_ro;
00105     double *m_io;
00106     void fft(unsigned int n, bool inverse,
00107              double *ri, double *ii, double *ro, double *io);
00108 #endif
00109 
00110     size_t makeBlockSizeAcceptable(size_t) const;
00111 };
00112 
00113 PluginInputDomainAdapter::PluginInputDomainAdapter(Plugin *plugin) :
00114     PluginWrapper(plugin)
00115 {
00116     m_impl = new Impl(plugin, m_inputSampleRate);
00117 }
00118 
00119 PluginInputDomainAdapter::~PluginInputDomainAdapter()
00120 {
00121     delete m_impl;
00122 }
00123   
00124 bool
00125 PluginInputDomainAdapter::initialise(size_t channels, size_t stepSize, size_t blockSize)
00126 {
00127     return m_impl->initialise(channels, stepSize, blockSize);
00128 }
00129 
00130 Plugin::InputDomain
00131 PluginInputDomainAdapter::getInputDomain() const
00132 {
00133     return TimeDomain;
00134 }
00135 
00136 size_t
00137 PluginInputDomainAdapter::getPreferredStepSize() const
00138 {
00139     return m_impl->getPreferredStepSize();
00140 }
00141 
00142 size_t
00143 PluginInputDomainAdapter::getPreferredBlockSize() const
00144 {
00145     return m_impl->getPreferredBlockSize();
00146 }
00147 
00148 Plugin::FeatureSet
00149 PluginInputDomainAdapter::process(const float *const *inputBuffers, RealTime timestamp)
00150 {
00151     return m_impl->process(inputBuffers, timestamp);
00152 }
00153 
00154 PluginInputDomainAdapter::Impl::Impl(Plugin *plugin, float inputSampleRate) :
00155     m_plugin(plugin),
00156     m_inputSampleRate(inputSampleRate),
00157     m_channels(0),
00158     m_blockSize(0),
00159     m_freqbuf(0),
00160     m_ri(0),
00161     m_window(0),
00162 #ifdef HAVE_FFTW3
00163     m_plan(0),
00164     m_cbuf(0)
00165 #else
00166     m_ro(0),
00167     m_io(0)
00168 #endif
00169 {
00170 }
00171 
00172 PluginInputDomainAdapter::Impl::~Impl()
00173 {
00174     // the adapter will delete the plugin
00175 
00176     if (m_channels > 0) {
00177         for (int c = 0; c < m_channels; ++c) {
00178             delete[] m_freqbuf[c];
00179         }
00180         delete[] m_freqbuf;
00181 #ifdef HAVE_FFTW3
00182         if (m_plan) {
00183             fftw_destroy_plan(m_plan);
00184             fftw_free(m_ri);
00185             fftw_free(m_cbuf);
00186             m_plan = 0;
00187         }
00188 #else
00189         delete[] m_ri;
00190         delete[] m_ro;
00191         delete[] m_io;
00192 #endif
00193         delete[] m_window;
00194     }
00195 }
00196 
00197 // for some visual studii apparently
00198 #ifndef M_PI
00199 #define M_PI 3.14159265358979232846
00200 #endif
00201     
00202 bool
00203 PluginInputDomainAdapter::Impl::initialise(size_t channels, size_t stepSize, size_t blockSize)
00204 {
00205     if (m_plugin->getInputDomain() == TimeDomain) {
00206 
00207         m_blockSize = int(blockSize);
00208         m_channels = int(channels);
00209 
00210         return m_plugin->initialise(channels, stepSize, blockSize);
00211     }
00212 
00213     if (blockSize < 2) {
00214         std::cerr << "ERROR: Vamp::HostExt::PluginInputDomainAdapter::Impl::initialise: blocksize < 2 not supported" << std::endl;
00215         return false;
00216     }                
00217         
00218     if (blockSize & (blockSize-1)) {
00219         std::cerr << "ERROR: Vamp::HostExt::PluginInputDomainAdapter::Impl::initialise: non-power-of-two\nblocksize " << blockSize << " not supported" << std::endl;
00220         return false;
00221     }
00222 
00223     if (m_channels > 0) {
00224         for (int c = 0; c < m_channels; ++c) {
00225             delete[] m_freqbuf[c];
00226         }
00227         delete[] m_freqbuf;
00228 #ifdef HAVE_FFTW3
00229         if (m_plan) {
00230             fftw_destroy_plan(m_plan);
00231             fftw_free(m_ri);
00232             fftw_free(m_cbuf);
00233             m_plan = 0;
00234         }
00235 #else
00236         delete[] m_ri;
00237         delete[] m_ro;
00238         delete[] m_io;
00239 #endif
00240         delete[] m_window;
00241     }
00242 
00243     m_blockSize = int(blockSize);
00244     m_channels = int(channels);
00245 
00246     m_freqbuf = new float *[m_channels];
00247     for (int c = 0; c < m_channels; ++c) {
00248         m_freqbuf[c] = new float[m_blockSize + 2];
00249     }
00250     m_window = new double[m_blockSize];
00251 
00252     for (int i = 0; i < m_blockSize; ++i) {
00253         // Hanning window
00254         m_window[i] = (0.50 - 0.50 * cos((2.0 * M_PI * i) / m_blockSize));
00255     }
00256 
00257 #ifdef HAVE_FFTW3
00258     m_ri = (double *)fftw_malloc(blockSize * sizeof(double));
00259     m_cbuf = (fftw_complex *)fftw_malloc((blockSize/2 + 1) * sizeof(fftw_complex));
00260     m_plan = fftw_plan_dft_r2c_1d(blockSize, m_ri, m_cbuf, FFTW_MEASURE);
00261 #else
00262     m_ri = new double[m_blockSize];
00263     m_ro = new double[m_blockSize];
00264     m_io = new double[m_blockSize];
00265 #endif
00266 
00267     return m_plugin->initialise(channels, stepSize, blockSize);
00268 }
00269 
00270 size_t
00271 PluginInputDomainAdapter::Impl::getPreferredStepSize() const
00272 {
00273     size_t step = m_plugin->getPreferredStepSize();
00274 
00275     if (step == 0 && (m_plugin->getInputDomain() == FrequencyDomain)) {
00276         step = getPreferredBlockSize() / 2;
00277     }
00278 
00279     return step;
00280 }
00281 
00282 size_t
00283 PluginInputDomainAdapter::Impl::getPreferredBlockSize() const
00284 {
00285     size_t block = m_plugin->getPreferredBlockSize();
00286 
00287     if (m_plugin->getInputDomain() == FrequencyDomain) {
00288         if (block == 0) {
00289             block = 1024;
00290         } else {
00291             block = makeBlockSizeAcceptable(block);
00292         }
00293     }
00294 
00295     return block;
00296 }
00297 
00298 size_t
00299 PluginInputDomainAdapter::Impl::makeBlockSizeAcceptable(size_t blockSize) const
00300 {
00301     if (blockSize < 2) {
00302 
00303         std::cerr << "WARNING: Vamp::HostExt::PluginInputDomainAdapter::Impl::initialise: blocksize < 2 not" << std::endl
00304                   << "supported, increasing from " << blockSize << " to 2" << std::endl;
00305         blockSize = 2;
00306         
00307     } else if (blockSize & (blockSize-1)) {
00308             
00309 #ifdef HAVE_FFTW3
00310         // not an issue with FFTW
00311 #else
00312 
00313         // not a power of two, can't handle that with our built-in FFT
00314         // implementation
00315 
00316         size_t nearest = blockSize;
00317         size_t power = 0;
00318         while (nearest > 1) {
00319             nearest >>= 1;
00320             ++power;
00321         }
00322         nearest = 1;
00323         while (power) {
00324             nearest <<= 1;
00325             --power;
00326         }
00327         
00328         if (blockSize - nearest > (nearest*2) - blockSize) {
00329             nearest = nearest*2;
00330         }
00331         
00332         std::cerr << "WARNING: Vamp::HostExt::PluginInputDomainAdapter::Impl::initialise: non-power-of-two\nblocksize " << blockSize << " not supported, using blocksize " << nearest << " instead" << std::endl;
00333         blockSize = nearest;
00334 
00335 #endif
00336     }
00337 
00338     return blockSize;
00339 }
00340 
00341 Plugin::FeatureSet
00342 PluginInputDomainAdapter::Impl::process(const float *const *inputBuffers,
00343                                         RealTime timestamp)
00344 {
00345     if (m_plugin->getInputDomain() == TimeDomain) {
00346         return m_plugin->process(inputBuffers, timestamp);
00347     }
00348 
00349     // The timestamp supplied should be (according to the Vamp::Plugin
00350     // spec) the time of the start of the time-domain input block.
00351     // However, we want to pass to the plugin an FFT output calculated
00352     // from the block of samples _centred_ on that timestamp.
00353     // 
00354     // We have two options:
00355     // 
00356     // 1. Buffer the input, calculating the fft of the values at the
00357     // passed-in block minus blockSize/2 rather than starting at the
00358     // passed-in block.  So each time we call process on the plugin,
00359     // we are passing in the same timestamp as was passed to our own
00360     // process plugin, but not (the frequency domain representation
00361     // of) the same set of samples.  Advantages: avoids confusion in
00362     // the host by ensuring the returned values have timestamps
00363     // comparable with that passed in to this function (in fact this
00364     // is pretty much essential for one-value-per-block outputs);
00365     // consistent with hosts such as SV that deal with the
00366     // frequency-domain transform themselves.  Disadvantages: means
00367     // making the not necessarily correct assumption that the samples
00368     // preceding the first official block are all zero (or some other
00369     // known value).
00370     //
00371     // 2. Increase the passed-in timestamps by half the blocksize.  So
00372     // when we call process, we are passing in the frequency domain
00373     // representation of the same set of samples as passed to us, but
00374     // with a different timestamp.  Advantages: simplicity; avoids
00375     // iffy assumption mentioned above.  Disadvantages: inconsistency
00376     // with SV in cases where stepSize != blockSize/2; potential
00377     // confusion arising from returned timestamps being calculated
00378     // from the adjusted input timestamps rather than the original
00379     // ones (and inaccuracy where the returned timestamp is implied,
00380     // as in one-value-per-block).
00381     //
00382     // Neither way is ideal, but I don't think either is strictly
00383     // incorrect either.  I think this is just a case where the same
00384     // plugin can legitimately produce differing results from the same
00385     // input data, depending on how that data is packaged.
00386     // 
00387     // We'll go for option 2, adjusting the timestamps.  Note in
00388     // particular that this means some results can differ from those
00389     // produced by SV.
00390 
00391 //    std::cerr << "PluginInputDomainAdapter: sampleRate " << m_inputSampleRate << ", blocksize " << m_blockSize << ", adjusting time from " << timestamp;
00392 
00393     timestamp = timestamp + RealTime::frame2RealTime
00394         (m_blockSize/2, int(m_inputSampleRate + 0.5));
00395 
00396 //    std::cerr << " to " << timestamp << std::endl;
00397 
00398     for (int c = 0; c < m_channels; ++c) {
00399 
00400         for (int i = 0; i < m_blockSize; ++i) {
00401             m_ri[i] = double(inputBuffers[c][i]) * m_window[i];
00402         }
00403 
00404         for (int i = 0; i < m_blockSize/2; ++i) {
00405             // FFT shift
00406             double value = m_ri[i];
00407             m_ri[i] = m_ri[i + m_blockSize/2];
00408             m_ri[i + m_blockSize/2] = value;
00409         }
00410 
00411 #ifdef HAVE_FFTW3
00412 
00413         fftw_execute(m_plan);
00414 
00415         for (int i = 0; i <= m_blockSize/2; ++i) {
00416             m_freqbuf[c][i * 2] = float(m_cbuf[i][0]);
00417             m_freqbuf[c][i * 2 + 1] = float(m_cbuf[i][1]);
00418         }
00419 
00420 #else
00421 
00422         fft(m_blockSize, false, m_ri, 0, m_ro, m_io);
00423 
00424         for (int i = 0; i <= m_blockSize/2; ++i) {
00425             m_freqbuf[c][i * 2] = float(m_ro[i]);
00426             m_freqbuf[c][i * 2 + 1] = float(m_io[i]);
00427         }
00428 
00429 #endif
00430     }
00431 
00432     return m_plugin->process(m_freqbuf, timestamp);
00433 }
00434 
00435 #ifndef HAVE_FFTW3
00436 
00437 void
00438 PluginInputDomainAdapter::Impl::fft(unsigned int n, bool inverse,
00439                                     double *ri, double *ii, double *ro, double *io)
00440 {
00441     if (!ri || !ro || !io) return;
00442 
00443     unsigned int bits;
00444     unsigned int i, j, k, m;
00445     unsigned int blockSize, blockEnd;
00446 
00447     double tr, ti;
00448 
00449     if (n < 2) return;
00450     if (n & (n-1)) return;
00451 
00452     double angle = 2.0 * M_PI;
00453     if (inverse) angle = -angle;
00454 
00455     for (i = 0; ; ++i) {
00456         if (n & (1 << i)) {
00457             bits = i;
00458             break;
00459         }
00460     }
00461 
00462     static unsigned int tableSize = 0;
00463     static int *table = 0;
00464 
00465     if (tableSize != n) {
00466 
00467         delete[] table;
00468 
00469         table = new int[n];
00470 
00471         for (i = 0; i < n; ++i) {
00472         
00473             m = i;
00474 
00475             for (j = k = 0; j < bits; ++j) {
00476                 k = (k << 1) | (m & 1);
00477                 m >>= 1;
00478             }
00479 
00480             table[i] = k;
00481         }
00482 
00483         tableSize = n;
00484     }
00485 
00486     if (ii) {
00487         for (i = 0; i < n; ++i) {
00488             ro[table[i]] = ri[i];
00489             io[table[i]] = ii[i];
00490         }
00491     } else {
00492         for (i = 0; i < n; ++i) {
00493             ro[table[i]] = ri[i];
00494             io[table[i]] = 0.0;
00495         }
00496     }
00497 
00498     blockEnd = 1;
00499 
00500     for (blockSize = 2; blockSize <= n; blockSize <<= 1) {
00501 
00502         double delta = angle / (double)blockSize;
00503         double sm2 = -sin(-2 * delta);
00504         double sm1 = -sin(-delta);
00505         double cm2 = cos(-2 * delta);
00506         double cm1 = cos(-delta);
00507         double w = 2 * cm1;
00508         double ar[3], ai[3];
00509 
00510         for (i = 0; i < n; i += blockSize) {
00511 
00512             ar[2] = cm2;
00513             ar[1] = cm1;
00514 
00515             ai[2] = sm2;
00516             ai[1] = sm1;
00517 
00518             for (j = i, m = 0; m < blockEnd; j++, m++) {
00519 
00520                 ar[0] = w * ar[1] - ar[2];
00521                 ar[2] = ar[1];
00522                 ar[1] = ar[0];
00523 
00524                 ai[0] = w * ai[1] - ai[2];
00525                 ai[2] = ai[1];
00526                 ai[1] = ai[0];
00527 
00528                 k = j + blockEnd;
00529                 tr = ar[0] * ro[k] - ai[0] * io[k];
00530                 ti = ar[0] * io[k] + ai[0] * ro[k];
00531 
00532                 ro[k] = ro[j] - tr;
00533                 io[k] = io[j] - ti;
00534 
00535                 ro[j] += tr;
00536                 io[j] += ti;
00537             }
00538         }
00539 
00540         blockEnd = blockSize;
00541     }
00542 
00543     if (inverse) {
00544 
00545         double denom = (double)n;
00546 
00547         for (i = 0; i < n; i++) {
00548             ro[i] /= denom;
00549             io[i] /= denom;
00550         }
00551     }
00552 }
00553 
00554 #endif
00555 
00556 }
00557         
00558 }
00559 

Generated on Fri Nov 7 13:10:34 2008 for VampPluginSDK by  doxygen 1.5.6