libavcodec/dcadsp.c

   1 /*
   2  * Copyright (C) 2016 foo86
   3  *
   4  * This file is part of FFmpeg.
   5  *
   6  * FFmpeg is free software; you can redistribute it and/or
   7  * modify it under the terms of the GNU Lesser General Public
   8  * License as published by the Free Software Foundation; either
   9  * version 2.1 of the License, or (at your option) any later version.
  10  *
  11  * FFmpeg 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 GNU
  14  * Lesser General Public License for more details.
  15  *
  16  * You should have received a copy of the GNU Lesser General Public
  17  * License along with FFmpeg; if not, write to the Free Software
  18  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
  19  */
  20
  21 #include "libavutil/mem_internal.h"
  22
  23 #include "dcadsp.h"
  24 #include "dcamath.h"
  25
  26 static void decode_hf_c(int32_t **dst,
  27                         const int32_t *vq_index,
  28                         const int8_t hf_vq[1024][32],
  29                         int32_t scale_factors[32][2],
  30                         ptrdiff_t sb_start, ptrdiff_t sb_end,
  31                         ptrdiff_t ofs, ptrdiff_t len)
  32 {
  33     int i, j;
  34
  35     for (i = sb_start; i < sb_end; i++) {
  36         const int8_t *coeff = hf_vq[vq_index[i]];
  37         int32_t scale = scale_factors[i][0];
  38         for (j = 0; j < len; j++)
  39             dst[i][j + ofs] = clip23(coeff[j] * scale + (1 << 3) >> 4);
  40     }
  41 }
  42
  43 static void decode_joint_c(int32_t **dst, int32_t **src,
  44                            const int32_t *scale_factors,
  45                            ptrdiff_t sb_start, ptrdiff_t sb_end,
  46                            ptrdiff_t ofs, ptrdiff_t len)
  47 {
  48     int i, j;
  49
  50     for (i = sb_start; i < sb_end; i++) {
  51         int32_t scale = scale_factors[i];
  52         for (j = 0; j < len; j++)
  53             dst[i][j + ofs] = clip23(mul17(src[i][j + ofs], scale));
  54     }
  55 }
  56
  57 static void lfe_fir_float_c(float *pcm_samples, int32_t *lfe_samples,
  58                             const float *filter_coeff, ptrdiff_t npcmblocks,
  59                             int dec_select)
  60 {
  61     // Select decimation factor
  62     int factor = 64 << dec_select;
  63     int ncoeffs = 8 >> dec_select;
  64     int nlfesamples = npcmblocks >> (dec_select + 1);
  65     int i, j, k;
  66
  67     for (i = 0; i < nlfesamples; i++) {
  68         // One decimated sample generates 64 or 128 interpolated ones
  69         for (j = 0; j < factor / 2; j++) {
  70             float a = 0;
  71             float b = 0;
  72
  73             for (k = 0; k < ncoeffs; k++) {
  74                 a += filter_coeff[      j * ncoeffs + k] * lfe_samples[-k];
  75                 b += filter_coeff[255 - j * ncoeffs - k] * lfe_samples[-k];
  76             }
  77
  78             pcm_samples[             j] = a;
  79             pcm_samples[factor / 2 + j] = b;
  80         }
  81
  82         lfe_samples++;
  83         pcm_samples += factor;
  84     }
  85 }
  86
  87 static void lfe_fir0_float_c(float *pcm_samples, int32_t *lfe_samples,
  88                              const float *filter_coeff, ptrdiff_t npcmblocks)
  89 {
  90     lfe_fir_float_c(pcm_samples, lfe_samples, filter_coeff, npcmblocks, 0);
  91 }
  92
  93 static void lfe_fir1_float_c(float *pcm_samples, int32_t *lfe_samples,
  94                              const float *filter_coeff, ptrdiff_t npcmblocks)
  95 {
  96     lfe_fir_float_c(pcm_samples, lfe_samples, filter_coeff, npcmblocks, 1);
  97 }
  98
  99 static void lfe_x96_float_c(float *dst, const float *src,
 100                             float *hist, ptrdiff_t len)
 101 {
 102     float prev = *hist;
 103     int i;
 104
 105     for (i = 0; i < len; i++) {
 106         float a = 0.25f * src[i] + 0.75f * prev;
 107         float b = 0.75f * src[i] + 0.25f * prev;
 108         prev = src[i];
 109         *dst++ = a;
 110         *dst++ = b;
 111     }
 112
 113     *hist = prev;
 114 }
 115
 116 static void sub_qmf32_float_c(SynthFilterContext *synth,
 117                               AVTXContext *imdct,
 118                               av_tx_fn imdct_fn,
 119                               float *pcm_samples,
 120                               int32_t **subband_samples_lo,
 121                               int32_t **subband_samples_hi,
 122                               float *hist1, int *offset, float *hist2,
 123                               const float *filter_coeff, ptrdiff_t npcmblocks,
 124                               float scale)
 125 {
 126     LOCAL_ALIGNED_32(float, input, [32]);
 127     int i, j;
 128
 129     for (j = 0; j < npcmblocks; j++) {
 130         // Load in one sample from each subband
 131         for (i = 0; i < 32; i++) {
 132             if ((i - 1) & 2)
 133                 input[i] = -subband_samples_lo[i][j];
 134             else
 135                 input[i] =  subband_samples_lo[i][j];
 136         }
 137
 138         // One subband sample generates 32 interpolated ones
 139         synth->synth_filter_float(imdct, hist1, offset,
 140                                   hist2, filter_coeff,
 141                                   pcm_samples, input, scale, imdct_fn);
 142         pcm_samples += 32;
 143     }
 144 }
 145
 146 static void sub_qmf64_float_c(SynthFilterContext *synth,
 147                               AVTXContext *imdct,
 148                               av_tx_fn imdct_fn,
 149                               float *pcm_samples,
 150                               int32_t **subband_samples_lo,
 151                               int32_t **subband_samples_hi,
 152                               float *hist1, int *offset, float *hist2,
 153                               const float *filter_coeff, ptrdiff_t npcmblocks,
 154                               float scale)
 155 {
 156     LOCAL_ALIGNED_32(float, input, [64]);
 157     int i, j;
 158
 159     if (!subband_samples_hi)
 160         memset(&input[32], 0, sizeof(input[0]) * 32);
 161
 162     for (j = 0; j < npcmblocks; j++) {
 163         // Load in one sample from each subband
 164         if (subband_samples_hi) {
 165             // Full 64 subbands, first 32 are residual coded
 166             for (i =  0; i < 32; i++) {
 167                 if ((i - 1) & 2)
 168                     input[i] = -subband_samples_lo[i][j] - subband_samples_hi[i][j];
 169                 else
 170                     input[i] =  subband_samples_lo[i][j] + subband_samples_hi[i][j];
 171             }
 172             for (i = 32; i < 64; i++) {
 173                 if ((i - 1) & 2)
 174                     input[i] = -subband_samples_hi[i][j];
 175                 else
 176                     input[i] =  subband_samples_hi[i][j];
 177             }
 178         } else {
 179             // Only first 32 subbands
 180             for (i =  0; i < 32; i++) {
 181                 if ((i - 1) & 2)
 182                     input[i] = -subband_samples_lo[i][j];
 183                 else
 184                     input[i] =  subband_samples_lo[i][j];
 185             }
 186         }
 187
 188         // One subband sample generates 64 interpolated ones
 189         synth->synth_filter_float_64(imdct, hist1, offset,
 190                                      hist2, filter_coeff,
 191                                      pcm_samples, input, scale, imdct_fn);
 192         pcm_samples += 64;
 193     }
 194 }
 195
 196 static void lfe_fir_fixed_c(int32_t *pcm_samples, int32_t *lfe_samples,
 197                             const int32_t *filter_coeff, ptrdiff_t npcmblocks)
 198 {
 199     // Select decimation factor
 200     int nlfesamples = npcmblocks >> 1;
 201     int i, j, k;
 202
 203     for (i = 0; i < nlfesamples; i++) {
 204         // One decimated sample generates 64 interpolated ones
 205         for (j = 0; j < 32; j++) {
 206             int64_t a = 0;
 207             int64_t b = 0;
 208
 209             for (k = 0; k < 8; k++) {
 210                 a += (int64_t)filter_coeff[      j * 8 + k] * lfe_samples[-k];
 211                 b += (int64_t)filter_coeff[255 - j * 8 - k] * lfe_samples[-k];
 212             }
 213
 214             pcm_samples[     j] = clip23(norm23(a));
 215             pcm_samples[32 + j] = clip23(norm23(b));
 216         }
 217
 218         lfe_samples++;
 219         pcm_samples += 64;
 220     }
 221 }
 222
 223 static void lfe_x96_fixed_c(int32_t *dst, const int32_t *src,
 224                             int32_t *hist, ptrdiff_t len)
 225 {
 226     int32_t prev = *hist;
 227     int i;
 228
 229     for (i = 0; i < len; i++) {
 230         int64_t a = INT64_C(2097471) * src[i] + INT64_C(6291137) * prev;
 231         int64_t b = INT64_C(6291137) * src[i] + INT64_C(2097471) * prev;
 232         prev = src[i];
 233         *dst++ = clip23(norm23(a));
 234         *dst++ = clip23(norm23(b));
 235     }
 236
 237     *hist = prev;
 238 }
 239
 240 static void sub_qmf32_fixed_c(SynthFilterContext *synth,
 241                               DCADCTContext *imdct,
 242                               int32_t *pcm_samples,
 243                               int32_t **subband_samples_lo,
 244                               int32_t **subband_samples_hi,
 245                               int32_t *hist1, int *offset, int32_t *hist2,
 246                               const int32_t *filter_coeff, ptrdiff_t npcmblocks)
 247 {
 248     LOCAL_ALIGNED_32(int32_t, input, [32]);
 249     int i, j;
 250
 251     for (j = 0; j < npcmblocks; j++) {
 252         // Load in one sample from each subband
 253         for (i = 0; i < 32; i++)
 254             input[i] = subband_samples_lo[i][j];
 255
 256         // One subband sample generates 32 interpolated ones
 257         synth->synth_filter_fixed(imdct, hist1, offset,
 258                                   hist2, filter_coeff,
 259                                   pcm_samples, input);
 260         pcm_samples += 32;
 261     }
 262 }
 263
 264 static void sub_qmf64_fixed_c(SynthFilterContext *synth,
 265                               DCADCTContext *imdct,
 266                               int32_t *pcm_samples,
 267                               int32_t **subband_samples_lo,
 268                               int32_t **subband_samples_hi,
 269                               int32_t *hist1, int *offset, int32_t *hist2,
 270                               const int32_t *filter_coeff, ptrdiff_t npcmblocks)
 271 {
 272     LOCAL_ALIGNED_32(int32_t, input, [64]);
 273     int i, j;
 274
 275     if (!subband_samples_hi)
 276         memset(&input[32], 0, sizeof(input[0]) * 32);
 277
 278     for (j = 0; j < npcmblocks; j++) {
 279         // Load in one sample from each subband
 280         if (subband_samples_hi) {
 281             // Full 64 subbands, first 32 are residual coded
 282             for (i =  0; i < 32; i++)
 283                 input[i] = subband_samples_lo[i][j] + subband_samples_hi[i][j];
 284             for (i = 32; i < 64; i++)
 285                 input[i] = subband_samples_hi[i][j];
 286         } else {
 287             // Only first 32 subbands
 288             for (i =  0; i < 32; i++)
 289                 input[i] = subband_samples_lo[i][j];
 290         }
 291
 292         // One subband sample generates 64 interpolated ones
 293         synth->synth_filter_fixed_64(imdct, hist1, offset,
 294                                      hist2, filter_coeff,
 295                                      pcm_samples, input);
 296         pcm_samples += 64;
 297     }
 298 }
 299
 300 static void decor_c(int32_t *dst, const int32_t *src, int coeff, ptrdiff_t len)
 301 {
 302     int i;
 303
 304     for (i = 0; i < len; i++)
 305         dst[i] += (SUINT)((int)(src[i] * (SUINT)coeff + (1 << 2)) >> 3);
 306 }
 307
 308 static void dmix_sub_xch_c(int32_t *dst1, int32_t *dst2,
 309                            const int32_t *src, ptrdiff_t len)
 310 {
 311     int i;
 312
 313     for (i = 0; i < len; i++) {
 314         int32_t cs = mul23(src[i], 5931520 /* M_SQRT1_2 * (1 << 23) */);
 315         dst1[i] -= cs;
 316         dst2[i] -= cs;
 317     }
 318 }
 319
 320 static void dmix_sub_c(int32_t *dst, const int32_t *src, int coeff, ptrdiff_t len)
 321 {
 322     int i;
 323
 324     for (i = 0; i < len; i++)
 325         dst[i] -= (unsigned)mul15(src[i], coeff);
 326 }
 327
 328 static void dmix_add_c(int32_t *dst, const int32_t *src, int coeff, ptrdiff_t len)
 329 {
 330     int i;
 331
 332     for (i = 0; i < len; i++)
 333         dst[i] += (unsigned)mul15(src[i], coeff);
 334 }
 335
 336 static void dmix_scale_c(int32_t *dst, int scale, ptrdiff_t len)
 337 {
 338     int i;
 339
 340     for (i = 0; i < len; i++)
 341         dst[i] = mul15(dst[i], scale);
 342 }
 343
 344 static void dmix_scale_inv_c(int32_t *dst, int scale_inv, ptrdiff_t len)
 345 {
 346     int i;
 347
 348     for (i = 0; i < len; i++)
 349         dst[i] = mul16(dst[i], scale_inv);
 350 }
 351
 352 static void filter0(SUINT32 *dst, const int32_t *src, int32_t coeff, ptrdiff_t len)
 353 {
 354     int i;
 355
 356     for (i = 0; i < len; i++)
 357         dst[i] -= mul22(src[i], coeff);
 358 }
 359
 360 static void filter1(SUINT32 *dst, const int32_t *src, int32_t coeff, ptrdiff_t len)
 361 {
 362     int i;
 363
 364     for (i = 0; i < len; i++)
 365         dst[i] -= mul23(src[i], coeff);
 366 }
 367
 368 static void assemble_freq_bands_c(int32_t *dst, int32_t *src0, int32_t *src1,
 369                                   const int32_t *coeff, ptrdiff_t len)
 370 {
 371     int i;
 372
 373     filter0(src0, src1, coeff[0], len);
 374     filter0(src1, src0, coeff[1], len);
 375     filter0(src0, src1, coeff[2], len);
 376     filter0(src1, src0, coeff[3], len);
 377
 378     for (i = 0; i < 8; i++, src0--) {
 379         filter1(src0, src1, coeff[i +  4], len);
 380         filter1(src1, src0, coeff[i + 12], len);
 381         filter1(src0, src1, coeff[i +  4], len);
 382     }
 383
 384     for (i = 0; i < len; i++) {
 385         *dst++ = *src1++;
 386         *dst++ = *++src0;
 387     }
 388 }
 389
 390 static void lbr_bank_c(float output[32][4], float **input,
 391                        const float *coeff, ptrdiff_t ofs, ptrdiff_t len)
 392 {
 393     float SW0 = coeff[0];
 394     float SW1 = coeff[1];
 395     float SW2 = coeff[2];
 396     float SW3 = coeff[3];
 397
 398     float C1  = coeff[4];
 399     float C2  = coeff[5];
 400     float C3  = coeff[6];
 401     float C4  = coeff[7];
 402
 403     float AL1 = coeff[8];
 404     float AL2 = coeff[9];
 405
 406     int i;
 407
 408     // Short window and 8 point forward MDCT
 409     for (i = 0; i < len; i++) {
 410         float *src = input[i] + ofs;
 411
 412         float a = src[-4] * SW0 - src[-1] * SW3;
 413         float b = src[-3] * SW1 - src[-2] * SW2;
 414         float c = src[ 2] * SW1 + src[ 1] * SW2;
 415         float d = src[ 3] * SW0 + src[ 0] * SW3;
 416
 417         output[i][0] = C1 * b - C2 * c + C4 * a - C3 * d;
 418         output[i][1] = C1 * d - C2 * a - C4 * b - C3 * c;
 419         output[i][2] = C3 * b + C2 * d - C4 * c + C1 * a;
 420         output[i][3] = C3 * a - C2 * b + C4 * d - C1 * c;
 421     }
 422
 423     // Aliasing cancellation for high frequencies
 424     for (i = 12; i < len - 1; i++) {
 425         float a = output[i  ][3] * AL1;
 426         float b = output[i+1][0] * AL1;
 427         output[i  ][3] += b - a;
 428         output[i+1][0] -= b + a;
 429         a = output[i  ][2] * AL2;
 430         b = output[i+1][1] * AL2;
 431         output[i  ][2] += b - a;
 432         output[i+1][1] -= b + a;
 433     }
 434 }
 435
 436 static void lfe_iir_c(float *output, const float *input,
 437                       const float iir[5][4], float hist[5][2],
 438                       ptrdiff_t factor)
 439 {
 440     float res, tmp;
 441     int i, j, k;
 442
 443     for (i = 0; i < 64; i++) {
 444         res = *input++;
 445
 446         for (j = 0; j < factor; j++) {
 447             for (k = 0; k < 5; k++) {
 448                 tmp = hist[k][0] * iir[k][0] + hist[k][1] * iir[k][1] + res;
 449                 res = hist[k][0] * iir[k][2] + hist[k][1] * iir[k][3] + tmp;
 450
 451                 hist[k][0] = hist[k][1];
 452                 hist[k][1] = tmp;
 453             }
 454
 455             *output++ = res;
 456             res = 0;
 457         }
 458     }
 459 }
 460
 461 av_cold void ff_dcadsp_init(DCADSPContext *s)
 462 {
 463     s->decode_hf     = decode_hf_c;
 464     s->decode_joint  = decode_joint_c;
 465
 466     s->lfe_fir_float[0] = lfe_fir0_float_c;
 467     s->lfe_fir_float[1] = lfe_fir1_float_c;
 468     s->lfe_x96_float    = lfe_x96_float_c;
 469     s->sub_qmf_float[0] = sub_qmf32_float_c;
 470     s->sub_qmf_float[1] = sub_qmf64_float_c;
 471
 472     s->lfe_fir_fixed    = lfe_fir_fixed_c;
 473     s->lfe_x96_fixed    = lfe_x96_fixed_c;
 474     s->sub_qmf_fixed[0] = sub_qmf32_fixed_c;
 475     s->sub_qmf_fixed[1] = sub_qmf64_fixed_c;
 476
 477     s->decor   = decor_c;
 478
 479     s->dmix_sub_xch   = dmix_sub_xch_c;
 480     s->dmix_sub       = dmix_sub_c;
 481     s->dmix_add       = dmix_add_c;
 482     s->dmix_scale     = dmix_scale_c;
 483     s->dmix_scale_inv = dmix_scale_inv_c;
 484
 485     s->assemble_freq_bands = assemble_freq_bands_c;
 486
 487     s->lbr_bank = lbr_bank_c;
 488     s->lfe_iir = lfe_iir_c;
 489
 490 #if ARCH_X86
 491     ff_dcadsp_init_x86(s);
 492 #endif
 493 }