Marsyas  0.6.0-alpha
/usr/src/RPM/BUILD/marsyas-0.6.0/src/marsyas/marsystems/Spectrum2Chroma.cpp
Go to the documentation of this file.
00001 /*
00002 ** Copyright (C) 1998-2010 George Tzanetakis <gtzan@cs.uvic.ca>
00003 **
00004 ** This program is free software; you can redistribute it and/or modify
00005 ** it under the terms of the GNU General Public License as published by
00006 ** the Free Software Foundation; either version 2 of the License, or
00007 ** (at your option) any later version.
00008 **
00009 ** This program is distributed in the hope that it will be useful,
00010 ** but WITHOUT ANY WARRANTY; without even the implied warranty of
00011 ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
00012 ** GNU General Public License for more details.
00013 **
00014 ** You should have received a copy of the GNU General Public License
00015 ** along with this program; if not, write to the Free Software
00016 ** Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA.
00017 */
00018 
00019 #include "Spectrum2Chroma.h"
00020 
00021 using namespace std;
00022 using namespace Marsyas;
00023 
00024 Spectrum2Chroma::Spectrum2Chroma(mrs_string name):MarSystem("Spectrum2Chroma", name)
00025 {
00026   addControls();
00027 
00028   noteNames_.push_back("A");
00029   noteNames_.push_back("A#");
00030   noteNames_.push_back("B");
00031   noteNames_.push_back("C");
00032   noteNames_.push_back("C#");
00033   noteNames_.push_back("D");
00034   noteNames_.push_back("D#");
00035   noteNames_.push_back("E");
00036   noteNames_.push_back("F");
00037   noteNames_.push_back("F#");
00038   noteNames_.push_back("G");
00039   noteNames_.push_back("G#");
00040 
00041   pnbins_ = 0;
00042   pmiddleAfreq_ = 0.0;
00043   pweightCenterFreq_ = 0.0;
00044   pweightStdDev_ = 0.0;
00045 }
00046 
00047 Spectrum2Chroma::Spectrum2Chroma(const Spectrum2Chroma& a) : MarSystem(a)
00048 {
00049   ctrl_nbins_ = getctrl("mrs_natural/nbins");
00050   ctrl_middleAfreq_ = getctrl("mrs_real/middleAfreq");
00051   ctrl_weightCenterFreq_ = getctrl("mrs_real/weightCenterFreq");
00052   ctrl_weightStdDev_ = getctrl("mrs_real/weightStdDev");
00053 
00054   noteNames_ = a.noteNames_;
00055   chromaMap_ = a.chromaMap_;
00056   pnbins_ = a.pnbins_;
00057   pmiddleAfreq_ = a.pmiddleAfreq_;
00058   pweightCenterFreq_ = a.pweightCenterFreq_;
00059   pweightStdDev_ = a.pweightStdDev_;
00060 }
00061 
00062 Spectrum2Chroma::~Spectrum2Chroma()
00063 {
00064 }
00065 
00066 MarSystem*
00067 Spectrum2Chroma::clone() const
00068 {
00069   return new Spectrum2Chroma(*this);
00070 }
00071 
00072 void
00073 Spectrum2Chroma::addControls()
00074 {
00075   addctrl("mrs_natural/nbins", 12, ctrl_nbins_); //diatonic chromatic scale by default
00076   addctrl("mrs_real/middleAfreq", 440.0, ctrl_middleAfreq_);
00077   mrs_real A0freq = 440.0/pow(2.0, 4.0);
00078   addctrl("mrs_real/weightCenterFreq", log(1000.0/A0freq)/log(2.0), ctrl_weightCenterFreq_);
00079   addctrl("mrs_real/weightStdDev", 0.0, ctrl_weightStdDev_);
00080 
00081   ctrl_nbins_->setState(true);
00082   ctrl_middleAfreq_->setState(true);
00083   ctrl_weightCenterFreq_->setState(true);
00084   ctrl_weightStdDev_->setState(true);
00085 }
00086 
00087 void
00088 Spectrum2Chroma::myUpdate(MarControlPtr sender)
00089 {
00090   mrs_natural t,o;
00091   (void) sender;  //suppress warning of unused parameter(s)
00092 
00093   ctrl_onSamples_->setValue(ctrl_inSamples_, NOUPDATE);
00094   ctrl_onObservations_->setValue(ctrl_nbins_, NOUPDATE);
00095   ctrl_osrate_->setValue(ctrl_israte_, NOUPDATE);
00096 
00097   //if (pnbins_ != ctrl_nbins_->to<mrs_natural>())
00098   // {
00099   pnbins_ = ctrl_nbins_->to<mrs_natural>();
00100   ostringstream oss;
00101   if(pnbins_ == 12)
00102   {
00103     for (mrs_natural n=0; n < pnbins_; n++)
00104     {
00105       oss << "Chroma_" << noteNames_[n] << "_" << ctrl_inObsNames_->to<mrs_string>();
00106     }
00107   }
00108   else
00109   {
00110     for (mrs_natural n=0; n < pnbins_; n++)
00111     {
00112       oss << "Chroma_" << n << "_" << ctrl_inObsNames_->to<mrs_string>();
00113     }
00114   }
00115   ctrl_onObsNames_->setValue(oss.str(), NOUPDATE);
00116 
00117 // }
00118 
00120   // calculate the Chroma map
00121   // based in the fft2chromamx.m MATLAB script by Dan Ellis
00122   // http://www.ee.columbia.edu/~dpwe/resources/matlab/chroma-ansyn/#1
00124   if(tinObservations_ != inObservations_ ||
00125       tonObservations_ != onObservations_ ||
00126       tisrate_ != israte_ ||
00127       pmiddleAfreq_ != ctrl_middleAfreq_->to<mrs_real>() ||
00128       pweightCenterFreq_ != ctrl_weightCenterFreq_->to<mrs_real>() ||
00129       pweightStdDev_ != ctrl_weightStdDev_->to<mrs_real>() )
00130   {
00131 
00132     pmiddleAfreq_ = ctrl_middleAfreq_->to<mrs_real>();
00133     pweightCenterFreq_ = ctrl_weightCenterFreq_->to<mrs_real>();
00134     pweightStdDev_ = ctrl_weightStdDev_->to<mrs_real>();
00135 
00136     mrs_natural nbins = ctrl_nbins_->to<mrs_natural>();
00137     mrs_natural nbins2 = (mrs_natural)floor(nbins/2.0+0.5); //equivalent to round()
00138     mrs_natural N2 = inObservations_; // we get N/2+1 spectrum points at the input...
00139     mrs_natural N = (N2-1)*2; //fft size
00140 
00141     //get the original audio sampling rate
00142     mrs_real srate = israte_*N;
00143 
00144     //calculate the frequencies (in octaves) for each spectrum bin
00145     //Note: make up value for 0Hz bin -> 1.5 octaves below bin 1
00146     //(so chroma is 50% rotated from bin 1, and bin width is broad)
00147     realvec fftfreqbins(N2);
00148     for(o=1; o < N2; ++o)
00149     {
00150       fftfreqbins(o) = nbins * hertz2octs(o * 1.0 / N2 * srate, ctrl_middleAfreq_->to<mrs_real>());
00151 
00152     }
00153     if (N2 > 1)
00154       fftfreqbins(0) = fftfreqbins(1)-1.5*nbins;
00155 
00156 
00157 
00158 
00159     //calculate the bin widths
00160     realvec binwidthbins(N2);
00161     for(o=0; o < N2-1; ++o)
00162     {
00163       binwidthbins(o) = fftfreqbins(o+1)-fftfreqbins(o);
00164       if(binwidthbins(o) < 1.0)
00165         binwidthbins(o) = 1.0;
00166     }
00167     binwidthbins(N2-1) = 1.0;
00168 
00169 
00170 
00171     //calculate chroma mapping
00172     chromaMap_.create(nbins, N2); //=wts in Dan Ellis MATLAB code
00173     realvec D(nbins, N2);
00174     for(o = 0; o < nbins; ++o)
00175     {
00176       for(t = 0; t < N2; ++t)
00177       {
00178         D(o,t) = fftfreqbins(t) - o;
00179 
00180         //project into range -nbins/2 .. nbins/2
00181         //(add on fixed offset of 10*nbins to ensure all values are positive
00182         //for the fmod remainder operation)
00183         D(o,t) = fmod(D(o,t) + nbins2 + 10*nbins, nbins) - nbins2;
00184 
00185         //Gaussian bumps (2*D(o,t) to make them narrower)
00186         chromaMap_(o,t) = exp(-0.5 * pow(2.0*D(o,t)/binwidthbins(t), 2.0) );
00187       }
00188     }
00189 
00190 
00191 
00192     //normalize each column by its RMS value
00193     mrs_real colRMS;
00194     for(t = 0; t < N2; ++t) //iterate over columns
00195     {
00196       colRMS = 0.0;
00197       //get RMS value for column t
00198       for(o = 0; o < nbins; ++o) //iterate over rows
00199       {
00200         colRMS += (chromaMap_(o,t)*chromaMap_(o,t));
00201       }
00202       //normalize column t
00203       if(colRMS != 0.0)
00204       {
00205         for(o = 0; o < nbins; ++o) //row
00206         {
00207           chromaMap_(o,t) = chromaMap_(o,t)/sqrt(colRMS);
00208         }
00209       }
00210     }
00211 
00212 
00213     //Maybe apply scaling for fft bins
00214     mrs_real ctroct = ctrl_weightCenterFreq_->to<mrs_real>();
00215     mrs_real octwidth = ctrl_weightStdDev_->to<mrs_real>();
00216 
00217     if(octwidth > 0.0)
00218     {
00219       for(o = 0; o < nbins; ++o)
00220       {
00221         for(t = 0; t < N2; ++t)
00222         {
00223           chromaMap_(o,t) = exp(-0.5 * pow((mrs_real) (fftfreqbins(t)/nbins - ctroct)/octwidth , (mrs_real)2.0) );
00224         }
00225       }
00226     }
00227 #ifdef MARSYAS_MATLAB
00228 //      MATLAB_PUT(chromaMap_, "marCM");
00229 //      MATLAB_EVAL("CM = fft2chromamx("<<N<<","<<nbins<<","<<srate<<","<<pmiddleAfreq_<<","<<pweightCenterFreq_<<","<<pweightStdDev_<<");");
00230 #endif //MARSYAS_MATLAB
00231   }
00232 }
00233 
00234 void
00235 Spectrum2Chroma::myProcess(realvec& in, realvec& out)
00236 {
00237   //input must contain spectral magnitude/power/density/etc
00238   //(e.g. output of PowerSpectrum MarSystem)
00239 
00240   mrs_natural o,t,i;
00241   mrs_real chroma_weight;
00242 
00243   out.setval(0.0);
00244 
00245 
00246   for(o=0; o< onObservations_; ++o)
00247   {
00248     for(i=0; i< inObservations_; ++i)
00249     {
00250       chroma_weight = chromaMap_(o,i);
00251 
00252       for(t=0; t< inSamples_; ++t)
00253       {
00254         out(o,t) += in(i,t)*chroma_weight;
00255       }
00256 
00257     }
00258   }
00259 
00260 #ifdef MARSYAS_MATLAB
00261 //  MATLAB_PUT(in, "spectrum");
00262 //  MATLAB_PUT(out,"marChromaVec");
00263 //  MATLAB_EVAL("chromaVec = CM*spectrum;");
00264 #endif //MARSYAS_MATLAB
00265 }