Marsyas
0.6.0-alpha
|
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 "PowerSpectrum.h" 00020 00021 using std::ostringstream; 00022 using namespace Marsyas; 00023 00024 #define PSD_POWER 1 00025 #define PSD_MAG 2 00026 #define PSD_DB 3 00027 #define PSD_WDB 4 00028 #define PSD_PD 5 00029 #define PSD_LOGMAG 6 00030 #define PSD_LOGMAG2 7 00031 00032 PowerSpectrum::PowerSpectrum(mrs_string name):MarSystem("PowerSpectrum",name) 00033 { 00034 ntype_ = PSD_POWER; 00035 N2_ = 0; 00036 re_ = 0.0; 00037 im_ = 0.0; 00038 dB_ = 0.0; 00039 pwr_ = 0.0; 00040 00041 addControls(); 00042 } 00043 00044 PowerSpectrum::PowerSpectrum(const PowerSpectrum& a):MarSystem(a) 00045 { 00046 ctrl_spectrumType_ = getctrl("mrs_string/spectrumType"); 00047 } 00048 00049 00050 PowerSpectrum::~PowerSpectrum() 00051 { 00052 } 00053 00054 void 00055 PowerSpectrum::addControls() 00056 { 00057 addctrl("mrs_string/spectrumType", "power", ctrl_spectrumType_); 00058 setctrlState("mrs_string/spectrumType", true); 00059 } 00060 00061 MarSystem* 00062 PowerSpectrum::clone() const 00063 { 00064 return new PowerSpectrum(*this); 00065 } 00066 00067 void 00068 PowerSpectrum::myUpdate(MarControlPtr sender) 00069 { 00070 (void) sender; //suppress warning of unused parameter(s) 00071 00072 //Spectrum outputs N values, corresponding to N/2+1 00073 //complex and unique spectrum points - see Spectrum.h documentation 00074 N2_ = ctrl_inObservations_->to<mrs_natural>()/2 + 1; 00075 00076 ctrl_onSamples_->setValue(ctrl_inSamples_, NOUPDATE); 00077 ctrl_onObservations_->setValue(N2_, NOUPDATE); //outputs N/2+1 real values 00078 ctrl_osrate_->setValue(ctrl_israte_->to<mrs_real>(), NOUPDATE); 00079 00080 stype_ = ctrl_spectrumType_->to<mrs_string>(); 00081 if (stype_ == "power") 00082 ntype_ = PSD_POWER; 00083 else if (stype_ == "magnitude") 00084 ntype_ = PSD_MAG; 00085 else if (stype_ == "decibels") 00086 ntype_ = PSD_DB; //PUT DB!! (DB->NEW) 00087 else if (stype_ == "wrongdBonsets") 00088 ntype_ = PSD_WDB; 00089 else if (stype_ == "powerdensity") 00090 ntype_ = PSD_PD; 00091 else if (stype_ == "logmagnitude") 00092 ntype_ = PSD_LOGMAG; 00093 else if (stype_ == "logmagnitude2") 00094 ntype_ = PSD_LOGMAG2; 00095 00096 // Add prefix to the observation names. 00097 mrs_string inObsNames = ctrl_inObsNames_->to<mrs_string>(); 00098 ctrl_onObsNames_->setValue("Power_" + stype_ + ctrl_inObsNames_->to<mrs_string>(), NOUPDATE); 00099 00100 // gtzan: removed cluttered output 00101 // ctrl_onObsNames_->setValue(obsNamesAddPrefix(inObsNames, "Power" + stype_ + "_"), NOUPDATE); 00102 } 00103 00104 void 00105 PowerSpectrum::myProcess(realvec& in, realvec& out) 00106 { 00107 mrs_natural t,o; 00108 for (t=0; t < inSamples_; ++t) 00109 { 00110 for (o=0; o < N2_; o++) 00111 { 00112 if (o==0) //DC bin (i.e. 0) 00113 { 00114 re_ = in(0,t); 00115 im_ = 0.0; 00116 } 00117 else if (o == N2_-1) //Nyquist bin (i.e. N/2) 00118 { 00119 re_ = in(1,t); 00120 im_ = 0.0; 00121 } 00122 else //all other bins 00123 { 00124 re_ = in(2*o, t); 00125 im_ = in(2*o+1, t); 00126 } 00127 00128 switch (ntype_) 00129 { 00130 case PSD_POWER: 00131 out(o, t) = re_*re_ + im_*im_; 00132 break; 00133 case PSD_MAG: 00134 out(o,t) = sqrt(re_ * re_ + im_ * im_); 00135 break; 00136 case PSD_DB: 00137 // TODO(sness) - Check validity of this magic number. Should probably be FLT_MIN 00138 dB_ = (mrs_real)(10*log10(re_ * re_ + im_ * im_ + 0.000000001)); 00139 out(o,t) = dB_; 00140 break; 00141 case PSD_WDB: 00142 // NOTE : FIXME!!! 00143 // 20*log10() seems to work better for toy_with_onsets 00144 // This is not good, and should be fixed inside the code that 00145 // looks for onsets. 00146 dB_ = (mrs_real)(20*log10(re_ * re_ + im_ * im_ + 0.000000001)); 00147 if (dB_ < -100) dB_ = -100; 00148 out(o,t) = dB_; 00149 break; 00150 case PSD_PD: 00151 pwr_ = re_ * re_ + im_ * im_; 00152 out(o,t) = (mrs_real)(2.0 * pwr_) / N2_; 00153 break; 00154 case PSD_LOGMAG2: 00155 out(o,t) = (mrs_real)log10(1+sqrt(re_ * re_ + im_ * im_)); 00156 break; 00157 case PSD_LOGMAG: 00158 out(o,t) = log(1+1000.0 * sqrt(re_ * re_ + im_ * im_)); 00159 // out(o,t) = pow(sqrt(re_ * re_ + im_ * im_), 0.5); 00160 // out(o,t) = asin(sqrt(re_ * re_ + im_ * im_)); 00161 00162 } 00163 } 00164 } 00165 00166 //MATLAB_PUT(out, "PowerSpectrum"); 00167 //MATLAB_EVAL("plot(PowerSpectrum)"); 00168 }