CLAM-Development  1.4.0
SynthSineSpectrum.cxx
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2004 MUSIC TECHNOLOGY GROUP (MTG)
3  * UNIVERSITAT POMPEU FABRA
4  *
5  *
6  * This program is free software; you can redistribute it and/or modify
7  * it under the terms of the GNU General Public License as published by
8  * the Free Software Foundation; either version 2 of the License, or
9  * (at your option) any later version.
10  *
11  * This program is distributed in the hope that it will be useful,
12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
14  * GNU General Public License for more details.
15  *
16  * You should have received a copy of the GNU General Public License
17  * along with this program; if not, write to the Free Software
18  * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19  *
20  */
21 
22 #include "SpectrumConfig.hxx"
23 #include "SpecTypeFlags.hxx"
24 #include "SynthSineSpectrum.hxx"
25 #include "BlackHarrisLobe.hxx"
26 
27 using namespace CLAM;
28 
29 
30 /******************************************************************/
31 /************************ SynthSineSpectrum ***********************/
32 
33 
34 
36  : mInput( "Input", this ),
37  mOutput("Output", this )
38 {
40 }
41 
43  : mInput( "Input", this ),
44  mOutput("Output", this )
45 {
46  Configure(cfg);
47 }
48 
49 
51 {
52  CopyAsConcreteConfig(mConfig, c);
54  wcfg.SetNormalize(EWindowNormalize::eNone);
55  wcfg.SetSize(MAINLOBE_TABLE_SIZE);
56  wcfg.SetType(EWindowType::eBlackmanHarris92);
57  mWndGen.Configure(wcfg);
58  //@TODO: there is some kind of bug in the window generator: when generating the
59  //BlackmannHarris transform it does not yield the same results as the above table!
60  //mBlackHarris92TransMainLobe.SetSize(MAINLOBE_TABLE_SIZE);
61  //mWndGen.Do(mBlackHarris92TransMainLobe);
62  return true;
63 }
64 
65 
67 {
68 }
69 
71 {
72  bool result = Do( mInput.GetData(), mOutput.GetData());
73  mInput.Consume();
74  mOutput.Produce();
75  return result;
76 }
77 
78 
79 //-------------------------------------------------------------------//
80 // note : for now we first fill the mSynthsineSpectrum. Then we convert this
81 // to a Spectrum. Its more cpu expensive, but as we dont know which processes
82 // has to be applied the mSynthSineSpectrum, we use it like this
83 // COULD BE OPTIMIZED LATER !!!! JO
84 bool SynthSineSpectrum::Do(const SpectralPeakArray& peakArray,Spectrum& residualSpectrumOut, double gain)
85 {
86  CLAM_DEBUG_ASSERT( AbleToExecute(), "SynthSineSpectrum::Do - processing is not running" );
87 
88  InitSynthSpec(mConfig.GetSpectrumSize()); // could be optimised with memset
89  FillSynthSineSpectrum(peakArray,gain);
90 
91  SpectrumConfig Scfg;
92  Scfg.SetScale(EScale::eLinear);
93  SpecTypeFlags sflags;
94  sflags.bComplex = 1;
95  sflags.bPolar = 0;
96  sflags.bMagPhase = 0;
97  sflags.bMagPhaseBPF = 0;
98  int spectralSize;
99  spectralSize = mSynthSineSpectrum.Size();
100  Scfg.SetType(sflags);
101  Scfg.SetSize(spectralSize);
102  Scfg.SetSpectralRange(residualSpectrumOut.GetSpectralRange());
103  residualSpectrumOut.Configure(Scfg);
104  residualSpectrumOut.SetComplexArray(mSynthSineSpectrum);
105  return true;
106 }
107 
108 void SynthSineSpectrum::FillSynthSineSpectrum ( const SpectralPeakArray& peakArray, double gain)
109 {
110  TSize mainLobeBins = TSize ( TData(8)*CLAM_pow(2.0,(double)mConfig.GetZeroPadding()) ) ; // Number of bins in the mainlope of a transformed BHarris92 window
111 
112  TIndex b,i,k,l,Index,numPeaks = peakArray.GetnPeaks();
113 
114  SpectralPeak currPeak;
115  double currMag,currBinPos,Sin,Cos,approxMag,fIndex,currFreq,phase;
116  Array< Complex > SpecPeakEnvelope(mainLobeBins);
117  SpecPeakEnvelope.SetSize(mainLobeBins);
118  DataArray& peakFreqBuffer=peakArray.GetFreqBuffer();
119  DataArray& peakMagBuffer=peakArray.GetMagBuffer();
120  DataArray& peakPhaseBuffer=peakArray.GetPhaseBuffer();
121 
122  TData samplingRate=mConfig.GetSamplingRate();
123  TSize spectrumSize=mConfig.GetSpectrumSize();
124  TSize zeroPadding=TSize( TData(mConfig.GetZeroPadding()) );
125 
126  TData binPosFactor=2*(spectrumSize-1)/samplingRate;
127  TData firstBinFactor= 3* CLAM_pow((TData)2,(TData)zeroPadding);
128 
129  int incr=int(CLAM_pow(2.0,9.0-mConfig.GetZeroPadding()));
130  int lastBin=4096-incr;
131  /* loop thru all the peaks ...*/
132  for(i=0;i<numPeaks;i++)
133  {
134  if(peakArray.GetScale()==EScale::eLog)//if it is in dB
135  {
136  currMag = CLAM_pow(10.0,(peakMagBuffer[i]/20.0));
137  }
138  else currMag = peakMagBuffer[i];
139  currFreq = peakFreqBuffer[i];
140  phase = peakPhaseBuffer[i];
141  Cos = CLAM_cos(phase); // we assume the phase is constant throughout the frame
142  Sin = CLAM_sin(phase);
143  // we use the frequency specify position of the peak, not the binpos !
144  currBinPos = currFreq*binPosFactor;
145 
146  double Loc = currBinPos;
147  TIndex firstBin = TIndex( int( Loc ) - firstBinFactor );
148  double BinRemainder = Loc - floor (Loc);
149  fIndex = (1.0 - BinRemainder);
150  Index=(int)(0.5 + fIndex*MAINLOBE_TABLE_SIZE/(mainLobeBins));
151 
152  /* SpecPeakEnvelope holds the 7*(2^ZeroPadding) Bin approximation of a Peak which is */
153  /* convoluted with the transform of a BlackmanHarris 92 window Main Lobe */
154  int n;
155  int currBin=Index;
156  double mainLobe = sBlackHarris92TransMainLobe[currBin];
157  for (n=0; currBin<=lastBin;n++)
158  {
159  approxMag = currMag * mainLobe;
160 
161  SpecPeakEnvelope[n].SetReal(TData(approxMag*Cos));
162  SpecPeakEnvelope[n].SetImag(TData(approxMag*Sin));
163  currBin+=incr;
164  mainLobe = sBlackHarris92TransMainLobe[currBin];
165  }
166 
167  if( (Index+4096-incr) >= MAINLOBE_TABLE_SIZE)
168  {
169  SpecPeakEnvelope[mainLobeBins-1].SetReal(0);
170  SpecPeakEnvelope[mainLobeBins-1].SetImag(0);
171  }
172  else
173  {
174  approxMag = currMag * sBlackHarris92TransMainLobe[Index+4096-incr];
175  SpecPeakEnvelope[mainLobeBins-1].SetReal(TData(approxMag*Cos));
176  SpecPeakEnvelope[mainLobeBins-1].SetImag(TData(approxMag*Sin));
177  }
178 
179  /* Add the current Peak to the spektrum and check for aliasing at the limits */
180  for (k = 0, l = firstBin; k < mainLobeBins ; k++,l++)
181  {
182  if (l > 0 && l < mConfig.GetSpectrumSize()-1) // we are inside the spectrum limits
183  {
184  mSynthSineSpectrum[l] += SpecPeakEnvelope[k];
185  }
186 
187  else if (l < 0) // part of peaklobe is below zero frequency, mirror it back
188  {
189  b = -l;
190  mSynthSineSpectrum[b].SetReal(mSynthSineSpectrum[b].Real()+SpecPeakEnvelope[k].Real());
191  mSynthSineSpectrum[b].SetImag(mSynthSineSpectrum[b].Imag()-SpecPeakEnvelope[k].Imag());
192  }
193 
194  else if (l == 0) // only real part, no phase
195  {
196  mSynthSineSpectrum[0].SetReal(mSynthSineSpectrum[0].Real()+2*(SpecPeakEnvelope[k].Real()));
197  }
198 
199  else if (l > mConfig.GetSpectrumSize()-1&&(l<mConfig.GetSpectrumSize()*2-1)) // part of peaklobe is over nyquist frequency, mirror it back
200  {
201  b = 2 * (spectrumSize-1) - l;
202  mSynthSineSpectrum[b].SetReal(mSynthSineSpectrum[b].Real()+SpecPeakEnvelope[k].Real());
203  mSynthSineSpectrum[b].SetImag(mSynthSineSpectrum[b].Imag()-SpecPeakEnvelope[k].Imag());
204  }
205 
206  else if (l == spectrumSize-1) // only real part, no phase
207  {
208  mSynthSineSpectrum[l].SetReal(mSynthSineSpectrum[l].Real()+2*(SpecPeakEnvelope[k].Real()));
209  }
210  }
211  }
212 }
213 //-------------------------------------------------------------------//
214 void SynthSineSpectrum::InitSynthSpec(TSize size)
215 {
216  TSize i;
217  mSynthSineSpectrum.Resize(size);
218  mSynthSineSpectrum.SetSize(size);
219  for(i=0;i<size;i++)
220  {
221  mSynthSineSpectrum[i] = Complex(0,0);
222  }
223 }
224