Marsyas  0.6.0-alpha
/usr/src/RPM/BUILD/marsyas-0.6.0/src/marsyas/marsystems/Spectrum2Mel.cpp
Go to the documentation of this file.
00001 /*
00002 ** Copyright (C) 1998-2006 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 "Spectrum2Mel.h"
00020 #include <algorithm>
00021 
00022 using namespace std;
00023 using namespace Marsyas;
00024 
00025 Spectrum2Mel::Spectrum2Mel(mrs_string name):MarSystem("Spectrum2Mel", name)
00026 {
00027   addControls();
00028 
00029   pmelBands_ = 0;
00030   pbandWidth_ = 0.0;
00031   pbandLowEdge_ = 0.0;
00032   pbandHighEdge_ = 0.0;
00033   phtkMel_ = false;
00034   pconstAmp_ = false;
00035 }
00036 
00037 Spectrum2Mel::Spectrum2Mel(const Spectrum2Mel& a) : MarSystem(a)
00038 {
00039   ctrl_melBands_ = getctrl("mrs_natural/melBands");
00040   ctrl_bandWidth_ = getctrl("mrs_real/bandWidth");
00041   ctrl_bandLowEdge_ = getctrl("mrs_real/bandLowEdge");
00042   ctrl_bandHighEdge_ = getctrl("mrs_real/bandHighEdge");
00043   ctrl_htkMel_ = getctrl("mrs_bool/htkMel");
00044   ctrl_constAmp_ = getctrl("mrs_bool/constAmp");
00045 
00046   melMap_ = a.melMap_;
00047   pmelBands_ = a.pmelBands_;
00048   pbandWidth_ = a.pbandWidth_;
00049   pbandLowEdge_ = a.pbandLowEdge_;
00050   pbandHighEdge_ = a.pbandHighEdge_;
00051   phtkMel_ = a.phtkMel_;
00052   pconstAmp_ = a.pconstAmp_;
00053 
00054 }
00055 
00056 Spectrum2Mel::~Spectrum2Mel()
00057 {
00058 }
00059 
00060 MarSystem*
00061 Spectrum2Mel::clone() const
00062 {
00063   return new Spectrum2Mel(*this);
00064 }
00065 
00066 void
00067 Spectrum2Mel::addControls()
00068 {
00069   addctrl("mrs_natural/melBands", 40, ctrl_melBands_);
00070   addctrl("mrs_real/bandWidth", 1.0, ctrl_bandWidth_);
00071   addctrl("mrs_real/bandLowEdge", 0.0, ctrl_bandLowEdge_);
00072   addctrl("mrs_real/bandHighEdge", -1.0, ctrl_bandHighEdge_);
00073   addctrl("mrs_bool/htkMel", false, ctrl_htkMel_);
00074   addctrl("mrs_bool/constAmp", false, ctrl_constAmp_);
00075 
00076   ctrl_melBands_->setState(true);
00077   ctrl_bandWidth_->setState(true);
00078   ctrl_bandLowEdge_->setState(true);
00079   ctrl_bandHighEdge_->setState(true);
00080   ctrl_htkMel_->setState(true);
00081   ctrl_constAmp_->setState(true);
00082 }
00083 
00084 void
00085 Spectrum2Mel::myUpdate(MarControlPtr sender)
00086 {
00087   mrs_natural t,o;
00088   (void) sender;  //suppress warning of unused parameter(s)
00089 
00090   ctrl_onSamples_->setValue(ctrl_inSamples_, NOUPDATE);
00091   ctrl_onObservations_->setValue(ctrl_melBands_, NOUPDATE);
00092   ctrl_osrate_->setValue(ctrl_israte_, NOUPDATE);
00093 
00094   // allocate memory for melMap_
00095   // I'm not familiar with this algorithm so this might be
00096   // wrong!  I'm just closing some memory faults!  -gp
00097   // size based on its usage in ::myProcess
00098   onObservations_ = ctrl_onObservations_->to<mrs_natural>();
00099 
00100 
00101   if (pmelBands_ != ctrl_melBands_->to<mrs_natural>())
00102   {
00103     pmelBands_ = ctrl_melBands_->to<mrs_natural>();
00104     ostringstream oss;
00105     for (mrs_natural n=0; n < pmelBands_; n++)
00106     {
00107       oss << "MelBand_" << n << ",";
00108     }
00109     ctrl_onObsNames_->setValue(oss.str(), NOUPDATE);
00110   }
00111 
00113   // calculate the Mel map
00114   // based in the fft2melmx.m MATLAB script by Dan Ellis
00115   // http://labrosa.ee.columbia.edu/projects/coversongs/
00117   if(tinObservations_ != inObservations_ ||
00118       tonObservations_ != onObservations_ ||
00119       tisrate_ != israte_ ||
00120       pbandWidth_ != ctrl_bandWidth_->to<mrs_real>() ||
00121       pbandLowEdge_ != ctrl_bandLowEdge_->to<mrs_real>() ||
00122       pbandHighEdge_ != ctrl_bandHighEdge_->to<mrs_real>() ||
00123       phtkMel_ != ctrl_htkMel_->to<mrs_bool>() ||
00124       pconstAmp_ != ctrl_constAmp_->to<mrs_bool>())
00125   {
00126     melMap_.allocate(onObservations_, inObservations_);
00127 
00128     mrs_natural nfilts = ctrl_melBands_->to<mrs_natural>();
00129     bool htkmel = ctrl_htkMel_->to<mrs_bool>();
00130     mrs_natural N2 = inObservations_;// we get N/2+1 spectrum points at the input...
00131     mrs_natural N = (N2-1)*2; //fft size
00132 
00133     //get the original audio sampling rate
00134     mrs_real srate = israte_*N;
00135 
00136     if(ctrl_bandHighEdge_->to<mrs_real>() == -1.0)
00137       ctrl_bandHighEdge_->setValue(srate/2.0, NOUPDATE);
00138 
00139     pbandWidth_ = ctrl_bandWidth_->to<mrs_real>();
00140     pbandLowEdge_ = ctrl_bandLowEdge_->to<mrs_real>();
00141     pbandHighEdge_ = ctrl_bandHighEdge_->to<mrs_real>();
00142     phtkMel_ = ctrl_htkMel_->to<mrs_bool>();
00143     pconstAmp_ = ctrl_constAmp_->to<mrs_bool>();
00144 
00145     //calculate the frequencies (in Hz) for each spectrum bin
00146     realvec fftfreqs(N2);
00147     for(o=0; o < N2; ++o) {
00148       fftfreqs(o) = o / (float)N * srate;
00149     }
00150 
00151     // 'Center freqs' of mel bands - uniformly spaced between limits
00152     mrs_real minmel = hertz2mel(ctrl_bandLowEdge_->to<mrs_real>(), htkmel);
00153     mrs_real maxmel = hertz2mel(ctrl_bandHighEdge_->to<mrs_real>(), htkmel);
00154     realvec binfrqs(nfilts+2);
00155     realvec binbin(nfilts+2);
00156     for(t=0; t < nfilts+2; ++t)
00157     {
00158       binfrqs(t) = mel2hertz(minmel + t/(nfilts+1.0)*(maxmel-minmel) , htkmel);
00159       binbin(t) = floor(binfrqs(t)/srate*(N-1) + 0.5);
00160     }
00161 
00162     realvec fs(3);
00163     mrs_real width = ctrl_bandWidth_->to<mrs_real>();
00164     realvec loslope(N2);
00165     realvec hislope(N2);
00166 
00167     for(mrs_natural i=0; i < nfilts; ++i)
00168     {
00169       for(t=0; t<3; ++t) {
00170         fs(t) = binfrqs(i+t);
00171       }
00172       //cout<<endl;
00173 
00174       //scale by width
00175       for(t=0; t<3; ++t)
00176         fs(t) = fs(1)+ width*(fs(t) - fs(1));
00177 
00178       //lower and upper slopes for all bins
00179       for(t=0; t < N2; ++t)
00180       {
00181         loslope(t) = (fftfreqs(t) - fs(0)) / (fs(1) - fs(0));
00182         hislope(t) = (fs(2) - fftfreqs(t)) / (fs(2) - fs(1));
00183       }
00184 
00185       // .. then intersect them with each other and zero
00186       for(t=0; t < N2; ++t)
00187         melMap_(i,t) = max((mrs_real)0.0, (mrs_real)min(loslope(t), hislope(t)));
00188     }
00189 
00190     if(ctrl_constAmp_->to<mrs_bool>() == false)
00191     {
00192       //Slaney-style mel is scaled to be approx constant E per channel
00193       mrs_real diagMatrix;
00194       for(o = 0; o < nfilts; ++o)
00195       {
00196         diagMatrix = 2.0 / (binfrqs(o+2) - binfrqs(o));
00197         for(t=0; t< N2; ++t)
00198           melMap_(o,t) *= diagMatrix;
00199       }
00200     }
00201   }
00202 }
00203 
00204 void
00205 Spectrum2Mel::myProcess(realvec& in, realvec& out)
00206 {
00207   //input must contain spectral magnitude/power/density/etc
00208   //(e.g. output of PowerSpectrum MarSystem)
00209   mrs_natural o,t;
00210   out.setval(0.0);
00211   for(t=0; t< inSamples_; ++t)
00212   {
00213     for(o=0; o< onObservations_; ++o)
00214     {
00215       for(mrs_natural i=0; i< inObservations_; ++i)
00216       {
00217         out(o,t)+= in(i,t)*melMap_(o,i);
00218       }
00219     }
00220   }
00221 }