MythTV  master
el_processor.cpp
Go to the documentation of this file.
1 /*
2 Copyright (C) 2007 Christian Kothe
3 
4 This program is free software; you can redistribute it and/or
5 modify it under the terms of the GNU General Public License
6 as published by the Free Software Foundation; either version 2
7 of the License, or (at your option) any later version.
8 
9 This program is distributed in the hope that it will be useful,
10 but WITHOUT ANY WARRANTY; without even the implied warranty of
11 MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
12 GNU General Public License for more details.
13 
14 You should have received a copy of the GNU General Public License
15 along with this program; if not, write to the Free Software
16 Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301, USA.
17 */
18 
19 #include "el_processor.h"
20 #include <cmath>
21 #include <complex>
22 #include <cstdlib>
23 #include <cstring>
24 #include <vector>
25 #ifdef USE_FFTW3
26 #include "fftw3.h"
27 #else
28 extern "C" {
29 #include "libavcodec/avfft.h"
30 #include "libavcodec/fft.h"
31 }
32 using FFTComplexArray = FFTSample[2];
33 #endif
34 
35 
36 #if defined(_WIN32) && defined(USE_FFTW3)
37 #pragma comment (lib,"libfftw3f-3.lib")
38 #endif
39 
40 using cfloat = std::complex<float>;
41 
42 static const float PI = 3.141592654;
43 static const float epsilon = 0.000001;
44 static const float center_level = 0.5*sqrt(0.5);
45 
46 // private implementation of the surround decoder
47 class decoder_impl {
48 public:
49  // create an instance of the decoder
50  // blocksize is fixed over the lifetime of this object for performance reasons
51  explicit decoder_impl(unsigned blocksize=8192): m_n(blocksize), m_halfN(blocksize/2) {
52 #ifdef USE_FFTW3
53  // create FFTW buffers
54  m_lt = (float*)fftwf_malloc(sizeof(float)*m_n);
55  m_rt = (float*)fftwf_malloc(sizeof(float)*m_n);
56  m_dst = (float*)fftwf_malloc(sizeof(float)*m_n);
57  m_dftL = (fftwf_complex*)fftwf_malloc(sizeof(fftwf_complex)*m_n);
58  m_dftR = (fftwf_complex*)fftwf_malloc(sizeof(fftwf_complex)*m_n);
59  m_src = (fftwf_complex*)fftwf_malloc(sizeof(fftwf_complex)*m_n);
60  m_loadL = fftwf_plan_dft_r2c_1d(m_n, m_lt, m_dftL,FFTW_MEASURE);
61  m_loadR = fftwf_plan_dft_r2c_1d(m_n, m_rt, m_dftR,FFTW_MEASURE);
62  m_store = fftwf_plan_dft_c2r_1d(m_n, m_src, m_dst,FFTW_MEASURE);
63 #else
64  // create lavc fft buffers
65  m_lt = (float*)av_malloc(sizeof(FFTSample)*m_n);
66  m_rt = (float*)av_malloc(sizeof(FFTSample)*m_n);
67  m_dftL = (FFTComplexArray*)av_malloc(sizeof(FFTComplex)*m_n*2);
68  m_dftR = (FFTComplexArray*)av_malloc(sizeof(FFTComplex)*m_n*2);
69  m_src = (FFTComplexArray*)av_malloc(sizeof(FFTComplex)*m_n*2);
70  m_fftContextForward = (FFTContext*)av_malloc(sizeof(FFTContext));
71  memset(m_fftContextForward, 0, sizeof(FFTContext));
72  m_fftContextReverse = (FFTContext*)av_malloc(sizeof(FFTContext));
73  memset(m_fftContextReverse, 0, sizeof(FFTContext));
74  ff_fft_init(m_fftContextForward, 13, 0);
75  ff_fft_init(m_fftContextReverse, 13, 1);
76 #endif
77  // resize our own buffers
78  m_frontR.resize(m_n);
79  m_frontL.resize(m_n);
80  m_avg.resize(m_n);
81  m_surR.resize(m_n);
82  m_surL.resize(m_n);
83  m_trueavg.resize(m_n);
84  m_xFs.resize(m_n);
85  m_yFs.resize(m_n);
86  m_inbuf[0].resize(m_n);
87  m_inbuf[1].resize(m_n);
88  for (unsigned c=0;c<6;c++) {
89  m_outbuf[c].resize(m_n);
90  m_filter[c].resize(m_n);
91  }
92  sample_rate(48000);
93  // generate the window function (square root of hann, b/c it is applied before and after the transform)
94  m_wnd.resize(m_n);
95  for (unsigned k=0;k<m_n;k++)
96  m_wnd[k] = sqrt(0.5*(1-std::cos(2*PI*k/m_n))/m_n);
97  m_currentBuf = 0;
98  memset(m_inbufs, 0, sizeof(m_inbufs));
99  memset(m_outbufs, 0, sizeof(m_outbufs));
100  // set the default coefficients
101  surround_coefficients(0.8165,0.5774);
102  phase_mode(0);
103  separation(1,1);
104  steering_mode(true);
105  }
106 
107  // destructor
109 #ifdef USE_FFTW3
110  // clean up the FFTW stuff
111  fftwf_destroy_plan(m_store);
112  fftwf_destroy_plan(m_loadR);
113  fftwf_destroy_plan(m_loadL);
114  fftwf_free(m_src);
115  fftwf_free(m_dftR);
116  fftwf_free(m_dftL);
117  fftwf_free(m_dst);
118  fftwf_free(m_rt);
119  fftwf_free(m_lt);
120 #else
121  ff_fft_end(m_fftContextForward);
122  ff_fft_end(m_fftContextReverse);
123  av_free(m_src);
124  av_free(m_dftR);
125  av_free(m_dftL);
126  av_free(m_rt);
127  av_free(m_lt);
130 #endif
131  }
132 
133  float ** getInputBuffers()
134  {
137  return m_inbufs;
138  }
139 
140  float ** getOutputBuffers()
141  {
148  return m_outbufs;
149  }
150 
151  // decode a chunk of stereo sound, has to contain exactly blocksize samples
152  // center_width [0..1] distributes the center information towards the front left/right channels, 1=full distribution, 0=no distribution
153  // dimension [0..1] moves the soundfield backwards, 0=front, 1=side
154  // adaption_rate [0..1] determines how fast the steering gets adapted, 1=instantaneous, 0.1 = very slow adaption
155  void decode(float center_width, float dimension, float adaption_rate) {
156  // process first part
157  int index = m_currentBuf*m_halfN;
158  float *in_second[2] = {&m_inbuf[0][index],&m_inbuf[1][index]};
159  m_currentBuf ^= 1;
160  index = m_currentBuf*m_halfN;
161  float *in_first[2] = {&m_inbuf[0][index],&m_inbuf[1][index]};
162  add_output(in_first,in_second,center_width,dimension,adaption_rate,true);
163  // shift last half of input buffer to the beginning
164  }
165 
166  // flush the internal buffers
167  void flush() {
168  for (unsigned k=0;k<m_n;k++) {
169  for (unsigned c=0;c<6;c++)
170  m_outbuf[c][k] = 0;
171  m_inbuf[0][k] = 0;
172  m_inbuf[1][k] = 0;
173  }
174  }
175 
176  // set lfe filter params
177  void sample_rate(unsigned int srate) {
178  // lfe filter is just straight through band limited
179  unsigned int cutoff = (30*m_n)/srate;
180  for (unsigned f=0;f<=m_halfN;f++) {
181  if (f<cutoff)
182  m_filter[5][f] = 0.5*sqrt(0.5);
183  else
184  m_filter[5][f] = 0.0;
185  }
186  }
187 
188  // set the assumed surround mixing coefficients
189  void surround_coefficients(float a, float b) {
190  // calc the simple coefficients
191  m_surroundHigh = a;
192  m_surroundLow = b;
193  m_surroundBalance = (a-b)/(a+b);
194  m_surroundLevel = 1/(a+b);
195  // calc the linear coefficients
196  cfloat i(0,1);
197  cfloat u((a+b)*i);
198  cfloat v((b-a)*i);
199  cfloat n(0.25,0);
200  cfloat o(1,0);
201  m_a = (v-o)*n; m_b = (o-u)*n; m_c = (-o-v)*n; m_d = (o+u)*n;
202  m_e = (o+v)*n; m_f = (o+u)*n; m_g = (o-v)*n; m_h = (o-u)*n;
203  }
204 
205  // set the phase shifting mode
206  void phase_mode(unsigned mode) {
207  const float modes[4][2] = {{0,0},{0,PI},{PI,0},{-PI/2,PI/2}};
208  m_phaseOffsetL = modes[mode][0];
209  m_phaseOffsetR = modes[mode][1];
210  }
211 
212  // what steering mode should be chosen
213  void steering_mode(bool mode) { m_linearSteering = mode; }
214 
215  // set front & rear separation controls
216  void separation(float front, float rear) {
217  m_frontSeparation = front;
218  m_rearSeparation = rear;
219  }
220 
221 private:
222  // polar <-> cartesian coodinates conversion
223  static inline float amplitude(const float cf[2]) { return std::sqrt(cf[0]*cf[0] + cf[1]*cf[1]); }
224  static inline float phase(const float cf[2]) { return std::atan2(cf[1],cf[0]); }
225  static inline cfloat polar(float a, float p) { return {static_cast<float>(a*std::cos(p)),static_cast<float>(a*std::sin(p))}; }
226  static inline float sqr(float x) { return x*x; }
227  // the dreaded min/max
228  static inline float min(float a, float b) { return a<b?a:b; }
229  static inline float max(float a, float b) { return a>b?a:b; }
230  static inline float clamp(float x) { return max(-1,min(1,x)); }
231 
232  // handle the output buffering for overlapped calls of block_decode
233  void add_output(float *input1[2], float *input2[2], float center_width, float dimension, float adaption_rate, bool /*result*/=false) {
234  // add the windowed data to the last 1/2 of the output buffer
235  float *out[6] = {&m_outbuf[0][0],&m_outbuf[1][0],&m_outbuf[2][0],&m_outbuf[3][0],&m_outbuf[4][0],&m_outbuf[5][0]};
236  block_decode(input1,input2,out,center_width,dimension,adaption_rate);
237  }
238 
239  // CORE FUNCTION: decode a block of data
240  void block_decode(float *input1[2], float *input2[2], float *output[6], float center_width, float dimension, float adaption_rate) {
241  // 1. scale the input by the window function; this serves a dual purpose:
242  // - first it improves the FFT resolution b/c boundary discontinuities (and their frequencies) get removed
243  // - second it allows for smooth blending of varying filters between the blocks
244  {
245  float* pWnd = &m_wnd[0];
246  float* pLt = &m_lt[0];
247  float* pRt = &m_rt[0];
248  float* pIn0 = input1[0];
249  float* pIn1 = input1[1];
250  for (unsigned k=0;k<m_halfN;k++) {
251  *pLt++ = *pIn0++ * *pWnd;
252  *pRt++ = *pIn1++ * *pWnd++;
253  }
254  pIn0 = input2[0];
255  pIn1 = input2[1];
256  for (unsigned k=0;k<m_halfN;k++) {
257  *pLt++ = *pIn0++ * *pWnd;
258  *pRt++ = *pIn1++ * *pWnd++;
259  }
260  }
261 
262 #ifdef USE_FFTW3
263  // ... and tranform it into the frequency domain
264  fftwf_execute(m_loadL);
265  fftwf_execute(m_loadR);
266 #else
267  ff_fft_permuteRC(m_fftContextForward, m_lt, (FFTComplex*)&m_dftL[0]);
268  av_fft_calc(m_fftContextForward, (FFTComplex*)&m_dftL[0]);
269 
270  ff_fft_permuteRC(m_fftContextForward, m_rt, (FFTComplex*)&m_dftR[0]);
271  av_fft_calc(m_fftContextForward, (FFTComplex*)&m_dftR[0]);
272 #endif
273 
274  // 2. compare amplitude and phase of each DFT bin and produce the X/Y coordinates in the sound field
275  // but dont do DC or N/2 component
276  for (unsigned f=0;f<m_halfN;f++) {
277  // get left/right amplitudes/phases
278  float ampL = amplitude(m_dftL[f]);
279  float ampR = amplitude(m_dftR[f]);
280  float phaseL = phase(m_dftL[f]);
281  float phaseR = phase(m_dftR[f]);
282 // if (ampL+ampR < epsilon)
283 // continue;
284 
285  // calculate the amplitude/phase difference
286  float ampDiff = clamp((ampL+ampR < epsilon) ? 0 : (ampR-ampL) / (ampR+ampL));
287  float phaseDiff = phaseL - phaseR;
288  if (phaseDiff < -PI) phaseDiff += 2*PI;
289  if (phaseDiff > PI) phaseDiff -= 2*PI;
290  phaseDiff = abs(phaseDiff);
291 
292  if (m_linearSteering) {
293  // --- this is the fancy new linear mode ---
294 
295  // get sound field x/y position
296  m_yFs[f] = get_yfs(ampDiff,phaseDiff);
297  m_xFs[f] = get_xfs(ampDiff,m_yFs[f]);
298 
299  // add dimension control
300  m_yFs[f] = clamp(m_yFs[f] - dimension);
301 
302  // add crossfeed control
303  m_xFs[f] = clamp(m_xFs[f] * (m_frontSeparation*(1+m_yFs[f])/2 + m_rearSeparation*(1-m_yFs[f])/2));
304 
305  // 3. generate frequency filters for each output channel
306  float left = (1-m_xFs[f])/2;
307  float right = (1+m_xFs[f])/2;
308  float front = (1+m_yFs[f])/2;
309  float back = (1-m_yFs[f])/2;
310  float volume[5] = {
311  front * (left * center_width + max(0,-m_xFs[f]) * (1-center_width)), // left
312  front * center_level*((1-abs(m_xFs[f])) * (1-center_width)), // center
313  front * (right * center_width + max(0, m_xFs[f]) * (1-center_width)), // right
314  back * m_surroundLevel * left, // left surround
315  back * m_surroundLevel * right // right surround
316  };
317 
318  // adapt the prior filter
319  for (unsigned c=0;c<5;c++)
320  m_filter[c][f] = (1-adaption_rate)*m_filter[c][f] + adaption_rate*volume[c];
321 
322  } else {
323  // --- this is the old & simple steering mode ---
324 
325  // calculate the amplitude/phase difference
326  float ampDiff = clamp((ampL+ampR < epsilon) ? 0 : (ampR-ampL) / (ampR+ampL));
327  float phaseDiff = phaseL - phaseR;
328  if (phaseDiff < -PI) phaseDiff += 2*PI;
329  if (phaseDiff > PI) phaseDiff -= 2*PI;
330  phaseDiff = abs(phaseDiff);
331 
332  // determine sound field x-position
333  m_xFs[f] = ampDiff;
334 
335  // determine preliminary sound field y-position from phase difference
336  m_yFs[f] = 1 - (phaseDiff/PI)*2;
337 
338  if (abs(m_xFs[f]) > m_surroundBalance) {
339  // blend linearly between the surrounds and the fronts if the balance exceeds the surround encoding balance
340  // this is necessary because the sound field is trapezoidal and will be stretched behind the listener
341  float frontness = (abs(m_xFs[f]) - m_surroundBalance)/(1-m_surroundBalance);
342  m_yFs[f] = (1-frontness) * m_yFs[f] + frontness * 1;
343  }
344 
345  // add dimension control
346  m_yFs[f] = clamp(m_yFs[f] - dimension);
347 
348  // add crossfeed control
349  m_xFs[f] = clamp(m_xFs[f] * (m_frontSeparation*(1+m_yFs[f])/2 + m_rearSeparation*(1-m_yFs[f])/2));
350 
351  // 3. generate frequency filters for each output channel, according to the signal position
352  // the sum of all channel volumes must be 1.0
353  float left = (1-m_xFs[f])/2;
354  float right = (1+m_xFs[f])/2;
355  float front = (1+m_yFs[f])/2;
356  float back = (1-m_yFs[f])/2;
357  float volume[5] = {
358  front * (left * center_width + max(0,-m_xFs[f]) * (1-center_width)), // left
359  front * center_level*((1-abs(m_xFs[f])) * (1-center_width)), // center
360  front * (right * center_width + max(0, m_xFs[f]) * (1-center_width)), // right
361  back * m_surroundLevel*max(0,min(1,((1-(m_xFs[f]/m_surroundBalance))/2))),// left surround
362  back * m_surroundLevel*max(0,min(1,((1+(m_xFs[f]/m_surroundBalance))/2))) // right surround
363  };
364 
365  // adapt the prior filter
366  for (unsigned c=0;c<5;c++)
367  m_filter[c][f] = (1-adaption_rate)*m_filter[c][f] + adaption_rate*volume[c];
368  }
369 
370  // ... and build the signal which we want to position
371  m_frontL[f] = polar(ampL+ampR,phaseL);
372  m_frontR[f] = polar(ampL+ampR,phaseR);
373  m_avg[f] = m_frontL[f] + m_frontR[f];
374  m_surL[f] = polar(ampL+ampR,phaseL+m_phaseOffsetL);
375  m_surR[f] = polar(ampL+ampR,phaseR+m_phaseOffsetR);
376  m_trueavg[f] = cfloat(m_dftL[f][0] + m_dftR[f][0], m_dftL[f][1] + m_dftR[f][1]);
377  }
378 
379  // 4. distribute the unfiltered reference signals over the channels
380  apply_filter(&m_frontL[0], &m_filter[0][0],&output[0][0]); // front left
381  apply_filter(&m_avg[0], &m_filter[1][0],&output[1][0]); // front center
382  apply_filter(&m_frontR[0], &m_filter[2][0],&output[2][0]); // front right
383  apply_filter(&m_surL[0], &m_filter[3][0],&output[3][0]); // surround left
384  apply_filter(&m_surR[0], &m_filter[4][0],&output[4][0]); // surround right
385  apply_filter(&m_trueavg[0],&m_filter[5][0],&output[5][0]); // lfe
386  }
387 
388 #define FASTER_CALC
389  // map from amplitude difference and phase difference to yfs
390  static inline double get_yfs(double ampDiff, double phaseDiff) {
391  double x = 1-(((1-sqr(ampDiff))*phaseDiff)/PI*2);
392 #ifdef FASTER_CALC
393  double tanX = tan(x);
394  return 0.16468622925824683 + 0.5009268347818189*x - 0.06462757726992101*x*x
395  + 0.09170680403453149*x*x*x + 0.2617754892323973*tanX - 0.04180413533856156*sqr(tanX);
396 #else
397  return 0.16468622925824683 + 0.5009268347818189*x - 0.06462757726992101*x*x
398  + 0.09170680403453149*x*x*x + 0.2617754892323973*tan(x) - 0.04180413533856156*sqr(tan(x));
399 #endif
400  }
401 
402  // map from amplitude difference and yfs to xfs
403  static inline double get_xfs(double ampDiff, double yfs) {
404  double x=ampDiff;
405  double y=yfs;
406 #ifdef FASTER_CALC
407  double tanX = tan(x);
408  double tanY = tan(y);
409  double asinX = asin(x);
410  double sinX = sin(x);
411  double sinY = sin(y);
412  double x3 = x*x*x;
413  double y2 = y*y;
414  double y3 = y*y2;
415  return 2.464833559224702*x - 423.52131153259404*x*y +
416  67.8557858606918*x3*y + 788.2429425544392*x*y2 -
417  79.97650354902909*x3*y2 - 513.8966153850349*x*y3 +
418  35.68117670186306*x3*y3 + 13867.406173420834*y*asinX -
419  2075.8237075786396*y2*asinX - 908.2722068360281*y3*asinX -
420  12934.654772878019*asinX*sinY - 13216.736529661162*y*tanX +
421  1288.6463247741938*y2*tanX + 1384.372969378453*y3*tanX +
422  12699.231471126128*sinY*tanX + 95.37131275594336*sinX*tanY -
423  91.21223198407546*tanX*tanY;
424 #else
425  return 2.464833559224702*x - 423.52131153259404*x*y +
426  67.8557858606918*x*x*x*y + 788.2429425544392*x*y*y -
427  79.97650354902909*x*x*x*y*y - 513.8966153850349*x*y*y*y +
428  35.68117670186306*x*x*x*y*y*y + 13867.406173420834*y*asin(x) -
429  2075.8237075786396*y*y*asin(x) - 908.2722068360281*y*y*y*asin(x) -
430  12934.654772878019*asin(x)*sin(y) - 13216.736529661162*y*tan(x) +
431  1288.6463247741938*y*y*tan(x) + 1384.372969378453*y*y*y*tan(x) +
432  12699.231471126128*sin(y)*tan(x) + 95.37131275594336*sin(x)*tan(y) -
433  91.21223198407546*tan(x)*tan(y);
434 #endif
435  }
436 
437  // filter the complex source signal and add it to target
438  void apply_filter(cfloat *signal, const float *flt, float *target) {
439  // filter the signal
440  for (unsigned f=0;f<=m_halfN;f++) {
441  m_src[f][0] = signal[f].real() * flt[f];
442  m_src[f][1] = signal[f].imag() * flt[f];
443  }
444 #ifdef USE_FFTW3
445  // transform into time domain
446  fftwf_execute(m_store);
447 
448  float* pT1 = &target[m_currentBuf*m_halfN];
449  float* pWnd1 = &m_wnd[0];
450  float* pDst1 = &m_dst[0];
451  float* pT2 = &target[(m_currentBuf^1)*m_halfN];
452  float* pWnd2 = &m_wnd[m_halfN];
453  float* pDst2 = &m_dst[m_halfN];
454  // add the result to target, windowed
455  for (unsigned int k=0;k<m_halfN;k++)
456  {
457  // 1st part is overlap add
458  *pT1++ += *pWnd1++ * *pDst1++;
459  // 2nd part is set as has no history
460  *pT2++ = *pWnd2++ * *pDst2++;
461  }
462 #else
463  // enforce odd symmetry
464  for (unsigned f=1;f<m_halfN;f++) {
465  m_src[m_n-f][0] = m_src[f][0];
466  m_src[m_n-f][1] = -m_src[f][1]; // complex conjugate
467  }
468  av_fft_permute(m_fftContextReverse, (FFTComplex*)&m_src[0]);
469  av_fft_calc(m_fftContextReverse, (FFTComplex*)&m_src[0]);
470 
471  float* pT1 = &target[m_currentBuf*m_halfN];
472  float* pWnd1 = &m_wnd[0];
473  float* pDst1 = &m_src[0][0];
474  float* pT2 = &target[(m_currentBuf^1)*m_halfN];
475  float* pWnd2 = &m_wnd[m_halfN];
476  float* pDst2 = &m_src[m_halfN][0];
477  // add the result to target, windowed
478  for (unsigned int k=0;k<m_halfN;k++)
479  {
480  // 1st part is overlap add
481  *pT1++ += *pWnd1++ * *pDst1; pDst1 += 2;
482  // 2nd part is set as has no history
483  *pT2++ = *pWnd2++ * *pDst2; pDst2 += 2;
484  }
485 #endif
486  }
487 
488 #ifndef USE_FFTW3
489 
493  void ff_fft_permuteRC(FFTContext *s, FFTSample *r, FFTComplex *z)
494  {
495  int j, k, np;
496  const uint16_t *revtab = s->revtab;
497 
498  /* reverse */
499  np = 1 << s->nbits;
500  for(j=0;j<np;j++) {
501  k = revtab[j];
502  z[k].re = r[j];
503  z[k].im = 0.0;
504  }
505  }
506 
512  void ff_fft_permuteCR(FFTContext *s, FFTComplex *z, FFTSample *r)
513  {
514  int j, k, np;
515  const uint16_t *revtab = s->revtab;
516 
517  /* reverse */
518  np = 1 << s->nbits;
519  for(j=0;j<np;j++) {
520  k = revtab[j];
521  if (k < j) {
522  r[k] = z[j].re;
523  r[j] = z[k].re;
524  }
525  }
526  }
527 #endif
528 
529  unsigned int m_n; // the block size
530  unsigned int m_halfN; // half block size precalculated
531 #ifdef USE_FFTW3
532  // FFTW data structures
533  float *m_lt,*m_rt,*m_dst; // left total, right total (source arrays), destination array
534  fftwf_complex *m_dftL,*m_dftR,*m_src; // intermediate arrays (FFTs of lt & rt, processing source)
535  fftwf_plan m_loadL,m_loadR,m_store; // plans for loading the data into the intermediate format and back
536 #else
538  FFTSample *m_lt,*m_rt; // left total, right total (source arrays), destination array
539  FFTComplexArray *m_dftL,*m_dftR,*m_src;// intermediate arrays (FFTs of lt & rt, processing source)
540 #endif
541  // buffers
542  std::vector<cfloat> m_frontL,m_frontR,m_avg,m_surL,m_surR; // the signal (phase-corrected) in the frequency domain
543  std::vector<cfloat> m_trueavg; // for lfe generation
544  std::vector<float> m_xFs,m_yFs; // the feature space positions for each frequency bin
545  std::vector<float> m_wnd; // the window function, precalculated
546  std::vector<float> m_filter[6]; // a frequency filter for each output channel
547  std::vector<float> m_inbuf[2]; // the sliding input buffers
548  std::vector<float> m_outbuf[6]; // the sliding output buffers
549  // coefficients
550  float m_surroundHigh,m_surroundLow; // high and low surround mixing coefficient (e.g. 0.8165/0.5774)
551  float m_surroundBalance; // the xfs balance that follows from the coeffs
552  float m_surroundLevel; // gain for the surround channels (follows from the coeffs
553  float m_phaseOffsetL, m_phaseOffsetR;// phase shifts to be applied to the rear channels
554  float m_frontSeparation; // front stereo separation
555  float m_rearSeparation; // rear stereo separation
556  bool m_linearSteering; // whether the steering should be linear or not
557  cfloat m_a,m_b,m_c,m_d,m_e,m_f,m_g,m_h; // coefficients for the linear steering
558  int m_currentBuf; // specifies which buffer is 2nd half of input sliding buffer
559  float * m_inbufs[2]; // for passing back to driver
560  float * m_outbufs[6]; // for passing back to driver
561 
562  friend class fsurround_decoder;
563 };
564 
565 
566 // implementation of the shell class
567 
568 fsurround_decoder::fsurround_decoder(unsigned blocksize): m_impl(new decoder_impl(blocksize)) { }
569 
571 
572 void fsurround_decoder::decode(float center_width, float dimension, float adaption_rate) {
573  m_impl->decode(center_width,dimension,adaption_rate);
574 }
575 
577 
579 
580 void fsurround_decoder::phase_mode(unsigned mode) { m_impl->phase_mode(mode); }
581 
583 
584 void fsurround_decoder::separation(float front, float rear) { m_impl->separation(front,rear); }
585 
587 {
588  return m_impl->getInputBuffers();
589 }
590 
592 {
593  return m_impl->getOutputBuffers();
594 }
595 
596 void fsurround_decoder::sample_rate(unsigned int samplerate)
597 {
598  m_impl->sample_rate(samplerate);
599 }
static float phase(const float cf[2])
void surround_coefficients(float a, float b)
void decode(float center_width, float dimension, float adaption_rate)
FFTSample[2] FFTComplexArray
float * m_inbufs[2]
static float clamp(float x)
unsigned int m_n
void separation(float front, float rear)
std::vector< float > m_xFs
void ff_fft_permuteCR(FFTContext *s, FFTComplex *z, FFTSample *r)
void steering_mode(bool mode)
std::vector< float > m_inbuf[2]
std::vector< float > m_filter[6]
std::vector< float > m_yFs
float * m_outbufs[6]
void sample_rate(unsigned int srate)
void apply_filter(cfloat *signal, const float *flt, float *target)
FFTComplexArray * m_src
float m_phaseOffsetL
std::vector< cfloat > m_frontL
void ff_fft_permuteRC(FFTContext *s, FFTSample *r, FFTComplex *z)
static float amplitude(const float cf[2])
std::vector< cfloat > m_surL
static cfloat polar(float a, float p)
static float max(float a, float b)
void add_output(float *input1[2], float *input2[2], float center_width, float dimension, float adaption_rate, bool=false)
float m_frontSeparation
std::vector< cfloat > m_frontR
void block_decode(float *input1[2], float *input2[2], float *output[6], float center_width, float dimension, float adaption_rate)
decoder_impl(unsigned blocksize=8192)
void steering_mode(bool mode)
float m_phaseOffsetR
void surround_coefficients(float a, float b)
float m_surroundLevel
unsigned short uint16_t
Definition: iso6937tables.h:1
float ** getInputBuffers()
std::vector< float > m_outbuf[6]
static double get_xfs(double ampDiff, double yfs)
FFTSample * m_lt
fsurround_decoder(unsigned blocksize=8192)
FFTComplexArray * m_dftL
void phase_mode(unsigned mode)
std::complex< float > cfloat
void * av_malloc(unsigned int size)
float ** getInputBuffers()
FFTContext * m_fftContextReverse
FFTContext * m_fftContextForward
float m_surroundHigh
FFTComplexArray * m_dftR
std::vector< cfloat > m_avg
static const float PI
static float sqr(float x)
static const float center_level
std::vector< float > m_wnd
static double get_yfs(double ampDiff, double phaseDiff)
void decode(float center_width=1, float dimension=0, float adaption_rate=1)
void separation(float front, float rear)
float ** getOutputBuffers()
std::vector< cfloat > m_trueavg
std::vector< cfloat > m_surR
void av_free(void *ptr)
void phase_mode(unsigned mode)
static float min(float a, float b)
static guint32 * back
Definition: goom_core.c:34
static const float epsilon
void sample_rate(unsigned int samplerate)
unsigned int m_halfN
static int x3
Definition: mythsocket.cpp:62
class decoder_impl * m_impl
Definition: el_processor.h:70
float m_surroundBalance
float ** getOutputBuffers()
#define output
float m_rearSeparation
FFTSample * m_rt