Marsyas
0.6.0-alpha
|
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 "Peaker.h" 00020 #include "../common_source.h" 00021 #include <algorithm> 00022 00023 using std::ostringstream; 00024 using std::min; 00025 using std::max; 00026 using std::cout; 00027 using std::endl; 00028 00029 00030 using namespace Marsyas; 00031 00032 Peaker::Peaker(mrs_string name):MarSystem("Peaker",name) 00033 { 00034 addControls(); 00035 } 00036 00037 Peaker::~Peaker() 00038 { 00039 } 00040 00041 MarSystem* 00042 Peaker::clone() const 00043 { 00044 return new Peaker(*this); 00045 } 00046 00047 00048 00049 void 00050 Peaker::addControls() 00051 { 00052 addctrl("mrs_real/peakSpacing", 0.0); 00053 addctrl("mrs_real/peakStrength", 0.0); 00054 addctrl("mrs_real/peakStrengthRelMax", 0.0); 00055 addctrl("mrs_real/peakStrengthRelThresh", 0.0); 00056 addctrl("mrs_real/peakStrengthAbs", 0.0); 00057 addctrl("mrs_real/peakStrengthThreshLpParam", 0.95); 00058 addctrl("mrs_natural/peakStart", (mrs_natural)0); 00059 addctrl("mrs_natural/peakEnd", (mrs_natural)0); 00060 addctrl("mrs_natural/interpolation", (mrs_natural)0); 00061 //addctrl("mrs_real/peakGain", 1.0); 00062 addctrl("mrs_bool/peakHarmonics", false); 00063 addctrl("mrs_bool/rmsNormalize", false); 00064 addctrl("mrs_natural/peakNeighbors", 2); 00065 } 00066 00067 00068 void 00069 Peaker::myUpdate(MarControlPtr sender) 00070 { 00071 (void) sender; //suppress warning of unused parameter(s) 00072 MRSDIAG("Peaker.cpp - Peaker:myUpdate"); 00073 00074 MarSystem::myUpdate (sender); 00075 00076 lpThresh_.stretch (getctrl ("mrs_natural/inSamples")->to<mrs_natural>()); 00077 } 00078 00079 void 00080 Peaker::myProcess(realvec& in, realvec& out) 00081 { 00082 mrs_natural t,o; 00083 00084 const mrs_natural peakSpacing = (mrs_natural)(inSamples_ *getctrl("mrs_real/peakSpacing")->to<mrs_real>() + .5); 00085 mrs_real peakStrengthRelRms, 00086 peakStrengthRelMax, 00087 peakStrengthRelThresh, 00088 peakStrengthAbs; 00089 //mrs_real peakGain; 00090 mrs_bool peakHarmonics; 00091 mrs_bool rmsNormalize; 00092 00093 mrs_natural peakStart; 00094 mrs_natural peakEnd; 00095 mrs_natural interpolationMode; 00096 mrs_natural peakNeighbors; 00097 00098 00099 peakStrengthRelRms = getctrl("mrs_real/peakStrength")->to<mrs_real>(); 00100 peakStrengthRelMax = getctrl("mrs_real/peakStrengthRelMax")->to<mrs_real>(); 00101 peakStrengthRelThresh = getctrl("mrs_real/peakStrengthRelThresh")->to<mrs_real>(); 00102 lpCoeff_ = getctrl("mrs_real/peakStrengthThreshLpParam")->to<mrs_real>(); 00103 peakStrengthAbs = getctrl("mrs_real/peakStrengthAbs")->to<mrs_real>(); 00104 peakStart = getctrl("mrs_natural/peakStart")->to<mrs_natural>(); 00105 peakEnd = getctrl("mrs_natural/peakEnd")->to<mrs_natural>(); 00106 interpolationMode = getctrl("mrs_natural/interpolation")->to<mrs_natural>(); 00107 //peakGain = getctrl("mrs_real/peakGain")->to<mrs_real>(); 00108 peakHarmonics = getctrl("mrs_bool/peakHarmonics")->to<mrs_bool>(); 00109 rmsNormalize = getctrl("mrs_bool/rmsNormalize")->to<mrs_bool>(); 00110 peakNeighbors = getctrl("mrs_natural/peakNeighbors")->to<mrs_natural>(); 00111 00112 00113 00114 if (peakEnd == 0) 00115 peakEnd = inSamples_; 00116 // FIXME This line defines an unused variable 00117 // mrs_real srate = getctrl("mrs_real/israte")->to<mrs_real>(); 00118 00119 out.setval(0.0); 00120 00121 00122 //peakStrengthRelRms = 0.0; 00123 00124 00125 00126 for (o = 0; o < inObservations_; o++) 00127 { 00128 rms_ = 0.0; 00129 max_ = -1e37; 00130 00131 for (t=peakStart; t < peakEnd; t++) 00132 { 00133 rms_ += in(o,t) * in(o,t); 00134 if (max_ < in(o,t)) 00135 max_ = in(o,t); 00136 } 00137 if (rms_ != 0.0) 00138 rms_ /= (peakEnd - peakStart); 00139 rms_ = sqrt(rms_); 00140 00141 mrs_real max; 00142 mrs_natural maxIndex; 00143 00144 bool peakFound = false; 00145 00146 if (peakStrengthRelThresh > .0) 00147 { 00148 in.getRow (o,lpThresh_); 00149 compLpThresh (lpThresh_, lpThresh_); // do it inplace to avoid another copy... 00150 } 00151 00152 for (t=peakStart; t < peakEnd; t++) 00153 { 00154 peakFound = true; 00155 00156 // peak has to be larger than neighbors 00157 for (int j = 1; j < peakNeighbors; j++) 00158 { 00159 mrs_natural index=t-j; 00160 if (index<0) index=0; 00161 if (in(o,index) >= in(o,t)) 00162 { 00163 peakFound = false; 00164 break; 00165 } 00166 index=t+j; 00167 if (index>=inSamples_) index=inSamples_-1; 00168 if (in(o,index) >= in(o,t)) 00169 { 00170 peakFound = false; 00171 break; 00172 } 00173 } 00174 00175 if (peakFound) 00176 { 00177 currThresh_ = lpThresh_(t); 00178 peakFound = doThresholding (in(o,t), peakStrengthRelRms, peakStrengthRelMax, peakStrengthRelThresh, peakStrengthAbs); 00179 } 00180 00181 if (peakFound) 00182 { 00183 // check for another peak in the peakSpacing area 00184 max = in(o,t); 00185 maxIndex = t; 00186 00187 00188 for (int j=0; j < peakSpacing; j++) 00189 { 00190 if (t+j < peakEnd-1) 00191 if (in(o,t+j) > max) 00192 { 00193 max = in(o,t+j); 00194 maxIndex = t+j; 00195 } 00196 } 00197 00198 t += peakSpacing; 00199 00200 if (rmsNormalize) 00201 { 00202 out(o,maxIndex) = in(o,maxIndex) / rms_; 00203 if(interpolationMode && maxIndex > 0 && maxIndex < inSamples_) 00204 { 00205 out(o,maxIndex-1) = in(o,maxIndex-1) /rms_; 00206 out(o,maxIndex+1) = in(o,maxIndex+1) / rms_; 00207 } 00208 } 00209 else 00210 { 00211 out(o,maxIndex) = in(o,maxIndex); 00212 if(interpolationMode && maxIndex > 0 && maxIndex < inSamples_) 00213 { 00214 out(o,maxIndex-1) = in(o,maxIndex-1); 00215 out(o,maxIndex+1) = in(o,maxIndex+1); 00216 } 00217 } 00218 00219 00220 // peakNeighbors = 0; 00221 // rms_ = 1.0; 00222 00223 00224 if (peakHarmonics) 00225 { 00226 twice_ = 2 * maxIndex; 00227 half_ = (mrs_natural) (0.5 * maxIndex + 0.5); 00228 triple_ = 3 * maxIndex; 00229 third_ = (mrs_natural) (0.33 * maxIndex + 0.5); 00230 00231 if (twice_ < (peakEnd - peakStart)) 00232 { 00233 peakFound = true; 00234 00235 // peak has to be larger than neighbors 00236 for (int j = 1; j < peakNeighbors; j++) 00237 { 00238 if (in(o,twice_-j) >= in(o,twice_)) 00239 { 00240 peakFound = false; 00241 break; 00242 } 00243 if (in(o,twice_+j) >= in(o,twice_)) 00244 { 00245 peakFound = false; 00246 break; 00247 } 00248 } 00249 00250 00251 if (peakFound) 00252 { 00253 out(o, maxIndex) += (in(o, twice_)/rms_); 00254 out(o, twice_) = in(o,twice_)/rms_ + 0.5 * out(o, maxIndex); 00255 } 00256 00257 } 00258 00259 if (half_ < (peakEnd - peakStart)) 00260 { 00261 peakFound = true; 00262 00263 // peak has to be larger than neighbors 00264 for (int j = 1; j < peakNeighbors; j++) 00265 { 00266 if (in(o,half_-j) >= in(o,half_)) 00267 { 00268 peakFound = false; 00269 break; 00270 } 00271 if (in(o,half_+j) >= in(o,half_)) 00272 { 00273 peakFound = false; 00274 break; 00275 } 00276 } 00277 00278 00279 if (peakFound) 00280 { 00281 out(o, maxIndex) += (in(o, half_)/rms_); 00282 out(o, half_) = in(o,half_)/rms_ + 0.5 * out(o, maxIndex); 00283 } 00284 00285 } 00286 00287 00288 if (triple_ < (peakEnd - peakStart)) 00289 { 00290 peakFound = true; 00291 00292 // peak has to be larger than neighbors 00293 for (int j = 1; j < peakNeighbors; j++) 00294 { 00295 if (in(o,triple_-j) >= in(o,triple_)) 00296 { 00297 peakFound = false; 00298 break; 00299 } 00300 if (in(o,triple_+j) >= in(o,triple_)) 00301 { 00302 peakFound = false; 00303 break; 00304 } 00305 } 00306 00307 00308 if (peakFound) 00309 { 00310 out(o, maxIndex) += (in(o, triple_)/rms_); 00311 out(o, triple_) = in(o,triple_)/rms_ + 0.5 * out(o, maxIndex); 00312 } 00313 00314 } 00315 00316 if (third_ < (peakEnd - peakStart)) 00317 { 00318 peakFound = true; 00319 00320 // peak has to be larger than neighbors 00321 00322 for (int j = 1; j < peakNeighbors; j++) 00323 { 00324 if (in(o,third_-j) >= in(o,third_)) 00325 { 00326 peakFound = false; 00327 break; 00328 } 00329 if (in(o,third_+j) >= in(o,triple_)) 00330 { 00331 peakFound = false; 00332 break; 00333 } 00334 } 00335 00336 if (peakFound) 00337 { 00338 out(o, maxIndex) += (in(o, third_)/rms_); 00339 out(o, third_) = in(o,third_)/rms_ + 0.5 * out(o, maxIndex); 00340 } 00341 00342 } 00343 } 00344 00345 peakFound = true; 00346 } 00347 00348 } 00349 } 00350 } 00351 00352 mrs_bool Peaker::doThresholding (mrs_real input, mrs_real peakStrengthRelRms, mrs_real peakStrengthRelMax, mrs_real peakStrengthRelThresh, mrs_real peakStrengthAbs) 00353 { 00354 mrs_real thresh = max (max (max (peakStrengthAbs, peakStrengthRelRms * rms_), peakStrengthRelMax * max_), peakStrengthRelThresh * currThresh_); 00355 00356 00357 if (input <= thresh) 00358 return false; 00359 else 00360 return true; 00361 } 00362 00363 void Peaker::compLpThresh (const mrs_realvec input, mrs_realvec &output) 00364 { 00365 mrs_natural i, 00366 len = input.getCols (); 00367 mrs_real buf = input(0); 00368 for (i = 0; i < len; ++i) 00369 { 00370 buf = input(i) + lpCoeff_ * (buf - input(i)); 00371 output(i) = buf; 00372 } 00373 // do reverse filtering to ensure zero group delay 00374 for (i = len-1; i >= 0; i--) 00375 { 00376 buf = output(i) + lpCoeff_ * (buf - output(i)); 00377 output(i) = buf; 00378 } 00379 }