libs/vamp-sdk/src/vamp-hostsdk/PluginInputDomainAdapter.cpp

   1 /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
   2
   3 /*
   4     Vamp
   5
   6     An API for audio analysis and feature extraction plugins.
   7
   8     Centre for Digital Music, Queen Mary, University of London.
   9     Copyright 2006-2009 Chris Cannam and QMUL.
  10
  11     This file is based in part on Don Cross's public domain FFT
  12     implementation.
  13
  14     Permission is hereby granted, free of charge, to any person
  15     obtaining a copy of this software and associated documentation
  16     files (the "Software"), to deal in the Software without
  17     restriction, including without limitation the rights to use, copy,
  18     modify, merge, publish, distribute, sublicense, and/or sell copies
  19     of the Software, and to permit persons to whom the Software is
  20     furnished to do so, subject to the following conditions:
  21
  22     The above copyright notice and this permission notice shall be
  23     included in all copies or substantial portions of the Software.
  24
  25     THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
  26     EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
  27     MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
  28     NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS BE LIABLE FOR
  29     ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF
  30     CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
  31     WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
  32
  33     Except as contained in this notice, the names of the Centre for
  34     Digital Music; Queen Mary, University of London; and Chris Cannam
  35     shall not be used in advertising or otherwise to promote the sale,
  36     use or other dealings in this Software without prior written
  37     authorization.
  38 */
  39
  40 #include "vamp-hostsdk/PluginInputDomainAdapter.h"
  41
  42 #include <cmath>
  43
  44
  45 /**
  46  * If you want to compile using FFTW instead of the built-in FFT
  47  * implementation for the PluginInputDomainAdapter, define HAVE_FFTW3
  48  * in the Makefile.
  49  *
  50  * Be aware that FFTW is licensed under the GPL -- unlike this SDK,
  51  * which is provided under a more liberal BSD license in order to
  52  * permit use in closed source applications.  The use of FFTW would
  53  * mean that your code would need to be licensed under the GPL as
  54  * well.  Do not define this symbol unless you understand and accept
  55  * the implications of this.
  56  *
  57  * Parties such as Linux distribution packagers who redistribute this
  58  * SDK for use in other programs should _not_ define this symbol, as
  59  * it would change the effective licensing terms under which the SDK
  60  * was available to third party developers.
  61  *
  62  * The default is not to use FFTW, and to use the built-in FFT instead.
  63  *
  64  * Note: The FFTW code uses FFTW_MEASURE, and so will perform badly on
  65  * its first invocation unless the host has saved and restored FFTW
  66  * wisdom (see the FFTW documentation).
  67  */
  68 #ifdef HAVE_FFTW3
  69 #include <fftw3.h>
  70 #endif
  71
  72
  73 _VAMP_SDK_HOSTSPACE_BEGIN(PluginInputDomainAdapter.cpp)
  74
  75 namespace Vamp {
  76
  77 namespace HostExt {
  78
  79 class PluginInputDomainAdapter::Impl
  80 {
  81 public:
  82     Impl(Plugin *plugin, float inputSampleRate);
  83     ~Impl();
  84
  85     bool initialise(size_t channels, size_t stepSize, size_t blockSize);
  86     void reset();
  87
  88     size_t getPreferredStepSize() const;
  89     size_t getPreferredBlockSize() const;
  90
  91     FeatureSet process(const float *const *inputBuffers, RealTime timestamp);
  92
  93     void setProcessTimestampMethod(ProcessTimestampMethod m);
  94     ProcessTimestampMethod getProcessTimestampMethod() const;
  95
  96     RealTime getTimestampAdjustment() const;
  97
  98 protected:
  99     Plugin *m_plugin;
 100     float m_inputSampleRate;
 101     int m_channels;
 102     int m_stepSize;
 103     int m_blockSize;
 104     float **m_freqbuf;
 105
 106     double *m_ri;
 107     double *m_window;
 108
 109     ProcessTimestampMethod m_method;
 110     int m_processCount;
 111     float **m_shiftBuffers;
 112
 113 #ifdef HAVE_FFTW3
 114     fftw_plan m_plan;
 115     fftw_complex *m_cbuf;
 116 #else
 117     double *m_ro;
 118     double *m_io;
 119     void fft(unsigned int n, bool inverse,
 120              double *ri, double *ii, double *ro, double *io);
 121 #endif
 122
 123     FeatureSet processShiftingTimestamp(const float *const *inputBuffers, RealTime timestamp);
 124     FeatureSet processShiftingData(const float *const *inputBuffers, RealTime timestamp);
 125
 126     size_t makeBlockSizeAcceptable(size_t) const;
 127 };
 128
 129 PluginInputDomainAdapter::PluginInputDomainAdapter(Plugin *plugin) :
 130     PluginWrapper(plugin)
 131 {
 132     m_impl = new Impl(plugin, m_inputSampleRate);
 133 }
 134
 135 PluginInputDomainAdapter::~PluginInputDomainAdapter()
 136 {
 137     delete m_impl;
 138 }
 139
 140 bool
 141 PluginInputDomainAdapter::initialise(size_t channels, size_t stepSize, size_t blockSize)
 142 {
 143     return m_impl->initialise(channels, stepSize, blockSize);
 144 }
 145
 146 void
 147 PluginInputDomainAdapter::reset()
 148 {
 149     m_impl->reset();
 150 }
 151
 152 Plugin::InputDomain
 153 PluginInputDomainAdapter::getInputDomain() const
 154 {
 155     return TimeDomain;
 156 }
 157
 158 size_t
 159 PluginInputDomainAdapter::getPreferredStepSize() const
 160 {
 161     return m_impl->getPreferredStepSize();
 162 }
 163
 164 size_t
 165 PluginInputDomainAdapter::getPreferredBlockSize() const
 166 {
 167     return m_impl->getPreferredBlockSize();
 168 }
 169
 170 Plugin::FeatureSet
 171 PluginInputDomainAdapter::process(const float *const *inputBuffers, RealTime timestamp)
 172 {
 173     return m_impl->process(inputBuffers, timestamp);
 174 }
 175
 176 void
 177 PluginInputDomainAdapter::setProcessTimestampMethod(ProcessTimestampMethod m)
 178 {
 179     m_impl->setProcessTimestampMethod(m);
 180 }
 181
 182 PluginInputDomainAdapter::ProcessTimestampMethod
 183 PluginInputDomainAdapter::getProcessTimestampMethod() const
 184 {
 185     return m_impl->getProcessTimestampMethod();
 186 }
 187
 188 RealTime
 189 PluginInputDomainAdapter::getTimestampAdjustment() const
 190 {
 191     return m_impl->getTimestampAdjustment();
 192 }
 193
 194
 195 PluginInputDomainAdapter::Impl::Impl(Plugin *plugin, float inputSampleRate) :
 196     m_plugin(plugin),
 197     m_inputSampleRate(inputSampleRate),
 198     m_channels(0),
 199     m_stepSize(0),
 200     m_blockSize(0),
 201     m_freqbuf(0),
 202     m_ri(0),
 203     m_window(0),
 204     m_method(ShiftTimestamp),
 205     m_processCount(0),
 206     m_shiftBuffers(0),
 207 #ifdef HAVE_FFTW3
 208     m_plan(0),
 209     m_cbuf(0)
 210 #else
 211     m_ro(0),
 212     m_io(0)
 213 #endif
 214 {
 215 }
 216
 217 PluginInputDomainAdapter::Impl::~Impl()
 218 {
 219     // the adapter will delete the plugin
 220
 221     if (m_shiftBuffers) {
 222         for (int c = 0; c < m_channels; ++c) {
 223             delete[] m_shiftBuffers[c];
 224         }
 225         delete[] m_shiftBuffers;
 226     }
 227
 228     if (m_channels > 0) {
 229         for (int c = 0; c < m_channels; ++c) {
 230             delete[] m_freqbuf[c];
 231         }
 232         delete[] m_freqbuf;
 233 #ifdef HAVE_FFTW3
 234         if (m_plan) {
 235             fftw_destroy_plan(m_plan);
 236             fftw_free(m_ri);
 237             fftw_free(m_cbuf);
 238             m_plan = 0;
 239         }
 240 #else
 241         delete[] m_ri;
 242         delete[] m_ro;
 243         delete[] m_io;
 244 #endif
 245         delete[] m_window;
 246     }
 247 }
 248
 249 // for some visual studii apparently
 250 #ifndef M_PI
 251 #define M_PI 3.14159265358979232846
 252 #endif
 253
 254 bool
 255 PluginInputDomainAdapter::Impl::initialise(size_t channels, size_t stepSize, size_t blockSize)
 256 {
 257     if (m_plugin->getInputDomain() == TimeDomain) {
 258
 259         m_stepSize = int(stepSize);
 260         m_blockSize = int(blockSize);
 261         m_channels = int(channels);
 262
 263         return m_plugin->initialise(channels, stepSize, blockSize);
 264     }
 265
 266     if (blockSize < 2) {
 267         std::cerr << "ERROR: PluginInputDomainAdapter::initialise: blocksize < 2 not supported" << std::endl;
 268         return false;
 269     }
 270
 271     if (blockSize & (blockSize-1)) {
 272         std::cerr << "ERROR: PluginInputDomainAdapter::initialise: non-power-of-two\nblocksize " << blockSize << " not supported" << std::endl;
 273         return false;
 274     }
 275
 276     if (m_channels > 0) {
 277         for (int c = 0; c < m_channels; ++c) {
 278             delete[] m_freqbuf[c];
 279         }
 280         delete[] m_freqbuf;
 281 #ifdef HAVE_FFTW3
 282         if (m_plan) {
 283             fftw_destroy_plan(m_plan);
 284             fftw_free(m_ri);
 285             fftw_free(m_cbuf);
 286             m_plan = 0;
 287         }
 288 #else
 289         delete[] m_ri;
 290         delete[] m_ro;
 291         delete[] m_io;
 292 #endif
 293         delete[] m_window;
 294     }
 295
 296     m_stepSize = int(stepSize);
 297     m_blockSize = int(blockSize);
 298     m_channels = int(channels);
 299
 300     m_freqbuf = new float *[m_channels];
 301     for (int c = 0; c < m_channels; ++c) {
 302         m_freqbuf[c] = new float[m_blockSize + 2];
 303     }
 304     m_window = new double[m_blockSize];
 305
 306     for (int i = 0; i < m_blockSize; ++i) {
 307         // Hanning window
 308         m_window[i] = (0.50 - 0.50 * cos((2.0 * M_PI * i) / m_blockSize));
 309     }
 310
 311 #ifdef HAVE_FFTW3
 312     m_ri = (double *)fftw_malloc(blockSize * sizeof(double));
 313     m_cbuf = (fftw_complex *)fftw_malloc((blockSize/2 + 1) * sizeof(fftw_complex));
 314     m_plan = fftw_plan_dft_r2c_1d(blockSize, m_ri, m_cbuf, FFTW_MEASURE);
 315 #else
 316     m_ri = new double[m_blockSize];
 317     m_ro = new double[m_blockSize];
 318     m_io = new double[m_blockSize];
 319 #endif
 320
 321     m_processCount = 0;
 322
 323     return m_plugin->initialise(channels, stepSize, blockSize);
 324 }
 325
 326 void
 327 PluginInputDomainAdapter::Impl::reset()
 328 {
 329     m_processCount = 0;
 330     m_plugin->reset();
 331 }
 332
 333 size_t
 334 PluginInputDomainAdapter::Impl::getPreferredStepSize() const
 335 {
 336     size_t step = m_plugin->getPreferredStepSize();
 337
 338     if (step == 0 && (m_plugin->getInputDomain() == FrequencyDomain)) {
 339         step = getPreferredBlockSize() / 2;
 340     }
 341
 342     return step;
 343 }
 344
 345 size_t
 346 PluginInputDomainAdapter::Impl::getPreferredBlockSize() const
 347 {
 348     size_t block = m_plugin->getPreferredBlockSize();
 349
 350     if (m_plugin->getInputDomain() == FrequencyDomain) {
 351         if (block == 0) {
 352             block = 1024;
 353         } else {
 354             block = makeBlockSizeAcceptable(block);
 355         }
 356     }
 357
 358     return block;
 359 }
 360
 361 size_t
 362 PluginInputDomainAdapter::Impl::makeBlockSizeAcceptable(size_t blockSize) const
 363 {
 364     if (blockSize < 2) {
 365
 366         std::cerr << "WARNING: PluginInputDomainAdapter::initialise: blocksize < 2 not" << std::endl
 367                   << "supported, increasing from " << blockSize << " to 2" << std::endl;
 368         blockSize = 2;
 369
 370     } else if (blockSize & (blockSize-1)) {
 371
 372 #ifdef HAVE_FFTW3
 373         // not an issue with FFTW
 374 #else
 375
 376         // not a power of two, can't handle that with our built-in FFT
 377         // implementation
 378
 379         size_t nearest = blockSize;
 380         size_t power = 0;
 381         while (nearest > 1) {
 382             nearest >>= 1;
 383             ++power;
 384         }
 385         nearest = 1;
 386         while (power) {
 387             nearest <<= 1;
 388             --power;
 389         }
 390
 391         if (blockSize - nearest > (nearest*2) - blockSize) {
 392             nearest = nearest*2;
 393         }
 394
 395         std::cerr << "WARNING: PluginInputDomainAdapter::initialise: non-power-of-two\nblocksize " << blockSize << " not supported, using blocksize " << nearest << " instead" << std::endl;
 396         blockSize = nearest;
 397
 398 #endif
 399     }
 400
 401     return blockSize;
 402 }
 403
 404 RealTime
 405 PluginInputDomainAdapter::Impl::getTimestampAdjustment() const
 406 {
 407     if (m_plugin->getInputDomain() == TimeDomain) {
 408         return RealTime::zeroTime;
 409     } else if (m_method == ShiftData || m_method == NoShift) {
 410         return RealTime::zeroTime;
 411     } else {
 412         return RealTime::frame2RealTime
 413             (m_blockSize/2, int(m_inputSampleRate + 0.5));
 414     }
 415 }
 416
 417 void
 418 PluginInputDomainAdapter::Impl::setProcessTimestampMethod(ProcessTimestampMethod m)
 419 {
 420     m_method = m;
 421 }
 422
 423 PluginInputDomainAdapter::ProcessTimestampMethod
 424 PluginInputDomainAdapter::Impl::getProcessTimestampMethod() const
 425 {
 426     return m_method;
 427 }
 428
 429 Plugin::FeatureSet
 430 PluginInputDomainAdapter::Impl::process(const float *const *inputBuffers,
 431                                         RealTime timestamp)
 432 {
 433     if (m_plugin->getInputDomain() == TimeDomain) {
 434         return m_plugin->process(inputBuffers, timestamp);
 435     }
 436
 437     if (m_method == ShiftTimestamp || m_method == NoShift) {
 438         return processShiftingTimestamp(inputBuffers, timestamp);
 439     } else {
 440         return processShiftingData(inputBuffers, timestamp);
 441     }
 442 }
 443
 444 Plugin::FeatureSet
 445 PluginInputDomainAdapter::Impl::processShiftingTimestamp(const float *const *inputBuffers,
 446                                                          RealTime timestamp)
 447 {
 448     if (m_method == ShiftTimestamp) {
 449         timestamp = timestamp + getTimestampAdjustment();
 450     }
 451
 452     for (int c = 0; c < m_channels; ++c) {
 453
 454         for (int i = 0; i < m_blockSize; ++i) {
 455             m_ri[i] = double(inputBuffers[c][i]) * m_window[i];
 456         }
 457
 458         for (int i = 0; i < m_blockSize/2; ++i) {
 459             // FFT shift
 460             double value = m_ri[i];
 461             m_ri[i] = m_ri[i + m_blockSize/2];
 462             m_ri[i + m_blockSize/2] = value;
 463         }
 464
 465 #ifdef HAVE_FFTW3
 466         fftw_execute(m_plan);
 467
 468         for (int i = 0; i <= m_blockSize/2; ++i) {
 469             m_freqbuf[c][i * 2] = float(m_cbuf[i][0]);
 470             m_freqbuf[c][i * 2 + 1] = float(m_cbuf[i][1]);
 471         }
 472 #else
 473         fft(m_blockSize, false, m_ri, 0, m_ro, m_io);
 474
 475         for (int i = 0; i <= m_blockSize/2; ++i) {
 476             m_freqbuf[c][i * 2] = float(m_ro[i]);
 477             m_freqbuf[c][i * 2 + 1] = float(m_io[i]);
 478         }
 479 #endif
 480     }
 481
 482     return m_plugin->process(m_freqbuf, timestamp);
 483 }
 484
 485 Plugin::FeatureSet
 486 PluginInputDomainAdapter::Impl::processShiftingData(const float *const *inputBuffers,
 487                                                     RealTime timestamp)
 488 {
 489     if (m_processCount == 0) {
 490         if (!m_shiftBuffers) {
 491             m_shiftBuffers = new float *[m_channels];
 492             for (int c = 0; c < m_channels; ++c) {
 493                 m_shiftBuffers[c] = new float[m_blockSize + m_blockSize/2];
 494             }
 495         }
 496         for (int c = 0; c < m_channels; ++c) {
 497             for (int i = 0; i < m_blockSize + m_blockSize/2; ++i) {
 498                 m_shiftBuffers[c][i] = 0.f;
 499             }
 500         }
 501     }
 502
 503     for (int c = 0; c < m_channels; ++c) {
 504         for (int i = m_stepSize; i < m_blockSize + m_blockSize/2; ++i) {
 505             m_shiftBuffers[c][i - m_stepSize] = m_shiftBuffers[c][i];
 506         }
 507         for (int i = 0; i < m_blockSize; ++i) {
 508             m_shiftBuffers[c][i + m_blockSize/2] = inputBuffers[c][i];
 509         }
 510     }
 511
 512     for (int c = 0; c < m_channels; ++c) {
 513
 514         for (int i = 0; i < m_blockSize; ++i) {
 515             m_ri[i] = double(m_shiftBuffers[c][i]) * m_window[i];
 516         }
 517
 518         for (int i = 0; i < m_blockSize/2; ++i) {
 519             // FFT shift
 520             double value = m_ri[i];
 521             m_ri[i] = m_ri[i + m_blockSize/2];
 522             m_ri[i + m_blockSize/2] = value;
 523         }
 524
 525 #ifdef HAVE_FFTW3
 526         fftw_execute(m_plan);
 527
 528         for (int i = 0; i <= m_blockSize/2; ++i) {
 529             m_freqbuf[c][i * 2] = float(m_cbuf[i][0]);
 530             m_freqbuf[c][i * 2 + 1] = float(m_cbuf[i][1]);
 531         }
 532 #else
 533         fft(m_blockSize, false, m_ri, 0, m_ro, m_io);
 534
 535         for (int i = 0; i <= m_blockSize/2; ++i) {
 536             m_freqbuf[c][i * 2] = float(m_ro[i]);
 537             m_freqbuf[c][i * 2 + 1] = float(m_io[i]);
 538         }
 539 #endif
 540     }
 541
 542     ++m_processCount;
 543
 544     return m_plugin->process(m_freqbuf, timestamp);
 545 }
 546
 547 #ifndef HAVE_FFTW3
 548
 549 void
 550 PluginInputDomainAdapter::Impl::fft(unsigned int n, bool inverse,
 551                                     double *ri, double *ii, double *ro, double *io)
 552 {
 553     if (!ri || !ro || !io) return;
 554
 555     unsigned int bits;
 556     unsigned int i, j, k, m;
 557     unsigned int blockSize, blockEnd;
 558
 559     double tr, ti;
 560
 561     if (n < 2) return;
 562     if (n & (n-1)) return;
 563
 564     double angle = 2.0 * M_PI;
 565     if (inverse) angle = -angle;
 566
 567     for (i = 0; ; ++i) {
 568         if (n & (1 << i)) {
 569             bits = i;
 570             break;
 571         }
 572     }
 573
 574     static unsigned int tableSize = 0;
 575     static int *table = 0;
 576
 577     if (tableSize != n) {
 578
 579         delete[] table;
 580
 581         table = new int[n];
 582
 583         for (i = 0; i < n; ++i) {
 584
 585             m = i;
 586
 587             for (j = k = 0; j < bits; ++j) {
 588                 k = (k << 1) | (m & 1);
 589                 m >>= 1;
 590             }
 591
 592             table[i] = k;
 593         }
 594
 595         tableSize = n;
 596     }
 597
 598     if (ii) {
 599         for (i = 0; i < n; ++i) {
 600             ro[table[i]] = ri[i];
 601             io[table[i]] = ii[i];
 602         }
 603     } else {
 604         for (i = 0; i < n; ++i) {
 605             ro[table[i]] = ri[i];
 606             io[table[i]] = 0.0;
 607         }
 608     }
 609
 610     blockEnd = 1;
 611
 612     for (blockSize = 2; blockSize <= n; blockSize <<= 1) {
 613
 614         double delta = angle / (double)blockSize;
 615         double sm2 = -sin(-2 * delta);
 616         double sm1 = -sin(-delta);
 617         double cm2 = cos(-2 * delta);
 618         double cm1 = cos(-delta);
 619         double w = 2 * cm1;
 620         double ar[3], ai[3];
 621
 622         for (i = 0; i < n; i += blockSize) {
 623
 624             ar[2] = cm2;
 625             ar[1] = cm1;
 626
 627             ai[2] = sm2;
 628             ai[1] = sm1;
 629
 630             for (j = i, m = 0; m < blockEnd; j++, m++) {
 631
 632                 ar[0] = w * ar[1] - ar[2];
 633                 ar[2] = ar[1];
 634                 ar[1] = ar[0];
 635
 636                 ai[0] = w * ai[1] - ai[2];
 637                 ai[2] = ai[1];
 638                 ai[1] = ai[0];
 639
 640                 k = j + blockEnd;
 641                 tr = ar[0] * ro[k] - ai[0] * io[k];
 642                 ti = ar[0] * io[k] + ai[0] * ro[k];
 643
 644                 ro[k] = ro[j] - tr;
 645                 io[k] = io[j] - ti;
 646
 647                 ro[j] += tr;
 648                 io[j] += ti;
 649             }
 650         }
 651
 652         blockEnd = blockSize;
 653     }
 654
 655     if (inverse) {
 656
 657         double denom = (double)n;
 658
 659         for (i = 0; i < n; i++) {
 660             ro[i] /= denom;
 661             io[i] /= denom;
 662         }
 663     }
 664 }
 665
 666 #endif
 667
 668 }
 669
 670 }
 671
 672 _VAMP_SDK_HOSTSPACE_END(PluginInputDomainAdapter.cpp)
 673