Marsyas  0.6.0-alpha
/usr/src/RPM/BUILD/marsyas-0.6.0/src/marsyas/marsystems/EnhADRess.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 "EnhADRess.h"
00020 #include "../common_source.h"
00021 
00022 using std::ostringstream;
00023 using std::abs;
00024 
00025 using namespace Marsyas;
00026 
00027 EnhADRess::EnhADRess(mrs_string name):MarSystem("EnhADRess", name)
00028 {
00029   addControls();
00030 }
00031 
00032 EnhADRess::EnhADRess(const EnhADRess& a) : MarSystem(a)
00033 {
00034 
00035 }
00036 
00037 EnhADRess::~EnhADRess()
00038 {
00039 }
00040 
00041 MarSystem*
00042 EnhADRess::clone() const
00043 {
00044   return new EnhADRess(*this);
00045 }
00046 
00047 void
00048 EnhADRess::addControls()
00049 {
00050 }
00051 
00052 void
00053 EnhADRess::myUpdate(MarControlPtr sender)
00054 {
00055   MRSDIAG("EnhADRess.cpp - EnhADRess:myUpdate");
00056   (void) sender;  //suppress warning of unused parameter(s)
00057 
00058   N2_ = inObservations_ / 2; //i.e. we get two vertically stacked spectrums at the input
00059   N4_ = N2_/2 + 1; //i.e. for each spectrum, we have N/2+1 spectral bins
00060 
00061   ctrl_onSamples_->setValue(ctrl_inSamples_, NOUPDATE);
00062   ctrl_onObservations_->setValue(N4_*3, NOUPDATE); //output data for N/2+1 spectral bins, stacked vertically for Mag, phase and pan
00063   ctrl_osrate_->setValue(ctrl_israte_, NOUPDATE);
00064 
00065   ostringstream oss;
00066   for(mrs_natural n=0; n< N4_; ++n)
00067   {
00068     oss << "EnhADRess_Mag_bin_" << n <<",";
00069   }
00070   for(mrs_natural n=0; n< N4_; ++n)
00071   {
00072     oss << "EnhADRess_Phase_bin_" << n <<",";
00073   }
00074   for(mrs_natural n=0; n< N4_; ++n)
00075   {
00076     oss << "EnhADRess_Pan_bin_" << n <<",";
00077   }
00078   ctrl_onObsNames_->setValue(oss.str(), NOUPDATE);
00079 }
00080 
00081 void
00082 EnhADRess::myProcess(realvec& in, realvec& out)
00083 {
00084   out.setval(0.0);
00085 
00086   for(mrs_natural t=0; t < inSamples_; ++t)
00087   {
00088     for (mrs_natural k=0; k < N4_; k++)
00089     {
00090       //get left channel spectrum bin
00091       if (k==0) //DC bin (i.e. 0)
00092       {
00093         rel_ = in(0,t);
00094         iml_ = 0.0;
00095       }
00096       else if (k == N4_-1) //Nyquist bin (i.e. N/2)
00097       {
00098         rel_ = in(1, t);
00099         iml_ = 0.0;
00100       }
00101       else //all other bins
00102       {
00103         rel_ = in(2*k, t);
00104         iml_ = in(2*k+1, t);
00105       }
00106 
00107       //get right channel spectrum bin
00108       if (k==0) //DC bin (i.e. 0)
00109       {
00110         rer_ = in(N2_,t);
00111         imr_ = 0.0;
00112       }
00113       else if (k == N4_-1) //Nyquist bin (i.e. N/2)
00114       {
00115         rer_ = in(N2_+1, t);
00116         imr_ = 0.0;
00117       }
00118       else //all other bins
00119       {
00120         rer_ = in(N2_ + 2*k, t);
00121         imr_ = in(N2_ + 2*k+1, t);
00122       }
00123 
00124       phaseL_ = atan2(iml_, rel_);      //left channel phases
00125       phaseR_ = atan2(imr_, rer_); //right channel phases
00126 
00127       deltaPhase_ = abs(phaseL_ - phaseR_);
00128 
00129       //wrap phase into the 0~2*PI range
00130       deltaPhase_ = (mrs_real)fmod((double)deltaPhase_, (double)2*PI);
00131 
00132       //left amplitude value
00133       Lk_ = sqrt(rel_*rel_ + iml_*iml_);
00134 
00135       //right amplitude value
00136       Rk_ = sqrt(rer_*rer_ + imr_*imr_);
00137 
00138       if(deltaPhase_ < PI/2)
00139       {
00140         minLk_ = Lk_ * sin(deltaPhase_);
00141         minRk_ = Rk_ * sin(deltaPhase_);
00142 
00143         if(Lk_ < Rk_) // i.e. minLk < minRk --> sound panned right
00144         {
00145           //store magnitude to output
00146           mrs_real mag = Rk_ - minLk_;
00147           // bins with amplitude inferior to -100dB are discarded
00148           if(20.0*log10(mag*mag+0.000000001) > -100.0)
00149           {
00150             out(k,t) = mag;
00151 
00152             //store phase to output
00153             out(k+N4_,t) = phaseR_;
00154 
00155             //store azimuth to output
00156             out(k+N4_*2,t) = 1.0 - Lk_ * cos(deltaPhase_) / Rk_ ; //-1.0-> L; 0.0-> C;
00157           }
00158         }
00159         else if(Lk_ > Rk_) // i.e. minLk > minRk --> sound panned left
00160         {
00161           //store magnitude to output
00162           mrs_real mag = Lk_ - minRk_;
00163           // bins with amplitude inferior to -100dB are discarded
00164           if(20.0*log10(mag*mag+0.000000001) > -100.0)
00165           {
00166             out(k,t) = mag;
00167 
00168             //store phase to output
00169             out(k+N4_,t) = phaseL_;
00170 
00171             //store azimuth to output
00172             out(k+N4_*2,t) = Rk_ * cos(deltaPhase_) / Lk_ -1.0; //0.0 -> C; 1.0-> R;
00173           }
00174         }
00175         else if(Lk_ == Rk_) //sound panned at the CENTER
00176         {
00177           //store magnitude to output
00178           mrs_real mag = Lk_ - minRk_;
00179           // bins with amplitude inferior to -100dB are discarded
00180           if(20.0*log10(mag*mag+0.000000001) > -100.0)
00181           {
00182             out(k,t) = mag;
00183 
00184             //store phase to output
00185             out(k+N4_,t) = phaseL_; //could have been phaseR_...
00186 
00187             //store azimuth to output
00188             out(k+N4_*2,t) = 0.0; //0.0 -> C;
00189           }
00190         }
00191       }
00192       else
00193       {
00194         if(20.0*log10(Lk_*Lk_+0.000000001) < -100.0)
00195           Lk_ = 0.0;
00196         if(20.0*log10(Rk_*Rk_+0.000000001) < -100.0)
00197           Rk_ = 0.0;
00198 
00199         if(Lk_ > Rk_)
00200         {
00201           out(k,t) = Lk_;
00202           out(k+N4_,t) = phaseL_;
00203           out(k+N4_*2,t) = 0.0; //-1.0;
00204         }
00205         else if(Rk_ > Lk_)
00206         {
00207           out(k,t) = Rk_;
00208           out(k+N4_,t) = phaseR_;
00209           out(k+N4_*2,t) = 0.0; //1.0;
00210         }
00211         else if(Lk_ == Rk_ && Lk_ != 0.0)
00212         {
00213           out(k,t) = Lk_;
00214           out(k+N4_,t) = phaseL_;
00215           out(k+N4_*2,t) = 0.0;
00216         }
00217       }
00218     }
00219   }
00220   //MATLAB_PUT(out, "out");
00221   //MATLAB_EVAL("figure(1);plot(out);figure(2)");
00222   //MATLAB_EVAL("plot(out(length(out)/3*2+1:end))")
00223 }