FFmpeg
opus_pvq.c
Go to the documentation of this file.
1 /*
2  * Copyright (c) 2007-2008 CSIRO
3  * Copyright (c) 2007-2009 Xiph.Org Foundation
4  * Copyright (c) 2008-2009 Gregory Maxwell
5  * Copyright (c) 2012 Andrew D'Addesio
6  * Copyright (c) 2013-2014 Mozilla Corporation
7  * Copyright (c) 2017 Rostislav Pehlivanov <atomnuker@gmail.com>
8  *
9  * This file is part of FFmpeg.
10  *
11  * FFmpeg is free software; you can redistribute it and/or
12  * modify it under the terms of the GNU Lesser General Public
13  * License as published by the Free Software Foundation; either
14  * version 2.1 of the License, or (at your option) any later version.
15  *
16  * FFmpeg is distributed in the hope that it will be useful,
17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
19  * Lesser General Public License for more details.
20  *
21  * You should have received a copy of the GNU Lesser General Public
22  * License along with FFmpeg; if not, write to the Free Software
23  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
24  */
25 
26 #include <float.h>
27 
28 #include "config_components.h"
29 
30 #include "libavutil/mem.h"
31 #include "mathops.h"
32 #include "opustab.h"
33 #include "opus_pvq.h"
34 
35 #define ROUND_MUL16(a,b) ((MUL16(a, b) + 16384) >> 15)
36 
37 #define CELT_PVQ_U(n, k) (ff_celt_pvq_u_row[FFMIN(n, k)][FFMAX(n, k)])
38 #define CELT_PVQ_V(n, k) (CELT_PVQ_U(n, k) + CELT_PVQ_U(n, (k) + 1))
39 
40 static inline int16_t celt_cos(int16_t x)
41 {
42  x = (MUL16(x, x) + 4096) >> 13;
43  x = (32767-x) + ROUND_MUL16(x, (-7651 + ROUND_MUL16(x, (8277 + ROUND_MUL16(-626, x)))));
44  return x + 1;
45 }
46 
47 static inline int celt_log2tan(int isin, int icos)
48 {
49  int lc, ls;
50  lc = opus_ilog(icos);
51  ls = opus_ilog(isin);
52  icos <<= 15 - lc;
53  isin <<= 15 - ls;
54  return (ls << 11) - (lc << 11) +
55  ROUND_MUL16(isin, ROUND_MUL16(isin, -2597) + 7932) -
56  ROUND_MUL16(icos, ROUND_MUL16(icos, -2597) + 7932);
57 }
58 
59 static inline int celt_bits2pulses(const uint8_t *cache, int bits)
60 {
61  // TODO: Find the size of cache and make it into an array in the parameters list
62  int i, low = 0, high;
63 
64  high = cache[0];
65  bits--;
66 
67  for (i = 0; i < 6; i++) {
68  int center = (low + high + 1) >> 1;
69  if (cache[center] >= bits)
70  high = center;
71  else
72  low = center;
73  }
74 
75  return (bits - (low == 0 ? -1 : cache[low]) <= cache[high] - bits) ? low : high;
76 }
77 
78 static inline int celt_pulses2bits(const uint8_t *cache, int pulses)
79 {
80  // TODO: Find the size of cache and make it into an array in the parameters list
81  return (pulses == 0) ? 0 : cache[pulses] + 1;
82 }
83 
84 static inline void celt_normalize_residual(const int * restrict iy, float * restrict X,
85  int N, float g)
86 {
87  int i;
88  for (i = 0; i < N; i++)
89  X[i] = g * iy[i];
90 }
91 
92 static void celt_exp_rotation_impl(float *X, uint32_t len, uint32_t stride,
93  float c, float s)
94 {
95  float *Xptr;
96  int i;
97 
98  Xptr = X;
99  for (i = 0; i < len - stride; i++) {
100  float x1 = Xptr[0];
101  float x2 = Xptr[stride];
102  Xptr[stride] = c * x2 + s * x1;
103  *Xptr++ = c * x1 - s * x2;
104  }
105 
106  Xptr = &X[len - 2 * stride - 1];
107  for (i = len - 2 * stride - 1; i >= 0; i--) {
108  float x1 = Xptr[0];
109  float x2 = Xptr[stride];
110  Xptr[stride] = c * x2 + s * x1;
111  *Xptr-- = c * x1 - s * x2;
112  }
113 }
114 
115 static inline void celt_exp_rotation(float *X, uint32_t len,
116  uint32_t stride, uint32_t K,
117  enum CeltSpread spread, const int encode)
118 {
119  uint32_t stride2 = 0;
120  float c, s;
121  float gain, theta;
122  int i;
123 
124  if (2*K >= len || spread == CELT_SPREAD_NONE)
125  return;
126 
127  gain = (float)len / (len + (20 - 5*spread) * K);
128  theta = M_PI * gain * gain / 4;
129 
130  c = cosf(theta);
131  s = sinf(theta);
132 
133  if (len >= stride << 3) {
134  stride2 = 1;
135  /* This is just a simple (equivalent) way of computing sqrt(len/stride) with rounding.
136  It's basically incrementing long as (stride2+0.5)^2 < len/stride. */
137  while ((stride2 * stride2 + stride2) * stride + (stride >> 2) < len)
138  stride2++;
139  }
140 
141  len /= stride;
142  for (i = 0; i < stride; i++) {
143  if (encode) {
144  celt_exp_rotation_impl(X + i * len, len, 1, c, -s);
145  if (stride2)
146  celt_exp_rotation_impl(X + i * len, len, stride2, s, -c);
147  } else {
148  if (stride2)
149  celt_exp_rotation_impl(X + i * len, len, stride2, s, c);
150  celt_exp_rotation_impl(X + i * len, len, 1, c, s);
151  }
152  }
153 }
154 
155 static inline uint32_t celt_extract_collapse_mask(const int *iy, uint32_t N, uint32_t B)
156 {
157  int i, j, N0 = N / B;
158  uint32_t collapse_mask = 0;
159 
160  if (B <= 1)
161  return 1;
162 
163  for (i = 0; i < B; i++)
164  for (j = 0; j < N0; j++)
165  collapse_mask |= (!!iy[i*N0+j]) << i;
166  return collapse_mask;
167 }
168 
169 static inline void celt_stereo_merge(float *X, float *Y, float mid, int N)
170 {
171  int i;
172  float xp = 0, side = 0;
173  float E[2];
174  float mid2;
175  float gain[2];
176 
177  /* Compute the norm of X+Y and X-Y as |X|^2 + |Y|^2 +/- sum(xy) */
178  for (i = 0; i < N; i++) {
179  xp += X[i] * Y[i];
180  side += Y[i] * Y[i];
181  }
182 
183  /* Compensating for the mid normalization */
184  xp *= mid;
185  mid2 = mid;
186  E[0] = mid2 * mid2 + side - 2 * xp;
187  E[1] = mid2 * mid2 + side + 2 * xp;
188  if (E[0] < 6e-4f || E[1] < 6e-4f) {
189  for (i = 0; i < N; i++)
190  Y[i] = X[i];
191  return;
192  }
193 
194  gain[0] = 1.0f / sqrtf(E[0]);
195  gain[1] = 1.0f / sqrtf(E[1]);
196 
197  for (i = 0; i < N; i++) {
198  float value[2];
199  /* Apply mid scaling (side is already scaled) */
200  value[0] = mid * X[i];
201  value[1] = Y[i];
202  X[i] = gain[0] * (value[0] - value[1]);
203  Y[i] = gain[1] * (value[0] + value[1]);
204  }
205 }
206 
207 static void celt_interleave_hadamard(float *tmp, float *X, int N0,
208  int stride, int hadamard)
209 {
210  int i, j, N = N0*stride;
211  const uint8_t *order = &ff_celt_hadamard_order[hadamard ? stride - 2 : 30];
212 
213  for (i = 0; i < stride; i++)
214  for (j = 0; j < N0; j++)
215  tmp[j*stride+i] = X[order[i]*N0+j];
216 
217  memcpy(X, tmp, N*sizeof(float));
218 }
219 
220 static void celt_deinterleave_hadamard(float *tmp, float *X, int N0,
221  int stride, int hadamard)
222 {
223  int i, j, N = N0*stride;
224  const uint8_t *order = &ff_celt_hadamard_order[hadamard ? stride - 2 : 30];
225 
226  for (i = 0; i < stride; i++)
227  for (j = 0; j < N0; j++)
228  tmp[order[i]*N0+j] = X[j*stride+i];
229 
230  memcpy(X, tmp, N*sizeof(float));
231 }
232 
233 static void celt_haar1(float *X, int N0, int stride)
234 {
235  int i, j;
236  N0 >>= 1;
237  for (i = 0; i < stride; i++) {
238  for (j = 0; j < N0; j++) {
239  float x0 = X[stride * (2 * j + 0) + i];
240  float x1 = X[stride * (2 * j + 1) + i];
241  X[stride * (2 * j + 0) + i] = (x0 + x1) * M_SQRT1_2;
242  X[stride * (2 * j + 1) + i] = (x0 - x1) * M_SQRT1_2;
243  }
244  }
245 }
246 
247 static inline int celt_compute_qn(int N, int b, int offset, int pulse_cap,
248  int stereo)
249 {
250  int qn, qb;
251  int N2 = 2 * N - 1;
252  if (stereo && N == 2)
253  N2--;
254 
255  /* The upper limit ensures that in a stereo split with itheta==16384, we'll
256  * always have enough bits left over to code at least one pulse in the
257  * side; otherwise it would collapse, since it doesn't get folded. */
258  qb = FFMIN3(b - pulse_cap - (4 << 3), (b + N2 * offset) / N2, 8 << 3);
259  qn = (qb < (1 << 3 >> 1)) ? 1 : ((ff_celt_qn_exp2[qb & 0x7] >> (14 - (qb >> 3))) + 1) >> 1 << 1;
260  return qn;
261 }
262 
263 /* Convert the quantized vector to an index */
264 static inline uint32_t celt_icwrsi(uint32_t N, uint32_t K, const int *y)
265 {
266  int i, idx = 0, sum = 0;
267  for (i = N - 1; i >= 0; i--) {
268  const uint32_t i_s = CELT_PVQ_U(N - i, sum + FFABS(y[i]) + 1);
269  idx += CELT_PVQ_U(N - i, sum) + (y[i] < 0)*i_s;
270  sum += FFABS(y[i]);
271  }
272  return idx;
273 }
274 
275 // this code was adapted from libopus
276 static inline uint64_t celt_cwrsi(uint32_t N, uint32_t K, uint32_t i, int *y)
277 {
278  uint64_t norm = 0;
279  uint32_t q, p;
280  int s, val;
281  int k0;
282 
283  while (N > 2) {
284  /*Lots of pulses case:*/
285  if (K >= N) {
286  const uint32_t *row = ff_celt_pvq_u_row[N];
287 
288  /* Are the pulses in this dimension negative? */
289  p = row[K + 1];
290  s = -(i >= p);
291  i -= p & s;
292 
293  /*Count how many pulses were placed in this dimension.*/
294  k0 = K;
295  q = row[N];
296  if (q > i) {
297  K = N;
298  do {
299  p = ff_celt_pvq_u_row[--K][N];
300  } while (p > i);
301  } else
302  for (p = row[K]; p > i; p = row[K])
303  K--;
304 
305  i -= p;
306  val = (k0 - K + s) ^ s;
307  norm += val * val;
308  *y++ = val;
309  } else { /*Lots of dimensions case:*/
310  /*Are there any pulses in this dimension at all?*/
311  p = ff_celt_pvq_u_row[K ][N];
312  q = ff_celt_pvq_u_row[K + 1][N];
313 
314  if (p <= i && i < q) {
315  i -= p;
316  *y++ = 0;
317  } else {
318  /*Are the pulses in this dimension negative?*/
319  s = -(i >= q);
320  i -= q & s;
321 
322  /*Count how many pulses were placed in this dimension.*/
323  k0 = K;
324  do p = ff_celt_pvq_u_row[--K][N];
325  while (p > i);
326 
327  i -= p;
328  val = (k0 - K + s) ^ s;
329  norm += val * val;
330  *y++ = val;
331  }
332  }
333  N--;
334  }
335 
336  /* N == 2 */
337  p = 2 * K + 1;
338  s = -(i >= p);
339  i -= p & s;
340  k0 = K;
341  K = (i + 1) / 2;
342 
343  if (K)
344  i -= 2 * K - 1;
345 
346  val = (k0 - K + s) ^ s;
347  norm += val * val;
348  *y++ = val;
349 
350  /* N==1 */
351  s = -i;
352  val = (K + s) ^ s;
353  norm += val * val;
354  *y = val;
355 
356  return norm;
357 }
358 
359 static inline void celt_encode_pulses(OpusRangeCoder *rc, int *y, uint32_t N, uint32_t K)
360 {
362 }
363 
364 static inline float celt_decode_pulses(OpusRangeCoder *rc, int *y, uint32_t N, uint32_t K)
365 {
366  const uint32_t idx = ff_opus_rc_dec_uint(rc, CELT_PVQ_V(N, K));
367  return celt_cwrsi(N, K, idx, y);
368 }
369 
370 #if CONFIG_OPUS_ENCODER
371 /*
372  * Faster than libopus's search, operates entirely in the signed domain.
373  * Slightly worse/better depending on N, K and the input vector.
374  */
375 static float ppp_pvq_search_c(float *X, int *y, int K, int N)
376 {
377  int i, y_norm = 0;
378  float res = 0.0f, xy_norm = 0.0f;
379 
380  for (i = 0; i < N; i++)
381  res += FFABS(X[i]);
382 
383  res = K/(res + FLT_EPSILON);
384 
385  for (i = 0; i < N; i++) {
386  y[i] = lrintf(res*X[i]);
387  y_norm += y[i]*y[i];
388  xy_norm += y[i]*X[i];
389  K -= FFABS(y[i]);
390  }
391 
392  while (K) {
393  int max_idx = 0, phase = FFSIGN(K);
394  float max_num = 0.0f;
395  float max_den = 1.0f;
396  y_norm += 1.0f;
397 
398  for (i = 0; i < N; i++) {
399  /* If the sum has been overshot and the best place has 0 pulses allocated
400  * to it, attempting to decrease it further will actually increase the
401  * sum. Prevent this by disregarding any 0 positions when decrementing. */
402  const int ca = 1 ^ ((y[i] == 0) & (phase < 0));
403  const int y_new = y_norm + 2*phase*FFABS(y[i]);
404  float xy_new = xy_norm + 1*phase*FFABS(X[i]);
405  xy_new = xy_new * xy_new;
406  if (ca && (max_den*xy_new) > (y_new*max_num)) {
407  max_den = y_new;
408  max_num = xy_new;
409  max_idx = i;
410  }
411  }
412 
413  K -= phase;
414 
415  phase *= FFSIGN(X[max_idx]);
416  xy_norm += 1*phase*X[max_idx];
417  y_norm += 2*phase*y[max_idx];
418  y[max_idx] += phase;
419  }
420 
421  return (float)y_norm;
422 }
423 #endif
424 
425 static uint32_t celt_alg_quant(OpusRangeCoder *rc, float *X, uint32_t N, uint32_t K,
426  enum CeltSpread spread, uint32_t blocks, float gain,
427  CeltPVQ *pvq)
428 {
429  int *y = pvq->qcoeff;
430 
431  celt_exp_rotation(X, N, blocks, K, spread, 1);
432  gain /= sqrtf(pvq->pvq_search(X, y, K, N));
433  celt_encode_pulses(rc, y, N, K);
434  celt_normalize_residual(y, X, N, gain);
435  celt_exp_rotation(X, N, blocks, K, spread, 0);
436  return celt_extract_collapse_mask(y, N, blocks);
437 }
438 
439 /** Decode pulse vector and combine the result with the pitch vector to produce
440  the final normalised signal in the current band. */
441 static uint32_t celt_alg_unquant(OpusRangeCoder *rc, float *X, uint32_t N, uint32_t K,
442  enum CeltSpread spread, uint32_t blocks, float gain,
443  CeltPVQ *pvq)
444 {
445  int *y = pvq->qcoeff;
446 
447  gain /= sqrtf(celt_decode_pulses(rc, y, N, K));
448  celt_normalize_residual(y, X, N, gain);
449  celt_exp_rotation(X, N, blocks, K, spread, 0);
450  return celt_extract_collapse_mask(y, N, blocks);
451 }
452 
453 static int celt_calc_theta(const float *X, const float *Y, int coupling, int N)
454 {
455  int i;
456  float e[2] = { 0.0f, 0.0f };
457  if (coupling) { /* Coupling case */
458  for (i = 0; i < N; i++) {
459  e[0] += (X[i] + Y[i])*(X[i] + Y[i]);
460  e[1] += (X[i] - Y[i])*(X[i] - Y[i]);
461  }
462  } else {
463  for (i = 0; i < N; i++) {
464  e[0] += X[i]*X[i];
465  e[1] += Y[i]*Y[i];
466  }
467  }
468  return lrintf(32768.0f*atan2f(sqrtf(e[1]), sqrtf(e[0]))/M_PI);
469 }
470 
471 static void celt_stereo_is_decouple(float *X, float *Y, float e_l, float e_r, int N)
472 {
473  int i;
474  const float energy_n = 1.0f/(sqrtf(e_l*e_l + e_r*e_r) + FLT_EPSILON);
475  e_l *= energy_n;
476  e_r *= energy_n;
477  for (i = 0; i < N; i++)
478  X[i] = e_l*X[i] + e_r*Y[i];
479 }
480 
481 static void celt_stereo_ms_decouple(float *X, float *Y, int N)
482 {
483  int i;
484  for (i = 0; i < N; i++) {
485  const float Xret = X[i];
486  X[i] = (X[i] + Y[i])*M_SQRT1_2;
487  Y[i] = (Y[i] - Xret)*M_SQRT1_2;
488  }
489 }
490 
492  OpusRangeCoder *rc,
493  const int band, float *X,
494  float *Y, int N, int b,
495  uint32_t blocks, float *lowband,
496  int duration, float *lowband_out,
497  int level, float gain,
498  float *lowband_scratch,
499  int fill, int quant)
500 {
501  int i;
502  const uint8_t *cache;
503  int stereo = !!Y, split = stereo;
504  int imid = 0, iside = 0;
505  uint32_t N0 = N;
506  int N_B = N / blocks;
507  int N_B0 = N_B;
508  int B0 = blocks;
509  int time_divide = 0;
510  int recombine = 0;
511  int inv = 0;
512  float mid = 0, side = 0;
513  int longblocks = (B0 == 1);
514  uint32_t cm = 0;
515 
516  if (N == 1) {
517  float *x = X;
518  for (i = 0; i <= stereo; i++) {
519  int sign = 0;
520  if (f->remaining2 >= 1 << 3) {
521  if (quant) {
522  sign = x[0] < 0;
523  ff_opus_rc_put_raw(rc, sign, 1);
524  } else {
525  sign = ff_opus_rc_get_raw(rc, 1);
526  }
527  f->remaining2 -= 1 << 3;
528  }
529  x[0] = 1.0f - 2.0f*sign;
530  x = Y;
531  }
532  if (lowband_out)
533  lowband_out[0] = X[0];
534  return 1;
535  }
536 
537  if (!stereo && level == 0) {
538  int tf_change = f->tf_change[band];
539  int k;
540  if (tf_change > 0)
541  recombine = tf_change;
542  /* Band recombining to increase frequency resolution */
543 
544  if (lowband &&
545  (recombine || ((N_B & 1) == 0 && tf_change < 0) || B0 > 1)) {
546  for (i = 0; i < N; i++)
547  lowband_scratch[i] = lowband[i];
548  lowband = lowband_scratch;
549  }
550 
551  for (k = 0; k < recombine; k++) {
552  if (quant || lowband)
553  celt_haar1(quant ? X : lowband, N >> k, 1 << k);
554  fill = ff_celt_bit_interleave[fill & 0xF] | ff_celt_bit_interleave[fill >> 4] << 2;
555  }
556  blocks >>= recombine;
557  N_B <<= recombine;
558 
559  /* Increasing the time resolution */
560  while ((N_B & 1) == 0 && tf_change < 0) {
561  if (quant || lowband)
562  celt_haar1(quant ? X : lowband, N_B, blocks);
563  fill |= fill << blocks;
564  blocks <<= 1;
565  N_B >>= 1;
566  time_divide++;
567  tf_change++;
568  }
569  B0 = blocks;
570  N_B0 = N_B;
571 
572  /* Reorganize the samples in time order instead of frequency order */
573  if (B0 > 1 && (quant || lowband))
575  N_B >> recombine, B0 << recombine,
576  longblocks);
577  }
578 
579  /* If we need 1.5 more bit than we can produce, split the band in two. */
580  cache = ff_celt_cache_bits +
582  if (!stereo && duration >= 0 && b > cache[cache[0]] + 12 && N > 2) {
583  N >>= 1;
584  Y = X + N;
585  split = 1;
586  duration -= 1;
587  if (blocks == 1)
588  fill = (fill & 1) | (fill << 1);
589  blocks = (blocks + 1) >> 1;
590  }
591 
592  if (split) {
593  int qn;
594  int itheta = quant ? celt_calc_theta(X, Y, stereo, N) : 0;
595  int mbits, sbits, delta;
596  int qalloc;
597  int pulse_cap;
598  int offset;
599  int orig_fill;
600  int tell;
601 
602  /* Decide on the resolution to give to the split parameter theta */
603  pulse_cap = ff_celt_log_freq_range[band] + duration * 8;
604  offset = (pulse_cap >> 1) - (stereo && N == 2 ? CELT_QTHETA_OFFSET_TWOPHASE :
606  qn = (stereo && band >= f->intensity_stereo) ? 1 :
607  celt_compute_qn(N, b, offset, pulse_cap, stereo);
608  tell = opus_rc_tell_frac(rc);
609  if (qn != 1) {
610  if (quant)
611  itheta = (itheta*qn + 8192) >> 14;
612  /* Entropy coding of the angle. We use a uniform pdf for the
613  * time split, a step for stereo, and a triangular one for the rest. */
614  if (quant) {
615  if (stereo && N > 2)
616  ff_opus_rc_enc_uint_step(rc, itheta, qn / 2);
617  else if (stereo || B0 > 1)
618  ff_opus_rc_enc_uint(rc, itheta, qn + 1);
619  else
620  ff_opus_rc_enc_uint_tri(rc, itheta, qn);
621  itheta = itheta * 16384 / qn;
622  if (stereo) {
623  if (itheta == 0)
624  celt_stereo_is_decouple(X, Y, f->block[0].lin_energy[band],
625  f->block[1].lin_energy[band], N);
626  else
628  }
629  } else {
630  if (stereo && N > 2)
631  itheta = ff_opus_rc_dec_uint_step(rc, qn / 2);
632  else if (stereo || B0 > 1)
633  itheta = ff_opus_rc_dec_uint(rc, qn+1);
634  else
635  itheta = ff_opus_rc_dec_uint_tri(rc, qn);
636  itheta = itheta * 16384 / qn;
637  }
638  } else if (stereo) {
639  if (quant) {
640  inv = f->apply_phase_inv ? itheta > 8192 : 0;
641  if (inv) {
642  for (i = 0; i < N; i++)
643  Y[i] *= -1;
644  }
645  celt_stereo_is_decouple(X, Y, f->block[0].lin_energy[band],
646  f->block[1].lin_energy[band], N);
647 
648  if (b > 2 << 3 && f->remaining2 > 2 << 3) {
649  ff_opus_rc_enc_log(rc, inv, 2);
650  } else {
651  inv = 0;
652  }
653  } else {
654  inv = (b > 2 << 3 && f->remaining2 > 2 << 3) ? ff_opus_rc_dec_log(rc, 2) : 0;
655  inv = f->apply_phase_inv ? inv : 0;
656  }
657  itheta = 0;
658  }
659  qalloc = opus_rc_tell_frac(rc) - tell;
660  b -= qalloc;
661 
662  orig_fill = fill;
663  if (itheta == 0) {
664  imid = 32767;
665  iside = 0;
666  fill = av_zero_extend(fill, blocks);
667  delta = -16384;
668  } else if (itheta == 16384) {
669  imid = 0;
670  iside = 32767;
671  fill &= ((1 << blocks) - 1) << blocks;
672  delta = 16384;
673  } else {
674  imid = celt_cos(itheta);
675  iside = celt_cos(16384-itheta);
676  /* This is the mid vs side allocation that minimizes squared error
677  in that band. */
678  delta = ROUND_MUL16((N - 1) << 7, celt_log2tan(iside, imid));
679  }
680 
681  mid = imid / 32768.0f;
682  side = iside / 32768.0f;
683 
684  /* This is a special case for N=2 that only works for stereo and takes
685  advantage of the fact that mid and side are orthogonal to encode
686  the side with just one bit. */
687  if (N == 2 && stereo) {
688  int c;
689  int sign = 0;
690  float tmp;
691  float *x2, *y2;
692  mbits = b;
693  /* Only need one bit for the side */
694  sbits = (itheta != 0 && itheta != 16384) ? 1 << 3 : 0;
695  mbits -= sbits;
696  c = (itheta > 8192);
697  f->remaining2 -= qalloc+sbits;
698 
699  x2 = c ? Y : X;
700  y2 = c ? X : Y;
701  if (sbits) {
702  if (quant) {
703  sign = x2[0]*y2[1] - x2[1]*y2[0] < 0;
704  ff_opus_rc_put_raw(rc, sign, 1);
705  } else {
706  sign = ff_opus_rc_get_raw(rc, 1);
707  }
708  }
709  sign = 1 - 2 * sign;
710  /* We use orig_fill here because we want to fold the side, but if
711  itheta==16384, we'll have cleared the low bits of fill. */
712  cm = pvq->quant_band(pvq, f, rc, band, x2, NULL, N, mbits, blocks, lowband, duration,
713  lowband_out, level, gain, lowband_scratch, orig_fill);
714  /* We don't split N=2 bands, so cm is either 1 or 0 (for a fold-collapse),
715  and there's no need to worry about mixing with the other channel. */
716  y2[0] = -sign * x2[1];
717  y2[1] = sign * x2[0];
718  X[0] *= mid;
719  X[1] *= mid;
720  Y[0] *= side;
721  Y[1] *= side;
722  tmp = X[0];
723  X[0] = tmp - Y[0];
724  Y[0] = tmp + Y[0];
725  tmp = X[1];
726  X[1] = tmp - Y[1];
727  Y[1] = tmp + Y[1];
728  } else {
729  /* "Normal" split code */
730  float *next_lowband2 = NULL;
731  float *next_lowband_out1 = NULL;
732  int next_level = 0;
733  int rebalance;
734  uint32_t cmt;
735 
736  /* Give more bits to low-energy MDCTs than they would
737  * otherwise deserve */
738  if (B0 > 1 && !stereo && (itheta & 0x3fff)) {
739  if (itheta > 8192)
740  /* Rough approximation for pre-echo masking */
741  delta -= delta >> (4 - duration);
742  else
743  /* Corresponds to a forward-masking slope of
744  * 1.5 dB per 10 ms */
745  delta = FFMIN(0, delta + (N << 3 >> (5 - duration)));
746  }
747  mbits = av_clip((b - delta) / 2, 0, b);
748  sbits = b - mbits;
749  f->remaining2 -= qalloc;
750 
751  if (lowband && !stereo)
752  next_lowband2 = lowband + N; /* >32-bit split case */
753 
754  /* Only stereo needs to pass on lowband_out.
755  * Otherwise, it's handled at the end */
756  if (stereo)
757  next_lowband_out1 = lowband_out;
758  else
759  next_level = level + 1;
760 
761  rebalance = f->remaining2;
762  if (mbits >= sbits) {
763  /* In stereo mode, we do not apply a scaling to the mid
764  * because we need the normalized mid for folding later */
765  cm = pvq->quant_band(pvq, f, rc, band, X, NULL, N, mbits, blocks,
766  lowband, duration, next_lowband_out1, next_level,
767  stereo ? 1.0f : (gain * mid), lowband_scratch, fill);
768  rebalance = mbits - (rebalance - f->remaining2);
769  if (rebalance > 3 << 3 && itheta != 0)
770  sbits += rebalance - (3 << 3);
771 
772  /* For a stereo split, the high bits of fill are always zero,
773  * so no folding will be done to the side. */
774  cmt = pvq->quant_band(pvq, f, rc, band, Y, NULL, N, sbits, blocks,
775  next_lowband2, duration, NULL, next_level,
776  gain * side, NULL, fill >> blocks);
777  cm |= cmt << ((B0 >> 1) & (stereo - 1));
778  } else {
779  /* For a stereo split, the high bits of fill are always zero,
780  * so no folding will be done to the side. */
781  cm = pvq->quant_band(pvq, f, rc, band, Y, NULL, N, sbits, blocks,
782  next_lowband2, duration, NULL, next_level,
783  gain * side, NULL, fill >> blocks);
784  cm <<= ((B0 >> 1) & (stereo - 1));
785  rebalance = sbits - (rebalance - f->remaining2);
786  if (rebalance > 3 << 3 && itheta != 16384)
787  mbits += rebalance - (3 << 3);
788 
789  /* In stereo mode, we do not apply a scaling to the mid because
790  * we need the normalized mid for folding later */
791  cm |= pvq->quant_band(pvq, f, rc, band, X, NULL, N, mbits, blocks,
792  lowband, duration, next_lowband_out1, next_level,
793  stereo ? 1.0f : (gain * mid), lowband_scratch, fill);
794  }
795  }
796  } else {
797  /* This is the basic no-split case */
798  uint32_t q = celt_bits2pulses(cache, b);
799  uint32_t curr_bits = celt_pulses2bits(cache, q);
800  f->remaining2 -= curr_bits;
801 
802  /* Ensures we can never bust the budget */
803  while (f->remaining2 < 0 && q > 0) {
804  f->remaining2 += curr_bits;
805  curr_bits = celt_pulses2bits(cache, --q);
806  f->remaining2 -= curr_bits;
807  }
808 
809  if (q != 0) {
810  /* Finally do the actual (de)quantization */
811  if (quant) {
812  cm = celt_alg_quant(rc, X, N, (q < 8) ? q : (8 + (q & 7)) << ((q >> 3) - 1),
813  f->spread, blocks, gain, pvq);
814  } else {
815  cm = celt_alg_unquant(rc, X, N, (q < 8) ? q : (8 + (q & 7)) << ((q >> 3) - 1),
816  f->spread, blocks, gain, pvq);
817  }
818  } else {
819  /* If there's no pulse, fill the band anyway */
820  uint32_t cm_mask = (1 << blocks) - 1;
821  fill &= cm_mask;
822  if (fill) {
823  if (!lowband) {
824  /* Noise */
825  for (i = 0; i < N; i++)
826  X[i] = (((int32_t)celt_rng(f)) >> 20);
827  cm = cm_mask;
828  } else {
829  /* Folded spectrum */
830  for (i = 0; i < N; i++) {
831  /* About 48 dB below the "normal" folding level */
832  X[i] = lowband[i] + (((celt_rng(f)) & 0x8000) ? 1.0f / 256 : -1.0f / 256);
833  }
834  cm = fill;
835  }
836  celt_renormalize_vector(X, N, gain);
837  } else {
838  memset(X, 0, N*sizeof(float));
839  }
840  }
841  }
842 
843  /* This code is used by the decoder and by the resynthesis-enabled encoder */
844  if (stereo) {
845  if (N > 2)
846  celt_stereo_merge(X, Y, mid, N);
847  if (inv) {
848  for (i = 0; i < N; i++)
849  Y[i] *= -1;
850  }
851  } else if (level == 0) {
852  int k;
853 
854  /* Undo the sample reorganization going from time order to frequency order */
855  if (B0 > 1)
856  celt_interleave_hadamard(pvq->hadamard_tmp, X, N_B >> recombine,
857  B0 << recombine, longblocks);
858 
859  /* Undo time-freq changes that we did earlier */
860  N_B = N_B0;
861  blocks = B0;
862  for (k = 0; k < time_divide; k++) {
863  blocks >>= 1;
864  N_B <<= 1;
865  cm |= cm >> blocks;
866  celt_haar1(X, N_B, blocks);
867  }
868 
869  for (k = 0; k < recombine; k++) {
871  celt_haar1(X, N0>>k, 1<<k);
872  }
873  blocks <<= recombine;
874 
875  /* Scale output for later folding */
876  if (lowband_out) {
877  float n = sqrtf(N0);
878  for (i = 0; i < N0; i++)
879  lowband_out[i] = n * X[i];
880  }
881  cm = av_zero_extend(cm, blocks);
882  }
883 
884  return cm;
885 }
886 
887 static QUANT_FN(pvq_decode_band)
888 {
889 #if CONFIG_OPUS_DECODER
890  return quant_band_template(pvq, f, rc, band, X, Y, N, b, blocks, lowband, duration,
891  lowband_out, level, gain, lowband_scratch, fill, 0);
892 #else
893  return 0;
894 #endif
895 }
896 
897 static QUANT_FN(pvq_encode_band)
898 {
899 #if CONFIG_OPUS_ENCODER
900  return quant_band_template(pvq, f, rc, band, X, Y, N, b, blocks, lowband, duration,
901  lowband_out, level, gain, lowband_scratch, fill, 1);
902 #else
903  return 0;
904 #endif
905 }
906 
908 {
909  CeltPVQ *s = av_malloc(sizeof(CeltPVQ));
910  if (!s)
911  return AVERROR(ENOMEM);
912 
913  s->quant_band = encode ? pvq_encode_band : pvq_decode_band;
914 
915 #if CONFIG_OPUS_ENCODER
916  s->pvq_search = ppp_pvq_search_c;
917 #if ARCH_X86
919 #endif
920 #endif
921 
922  *pvq = s;
923 
924  return 0;
925 }
926 
928 {
929  av_freep(pvq);
930 }
ff_celt_cache_bits
const uint8_t ff_celt_cache_bits[392]
Definition: opustab.c:890
ff_celt_bit_deinterleave
const uint8_t ff_celt_bit_deinterleave[]
Definition: opustab.c:938
level
uint8_t level
Definition: svq3.c:205
av_clip
#define av_clip
Definition: common.h:100
AVERROR
Filter the word “frame” indicates either a video frame or a group of audio as stored in an AVFrame structure Format for each input and each output the list of supported formats For video that means pixel format For audio that means channel sample they are references to shared objects When the negotiation mechanism computes the intersection of the formats supported at each end of a all references to both lists are replaced with a reference to the intersection And when a single format is eventually chosen for a link amongst the remaining all references to the list are updated That means that if a filter requires that its input and output have the same format amongst a supported all it has to do is use a reference to the same list of formats query_formats can leave some formats unset and return AVERROR(EAGAIN) to cause the negotiation mechanism toagain later. That can be used by filters with complex requirements to use the format negotiated on one link to set the formats supported on another. Frame references ownership and permissions
ff_celt_qn_exp2
const uint16_t ff_celt_qn_exp2[]
Definition: opustab.c:951
celt_alg_unquant
static uint32_t celt_alg_unquant(OpusRangeCoder *rc, float *X, uint32_t N, uint32_t K, enum CeltSpread spread, uint32_t blocks, float gain, CeltPVQ *pvq)
Decode pulse vector and combine the result with the pitch vector to produce the final normalised sign...
Definition: opus_pvq.c:441
opus_ilog
#define opus_ilog(i)
Definition: opus_rc.h:30
ff_opus_rc_get_raw
uint32_t ff_opus_rc_get_raw(OpusRangeCoder *rc, uint32_t count)
CELT: read 1-25 raw bits at the end of the frame, backwards byte-wise.
Definition: opus_rc.c:140
quant_band_template
static av_always_inline uint32_t quant_band_template(CeltPVQ *pvq, CeltFrame *f, OpusRangeCoder *rc, const int band, float *X, float *Y, int N, int b, uint32_t blocks, float *lowband, int duration, float *lowband_out, int level, float gain, float *lowband_scratch, int fill, int quant)
Definition: opus_pvq.c:491
celt_cwrsi
static uint64_t celt_cwrsi(uint32_t N, uint32_t K, uint32_t i, int *y)
Definition: opus_pvq.c:276
ff_opus_rc_enc_uint
void ff_opus_rc_enc_uint(OpusRangeCoder *rc, uint32_t val, uint32_t size)
CELT: write a uniformly distributed integer.
Definition: opus_rc.c:204
atan2f
#define atan2f(y, x)
Definition: libm.h:45
tmp
static uint8_t tmp[11]
Definition: aes_ctr.c:28
celt_haar1
static void celt_haar1(float *X, int N0, int stride)
Definition: opus_pvq.c:233
CELT_SPREAD_NONE
@ CELT_SPREAD_NONE
Definition: opus_celt.h:55
b
#define b
Definition: input.c:41
ff_celt_log_freq_range
const uint8_t ff_celt_log_freq_range[]
Definition: opustab.c:776
celt_stereo_is_decouple
static void celt_stereo_is_decouple(float *X, float *Y, float e_l, float e_r, int N)
Definition: opus_pvq.c:471
high
int high
Definition: dovi_rpuenc.c:38
celt_exp_rotation_impl
static void celt_exp_rotation_impl(float *X, uint32_t len, uint32_t stride, float c, float s)
Definition: opus_pvq.c:92
ff_celt_hadamard_order
const uint8_t ff_celt_hadamard_order[]
Definition: opustab.c:943
float.h
celt_calc_theta
static int celt_calc_theta(const float *X, const float *Y, int coupling, int N)
Definition: opus_pvq.c:453
celt_decode_pulses
static float celt_decode_pulses(OpusRangeCoder *rc, int *y, uint32_t N, uint32_t K)
Definition: opus_pvq.c:364
N2
#define N2
Definition: vf_pp7.c:74
QUANT_FN
static QUANT_FN(pvq_decode_band)
Definition: opus_pvq.c:887
ff_celt_pvq_init
int av_cold ff_celt_pvq_init(CeltPVQ **pvq, int encode)
Definition: opus_pvq.c:907
ff_celt_pvq_u_row
const uint32_t *const ff_celt_pvq_u_row[15]
Definition: opustab.c:1155
av_malloc
#define av_malloc(s)
Definition: tableprint_vlc.h:30
celt_exp_rotation
static void celt_exp_rotation(float *X, uint32_t len, uint32_t stride, uint32_t K, enum CeltSpread spread, const int encode)
Definition: opus_pvq.c:115
CeltPVQ
Definition: opus_pvq.h:37
celt_log2tan
static int celt_log2tan(int isin, int icos)
Definition: opus_pvq.c:47
cosf
#define cosf(x)
Definition: libm.h:78
FFSIGN
#define FFSIGN(a)
Definition: common.h:75
celt_pulses2bits
static int celt_pulses2bits(const uint8_t *cache, int pulses)
Definition: opus_pvq.c:78
val
static double val(void *priv, double ch)
Definition: aeval.c:77
celt_bits2pulses
static int celt_bits2pulses(const uint8_t *cache, int bits)
Definition: opus_pvq.c:59
quant
static const uint8_t quant[64]
Definition: vmixdec.c:71
celt_renormalize_vector
static av_always_inline void celt_renormalize_vector(float *X, int N, float gain)
Definition: opus_celt.h:155
MUL16
#define MUL16(ra, rb)
Definition: mathops.h:87
av_cold
#define av_cold
Definition: attributes.h:90
CeltPVQ::quant_band
QUANT_FN * quant_band
Definition: opus_pvq.h:42
duration
int64_t duration
Definition: movenc.c:65
float
float
Definition: af_crystalizer.c:122
CELT_PVQ_U
#define CELT_PVQ_U(n, k)
Definition: opus_pvq.c:37
ff_opus_rc_dec_uint
uint32_t ff_opus_rc_dec_uint(OpusRangeCoder *rc, uint32_t size)
CELT: read a uniform distribution.
Definition: opus_rc.c:182
s
#define s(width, name)
Definition: cbs_vp9.c:198
celt_alg_quant
static uint32_t celt_alg_quant(OpusRangeCoder *rc, float *X, uint32_t N, uint32_t K, enum CeltSpread spread, uint32_t blocks, float gain, CeltPVQ *pvq)
Definition: opus_pvq.c:425
g
const char * g
Definition: vf_curves.c:128
bits
uint8_t bits
Definition: vp3data.h:128
CELT_QTHETA_OFFSET_TWOPHASE
#define CELT_QTHETA_OFFSET_TWOPHASE
Definition: opus_celt.h:50
B
#define B
Definition: huffyuv.h:42
CeltPVQ::hadamard_tmp
float hadamard_tmp[256]
Definition: opus_pvq.h:39
ff_opus_rc_enc_uint_tri
void ff_opus_rc_enc_uint_tri(OpusRangeCoder *rc, uint32_t k, int qn)
Definition: opus_rc.c:258
celt_icwrsi
static uint32_t celt_icwrsi(uint32_t N, uint32_t K, const int *y)
Definition: opus_pvq.c:264
E
#define E
Definition: avdct.c:33
FFABS
#define FFABS(a)
Absolute value, Note, INT_MIN / INT64_MIN result in undefined behavior as they are not representable ...
Definition: common.h:74
CeltPVQ::qcoeff
int qcoeff[256]
Definition: opus_pvq.h:38
ff_celt_cache_index
const int16_t ff_celt_cache_index[105]
Definition: opustab.c:920
NULL
#define NULL
Definition: coverity.c:32
ff_celt_pvq_uninit
void av_cold ff_celt_pvq_uninit(CeltPVQ **pvq)
Definition: opus_pvq.c:927
tell
static int BS_FUNC() tell(const BSCTX *bc)
Return number of bits already read.
Definition: bitstream_template.h:146
mathops.h
CELT_QTHETA_OFFSET
#define CELT_QTHETA_OFFSET
Definition: opus_celt.h:49
sqrtf
static __device__ float sqrtf(float a)
Definition: cuda_runtime.h:184
sinf
#define sinf(x)
Definition: libm.h:419
ff_opus_rc_dec_uint_tri
uint32_t ff_opus_rc_dec_uint_tri(OpusRangeCoder *rc, int qn)
Definition: opus_rc.c:234
celt_deinterleave_hadamard
static void celt_deinterleave_hadamard(float *tmp, float *X, int N0, int stride, int hadamard)
Definition: opus_pvq.c:220
c
Undefined Behavior In the C some operations are like signed integer dereferencing freed accessing outside allocated Undefined Behavior must not occur in a C it is not safe even if the output of undefined operations is unused The unsafety may seem nit picking but Optimizing compilers have in fact optimized code on the assumption that no undefined Behavior occurs Optimizing code based on wrong assumptions can and has in some cases lead to effects beyond the output of computations The signed integer overflow problem in speed critical code Code which is highly optimized and works with signed integers sometimes has the problem that often the output of the computation does not c
Definition: undefined.txt:32
opustab.h
CELT_MAX_BANDS
#define CELT_MAX_BANDS
Definition: opus_celt.h:42
f
f
Definition: af_crystalizer.c:122
opus_pvq.h
N0
#define N0
Definition: vf_pp7.c:72
celt_cos
static int16_t celt_cos(int16_t x)
Definition: opus_pvq.c:40
OpusRangeCoder
Definition: opus_rc.h:39
split
static char * split(char *message, char delim)
Definition: af_channelmap.c:81
celt_rng
static av_always_inline uint32_t celt_rng(CeltFrame *f)
Definition: opus_celt.h:149
encode
static void encode(AVCodecContext *ctx, AVFrame *frame, AVPacket *pkt, FILE *output)
Definition: encode_audio.c:94
offset
it s the only field you need to keep assuming you have a context There is some magic you don t need to care about around this just let it vf offset
Definition: writing_filters.txt:86
N
#define N
Definition: af_mcompand.c:54
av_zero_extend
#define av_zero_extend
Definition: common.h:151
M_PI
#define M_PI
Definition: mathematics.h:67
Y
#define Y
Definition: boxblur.h:37
ROUND_MUL16
#define ROUND_MUL16(a, b)
Definition: opus_pvq.c:35
celt_normalize_residual
static void celt_normalize_residual(const int *restrict iy, float *restrict X, int N, float g)
Definition: opus_pvq.c:84
ff_opus_rc_enc_uint_step
void ff_opus_rc_enc_uint_step(OpusRangeCoder *rc, uint32_t val, int k0)
Definition: opus_rc.c:226
lrintf
#define lrintf(x)
Definition: libm_mips.h:72
i
#define i(width, name, range_min, range_max)
Definition: cbs_h2645.c:256
FFMIN3
#define FFMIN3(a, b, c)
Definition: macros.h:50
ff_opus_rc_put_raw
void ff_opus_rc_put_raw(OpusRangeCoder *rc, uint32_t val, uint32_t count)
CELT: write 0 - 31 bits to the rawbits buffer.
Definition: opus_rc.c:161
delta
float delta
Definition: vorbis_enc_data.h:430
av_always_inline
#define av_always_inline
Definition: attributes.h:49
value
it s the only field you need to keep assuming you have a context There is some magic you don t need to care about around this just let it vf default value
Definition: writing_filters.txt:86
celt_encode_pulses
static void celt_encode_pulses(OpusRangeCoder *rc, int *y, uint32_t N, uint32_t K)
Definition: opus_pvq.c:359
FFMIN
#define FFMIN(a, b)
Definition: macros.h:49
len
int len
Definition: vorbis_enc_data.h:426
opus_rc_tell_frac
static av_always_inline uint32_t opus_rc_tell_frac(const OpusRangeCoder *rc)
Definition: opus_rc.h:65
B0
@ B0
Definition: mvs.c:526
stride
#define stride
Definition: h264pred_template.c:537
ff_opus_rc_dec_uint_step
uint32_t ff_opus_rc_dec_uint_step(OpusRangeCoder *rc, int k0)
Definition: opus_rc.c:211
celt_stereo_merge
static void celt_stereo_merge(float *X, float *Y, float mid, int N)
Definition: opus_pvq.c:169
M_SQRT1_2
#define M_SQRT1_2
Definition: mathematics.h:103
CeltSpread
CeltSpread
Definition: opus_celt.h:54
cm
#define cm
Definition: dvbsubdec.c:40
ff_opus_rc_enc_log
void ff_opus_rc_enc_log(OpusRangeCoder *rc, int val, uint32_t bits)
Definition: opus_rc.c:131
celt_compute_qn
static int celt_compute_qn(int N, int b, int offset, int pulse_cap, int stereo)
Definition: opus_pvq.c:247
ff_celt_pvq_init_x86
void ff_celt_pvq_init_x86(struct CeltPVQ *s)
Definition: celt_pvq_init.c:32
celt_interleave_hadamard
static void celt_interleave_hadamard(float *tmp, float *X, int N0, int stride, int hadamard)
Definition: opus_pvq.c:207
mem.h
av_freep
#define av_freep(p)
Definition: tableprint_vlc.h:34
K
#define K
Definition: palette.c:25
int32_t
int32_t
Definition: audioconvert.c:56
ff_celt_bit_interleave
const uint8_t ff_celt_bit_interleave[]
Definition: opustab.c:934
X
@ X
Definition: vf_addroi.c:27
ff_opus_rc_dec_log
uint32_t ff_opus_rc_dec_log(OpusRangeCoder *rc, uint32_t bits)
Definition: opus_rc.c:114
CeltPVQ::pvq_search
float(* pvq_search)(float *X, int *y, int K, int N)
Definition: opus_pvq.h:41
pulses
static const int8_t pulses[4]
Number of non-zero pulses in the MP-MLQ excitation.
Definition: g723_1.h:260
celt_extract_collapse_mask
static uint32_t celt_extract_collapse_mask(const int *iy, uint32_t N, uint32_t B)
Definition: opus_pvq.c:155
celt_stereo_ms_decouple
static void celt_stereo_ms_decouple(float *X, float *Y, int N)
Definition: opus_pvq.c:481
CeltFrame
Definition: opus_celt.h:97
CELT_PVQ_V
#define CELT_PVQ_V(n, k)
Definition: opus_pvq.c:38