Marsyas  0.5.0-beta1
/Users/jleben/code/marsyas/src/marsyas/Conversions.cpp
Go to the documentation of this file.
00001 /*
00002 ** Copyright (C) 1998-2005 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 <marsyas/Conversions.h>
00020 #include <marsyas/realvec.h>
00021 
00022 using std::ostringstream;
00023 using namespace Marsyas;
00024 
00025 mrs_real
00026 Marsyas::pitch2hertz(mrs_real pitch) {
00027   return mrs_real(440.0 * pow(2.0, ((pitch - 69.0) / 12.0)));
00028 }
00029 
00030 mrs_real
00031 Marsyas::hertz2pitch(mrs_real hz) {
00032   return (hz == 0.0) ? (mrs_real)0.0 : (mrs_real)(69.0 + 12.0 * (log(hz/440.0)/log(2.0)));
00033 }
00034 
00035 mrs_real
00036 Marsyas::samples2hertz(mrs_natural samples, mrs_real srate) {
00037   return (samples <= 0)  ? (mrs_real)0.0 : (mrs_real) (srate * 1.0) / (samples);
00038 }
00039 
00040 mrs_real
00041 Marsyas::samples2hertz(mrs_real samples, mrs_real srate) {
00042   return (samples <= 0.001)  ? (mrs_real)0.0 : (mrs_real) (srate * 1.0) / (samples);
00043 }
00044 
00045 
00046 mrs_natural
00047 Marsyas::hertz2samples(mrs_real hz, mrs_real srate) {
00048   return (hz == 0.0) ? (mrs_natural)0 : (mrs_natural) (srate / hz);
00049 }
00050 
00051 /* convert a string representing time to number of samples base on the
00052 given sample rate. Format "123.456#" where # is the time division.
00053 Valid time divisions: { h, m, s, ms, us }.
00054 On a format error,
00055 Errors: -1 is returned. ie more than 1 decimal point, invalid time
00056 division.
00057 */
00058 mrs_natural
00059 Marsyas::time2samples(mrs_string time, mrs_real srate) {
00060   //example times: { "10us", "10ms", "10s", "10m", "10h" }
00061   if (time=="") { return 0; }
00062   // calculate time value
00063   mrs_real samples=0;
00064   int i=0;
00065   int len=(int)time.length();
00066   bool decimal_point=false;
00067   mrs_real divisor = 10.0;
00068   for (i=0; i<len && (time[i]=='.' || (time[i]>='0' && time[i]<='9')); ++i) {
00069     if (decimal_point) {
00070       if (time[i]=='.') { return -1; }
00071       samples = samples + ((mrs_real)(time[i]-'0'))/divisor;
00072       divisor = divisor * 10.0;
00073     } else if (time[i]=='.') {
00074       decimal_point=true;
00075     } else {
00076       samples = samples * 10.0 + (time[i]-'0');
00077     }
00078   }
00079   //
00080   if (i<len) {
00081     char a=time[++i];
00082     if (i>=len) {
00083       if (a=='h') { // hours
00084         samples= 120.0*samples*srate;
00085       } else if (a=='m') { // minutes
00086         samples=  60.0*samples*srate;
00087       } else if (a=='s') { // seconds
00088         samples=       samples*srate;
00089       } else {
00090         return -1;
00091       }
00092     } else {
00093       char b=time[i];
00094       if ((i+1)>=len) {
00095         if (a=='u' && b=='s') { // micro-seconds
00096           samples= samples/1000000.0*srate;
00097         } else if (a=='m' && b=='s') { // milli-seconds
00098           samples= samples/1000.0*srate;
00099         } else {
00100           return -1;
00101         }
00102       }
00103     }
00104   }
00105   return (mrs_natural)samples;
00106 }
00107 mrs_natural
00108 Marsyas::time2usecs(mrs_string time) {
00109   //example times: { "10us", "10ms", "10s", "10m", "10h" }
00110   if (time=="") { return 0; }
00111   // calculate time value
00112   mrs_real usecs=0;
00113   int i=0;
00114   int len=(int)time.length();
00115   bool decimal_point=false;
00116   mrs_real divisor = 10.0;
00117   for (i=0; i<len && (time[i]=='.' || (time[i]>='0' && time[i]<='9')); ++i) {
00118     if (decimal_point) {
00119       if (time[i]=='.') { return -1; }
00120       usecs = usecs + ((mrs_real)(time[i]-'0'))/divisor;
00121       divisor = divisor * 10.0;
00122     } else if (time[i]=='.') {
00123       decimal_point=true;
00124     } else {
00125       usecs = usecs * 10.0 + (time[i]-'0');
00126     }
00127   }
00128   //
00129   if (i<len) {
00130     char a=time[++i];
00131     if (i>=len) {
00132       if (a=='h') { // hours
00133         usecs *= 1000.0 * 1000.0 * 60.0 * 60.0;
00134       } else if (a=='m') { // minutes
00135         usecs *= 1000.0 * 1000.0 * 60.0;
00136       } else if (a=='s') { // seconds
00137         usecs *= 1000.0 * 1000.0;
00138       } else {
00139         return -1;
00140       }
00141     } else {
00142       char b=time[i];
00143       if ((i+1)>=len) {
00144         if (a=='u' && b=='s') { // micro-seconds
00145           ;
00146         } else if (a=='m' && b=='s') { // milli-seconds
00147           usecs *= 1000.0;
00148         } else {
00149           return -1;
00150         }
00151       }
00152     }
00153   }
00154   return (mrs_natural)usecs;
00155 }
00156 
00157 mrs_real Marsyas::amplitude2dB(mrs_real a)
00158 {
00159   MRSASSERT (a > 0);
00160   return 20*log10(a);
00161 }
00162 
00163 mrs_real Marsyas::dB2amplitude(mrs_real a)
00164 {
00165   return pow(10.0, a/20);
00166 }
00167 
00168 mrs_real
00169 Marsyas::hertz2octs(mrs_real f, mrs_real middleAfreq)
00170 {
00171   //adapted from Dan Ellis fft2chromamx.m MATLAB routine
00172   //
00173   // octs = hz2octs(freq, A440)
00174   // Convert a frequency in Hz into a real number counting
00175   // the octaves above A0. So hz2octs(440) = 4.0
00176   // Optional A440 specifies the Hz to be treated as middle A (default 440).
00177   // 2006-06-29 dpwe@ee.columbia.edu for fft2chromamx
00178   //
00179   // A4 = 440 Hz, so A0 = 440/16 Hz
00180   // octs = log(freq./(A440/16))./log(2);
00181   //
00182   return log(f/(middleAfreq/16.0))/log(2.0);
00183 }
00184 
00185 mrs_real Marsyas::hertz2bark(mrs_real f, mrs_natural mode)
00186 {
00187   switch (mode)
00188   {
00189   case 0:
00190   default:
00191     return  6 * log(f/600 + sqrt(1+ (pow(f/600,2)))); // 6*asinh(f/600);
00192 
00193   case  1:  //  zwicker
00194     return  13 * atan  (0.00076*f)  +  3.5  *  atan  ((f*0.000133333333333333)*(f*0.000133333333333333));
00195 
00196   case  2:  //  terhardt
00197     return  13.3  *  atan  (0.00075*f);
00198 
00199   case  3:  //  schroeder
00200     return  7.0  *  log  (f*0.00153846153846154  +  sqrt  ((f*0.00153846153846154)*(f*0.00153846153846154)  +  1));  //  7*asinh  (fFrequency/650);
00201   }
00202 }
00203 
00204 mrs_real Marsyas::bark2hertz(mrs_real f, mrs_natural mode)
00205 {
00206   switch  (mode)
00207   {
00208   case 0:
00209   case 1:
00210   default:
00211     return 600*sinh(f/6.0);
00212 
00213   case  2:  //  terhardt
00214     return  4*1000.0/3.0  *  tan  (f*0.075187969924812);
00215   case  3:  //  schroeder
00216     return  650.0  *  sinh  (f*0.142857142857143);
00217   }
00218 }
00219 
00220 mrs_real Marsyas::hertz2erb(mrs_real hz)
00221 {
00222   return 21.4 * log10(1+ (hz / 229.0));
00223 }
00224 
00225 mrs_real Marsyas::erb2hertz(mrs_real erb)
00226 {
00227   return (pow(10, (erb/21.4)) - 1.0) * 229.0;
00228 }
00229 
00230 
00231 mrs_real
00232 Marsyas::hertz2mel(mrs_real f, bool htk)
00233 {
00234   //  z = hertz2mel(f,htk)
00235   //  Convert frequencies f (in Hz) to mel 'scale'.
00236   //  Optional htk = 1 uses the mel axis defined in the HTKBook
00237   //  otherwise use Slaney's formula.
00238   //
00239   //  adapted from Dan Ellis fft2melmx.m MATLAB code
00240 
00241   if(htk)
00242   {
00243     return 2595.0 * log10(1.0 + f / 700.0);
00244   }
00245   else
00246   {
00247     // Mel fn to match Slaney's Auditory Toolbox mfcc.m
00248     mrs_real f_0 = 0.0; //133.33333;
00249     mrs_real f_sp = 200.0/3.0; //66.66667;
00250     mrs_real brkfrq = 1000.0;
00251     mrs_real brkpt  = (brkfrq - f_0)/f_sp;  //starting mel value for log region
00252 
00253     //the magic 1.0711703 which is the ratio needed to get from
00254     //1000 Hz to 6400 Hz in 27 steps, and is *almost* the ratio between
00255     //1000 Hz and the preceding linear filter center at 933.33333 Hz
00256     //(actually 1000/933.33333 = 1.07142857142857 and
00257     //exp(log(6.4)/27) = 1.07117028749447)
00258     mrs_real logstep = exp(log(6.4)/27.0);
00259 
00260     if(f < brkfrq)
00261       return (f - f_0) / f_sp; //linear
00262     else
00263       return brkpt + log(f / brkfrq) / log(logstep); //non-linear
00264   }
00265 }
00266 
00267 mrs_real
00268 Marsyas::mel2hertz(mrs_real z, bool htk)
00269 {
00270   //   f = mel2hz(z, htk)
00271   //   Convert 'mel scale' frequencies into Hz
00272   //   Optional htk = 1 means use the HTK formula
00273   //   else use the formula from Slaney's mfcc.m
00274   //
00275   //     Adapted from Dan Ellis fft2melmx.m MATLAB code
00276 
00277   if(htk)
00278   {
00279     return 700.0 * (pow(10.0, z/2595.0) - 1.0);
00280   }
00281   else
00282   {
00283     mrs_real f_0 = 0.0; //133.33333;
00284     mrs_real f_sp = 200.0 / 3.0; //66.66667;
00285     mrs_real brkfrq = 1000.0;
00286     mrs_real brkpt  = (brkfrq - f_0)/f_sp; //starting mel value for log region
00287 
00288     //the magic 1.0711703 which is the ratio needed to get from 1000 Hz
00289     //to 6400 Hz in 27 steps, and is *almost* the ratio between 1000 Hz
00290     //and the preceding linear filter center at 933.33333 Hz
00291     //(actually 1000/933.33333 = 1.07142857142857 and
00292     //exp(log(6.4)/27) = 1.07117028749447)
00293     mrs_real logstep = exp(log(6.4)/27.0);
00294 
00295     if(z < brkpt)
00296       return f_0 + f_sp * z;
00297     else
00298       return brkfrq * exp(log(logstep)*(z-brkpt));
00299   }
00300 }
00301 
00302 // Returns the frequency in Hz, represented by a PowerSpectrum bin number
00303 // (as returned from the PowerSpectrum Marsystem).
00304 // * bin is the bin number
00305 // * nyquistFreq is the Nyquist frequency (the sampling rate / 2)
00306 // * framerate is the samplerate / total # bins (usually we can pass in israte_)
00307 //
00308 // For example, if the nyquistFreq = 4000, the original sample rate was 8000Hz,
00309 // and we are using a framesize of 8 samples (8 FFT bins), then we should pass
00310 // in framerate=1000(Hz). Then there are 5 PowerSpect bins ((8/2)+1)
00311 // corresponding to frequencies as follows:
00312 // bin=0: returns 0
00313 // bin=1: returns 4000   (+/- 4000 Hz)
00314 // bin=2: returns 1000   (+/- 1000 Hz)
00315 // bin=3: returns 2000   (+/- 2000 Hz)
00316 // bin=4: returns 3000   (+/- 3000 Hz)
00317 mrs_real
00318 Marsyas::bin2hertz(mrs_natural bin, mrs_real nyquistFreq, mrs_real framerate)
00319 {
00320   if (bin == 0)
00321     return 0;
00322   if (bin == 1)
00323     return nyquistFreq;
00324   else
00325     return (bin-1) * framerate;
00326 }
00327 
00328 mrs_natural Marsyas::powerOfTwo(mrs_real v)
00329 {
00330   mrs_natural n=1, res=0;
00331   while(res < v)
00332   {
00333     res = (mrs_natural) pow(2.0, n+.0);
00334     n++;
00335   }
00336   return res;
00337 }
00338 
00339 void
00340 Marsyas::string2parameters(mrs_string s, realvec &v, char d)
00341 {
00342   mrs_natural i =0, pos=0, newPos=0;
00343   mrs_string tmp;
00344   while(newPos != -1 )
00345   {
00346     newPos = (mrs_natural) s.find_first_of(&d, pos, 1);
00347     tmp = s.substr(pos, newPos);
00348     v(i++) = atof(tmp.c_str());
00349     pos = newPos+1;
00350   }
00351 }