Marsyas
0.6.0-alpha
|
00001 00002 /* 00003 ** Copyright (C) 1998-2010 George Tzanetakis <gtzan@cs.uvic.ca> 00004 ** 00005 ** This program is free software; you can redistribute it and/or modify 00006 ** it under the terms of the GNU General Public License as published by 00007 ** the Free Software Foundation; either version 2 of the License, or 00008 ** (at your option) any later version. 00009 ** 00010 ** This program is distributed in the hope that it will be useful, 00011 ** but WITHOUT ANY WARRANTY; without even the implied warranty of 00012 ** MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 00013 ** GNU General Public License for more details. 00014 ** 00015 ** You should have received a copy of the GNU General Public License 00016 ** along with this program; if not, write to the Free Software 00017 ** Foundation, Inc., 59 Temple Place - Suite 330, Boston, MA 02111-1307, USA. 00018 */ 00019 00020 #include "LSP.h" 00021 #include "../common_source.h" 00022 00023 #include <marsyas/NumericLib.h> 00024 00025 #include <algorithm> 00026 #include <vector> 00027 00028 // #define _MATLAB_LSP_ 00029 00030 00031 using std::ostringstream; 00032 using std::vector; 00033 using std::polar; 00034 00035 using namespace Marsyas; 00036 00037 LSP::LSP(mrs_string name):MarSystem("LSP",name) 00038 { 00039 addControls(); 00040 } 00041 00042 LSP::~LSP() 00043 { 00044 } 00045 00046 MarSystem* 00047 LSP::clone() const 00048 { 00049 return new LSP(*this); 00050 } 00051 00052 void 00053 LSP::addControls() 00054 { 00055 addctrl("mrs_natural/order", (mrs_natural)10);//read-only control 00056 addctrl("mrs_real/gamma", (mrs_real)1.0); 00057 } 00058 00059 void 00060 LSP::myUpdate(MarControlPtr sender) 00061 { 00062 (void) sender; //suppress warning of unused parameter(s) 00063 MRSDIAG("LSP.cpp - LSP:myUpdate"); 00064 00065 order_ = getctrl("mrs_natural/inObservations")->to<mrs_natural>() - 2; 00066 setctrl("mrs_natural/order", order_);//read-only control 00067 00068 setctrl("mrs_natural/onObservations", order_); 00069 setctrl("mrs_natural/onSamples", getctrl("mrs_natural/inSamples")); 00070 setctrl("mrs_real/osrate", getctrl("mrs_real/israte")); 00071 00072 //LSP features names 00073 ostringstream oss; 00074 for (mrs_natural i = 0; i < order_; ++i) 00075 oss << "LSP_" << i+1 << ","; 00076 setctrl("mrs_string/onObsNames", oss.str()); 00077 } 00078 00079 void 00080 LSP::myProcess(realvec& in, realvec& out) 00081 { 00082 NumericLib numLib; 00083 00084 mrs_real gamma = getctrl("mrs_real/gamma")->to<mrs_real>(); 00085 vector<mrs_real> ak(order_); 00086 00087 if( gamma != 1.0) 00088 for(mrs_natural j = 0; j < order_ ; j++) 00089 { 00090 ak[j] = in(j)*pow((mrs_real)gamma, (mrs_real)j+1);//*(-1.0);//apply pole-shifting 00091 } 00092 else 00093 for(mrs_natural j = 0; j < order_ ; j++) 00094 { 00095 ak[j] = in(j);//*(-1.0); //no pole-shifting applied 00096 } 00097 00098 vector<mrs_complex> P(order_+2); 00099 vector<mrs_complex> Q(order_+2); 00100 vector<mrs_complex> Proots(order_+1); 00101 vector<mrs_complex> Qroots(order_+1); 00102 00103 P[order_+1] = polar(1.0, 0.0); 00104 Q[order_+1] = polar(1.0, 0.0); 00105 for(mrs_natural k = 0; k < order_; k++) 00106 { 00107 P[order_-k] = polar((double)(ak[k] + ak[order_-1-k]), 0.0); 00108 Q[order_-k] = polar((double)(ak[k] - ak[order_-1-k]), 0.0); 00109 } 00110 P[0] = polar(1.0, 0.0); 00111 Q[0] = polar(-1.0, 0.0); 00112 00113 if (!numLib.polyRoots(P, false, order_+1, Proots)) { 00114 //P has only real coefs => complexCoefs = false 00115 MRSERR("LSP::myProcess() - numerical error in polynomial root calculation!"); 00116 } 00117 if(!numLib.polyRoots(Q, false, order_+1, Qroots)) { 00118 //Q has only real coefs => complexCoefs = false 00119 MRSERR("LSP::myProcess() - numerical error in polynomial root calculation!"); 00120 } 00121 00122 mrs_real phase; 00123 vector<mrs_real> out_vec; 00124 for(mrs_natural k = 0; k <= order_; k++) 00125 { 00126 phase = arg(Proots[k]); 00127 if((phase > 0) && (phase < PI)) 00128 { 00129 out_vec.push_back(phase); 00130 } 00131 } 00132 for(mrs_natural k = 0; k <= order_; k++) 00133 { 00134 phase = arg(Qroots[k]); 00135 if((phase > 0) && (phase < PI)) 00136 { 00137 out_vec.push_back(phase); 00138 } 00139 } 00140 sort(out_vec.begin(), out_vec.end()); //sorts LSP freqs into ascending order 00141 00142 //output sorted LSP frequencies 00143 for(mrs_natural i = 0; i < order_; ++i) 00144 out(i) = out_vec[i]; 00145 00146 #ifdef _MATLAB_LSP_ 00147 MATLAB_PUT(order_, "LSP_order"); 00148 MATLAB_PUT(in, "LSP_in"); 00149 MATLAB_PUT(P, "LSP_P"); 00150 MATLAB_PUT(Q, "LSP_Q"); 00151 MATLAB_PUT(Proots, "LSP_Proots"); 00152 MATLAB_PUT(Qroots, "LSP_Qroots"); 00153 MATLAB_PUT(out_vec, "LSP_out1"); 00154 MATLAB_PUT(out, "LSP_out2"); 00155 MATLAB_EVAL("LSP_test(LSP_order, LSP_in, LSP_P, LSP_Q, LSP_Proots, LSP_Qroots, LSP_out1, LSP_out2);"); 00156 #endif 00157 }