cannam@1: /* -*- c-basic-offset: 4 indent-tabs-mode: nil -*-  vi:set ts=8 sts=4 sw=4: */
cannam@1: 
cannam@1: /*
cannam@1:     Vamp
cannam@1: 
cannam@1:     An API for audio analysis and feature extraction plugins.
cannam@1: 
cannam@1:     Centre for Digital Music, Queen Mary, University of London.
cannam@1:     Copyright 2006 Chris Cannam.
cannam@16:     FFT code from Don Cross's public domain FFT implementation.
cannam@1:   
cannam@1:     Permission is hereby granted, free of charge, to any person
cannam@1:     obtaining a copy of this software and associated documentation
cannam@1:     files (the "Software"), to deal in the Software without
cannam@1:     restriction, including without limitation the rights to use, copy,
cannam@1:     modify, merge, publish, distribute, sublicense, and/or sell copies
cannam@1:     of the Software, and to permit persons to whom the Software is
cannam@1:     furnished to do so, subject to the following conditions:
cannam@1: 
cannam@1:     The above copyright notice and this permission notice shall be
cannam@1:     included in all copies or substantial portions of the Software.
cannam@1: 
cannam@1:     THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
cannam@1:     EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
cannam@1:     MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
cannam@6:     NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS BE LIABLE FOR
cannam@1:     ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF
cannam@1:     CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
cannam@1:     WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
cannam@1: 
cannam@1:     Except as contained in this notice, the names of the Centre for
cannam@1:     Digital Music; Queen Mary, University of London; and Chris Cannam
cannam@1:     shall not be used in advertising or otherwise to promote the sale,
cannam@1:     use or other dealings in this Software without prior written
cannam@1:     authorization.
cannam@1: */
cannam@1: 
cannam@16: #include "PluginHostAdapter.h"
cannam@1: #include "vamp.h"
cannam@1: 
cannam@16: #include <iostream>
cannam@16: #include <sndfile.h>
cannam@40: #include <dirent.h> // POSIX directory open and read
cannam@1: 
cannam@1: #include "system.h"
cannam@1: 
cannam@19: #include <cmath>
cannam@19: 
cannam@16: using std::cout;
cannam@16: using std::cerr;
cannam@16: using std::endl;
cannam@16: using std::string;
cannam@32: using std::vector;
cannam@16: 
cannam@43: #define HOST_VERSION "1.0"
cannam@40: 
cannam@16: void printFeatures(int, int, int, Vamp::Plugin::FeatureSet);
cannam@16: void transformInput(float *, size_t);
cannam@16: void fft(unsigned int, bool, double *, double *, double *, double *);
cannam@40: void printPluginPath();
cannam@40: 
cannam@40: #ifdef HAVE_OPENDIR
cannam@40: void enumeratePlugins();
cannam@40: #endif
cannam@16: 
cannam@1: /*
cannam@16:     A very simple Vamp plugin host.  Given the name of a plugin
cannam@16:     library and the name of a sound file on the command line, it loads
cannam@16:     the first plugin in the library and runs it on the sound file,
cannam@16:     dumping the plugin's first output to stdout.
cannam@1: */
cannam@1: 
cannam@1: int main(int argc, char **argv)
cannam@1: {
cannam@16:     if (argc < 2 || argc > 4) {
cannam@40:         char *scooter = argv[0];
cannam@40:         char *name = 0;
cannam@40:         while (scooter && *scooter) {
cannam@40:             if (*scooter == '/' || *scooter == '\\') name = ++scooter;
cannam@40:             else ++scooter;
cannam@40:         }
cannam@40:         if (!name || !*name) name = argv[0];
cannam@40:         cerr << "\n"
cannam@40:              << name << ": A simple Vamp plugin host.\n\n"
cannam@40:             "Centre for Digital Music, Queen Mary, University of London.\n"
cannam@40:             "Copyright 2006 Chris Cannam and QMUL.\n"
cannam@40:             "Freely redistributable; published under a BSD-style license.\n\n"
cannam@40:             "Usage:\n\n"
cannam@40:             "  " << name << " pluginlibrary." << PLUGIN_SUFFIX << "\n\n"
cannam@40:             "    -- Load \"pluginlibrary\" and list the Vamp plugins it contains.\n\n"
cannam@40:             "  " << name << " pluginlibrary." << PLUGIN_SUFFIX << ":plugin file.wav [outputno]\n\n"
cannam@40:             "    -- Load plugin id \"plugin\" from \"pluginlibrary\" and run it on the\n"
cannam@40:             "       audio data in \"file.wav\", dumping the output from \"outputno\"\n"
cannam@40:             "       (default 0) to standard output.\n\n"
cannam@40: #ifdef HAVE_OPENDIR
cannam@40:             "  " << name << " -l\n\n"
cannam@40:             "    -- List the plugin libraries and Vamp plugins in the plugin search path.\n\n"
cannam@40: #endif
cannam@40:             "  " << name << " -p\n\n"
cannam@40:             "    -- Print out the Vamp plugin search path.\n\n"
cannam@43: 	    "  " << name << " -v\n\n"
cannam@43: 	    "    -- Display version information only.\n\n"
cannam@43: 	    "Note that this host does not use the plugin search path when loadinga plugin.\nIf a plugin library is specified, it should be with a full file path.\n"
cannam@40:              << endl;
cannam@1:         return 2;
cannam@1:     }
cannam@43:     
cannam@43:     if (argc == 2 && !strcmp(argv[1], "-v")) {
cannam@43: 	cout << "Simple Vamp plugin host version: " << HOST_VERSION << endl
cannam@43: 	     << "Vamp API version: " << VAMP_API_VERSION << endl
cannam@43: 	     << "Vamp SDK version: " << VAMP_SDK_VERSION << endl;
cannam@43: 	return 0;
cannam@43:     }
cannam@43:     
cannam@40:     if (argc == 2 && !strcmp(argv[1], "-l")) {
cannam@40: #ifdef HAVE_OPENDIR
cannam@40:         enumeratePlugins();
cannam@40: #endif
cannam@40:         return 0;
cannam@40:     }
cannam@40:     if (argc == 2 && !strcmp(argv[1], "-p")) {
cannam@40:         printPluginPath();
cannam@40:         return 0;
cannam@40:     }
cannam@40: 
cannam@16:     cerr << endl << argv[0] << ": Running..." << endl;
cannam@1: 
cannam@16:     string soname = argv[1];
cannam@49:     string plugid = "";
cannam@16:     string wavname;
cannam@16:     if (argc >= 3) wavname = argv[2];
cannam@16: 
cannam@20:     int sep = soname.find(":");
cannam@40:     if (sep >= 0 && sep < int(soname.length())) {
cannam@49:         plugid = soname.substr(sep + 1);
cannam@20:         soname = soname.substr(0, sep);
cannam@16:     }
cannam@1: 
cannam@1:     void *libraryHandle = DLOPEN(soname, RTLD_LAZY);
cannam@1: 
cannam@1:     if (!libraryHandle) {
cannam@16:         cerr << argv[0] << ": Failed to open plugin library " 
cannam@16:                   << soname << ": " << DLERROR() << endl;
cannam@1:         return 1;
cannam@1:     }
cannam@1: 
cannam@16:     cerr << argv[0] << ": Opened plugin library " << soname << endl;
cannam@1: 
cannam@1:     VampGetPluginDescriptorFunction fn = (VampGetPluginDescriptorFunction)
cannam@1:         DLSYM(libraryHandle, "vampGetPluginDescriptor");
cannam@1:     
cannam@1:     if (!fn) {
cannam@16:         cerr << argv[0] << ": No Vamp descriptor function in library "
cannam@16:                   << soname << endl;
cannam@1:         DLCLOSE(libraryHandle);
cannam@1:         return 1;
cannam@1:     }
cannam@1: 
cannam@16:     cerr << argv[0] << ": Found plugin descriptor function" << endl;
cannam@1: 
cannam@1:     int index = 0;
cannam@16:     int plugnumber = -1;
cannam@1:     const VampPluginDescriptor *descriptor = 0;
cannam@1: 
cannam@50:     while ((descriptor = fn(VAMP_API_VERSION, index))) {
cannam@1: 
cannam@16:         Vamp::PluginHostAdapter plugin(descriptor, 48000);
cannam@16:         cerr << argv[0] << ": Plugin " << (index+1)
cannam@49:                   << " is \"" << plugin.getIdentifier() << "\"" << endl;
cannam@16: 
cannam@49:         if (plugin.getIdentifier() == plugid) plugnumber = index;
cannam@1:         
cannam@1:         ++index;
cannam@1:     }
cannam@1: 
cannam@16:     cerr << argv[0] << ": Done\n" << endl;
cannam@16: 
cannam@16:     if (wavname == "") {
cannam@16:         DLCLOSE(libraryHandle);
cannam@16:         return 0;
cannam@16:     }
cannam@16: 
cannam@16:     if (plugnumber < 0) {
cannam@49:         if (plugid != "") {
cannam@49:             cerr << "ERROR: No such plugin as " << plugid << " in library"
cannam@16:                  << endl;
cannam@16:             DLCLOSE(libraryHandle);
cannam@16:             return 0;
cannam@16:         } else {
cannam@16:             plugnumber = 0;
cannam@16:         }
cannam@16:     }
cannam@16: 
cannam@50:     descriptor = fn(VAMP_API_VERSION, plugnumber);
cannam@16:     if (!descriptor) {
cannam@16:         DLCLOSE(libraryHandle);
cannam@16:         return 0;
cannam@16:     }
cannam@16:     
cannam@16:     SNDFILE *sndfile;
cannam@16:     SF_INFO sfinfo;
cannam@16:     memset(&sfinfo, 0, sizeof(SF_INFO));
cannam@16: 
cannam@16:     sndfile = sf_open(wavname.c_str(), SFM_READ, &sfinfo);
cannam@16:     if (!sndfile) {
cannam@16: 	cerr << "ERROR: Failed to open input file \"" << wavname << "\": "
cannam@16: 	     << sf_strerror(sndfile) << endl;
cannam@16:         DLCLOSE(libraryHandle);
cannam@16: 	return 1;
cannam@16:     }
cannam@16: 
cannam@16:     Vamp::PluginHostAdapter *plugin =
cannam@16:         new Vamp::PluginHostAdapter(descriptor, sfinfo.samplerate);
cannam@16: 
cannam@49:     cerr << "Running " << plugin->getIdentifier() << "..." << endl;
cannam@16: 
cannam@16:     int blockSize = plugin->getPreferredBlockSize();
cannam@16:     int stepSize = plugin->getPreferredStepSize();
cannam@16: 
cannam@16:     cerr << "Preferred block size = " << blockSize << ", step size = "
cannam@29:          << stepSize << endl;
cannam@16: 
cannam@16:     if (blockSize == 0) blockSize = 1024;
cannam@16: 
cannam@29:     bool rightBlockSize = true;
cannam@42: 
cannam@29:     if (plugin->getInputDomain() == Vamp::Plugin::FrequencyDomain) {
cannam@42: 
cannam@29:         int p = 1, b = blockSize;
cannam@29:         while (b) {
cannam@29:             p <<= 1;
cannam@29:             b >>= 1;
cannam@29:         }
cannam@29:         if (p != blockSize * 2) {
cannam@29:             cerr << "WARNING: Plugin requested non-power-of-two block size of "
cannam@29:                  << blockSize << ",\nwhich is not supported by this host.  ";
cannam@29:             blockSize = p;
cannam@29:             cerr << "Rounding up to " << blockSize << "." << endl;
cannam@29:             rightBlockSize = false;
cannam@29:         }
cannam@42:         if (stepSize == 0) stepSize = blockSize / 2;
cannam@42: 
cannam@42:     } else {
cannam@42: 
cannam@42:         if (stepSize == 0) stepSize = blockSize;
cannam@29:     }
cannam@29: 
cannam@16:     int channels = sfinfo.channels;
cannam@16: 
cannam@16:     float *filebuf = new float[blockSize * channels];
cannam@16:     float **plugbuf = new float*[channels];
cannam@47:     for (int c = 0; c < channels; ++c) plugbuf[c] = new float[blockSize + 2];
cannam@16: 
cannam@16:     cerr << "Using block size = " << blockSize << ", step size = "
cannam@16:               << stepSize << endl;
cannam@16: 
cannam@16:     int minch = plugin->getMinChannelCount();
cannam@16:     int maxch = plugin->getMaxChannelCount();
cannam@16:     cerr << "Plugin accepts " << minch << " -> " << maxch << " channel(s)" << endl;
cannam@16: 
cannam@16:     Vamp::Plugin::OutputList outputs = plugin->getOutputDescriptors();
cannam@16:     Vamp::Plugin::OutputDescriptor od;
cannam@16: 
cannam@29:     int returnValue = 1;
cannam@29: 
cannam@16:     int output = 0;
cannam@16:     if (argc == 4) output = atoi(argv[3]);
cannam@16: 
cannam@16:     bool mix = false;
cannam@16: 
cannam@16:     if (minch > channels || maxch < channels) {
cannam@16:         if (minch == 1) {
cannam@16:             cerr << "WARNING: Sound file has " << channels << " channels, mixing down to 1" << endl;
cannam@16:             mix = true;
cannam@16:             channels = 1;
cannam@16:         } else {
cannam@16:             cerr << "ERROR: Sound file has " << channels << " channels, out of range for plugin" << endl;
cannam@16:             goto done;
cannam@16:         }
cannam@16:     }
cannam@16: 
cannam@16:     if (outputs.empty()) {
cannam@16: 	cerr << "Plugin has no outputs!" << endl;
cannam@16:         goto done;
cannam@16:     }
cannam@16: 
cannam@16:     if (int(outputs.size()) <= output) {
cannam@16: 	cerr << "Output " << output << " requested, but plugin has only " << outputs.size() << " output(s)" << endl;
cannam@16:         goto done;
cannam@16:     }        
cannam@16: 
cannam@16:     od = outputs[output];
cannam@49:     cerr << "Output is " << od.identifier << endl;
cannam@16: 
cannam@29:     if (!plugin->initialise(channels, stepSize, blockSize)) {
cannam@29:         cerr << "ERROR: Plugin initialise (channels = " << channels
cannam@29:              << ", stepSize = " << stepSize << ", blockSize = "
cannam@29:              << blockSize << ") failed." << endl;
cannam@29:         if (!rightBlockSize) {
cannam@29:             cerr << "(Probably because I couldn't provide the plugin's preferred block size.)" << endl;
cannam@29:         }
cannam@29:         goto done;
cannam@29:     }
cannam@16: 
cannam@16:     for (size_t i = 0; i < sfinfo.frames; i += stepSize) {
cannam@16: 
cannam@16:         int count;
cannam@16: 
cannam@16:         if (sf_seek(sndfile, i, SEEK_SET) < 0) {
cannam@16:             cerr << "ERROR: sf_seek failed: " << sf_strerror(sndfile) << endl;
cannam@16:             break;
cannam@16:         }
cannam@16:         
cannam@16:         if ((count = sf_readf_float(sndfile, filebuf, blockSize)) < 0) {
cannam@16:             cerr << "ERROR: sf_readf_float failed: " << sf_strerror(sndfile) << endl;
cannam@16:             break;
cannam@16:         }
cannam@16: 
cannam@16:         for (int c = 0; c < channels; ++c) {
cannam@16:             for (int j = 0; j < blockSize; ++j) {
cannam@16:                 plugbuf[c][j] = 0.0f;
cannam@16:             }
cannam@16:         }
cannam@16: 
cannam@52:         for (int j = 0; j < blockSize && j < count; ++j) {
cannam@52:             int tc = 0;
cannam@52:             for (int c = 0; c < sfinfo.channels; ++c) {
cannam@52:                 tc = c;
cannam@52:                 if (mix) tc = 0;
cannam@52:                 plugbuf[tc][j] += filebuf[j * sfinfo.channels + c];
cannam@16:             }
cannam@52:             if (mix) {
cannam@52:                 plugbuf[0][j] /= sfinfo.channels;
cannam@52:             }
cannam@52:         }
cannam@16: 
cannam@52:         if (plugin->getInputDomain() == Vamp::Plugin::FrequencyDomain) {
cannam@52:             for (int c = 0; c < sfinfo.channels; ++c) {
cannam@52:                 transformInput(plugbuf[c], blockSize);
cannam@52:                 if (mix) break;
cannam@16:             }
cannam@16:         }
cannam@16: 
cannam@16:         printFeatures
cannam@16:             (i, sfinfo.samplerate, output, plugin->process
cannam@16:              (plugbuf, Vamp::RealTime::frame2RealTime(i, sfinfo.samplerate)));
cannam@16:     }
cannam@16: 
cannam@16:     printFeatures(sfinfo.frames, sfinfo.samplerate, output,
cannam@16:                   plugin->getRemainingFeatures());
cannam@16: 
cannam@29:     returnValue = 0;
cannam@29: 
cannam@16: done:
cannam@16:     delete plugin;
cannam@1: 
cannam@1:     DLCLOSE(libraryHandle);
cannam@16:     sf_close(sndfile);
cannam@29:     return returnValue;
cannam@1: }
cannam@1: 
cannam@16: void
cannam@40: printPluginPath()
cannam@40: {
cannam@40:     vector<string> path = Vamp::PluginHostAdapter::getPluginPath();
cannam@40:     for (size_t i = 0; i < path.size(); ++i) {
cannam@40:         cerr << path[i] << endl;
cannam@40:     }
cannam@40: }
cannam@40: 
cannam@40: #ifdef HAVE_OPENDIR
cannam@40: 
cannam@40: void
cannam@40: enumeratePlugins()
cannam@40: {
cannam@40:     cerr << endl << "Vamp plugin libraries found in search path:" << endl;
cannam@40:     vector<string> path = Vamp::PluginHostAdapter::getPluginPath();
cannam@40:     for (size_t i = 0; i < path.size(); ++i) {
cannam@40:         cerr << "\n" << path[i] << ":" << endl;
cannam@40:         DIR *d = opendir(path[i].c_str());
cannam@40:         if (!d) {
cannam@40:             perror("Failed to open directory");
cannam@40:             continue;
cannam@40:         }
cannam@40:         struct dirent *e = 0;
cannam@40:         while ((e = readdir(d))) {
cannam@47: //            cerr << "reading: " << e->d_name << endl;
cannam@47:             if (!(e->d_type & DT_REG)) {
cannam@47: //                cerr << e->d_name << ": not a regular file" << endl;
cannam@47:                 continue;
cannam@47:             }
cannam@40:             int len = strlen(e->d_name);
cannam@40:             if (len < int(strlen(PLUGIN_SUFFIX) + 2) ||
cannam@40:                 e->d_name[len - strlen(PLUGIN_SUFFIX) - 1] != '.' ||
cannam@40:                 strcmp(e->d_name + len - strlen(PLUGIN_SUFFIX), PLUGIN_SUFFIX)) {
cannam@47: //                cerr << e->d_name << ": not a library file" << endl;
cannam@40:                 continue;
cannam@40:             }
cannam@40:             char *fp = new char[path[i].length() + len + 3];
cannam@40:             sprintf(fp, "%s/%s", path[i].c_str(), e->d_name);
cannam@40:             void *handle = DLOPEN(string(fp), RTLD_LAZY);
cannam@40:             if (handle) {
cannam@40:                 VampGetPluginDescriptorFunction fn =
cannam@40:                     (VampGetPluginDescriptorFunction)DLSYM
cannam@40:                     (handle, "vampGetPluginDescriptor");
cannam@40:                 if (fn) {
cannam@40:                     cerr << "\n  " << e->d_name << ":" << endl;
cannam@40:                     int index = 0;
cannam@40:                     const VampPluginDescriptor *descriptor = 0;
cannam@50:                     while ((descriptor = fn(VAMP_API_VERSION, index))) {
cannam@40:                         Vamp::PluginHostAdapter plugin(descriptor, 48000);
cannam@47:                         char c = char('A' + index);
cannam@47:                         if (c > 'Z') c = char('a' + (index - 26));
cannam@50:                         cerr << "    [" << c << "] [v"
cannam@50:                              << plugin.getVampApiVersion() << "] "
cannam@49:                              << plugin.getName()
cannam@49:                              << ", \"" << plugin.getIdentifier() << "\""
cannam@40:                              << " [" << plugin.getMaker()
cannam@49:                              << "]" << endl;
cannam@49:                         if (plugin.getDescription() != "") {
cannam@49:                             cerr << "        - " << plugin.getDescription() << endl;
cannam@49:                         }
cannam@40:                         Vamp::Plugin::OutputList outputs =
cannam@40:                             plugin.getOutputDescriptors();
cannam@40:                         if (outputs.size() > 1) {
cannam@40:                             for (size_t j = 0; j < outputs.size(); ++j) {
cannam@40:                                 cerr << "         (" << j << ") "
cannam@49:                                      << outputs[j].name
cannam@49:                                      << ", \"" << outputs[j].identifier << "\""
cannam@45:                                      << endl;
cannam@49:                                 if (outputs[j].description != "") {
cannam@49:                                     cerr << "             - " 
cannam@49:                                          << outputs[j].description << endl;
cannam@49:                                 }
cannam@40:                             }
cannam@40:                         }
cannam@40:                         ++index;
cannam@40:                     }
cannam@47:                 } else {
cannam@47: //                    cerr << e->d_name << ": no Vamp descriptor function" << endl;
cannam@40:                 }
cannam@40:                 DLCLOSE(handle);
cannam@47:             } else {
cannam@47:                 cerr << "\n" << e->d_name << ": unable to load library (" << DLERROR() << ")" << endl;
cannam@47:             }                
cannam@40:         }
cannam@40:         closedir(d);
cannam@40:     }
cannam@40:     cerr << endl;
cannam@40: }
cannam@40: 
cannam@40: #endif
cannam@40: 
cannam@40: 
cannam@40: void
cannam@16: printFeatures(int frame, int sr, int output, Vamp::Plugin::FeatureSet features)
cannam@16: {
cannam@16:     for (unsigned int i = 0; i < features[output].size(); ++i) {
cannam@16:         Vamp::RealTime rt = Vamp::RealTime::frame2RealTime(frame, sr);
cannam@16:         if (features[output][i].hasTimestamp) {
cannam@16:             rt = features[output][i].timestamp;
cannam@16:         }
cannam@16:         cout << rt.toString() << ":";
cannam@16:         for (unsigned int j = 0; j < features[output][i].values.size(); ++j) {
cannam@16:             cout << " " << features[output][i].values[j];
cannam@16:         }
cannam@16:         cout << endl;
cannam@16:     }
cannam@16: }
cannam@16: 
cannam@16: void
cannam@16: transformInput(float *buffer, size_t size)
cannam@16: {
cannam@16:     double *inbuf = new double[size * 2];
cannam@16:     double *outbuf = new double[size * 2];
cannam@16: 
cannam@16:     // Copy across with Hanning window
cannam@16:     for (size_t i = 0; i < size; ++i) {
cannam@16:         inbuf[i] = double(buffer[i]) * (0.50 - 0.50 * cos(2 * M_PI * i / size));
cannam@16:         inbuf[i + size] = 0.0;
cannam@16:     }
cannam@16:     
cannam@16:     for (size_t i = 0; i < size/2; ++i) {
cannam@16:         double temp = inbuf[i];
cannam@16:         inbuf[i] = inbuf[i + size/2];
cannam@16:         inbuf[i + size/2] = temp;
cannam@16:     }
cannam@16: 
cannam@16:     fft(size, false, inbuf, inbuf + size, outbuf, outbuf + size);
cannam@16: 
cannam@47:     for (size_t i = 0; i <= size/2; ++i) {
cannam@16:         buffer[i * 2] = outbuf[i];
cannam@16:         buffer[i * 2 + 1] = outbuf[i + size];
cannam@16:     }
cannam@16:     
cannam@46:     delete[] inbuf;
cannam@46:     delete[] outbuf;
cannam@16: }
cannam@16: 
cannam@16: void
cannam@16: fft(unsigned int n, bool inverse, double *ri, double *ii, double *ro, double *io)
cannam@16: {
cannam@16:     if (!ri || !ro || !io) return;
cannam@16: 
cannam@16:     unsigned int bits;
cannam@16:     unsigned int i, j, k, m;
cannam@16:     unsigned int blockSize, blockEnd;
cannam@16: 
cannam@16:     double tr, ti;
cannam@16: 
cannam@16:     if (n < 2) return;
cannam@16:     if (n & (n-1)) return;
cannam@16: 
cannam@16:     double angle = 2.0 * M_PI;
cannam@16:     if (inverse) angle = -angle;
cannam@16: 
cannam@16:     for (i = 0; ; ++i) {
cannam@16: 	if (n & (1 << i)) {
cannam@16: 	    bits = i;
cannam@16: 	    break;
cannam@16: 	}
cannam@16:     }
cannam@16: 
cannam@16:     static unsigned int tableSize = 0;
cannam@16:     static int *table = 0;
cannam@16: 
cannam@16:     if (tableSize != n) {
cannam@16: 
cannam@16: 	delete[] table;
cannam@16: 
cannam@16: 	table = new int[n];
cannam@16: 
cannam@16: 	for (i = 0; i < n; ++i) {
cannam@16: 	
cannam@16: 	    m = i;
cannam@16: 
cannam@16: 	    for (j = k = 0; j < bits; ++j) {
cannam@16: 		k = (k << 1) | (m & 1);
cannam@16: 		m >>= 1;
cannam@16: 	    }
cannam@16: 
cannam@16: 	    table[i] = k;
cannam@16: 	}
cannam@16: 
cannam@16: 	tableSize = n;
cannam@16:     }
cannam@16: 
cannam@16:     if (ii) {
cannam@16: 	for (i = 0; i < n; ++i) {
cannam@16: 	    ro[table[i]] = ri[i];
cannam@16: 	    io[table[i]] = ii[i];
cannam@16: 	}
cannam@16:     } else {
cannam@16: 	for (i = 0; i < n; ++i) {
cannam@16: 	    ro[table[i]] = ri[i];
cannam@16: 	    io[table[i]] = 0.0;
cannam@16: 	}
cannam@16:     }
cannam@16: 
cannam@16:     blockEnd = 1;
cannam@16: 
cannam@16:     for (blockSize = 2; blockSize <= n; blockSize <<= 1) {
cannam@16: 
cannam@16: 	double delta = angle / (double)blockSize;
cannam@16: 	double sm2 = -sin(-2 * delta);
cannam@16: 	double sm1 = -sin(-delta);
cannam@16: 	double cm2 = cos(-2 * delta);
cannam@16: 	double cm1 = cos(-delta);
cannam@16: 	double w = 2 * cm1;
cannam@16: 	double ar[3], ai[3];
cannam@16: 
cannam@16: 	for (i = 0; i < n; i += blockSize) {
cannam@16: 
cannam@16: 	    ar[2] = cm2;
cannam@16: 	    ar[1] = cm1;
cannam@16: 
cannam@16: 	    ai[2] = sm2;
cannam@16: 	    ai[1] = sm1;
cannam@16: 
cannam@16: 	    for (j = i, m = 0; m < blockEnd; j++, m++) {
cannam@16: 
cannam@16: 		ar[0] = w * ar[1] - ar[2];
cannam@16: 		ar[2] = ar[1];
cannam@16: 		ar[1] = ar[0];
cannam@16: 
cannam@16: 		ai[0] = w * ai[1] - ai[2];
cannam@16: 		ai[2] = ai[1];
cannam@16: 		ai[1] = ai[0];
cannam@16: 
cannam@16: 		k = j + blockEnd;
cannam@16: 		tr = ar[0] * ro[k] - ai[0] * io[k];
cannam@16: 		ti = ar[0] * io[k] + ai[0] * ro[k];
cannam@16: 
cannam@16: 		ro[k] = ro[j] - tr;
cannam@16: 		io[k] = io[j] - ti;
cannam@16: 
cannam@16: 		ro[j] += tr;
cannam@16: 		io[j] += ti;
cannam@16: 	    }
cannam@16: 	}
cannam@16: 
cannam@16: 	blockEnd = blockSize;
cannam@16:     }
cannam@16: 
cannam@16:     if (inverse) {
cannam@16: 
cannam@16: 	double denom = (double)n;
cannam@16: 
cannam@16: 	for (i = 0; i < n; i++) {
cannam@16: 	    ro[i] /= denom;
cannam@16: 	    io[i] /= denom;
cannam@16: 	}
cannam@16:     }
cannam@16: }
cannam@16: 
cannam@16: