git.sesse.net Git - ffmpeg/blob - libavcodec/ac3enc.c

   1 /*
   2  * The simplest AC3 encoder
   3  * Copyright (c) 2000 Fabrice Bellard.
   4  *
   5  * This file is part of FFmpeg.
   6  *
   7  * FFmpeg is free software; you can redistribute it and/or
   8  * modify it under the terms of the GNU Lesser General Public
   9  * License as published by the Free Software Foundation; either
  10  * version 2.1 of the License, or (at your option) any later version.
  11  *
  12  * FFmpeg is distributed in the hope that it will be useful,
  13  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  14  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  15  * Lesser General Public License for more details.
  16  *
  17  * You should have received a copy of the GNU Lesser General Public
  18  * License along with FFmpeg; if not, write to the Free Software
  19  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
  20  */
  21
  22 /**
  23  * @file ac3enc.c
  24  * The simplest AC3 encoder.
  25  */
  26 //#define DEBUG
  27 //#define DEBUG_BITALLOC
  28 #include "avcodec.h"
  29 #include "bitstream.h"
  30 #include "crc.h"
  31 #include "ac3.h"
  32
  33 typedef struct AC3EncodeContext {
  34     PutBitContext pb;
  35     int nb_channels;
  36     int nb_all_channels;
  37     int lfe_channel;
  38     int bit_rate;
  39     unsigned int sample_rate;
  40     unsigned int bitstream_id;
  41     unsigned int frame_size_min; /* minimum frame size in case rounding is necessary */
  42     unsigned int frame_size; /* current frame size in words */
  43     unsigned int bits_written;
  44     unsigned int samples_written;
  45     int sr_shift;
  46     unsigned int frame_size_code;
  47     unsigned int sr_code; /* frequency */
  48     unsigned int channel_mode;
  49     int lfe;
  50     unsigned int bitstream_mode;
  51     short last_samples[AC3_MAX_CHANNELS][256];
  52     unsigned int chbwcod[AC3_MAX_CHANNELS];
  53     int nb_coefs[AC3_MAX_CHANNELS];
  54
  55     /* bitrate allocation control */
  56     int slow_gain_code, slow_decay_code, fast_decay_code, db_per_bit_code, floor_code;
  57     AC3BitAllocParameters bit_alloc;
  58     int coarse_snr_offset;
  59     int fast_gain_code[AC3_MAX_CHANNELS];
  60     int fine_snr_offset[AC3_MAX_CHANNELS];
  61     /* mantissa encoding */
  62     int mant1_cnt, mant2_cnt, mant4_cnt;
  63 } AC3EncodeContext;
  64
  65 static int16_t costab[64];
  66 static int16_t sintab[64];
  67 static int16_t fft_rev[512];
  68 static int16_t xcos1[128];
  69 static int16_t xsin1[128];
  70
  71 #define MDCT_NBITS 9
  72 #define N         (1 << MDCT_NBITS)
  73
  74 /* new exponents are sent if their Norm 1 exceed this number */
  75 #define EXP_DIFF_THRESHOLD 1000
  76
  77 static void fft_init(int ln);
  78
  79 static inline int16_t fix15(float a)
  80 {
  81     int v;
  82     v = (int)(a * (float)(1 << 15));
  83     if (v < -32767)
  84         v = -32767;
  85     else if (v > 32767)
  86         v = 32767;
  87     return v;
  88 }
  89
  90 typedef struct IComplex {
  91     short re,im;
  92 } IComplex;
  93
  94 static void fft_init(int ln)
  95 {
  96     int i, j, m, n;
  97     float alpha;
  98
  99     n = 1 << ln;
 100
 101     for(i=0;i<(n/2);i++) {
 102         alpha = 2 * M_PI * (float)i / (float)n;
 103         costab[i] = fix15(cos(alpha));
 104         sintab[i] = fix15(sin(alpha));
 105     }
 106
 107     for(i=0;i<n;i++) {
 108         m=0;
 109         for(j=0;j<ln;j++) {
 110             m |= ((i >> j) & 1) << (ln-j-1);
 111         }
 112         fft_rev[i]=m;
 113     }
 114 }
 115
 116 /* butter fly op */
 117 #define BF(pre, pim, qre, qim, pre1, pim1, qre1, qim1) \
 118 {\
 119   int ax, ay, bx, by;\
 120   bx=pre1;\
 121   by=pim1;\
 122   ax=qre1;\
 123   ay=qim1;\
 124   pre = (bx + ax) >> 1;\
 125   pim = (by + ay) >> 1;\
 126   qre = (bx - ax) >> 1;\
 127   qim = (by - ay) >> 1;\
 128 }
 129
 130 #define MUL16(a,b) ((a) * (b))
 131
 132 #define CMUL(pre, pim, are, aim, bre, bim) \
 133 {\
 134    pre = (MUL16(are, bre) - MUL16(aim, bim)) >> 15;\
 135    pim = (MUL16(are, bim) + MUL16(bre, aim)) >> 15;\
 136 }
 137
 138
 139 /* do a 2^n point complex fft on 2^ln points. */
 140 static void fft(IComplex *z, int ln)
 141 {
 142     int        j, l, np, np2;
 143     int        nblocks, nloops;
 144     register IComplex *p,*q;
 145     int tmp_re, tmp_im;
 146
 147     np = 1 << ln;
 148
 149     /* reverse */
 150     for(j=0;j<np;j++) {
 151         int k;
 152         IComplex tmp;
 153         k = fft_rev[j];
 154         if (k < j) {
 155             tmp = z[k];
 156             z[k] = z[j];
 157             z[j] = tmp;
 158         }
 159     }
 160
 161     /* pass 0 */
 162
 163     p=&z[0];
 164     j=(np >> 1);
 165     do {
 166         BF(p[0].re, p[0].im, p[1].re, p[1].im,
 167            p[0].re, p[0].im, p[1].re, p[1].im);
 168         p+=2;
 169     } while (--j != 0);
 170
 171     /* pass 1 */
 172
 173     p=&z[0];
 174     j=np >> 2;
 175     do {
 176         BF(p[0].re, p[0].im, p[2].re, p[2].im,
 177            p[0].re, p[0].im, p[2].re, p[2].im);
 178         BF(p[1].re, p[1].im, p[3].re, p[3].im,
 179            p[1].re, p[1].im, p[3].im, -p[3].re);
 180         p+=4;
 181     } while (--j != 0);
 182
 183     /* pass 2 .. ln-1 */
 184
 185     nblocks = np >> 3;
 186     nloops = 1 << 2;
 187     np2 = np >> 1;
 188     do {
 189         p = z;
 190         q = z + nloops;
 191         for (j = 0; j < nblocks; ++j) {
 192
 193             BF(p->re, p->im, q->re, q->im,
 194                p->re, p->im, q->re, q->im);
 195
 196             p++;
 197             q++;
 198             for(l = nblocks; l < np2; l += nblocks) {
 199                 CMUL(tmp_re, tmp_im, costab[l], -sintab[l], q->re, q->im);
 200                 BF(p->re, p->im, q->re, q->im,
 201                    p->re, p->im, tmp_re, tmp_im);
 202                 p++;
 203                 q++;
 204             }
 205             p += nloops;
 206             q += nloops;
 207         }
 208         nblocks = nblocks >> 1;
 209         nloops = nloops << 1;
 210     } while (nblocks != 0);
 211 }
 212
 213 /* do a 512 point mdct */
 214 static void mdct512(int32_t *out, int16_t *in)
 215 {
 216     int i, re, im, re1, im1;
 217     int16_t rot[N];
 218     IComplex x[N/4];
 219
 220     /* shift to simplify computations */
 221     for(i=0;i<N/4;i++)
 222         rot[i] = -in[i + 3*N/4];
 223     for(i=N/4;i<N;i++)
 224         rot[i] = in[i - N/4];
 225
 226     /* pre rotation */
 227     for(i=0;i<N/4;i++) {
 228         re = ((int)rot[2*i] - (int)rot[N-1-2*i]) >> 1;
 229         im = -((int)rot[N/2+2*i] - (int)rot[N/2-1-2*i]) >> 1;
 230         CMUL(x[i].re, x[i].im, re, im, -xcos1[i], xsin1[i]);
 231     }
 232
 233     fft(x, MDCT_NBITS - 2);
 234
 235     /* post rotation */
 236     for(i=0;i<N/4;i++) {
 237         re = x[i].re;
 238         im = x[i].im;
 239         CMUL(re1, im1, re, im, xsin1[i], xcos1[i]);
 240         out[2*i] = im1;
 241         out[N/2-1-2*i] = re1;
 242     }
 243 }
 244
 245 /* XXX: use another norm ? */
 246 static int calc_exp_diff(uint8_t *exp1, uint8_t *exp2, int n)
 247 {
 248     int sum, i;
 249     sum = 0;
 250     for(i=0;i<n;i++) {
 251         sum += abs(exp1[i] - exp2[i]);
 252     }
 253     return sum;
 254 }
 255
 256 static void compute_exp_strategy(uint8_t exp_strategy[NB_BLOCKS][AC3_MAX_CHANNELS],
 257                                  uint8_t exp[NB_BLOCKS][AC3_MAX_CHANNELS][N/2],
 258                                  int ch, int is_lfe)
 259 {
 260     int i, j;
 261     int exp_diff;
 262
 263     /* estimate if the exponent variation & decide if they should be
 264        reused in the next frame */
 265     exp_strategy[0][ch] = EXP_NEW;
 266     for(i=1;i<NB_BLOCKS;i++) {
 267         exp_diff = calc_exp_diff(exp[i][ch], exp[i-1][ch], N/2);
 268 #ifdef DEBUG
 269         av_log(NULL, AV_LOG_DEBUG, "exp_diff=%d\n", exp_diff);
 270 #endif
 271         if (exp_diff > EXP_DIFF_THRESHOLD)
 272             exp_strategy[i][ch] = EXP_NEW;
 273         else
 274             exp_strategy[i][ch] = EXP_REUSE;
 275     }
 276     if (is_lfe)
 277         return;
 278
 279     /* now select the encoding strategy type : if exponents are often
 280        recoded, we use a coarse encoding */
 281     i = 0;
 282     while (i < NB_BLOCKS) {
 283         j = i + 1;
 284         while (j < NB_BLOCKS && exp_strategy[j][ch] == EXP_REUSE)
 285             j++;
 286         switch(j - i) {
 287         case 1:
 288             exp_strategy[i][ch] = EXP_D45;
 289             break;
 290         case 2:
 291         case 3:
 292             exp_strategy[i][ch] = EXP_D25;
 293             break;
 294         default:
 295             exp_strategy[i][ch] = EXP_D15;
 296             break;
 297         }
 298         i = j;
 299     }
 300 }
 301
 302 /* set exp[i] to min(exp[i], exp1[i]) */
 303 static void exponent_min(uint8_t exp[N/2], uint8_t exp1[N/2], int n)
 304 {
 305     int i;
 306
 307     for(i=0;i<n;i++) {
 308         if (exp1[i] < exp[i])
 309             exp[i] = exp1[i];
 310     }
 311 }
 312
 313 /* update the exponents so that they are the ones the decoder will
 314    decode. Return the number of bits used to code the exponents */
 315 static int encode_exp(uint8_t encoded_exp[N/2],
 316                       uint8_t exp[N/2],
 317                       int nb_exps,
 318                       int exp_strategy)
 319 {
 320     int group_size, nb_groups, i, j, k, exp_min;
 321     uint8_t exp1[N/2];
 322
 323     switch(exp_strategy) {
 324     case EXP_D15:
 325         group_size = 1;
 326         break;
 327     case EXP_D25:
 328         group_size = 2;
 329         break;
 330     default:
 331     case EXP_D45:
 332         group_size = 4;
 333         break;
 334     }
 335     nb_groups = ((nb_exps + (group_size * 3) - 4) / (3 * group_size)) * 3;
 336
 337     /* for each group, compute the minimum exponent */
 338     exp1[0] = exp[0]; /* DC exponent is handled separately */
 339     k = 1;
 340     for(i=1;i<=nb_groups;i++) {
 341         exp_min = exp[k];
 342         assert(exp_min >= 0 && exp_min <= 24);
 343         for(j=1;j<group_size;j++) {
 344             if (exp[k+j] < exp_min)
 345                 exp_min = exp[k+j];
 346         }
 347         exp1[i] = exp_min;
 348         k += group_size;
 349     }
 350
 351     /* constraint for DC exponent */
 352     if (exp1[0] > 15)
 353         exp1[0] = 15;
 354
 355     /* Decrease the delta between each groups to within 2
 356      * so that they can be differentially encoded */
 357     for (i=1;i<=nb_groups;i++)
 358         exp1[i] = FFMIN(exp1[i], exp1[i-1] + 2);
 359     for (i=nb_groups-1;i>=0;i--)
 360         exp1[i] = FFMIN(exp1[i], exp1[i+1] + 2);
 361
 362     /* now we have the exponent values the decoder will see */
 363     encoded_exp[0] = exp1[0];
 364     k = 1;
 365     for(i=1;i<=nb_groups;i++) {
 366         for(j=0;j<group_size;j++) {
 367             encoded_exp[k+j] = exp1[i];
 368         }
 369         k += group_size;
 370     }
 371
 372 #if defined(DEBUG)
 373     av_log(NULL, AV_LOG_DEBUG, "exponents: strategy=%d\n", exp_strategy);
 374     for(i=0;i<=nb_groups * group_size;i++) {
 375         av_log(NULL, AV_LOG_DEBUG, "%d ", encoded_exp[i]);
 376     }
 377     av_log(NULL, AV_LOG_DEBUG, "\n");
 378 #endif
 379
 380     return 4 + (nb_groups / 3) * 7;
 381 }
 382
 383 /* return the size in bits taken by the mantissa */
 384 static int compute_mantissa_size(AC3EncodeContext *s, uint8_t *m, int nb_coefs)
 385 {
 386     int bits, mant, i;
 387
 388     bits = 0;
 389     for(i=0;i<nb_coefs;i++) {
 390         mant = m[i];
 391         switch(mant) {
 392         case 0:
 393             /* nothing */
 394             break;
 395         case 1:
 396             /* 3 mantissa in 5 bits */
 397             if (s->mant1_cnt == 0)
 398                 bits += 5;
 399             if (++s->mant1_cnt == 3)
 400                 s->mant1_cnt = 0;
 401             break;
 402         case 2:
 403             /* 3 mantissa in 7 bits */
 404             if (s->mant2_cnt == 0)
 405                 bits += 7;
 406             if (++s->mant2_cnt == 3)
 407                 s->mant2_cnt = 0;
 408             break;
 409         case 3:
 410             bits += 3;
 411             break;
 412         case 4:
 413             /* 2 mantissa in 7 bits */
 414             if (s->mant4_cnt == 0)
 415                 bits += 7;
 416             if (++s->mant4_cnt == 2)
 417                 s->mant4_cnt = 0;
 418             break;
 419         case 14:
 420             bits += 14;
 421             break;
 422         case 15:
 423             bits += 16;
 424             break;
 425         default:
 426             bits += mant - 1;
 427             break;
 428         }
 429     }
 430     return bits;
 431 }
 432
 433
 434 static void bit_alloc_masking(AC3EncodeContext *s,
 435                               uint8_t encoded_exp[NB_BLOCKS][AC3_MAX_CHANNELS][N/2],
 436                               uint8_t exp_strategy[NB_BLOCKS][AC3_MAX_CHANNELS],
 437                               int16_t psd[NB_BLOCKS][AC3_MAX_CHANNELS][N/2],
 438                               int16_t mask[NB_BLOCKS][AC3_MAX_CHANNELS][50])
 439 {
 440     int blk, ch;
 441     int16_t band_psd[NB_BLOCKS][AC3_MAX_CHANNELS][50];
 442
 443     for(blk=0; blk<NB_BLOCKS; blk++) {
 444         for(ch=0;ch<s->nb_all_channels;ch++) {
 445             if(exp_strategy[blk][ch] == EXP_REUSE) {
 446                 memcpy(psd[blk][ch], psd[blk-1][ch], (N/2)*sizeof(int16_t));
 447                 memcpy(mask[blk][ch], mask[blk-1][ch], 50*sizeof(int16_t));
 448             } else {
 449                 ff_ac3_bit_alloc_calc_psd(encoded_exp[blk][ch], 0,
 450                                           s->nb_coefs[ch],
 451                                           psd[blk][ch], band_psd[blk][ch]);
 452                 ff_ac3_bit_alloc_calc_mask(&s->bit_alloc, band_psd[blk][ch],
 453                                            0, s->nb_coefs[ch],
 454                                            ff_ac3_fast_gain_tab[s->fast_gain_code[ch]],
 455                                            ch == s->lfe_channel,
 456                                            DBA_NONE, 0, NULL, NULL, NULL,
 457                                            mask[blk][ch]);
 458             }
 459         }
 460     }
 461 }
 462
 463 static int bit_alloc(AC3EncodeContext *s,
 464                      int16_t mask[NB_BLOCKS][AC3_MAX_CHANNELS][50],
 465                      int16_t psd[NB_BLOCKS][AC3_MAX_CHANNELS][N/2],
 466                      uint8_t bap[NB_BLOCKS][AC3_MAX_CHANNELS][N/2],
 467                      int frame_bits, int coarse_snr_offset, int fine_snr_offset)
 468 {
 469     int i, ch;
 470     int snr_offset;
 471
 472     snr_offset = (((coarse_snr_offset - 15) << 4) + fine_snr_offset) << 2;
 473
 474     /* compute size */
 475     for(i=0;i<NB_BLOCKS;i++) {
 476         s->mant1_cnt = 0;
 477         s->mant2_cnt = 0;
 478         s->mant4_cnt = 0;
 479         for(ch=0;ch<s->nb_all_channels;ch++) {
 480             ff_ac3_bit_alloc_calc_bap(mask[i][ch], psd[i][ch], 0,
 481                                       s->nb_coefs[ch], snr_offset,
 482                                       s->bit_alloc.floor, bap[i][ch]);
 483             frame_bits += compute_mantissa_size(s, bap[i][ch],
 484                                                  s->nb_coefs[ch]);
 485         }
 486     }
 487 #if 0
 488     printf("csnr=%d fsnr=%d frame_bits=%d diff=%d\n",
 489            coarse_snr_offset, fine_snr_offset, frame_bits,
 490            16 * s->frame_size - ((frame_bits + 7) & ~7));
 491 #endif
 492     return 16 * s->frame_size - frame_bits;
 493 }
 494
 495 #define SNR_INC1 4
 496
 497 static int compute_bit_allocation(AC3EncodeContext *s,
 498                                   uint8_t bap[NB_BLOCKS][AC3_MAX_CHANNELS][N/2],
 499                                   uint8_t encoded_exp[NB_BLOCKS][AC3_MAX_CHANNELS][N/2],
 500                                   uint8_t exp_strategy[NB_BLOCKS][AC3_MAX_CHANNELS],
 501                                   int frame_bits)
 502 {
 503     int i, ch;
 504     int coarse_snr_offset, fine_snr_offset;
 505     uint8_t bap1[NB_BLOCKS][AC3_MAX_CHANNELS][N/2];
 506     int16_t psd[NB_BLOCKS][AC3_MAX_CHANNELS][N/2];
 507     int16_t mask[NB_BLOCKS][AC3_MAX_CHANNELS][50];
 508     static int frame_bits_inc[8] = { 0, 0, 2, 2, 2, 4, 2, 4 };
 509
 510     /* init default parameters */
 511     s->slow_decay_code = 2;
 512     s->fast_decay_code = 1;
 513     s->slow_gain_code = 1;
 514     s->db_per_bit_code = 2;
 515     s->floor_code = 4;
 516     for(ch=0;ch<s->nb_all_channels;ch++)
 517         s->fast_gain_code[ch] = 4;
 518
 519     /* compute real values */
 520     s->bit_alloc.sr_code = s->sr_code;
 521     s->bit_alloc.sr_shift = s->sr_shift;
 522     s->bit_alloc.slow_decay = ff_ac3_slow_decay_tab[s->slow_decay_code] >> s->sr_shift;
 523     s->bit_alloc.fast_decay = ff_ac3_fast_decay_tab[s->fast_decay_code] >> s->sr_shift;
 524     s->bit_alloc.slow_gain = ff_ac3_slow_gain_tab[s->slow_gain_code];
 525     s->bit_alloc.db_per_bit = ff_ac3_db_per_bit_tab[s->db_per_bit_code];
 526     s->bit_alloc.floor = ff_ac3_floor_tab[s->floor_code];
 527
 528     /* header size */
 529     frame_bits += 65;
 530     // if (s->channel_mode == 2)
 531     //    frame_bits += 2;
 532     frame_bits += frame_bits_inc[s->channel_mode];
 533
 534     /* audio blocks */
 535     for(i=0;i<NB_BLOCKS;i++) {
 536         frame_bits += s->nb_channels * 2 + 2; /* blksw * c, dithflag * c, dynrnge, cplstre */
 537         if (s->channel_mode == AC3_CHMODE_STEREO) {
 538             frame_bits++; /* rematstr */
 539             if(i==0) frame_bits += 4;
 540         }
 541         frame_bits += 2 * s->nb_channels; /* chexpstr[2] * c */
 542         if (s->lfe)
 543             frame_bits++; /* lfeexpstr */
 544         for(ch=0;ch<s->nb_channels;ch++) {
 545             if (exp_strategy[i][ch] != EXP_REUSE)
 546                 frame_bits += 6 + 2; /* chbwcod[6], gainrng[2] */
 547         }
 548         frame_bits++; /* baie */
 549         frame_bits++; /* snr */
 550         frame_bits += 2; /* delta / skip */
 551     }
 552     frame_bits++; /* cplinu for block 0 */
 553     /* bit alloc info */
 554     /* sdcycod[2], fdcycod[2], sgaincod[2], dbpbcod[2], floorcod[3] */
 555     /* csnroffset[6] */
 556     /* (fsnoffset[4] + fgaincod[4]) * c */
 557     frame_bits += 2*4 + 3 + 6 + s->nb_all_channels * (4 + 3);
 558
 559     /* auxdatae, crcrsv */
 560     frame_bits += 2;
 561
 562     /* CRC */
 563     frame_bits += 16;
 564
 565     /* calculate psd and masking curve before doing bit allocation */
 566     bit_alloc_masking(s, encoded_exp, exp_strategy, psd, mask);
 567
 568     /* now the big work begins : do the bit allocation. Modify the snr
 569        offset until we can pack everything in the requested frame size */
 570
 571     coarse_snr_offset = s->coarse_snr_offset;
 572     while (coarse_snr_offset >= 0 &&
 573            bit_alloc(s, mask, psd, bap, frame_bits, coarse_snr_offset, 0) < 0)
 574         coarse_snr_offset -= SNR_INC1;
 575     if (coarse_snr_offset < 0) {
 576         av_log(NULL, AV_LOG_ERROR, "Bit allocation failed. Try increasing the bitrate.\n");
 577         return -1;
 578     }
 579     while ((coarse_snr_offset + SNR_INC1) <= 63 &&
 580            bit_alloc(s, mask, psd, bap1, frame_bits,
 581                      coarse_snr_offset + SNR_INC1, 0) >= 0) {
 582         coarse_snr_offset += SNR_INC1;
 583         memcpy(bap, bap1, sizeof(bap1));
 584     }
 585     while ((coarse_snr_offset + 1) <= 63 &&
 586            bit_alloc(s, mask, psd, bap1, frame_bits, coarse_snr_offset + 1, 0) >= 0) {
 587         coarse_snr_offset++;
 588         memcpy(bap, bap1, sizeof(bap1));
 589     }
 590
 591     fine_snr_offset = 0;
 592     while ((fine_snr_offset + SNR_INC1) <= 15 &&
 593            bit_alloc(s, mask, psd, bap1, frame_bits,
 594                      coarse_snr_offset, fine_snr_offset + SNR_INC1) >= 0) {
 595         fine_snr_offset += SNR_INC1;
 596         memcpy(bap, bap1, sizeof(bap1));
 597     }
 598     while ((fine_snr_offset + 1) <= 15 &&
 599            bit_alloc(s, mask, psd, bap1, frame_bits,
 600                      coarse_snr_offset, fine_snr_offset + 1) >= 0) {
 601         fine_snr_offset++;
 602         memcpy(bap, bap1, sizeof(bap1));
 603     }
 604
 605     s->coarse_snr_offset = coarse_snr_offset;
 606     for(ch=0;ch<s->nb_all_channels;ch++)
 607         s->fine_snr_offset[ch] = fine_snr_offset;
 608 #if defined(DEBUG_BITALLOC)
 609     {
 610         int j;
 611
 612         for(i=0;i<6;i++) {
 613             for(ch=0;ch<s->nb_all_channels;ch++) {
 614                 printf("Block #%d Ch%d:\n", i, ch);
 615                 printf("bap=");
 616                 for(j=0;j<s->nb_coefs[ch];j++) {
 617                     printf("%d ",bap[i][ch][j]);
 618                 }
 619                 printf("\n");
 620             }
 621         }
 622     }
 623 #endif
 624     return 0;
 625 }
 626
 627 static int AC3_encode_init(AVCodecContext *avctx)
 628 {
 629     int freq = avctx->sample_rate;
 630     int bitrate = avctx->bit_rate;
 631     int channels = avctx->channels;
 632     AC3EncodeContext *s = avctx->priv_data;
 633     int i, j, ch;
 634     float alpha;
 635     static const uint8_t channel_mode_defs[6] = {
 636         0x01, /* C */
 637         0x02, /* L R */
 638         0x03, /* L C R */
 639         0x06, /* L R SL SR */
 640         0x07, /* L C R SL SR */
 641         0x07, /* L C R SL SR (+LFE) */
 642     };
 643
 644     avctx->frame_size = AC3_FRAME_SIZE;
 645
 646     ac3_common_init();
 647
 648     /* number of channels */
 649     if (channels < 1 || channels > 6)
 650         return -1;
 651     s->channel_mode = channel_mode_defs[channels - 1];
 652     s->lfe = (channels == 6) ? 1 : 0;
 653     s->nb_all_channels = channels;
 654     s->nb_channels = channels > 5 ? 5 : channels;
 655     s->lfe_channel = s->lfe ? 5 : -1;
 656
 657     /* frequency */
 658     for(i=0;i<3;i++) {
 659         for(j=0;j<3;j++)
 660             if ((ff_ac3_sample_rate_tab[j] >> i) == freq)
 661                 goto found;
 662     }
 663     return -1;
 664  found:
 665     s->sample_rate = freq;
 666     s->sr_shift = i;
 667     s->sr_code = j;
 668     s->bitstream_id = 8 + s->sr_shift;
 669     s->bitstream_mode = 0; /* complete main audio service */
 670
 671     /* bitrate & frame size */
 672     for(i=0;i<19;i++) {
 673         if ((ff_ac3_bitrate_tab[i] >> s->sr_shift)*1000 == bitrate)
 674             break;
 675     }
 676     if (i == 19)
 677         return -1;
 678     s->bit_rate = bitrate;
 679     s->frame_size_code = i << 1;
 680     s->frame_size_min = ff_ac3_frame_size_tab[s->frame_size_code][s->sr_code];
 681     s->bits_written = 0;
 682     s->samples_written = 0;
 683     s->frame_size = s->frame_size_min;
 684
 685     /* bit allocation init */
 686     for(ch=0;ch<s->nb_channels;ch++) {
 687         /* bandwidth for each channel */
 688         /* XXX: should compute the bandwidth according to the frame
 689            size, so that we avoid anoying high freq artefacts */
 690         s->chbwcod[ch] = 50; /* sample bandwidth as mpeg audio layer 2 table 0 */
 691         s->nb_coefs[ch] = ((s->chbwcod[ch] + 12) * 3) + 37;
 692     }
 693     if (s->lfe) {
 694         s->nb_coefs[s->lfe_channel] = 7; /* fixed */
 695     }
 696     /* initial snr offset */
 697     s->coarse_snr_offset = 40;
 698
 699     /* mdct init */
 700     fft_init(MDCT_NBITS - 2);
 701     for(i=0;i<N/4;i++) {
 702         alpha = 2 * M_PI * (i + 1.0 / 8.0) / (float)N;
 703         xcos1[i] = fix15(-cos(alpha));
 704         xsin1[i] = fix15(-sin(alpha));
 705     }
 706
 707     avctx->coded_frame= avcodec_alloc_frame();
 708     avctx->coded_frame->key_frame= 1;
 709
 710     return 0;
 711 }
 712
 713 /* output the AC3 frame header */
 714 static void output_frame_header(AC3EncodeContext *s, unsigned char *frame)
 715 {
 716     init_put_bits(&s->pb, frame, AC3_MAX_CODED_FRAME_SIZE);
 717
 718     put_bits(&s->pb, 16, 0x0b77); /* frame header */
 719     put_bits(&s->pb, 16, 0); /* crc1: will be filled later */
 720     put_bits(&s->pb, 2, s->sr_code);
 721     put_bits(&s->pb, 6, s->frame_size_code + (s->frame_size - s->frame_size_min));
 722     put_bits(&s->pb, 5, s->bitstream_id);
 723     put_bits(&s->pb, 3, s->bitstream_mode);
 724     put_bits(&s->pb, 3, s->channel_mode);
 725     if ((s->channel_mode & 0x01) && s->channel_mode != AC3_CHMODE_MONO)
 726         put_bits(&s->pb, 2, 1); /* XXX -4.5 dB */
 727     if (s->channel_mode & 0x04)
 728         put_bits(&s->pb, 2, 1); /* XXX -6 dB */
 729     if (s->channel_mode == AC3_CHMODE_STEREO)
 730         put_bits(&s->pb, 2, 0); /* surround not indicated */
 731     put_bits(&s->pb, 1, s->lfe); /* LFE */
 732     put_bits(&s->pb, 5, 31); /* dialog norm: -31 db */
 733     put_bits(&s->pb, 1, 0); /* no compression control word */
 734     put_bits(&s->pb, 1, 0); /* no lang code */
 735     put_bits(&s->pb, 1, 0); /* no audio production info */
 736     put_bits(&s->pb, 1, 0); /* no copyright */
 737     put_bits(&s->pb, 1, 1); /* original bitstream */
 738     put_bits(&s->pb, 1, 0); /* no time code 1 */
 739     put_bits(&s->pb, 1, 0); /* no time code 2 */
 740     put_bits(&s->pb, 1, 0); /* no additional bit stream info */
 741 }
 742
 743 /* symetric quantization on 'levels' levels */
 744 static inline int sym_quant(int c, int e, int levels)
 745 {
 746     int v;
 747
 748     if (c >= 0) {
 749         v = (levels * (c << e)) >> 24;
 750         v = (v + 1) >> 1;
 751         v = (levels >> 1) + v;
 752     } else {
 753         v = (levels * ((-c) << e)) >> 24;
 754         v = (v + 1) >> 1;
 755         v = (levels >> 1) - v;
 756     }
 757     assert (v >= 0 && v < levels);
 758     return v;
 759 }
 760
 761 /* asymetric quantization on 2^qbits levels */
 762 static inline int asym_quant(int c, int e, int qbits)
 763 {
 764     int lshift, m, v;
 765
 766     lshift = e + qbits - 24;
 767     if (lshift >= 0)
 768         v = c << lshift;
 769     else
 770         v = c >> (-lshift);
 771     /* rounding */
 772     v = (v + 1) >> 1;
 773     m = (1 << (qbits-1));
 774     if (v >= m)
 775         v = m - 1;
 776     assert(v >= -m);
 777     return v & ((1 << qbits)-1);
 778 }
 779
 780 /* Output one audio block. There are NB_BLOCKS audio blocks in one AC3
 781    frame */
 782 static void output_audio_block(AC3EncodeContext *s,
 783                                uint8_t exp_strategy[AC3_MAX_CHANNELS],
 784                                uint8_t encoded_exp[AC3_MAX_CHANNELS][N/2],
 785                                uint8_t bap[AC3_MAX_CHANNELS][N/2],
 786                                int32_t mdct_coefs[AC3_MAX_CHANNELS][N/2],
 787                                int8_t global_exp[AC3_MAX_CHANNELS],
 788                                int block_num)
 789 {
 790     int ch, nb_groups, group_size, i, baie, rbnd;
 791     uint8_t *p;
 792     uint16_t qmant[AC3_MAX_CHANNELS][N/2];
 793     int exp0, exp1;
 794     int mant1_cnt, mant2_cnt, mant4_cnt;
 795     uint16_t *qmant1_ptr, *qmant2_ptr, *qmant4_ptr;
 796     int delta0, delta1, delta2;
 797
 798     for(ch=0;ch<s->nb_channels;ch++)
 799         put_bits(&s->pb, 1, 0); /* 512 point MDCT */
 800     for(ch=0;ch<s->nb_channels;ch++)
 801         put_bits(&s->pb, 1, 1); /* no dither */
 802     put_bits(&s->pb, 1, 0); /* no dynamic range */
 803     if (block_num == 0) {
 804         /* for block 0, even if no coupling, we must say it. This is a
 805            waste of bit :-) */
 806         put_bits(&s->pb, 1, 1); /* coupling strategy present */
 807         put_bits(&s->pb, 1, 0); /* no coupling strategy */
 808     } else {
 809         put_bits(&s->pb, 1, 0); /* no new coupling strategy */
 810     }
 811
 812     if (s->channel_mode == AC3_CHMODE_STEREO)
 813       {
 814         if(block_num==0)
 815           {
 816             /* first block must define rematrixing (rematstr)  */
 817             put_bits(&s->pb, 1, 1);
 818
 819             /* dummy rematrixing rematflg(1:4)=0 */
 820             for (rbnd=0;rbnd<4;rbnd++)
 821               put_bits(&s->pb, 1, 0);
 822           }
 823         else
 824           {
 825             /* no matrixing (but should be used in the future) */
 826             put_bits(&s->pb, 1, 0);
 827           }
 828       }
 829
 830 #if defined(DEBUG)
 831     {
 832       static int count = 0;
 833       av_log(NULL, AV_LOG_DEBUG, "Block #%d (%d)\n", block_num, count++);
 834     }
 835 #endif
 836     /* exponent strategy */
 837     for(ch=0;ch<s->nb_channels;ch++) {
 838         put_bits(&s->pb, 2, exp_strategy[ch]);
 839     }
 840
 841     if (s->lfe) {
 842         put_bits(&s->pb, 1, exp_strategy[s->lfe_channel]);
 843     }
 844
 845     for(ch=0;ch<s->nb_channels;ch++) {
 846         if (exp_strategy[ch] != EXP_REUSE)
 847             put_bits(&s->pb, 6, s->chbwcod[ch]);
 848     }
 849
 850     /* exponents */
 851     for (ch = 0; ch < s->nb_all_channels; ch++) {
 852         switch(exp_strategy[ch]) {
 853         case EXP_REUSE:
 854             continue;
 855         case EXP_D15:
 856             group_size = 1;
 857             break;
 858         case EXP_D25:
 859             group_size = 2;
 860             break;
 861         default:
 862         case EXP_D45:
 863             group_size = 4;
 864             break;
 865         }
 866         nb_groups = (s->nb_coefs[ch] + (group_size * 3) - 4) / (3 * group_size);
 867         p = encoded_exp[ch];
 868
 869         /* first exponent */
 870         exp1 = *p++;
 871         put_bits(&s->pb, 4, exp1);
 872
 873         /* next ones are delta encoded */
 874         for(i=0;i<nb_groups;i++) {
 875             /* merge three delta in one code */
 876             exp0 = exp1;
 877             exp1 = p[0];
 878             p += group_size;
 879             delta0 = exp1 - exp0 + 2;
 880
 881             exp0 = exp1;
 882             exp1 = p[0];
 883             p += group_size;
 884             delta1 = exp1 - exp0 + 2;
 885
 886             exp0 = exp1;
 887             exp1 = p[0];
 888             p += group_size;
 889             delta2 = exp1 - exp0 + 2;
 890
 891             put_bits(&s->pb, 7, ((delta0 * 5 + delta1) * 5) + delta2);
 892         }
 893
 894         if (ch != s->lfe_channel)
 895             put_bits(&s->pb, 2, 0); /* no gain range info */
 896     }
 897
 898     /* bit allocation info */
 899     baie = (block_num == 0);
 900     put_bits(&s->pb, 1, baie);
 901     if (baie) {
 902         put_bits(&s->pb, 2, s->slow_decay_code);
 903         put_bits(&s->pb, 2, s->fast_decay_code);
 904         put_bits(&s->pb, 2, s->slow_gain_code);
 905         put_bits(&s->pb, 2, s->db_per_bit_code);
 906         put_bits(&s->pb, 3, s->floor_code);
 907     }
 908
 909     /* snr offset */
 910     put_bits(&s->pb, 1, baie); /* always present with bai */
 911     if (baie) {
 912         put_bits(&s->pb, 6, s->coarse_snr_offset);
 913         for(ch=0;ch<s->nb_all_channels;ch++) {
 914             put_bits(&s->pb, 4, s->fine_snr_offset[ch]);
 915             put_bits(&s->pb, 3, s->fast_gain_code[ch]);
 916         }
 917     }
 918
 919     put_bits(&s->pb, 1, 0); /* no delta bit allocation */
 920     put_bits(&s->pb, 1, 0); /* no data to skip */
 921
 922     /* mantissa encoding : we use two passes to handle the grouping. A
 923        one pass method may be faster, but it would necessitate to
 924        modify the output stream. */
 925
 926     /* first pass: quantize */
 927     mant1_cnt = mant2_cnt = mant4_cnt = 0;
 928     qmant1_ptr = qmant2_ptr = qmant4_ptr = NULL;
 929
 930     for (ch = 0; ch < s->nb_all_channels; ch++) {
 931         int b, c, e, v;
 932
 933         for(i=0;i<s->nb_coefs[ch];i++) {
 934             c = mdct_coefs[ch][i];
 935             e = encoded_exp[ch][i] - global_exp[ch];
 936             b = bap[ch][i];
 937             switch(b) {
 938             case 0:
 939                 v = 0;
 940                 break;
 941             case 1:
 942                 v = sym_quant(c, e, 3);
 943                 switch(mant1_cnt) {
 944                 case 0:
 945                     qmant1_ptr = &qmant[ch][i];
 946                     v = 9 * v;
 947                     mant1_cnt = 1;
 948                     break;
 949                 case 1:
 950                     *qmant1_ptr += 3 * v;
 951                     mant1_cnt = 2;
 952                     v = 128;
 953                     break;
 954                 default:
 955                     *qmant1_ptr += v;
 956                     mant1_cnt = 0;
 957                     v = 128;
 958                     break;
 959                 }
 960                 break;
 961             case 2:
 962                 v = sym_quant(c, e, 5);
 963                 switch(mant2_cnt) {
 964                 case 0:
 965                     qmant2_ptr = &qmant[ch][i];
 966                     v = 25 * v;
 967                     mant2_cnt = 1;
 968                     break;
 969                 case 1:
 970                     *qmant2_ptr += 5 * v;
 971                     mant2_cnt = 2;
 972                     v = 128;
 973                     break;
 974                 default:
 975                     *qmant2_ptr += v;
 976                     mant2_cnt = 0;
 977                     v = 128;
 978                     break;
 979                 }
 980                 break;
 981             case 3:
 982                 v = sym_quant(c, e, 7);
 983                 break;
 984             case 4:
 985                 v = sym_quant(c, e, 11);
 986                 switch(mant4_cnt) {
 987                 case 0:
 988                     qmant4_ptr = &qmant[ch][i];
 989                     v = 11 * v;
 990                     mant4_cnt = 1;
 991                     break;
 992                 default:
 993                     *qmant4_ptr += v;
 994                     mant4_cnt = 0;
 995                     v = 128;
 996                     break;
 997                 }
 998                 break;
 999             case 5:
1000                 v = sym_quant(c, e, 15);
1001                 break;
1002             case 14:
1003                 v = asym_quant(c, e, 14);
1004                 break;
1005             case 15:
1006                 v = asym_quant(c, e, 16);
1007                 break;
1008             default:
1009                 v = asym_quant(c, e, b - 1);
1010                 break;
1011             }
1012             qmant[ch][i] = v;
1013         }
1014     }
1015
1016     /* second pass : output the values */
1017     for (ch = 0; ch < s->nb_all_channels; ch++) {
1018         int b, q;
1019
1020         for(i=0;i<s->nb_coefs[ch];i++) {
1021             q = qmant[ch][i];
1022             b = bap[ch][i];
1023             switch(b) {
1024             case 0:
1025                 break;
1026             case 1:
1027                 if (q != 128)
1028                     put_bits(&s->pb, 5, q);
1029                 break;
1030             case 2:
1031                 if (q != 128)
1032                     put_bits(&s->pb, 7, q);
1033                 break;
1034             case 3:
1035                 put_bits(&s->pb, 3, q);
1036                 break;
1037             case 4:
1038                 if (q != 128)
1039                     put_bits(&s->pb, 7, q);
1040                 break;
1041             case 14:
1042                 put_bits(&s->pb, 14, q);
1043                 break;
1044             case 15:
1045                 put_bits(&s->pb, 16, q);
1046                 break;
1047             default:
1048                 put_bits(&s->pb, b - 1, q);
1049                 break;
1050             }
1051         }
1052     }
1053 }
1054
1055 #define CRC16_POLY ((1 << 0) | (1 << 2) | (1 << 15) | (1 << 16))
1056
1057 static unsigned int mul_poly(unsigned int a, unsigned int b, unsigned int poly)
1058 {
1059     unsigned int c;
1060
1061     c = 0;
1062     while (a) {
1063         if (a & 1)
1064             c ^= b;
1065         a = a >> 1;
1066         b = b << 1;
1067         if (b & (1 << 16))
1068             b ^= poly;
1069     }
1070     return c;
1071 }
1072
1073 static unsigned int pow_poly(unsigned int a, unsigned int n, unsigned int poly)
1074 {
1075     unsigned int r;
1076     r = 1;
1077     while (n) {
1078         if (n & 1)
1079             r = mul_poly(r, a, poly);
1080         a = mul_poly(a, a, poly);
1081         n >>= 1;
1082     }
1083     return r;
1084 }
1085
1086
1087 /* compute log2(max(abs(tab[]))) */
1088 static int log2_tab(int16_t *tab, int n)
1089 {
1090     int i, v;
1091
1092     v = 0;
1093     for(i=0;i<n;i++) {
1094         v |= abs(tab[i]);
1095     }
1096     return av_log2(v);
1097 }
1098
1099 static void lshift_tab(int16_t *tab, int n, int lshift)
1100 {
1101     int i;
1102
1103     if (lshift > 0) {
1104         for(i=0;i<n;i++) {
1105             tab[i] <<= lshift;
1106         }
1107     } else if (lshift < 0) {
1108         lshift = -lshift;
1109         for(i=0;i<n;i++) {
1110             tab[i] >>= lshift;
1111         }
1112     }
1113 }
1114
1115 /* fill the end of the frame and compute the two crcs */
1116 static int output_frame_end(AC3EncodeContext *s)
1117 {
1118     int frame_size, frame_size_58, n, crc1, crc2, crc_inv;
1119     uint8_t *frame;
1120
1121     frame_size = s->frame_size; /* frame size in words */
1122     /* align to 8 bits */
1123     flush_put_bits(&s->pb);
1124     /* add zero bytes to reach the frame size */
1125     frame = s->pb.buf;
1126     n = 2 * s->frame_size - (pbBufPtr(&s->pb) - frame) - 2;
1127     assert(n >= 0);
1128     if(n>0)
1129       memset(pbBufPtr(&s->pb), 0, n);
1130
1131     /* Now we must compute both crcs : this is not so easy for crc1
1132        because it is at the beginning of the data... */
1133     frame_size_58 = (frame_size >> 1) + (frame_size >> 3);
1134     crc1 = bswap_16(av_crc(av_crc8005, 0, frame + 4, 2 * frame_size_58 - 4));
1135     /* XXX: could precompute crc_inv */
1136     crc_inv = pow_poly((CRC16_POLY >> 1), (16 * frame_size_58) - 16, CRC16_POLY);
1137     crc1 = mul_poly(crc_inv, crc1, CRC16_POLY);
1138     AV_WB16(frame+2,crc1);
1139
1140     crc2 = bswap_16(av_crc(av_crc8005, 0, frame + 2 * frame_size_58, (frame_size - frame_size_58) * 2 - 2));
1141     AV_WB16(frame+2*frame_size-2,crc2);
1142
1143     //    printf("n=%d frame_size=%d\n", n, frame_size);
1144     return frame_size * 2;
1145 }
1146
1147 static int AC3_encode_frame(AVCodecContext *avctx,
1148                             unsigned char *frame, int buf_size, void *data)
1149 {
1150     AC3EncodeContext *s = avctx->priv_data;
1151     int16_t *samples = data;
1152     int i, j, k, v, ch;
1153     int16_t input_samples[N];
1154     int32_t mdct_coef[NB_BLOCKS][AC3_MAX_CHANNELS][N/2];
1155     uint8_t exp[NB_BLOCKS][AC3_MAX_CHANNELS][N/2];
1156     uint8_t exp_strategy[NB_BLOCKS][AC3_MAX_CHANNELS];
1157     uint8_t encoded_exp[NB_BLOCKS][AC3_MAX_CHANNELS][N/2];
1158     uint8_t bap[NB_BLOCKS][AC3_MAX_CHANNELS][N/2];
1159     int8_t exp_samples[NB_BLOCKS][AC3_MAX_CHANNELS];
1160     int frame_bits;
1161
1162     frame_bits = 0;
1163     for(ch=0;ch<s->nb_all_channels;ch++) {
1164         /* fixed mdct to the six sub blocks & exponent computation */
1165         for(i=0;i<NB_BLOCKS;i++) {
1166             int16_t *sptr;
1167             int sinc;
1168
1169             /* compute input samples */
1170             memcpy(input_samples, s->last_samples[ch], N/2 * sizeof(int16_t));
1171             sinc = s->nb_all_channels;
1172             sptr = samples + (sinc * (N/2) * i) + ch;
1173             for(j=0;j<N/2;j++) {
1174                 v = *sptr;
1175                 input_samples[j + N/2] = v;
1176                 s->last_samples[ch][j] = v;
1177                 sptr += sinc;
1178             }
1179
1180             /* apply the MDCT window */
1181             for(j=0;j<N/2;j++) {
1182                 input_samples[j] = MUL16(input_samples[j],
1183                                          ff_ac3_window[j]) >> 15;
1184                 input_samples[N-j-1] = MUL16(input_samples[N-j-1],
1185                                              ff_ac3_window[j]) >> 15;
1186             }
1187
1188             /* Normalize the samples to use the maximum available
1189                precision */
1190             v = 14 - log2_tab(input_samples, N);
1191             if (v < 0)
1192                 v = 0;
1193             exp_samples[i][ch] = v - 9;
1194             lshift_tab(input_samples, N, v);
1195
1196             /* do the MDCT */
1197             mdct512(mdct_coef[i][ch], input_samples);
1198
1199             /* compute "exponents". We take into account the
1200                normalization there */
1201             for(j=0;j<N/2;j++) {
1202                 int e;
1203                 v = abs(mdct_coef[i][ch][j]);
1204                 if (v == 0)
1205                     e = 24;
1206                 else {
1207                     e = 23 - av_log2(v) + exp_samples[i][ch];
1208                     if (e >= 24) {
1209                         e = 24;
1210                         mdct_coef[i][ch][j] = 0;
1211                     }
1212                 }
1213                 exp[i][ch][j] = e;
1214             }
1215         }
1216
1217         compute_exp_strategy(exp_strategy, exp, ch, ch == s->lfe_channel);
1218
1219         /* compute the exponents as the decoder will see them. The
1220            EXP_REUSE case must be handled carefully : we select the
1221            min of the exponents */
1222         i = 0;
1223         while (i < NB_BLOCKS) {
1224             j = i + 1;
1225             while (j < NB_BLOCKS && exp_strategy[j][ch] == EXP_REUSE) {
1226                 exponent_min(exp[i][ch], exp[j][ch], s->nb_coefs[ch]);
1227                 j++;
1228             }
1229             frame_bits += encode_exp(encoded_exp[i][ch],
1230                                      exp[i][ch], s->nb_coefs[ch],
1231                                      exp_strategy[i][ch]);
1232             /* copy encoded exponents for reuse case */
1233             for(k=i+1;k<j;k++) {
1234                 memcpy(encoded_exp[k][ch], encoded_exp[i][ch],
1235                        s->nb_coefs[ch] * sizeof(uint8_t));
1236             }
1237             i = j;
1238         }
1239     }
1240
1241     /* adjust for fractional frame sizes */
1242     while(s->bits_written >= s->bit_rate && s->samples_written >= s->sample_rate) {
1243         s->bits_written -= s->bit_rate;
1244         s->samples_written -= s->sample_rate;
1245     }
1246     s->frame_size = s->frame_size_min + (s->bits_written * s->sample_rate < s->samples_written * s->bit_rate);
1247     s->bits_written += s->frame_size * 16;
1248     s->samples_written += AC3_FRAME_SIZE;
1249
1250     compute_bit_allocation(s, bap, encoded_exp, exp_strategy, frame_bits);
1251     /* everything is known... let's output the frame */
1252     output_frame_header(s, frame);
1253
1254     for(i=0;i<NB_BLOCKS;i++) {
1255         output_audio_block(s, exp_strategy[i], encoded_exp[i],
1256                            bap[i], mdct_coef[i], exp_samples[i], i);
1257     }
1258     return output_frame_end(s);
1259 }
1260
1261 static int AC3_encode_close(AVCodecContext *avctx)
1262 {
1263     av_freep(&avctx->coded_frame);
1264     return 0;
1265 }
1266
1267 #if 0
1268 /*************************************************************************/
1269 /* TEST */
1270
1271 #undef random
1272 #define FN (N/4)
1273
1274 void fft_test(void)
1275 {
1276     IComplex in[FN], in1[FN];
1277     int k, n, i;
1278     float sum_re, sum_im, a;
1279
1280     /* FFT test */
1281
1282     for(i=0;i<FN;i++) {
1283         in[i].re = random() % 65535 - 32767;
1284         in[i].im = random() % 65535 - 32767;
1285         in1[i] = in[i];
1286     }
1287     fft(in, 7);
1288
1289     /* do it by hand */
1290     for(k=0;k<FN;k++) {
1291         sum_re = 0;
1292         sum_im = 0;
1293         for(n=0;n<FN;n++) {
1294             a = -2 * M_PI * (n * k) / FN;
1295             sum_re += in1[n].re * cos(a) - in1[n].im * sin(a);
1296             sum_im += in1[n].re * sin(a) + in1[n].im * cos(a);
1297         }
1298         printf("%3d: %6d,%6d %6.0f,%6.0f\n",
1299                k, in[k].re, in[k].im, sum_re / FN, sum_im / FN);
1300     }
1301 }
1302
1303 void mdct_test(void)
1304 {
1305     int16_t input[N];
1306     int32_t output[N/2];
1307     float input1[N];
1308     float output1[N/2];
1309     float s, a, err, e, emax;
1310     int i, k, n;
1311
1312     for(i=0;i<N;i++) {
1313         input[i] = (random() % 65535 - 32767) * 9 / 10;
1314         input1[i] = input[i];
1315     }
1316
1317     mdct512(output, input);
1318
1319     /* do it by hand */
1320     for(k=0;k<N/2;k++) {
1321         s = 0;
1322         for(n=0;n<N;n++) {
1323             a = (2*M_PI*(2*n+1+N/2)*(2*k+1) / (4 * N));
1324             s += input1[n] * cos(a);
1325         }
1326         output1[k] = -2 * s / N;
1327     }
1328
1329     err = 0;
1330     emax = 0;
1331     for(i=0;i<N/2;i++) {
1332         printf("%3d: %7d %7.0f\n", i, output[i], output1[i]);
1333         e = output[i] - output1[i];
1334         if (e > emax)
1335             emax = e;
1336         err += e * e;
1337     }
1338     printf("err2=%f emax=%f\n", err / (N/2), emax);
1339 }
1340
1341 void test_ac3(void)
1342 {
1343     AC3EncodeContext ctx;
1344     unsigned char frame[AC3_MAX_CODED_FRAME_SIZE];
1345     short samples[AC3_FRAME_SIZE];
1346     int ret, i;
1347
1348     AC3_encode_init(&ctx, 44100, 64000, 1);
1349
1350     fft_test();
1351     mdct_test();
1352
1353     for(i=0;i<AC3_FRAME_SIZE;i++)
1354         samples[i] = (int)(sin(2*M_PI*i*1000.0/44100) * 10000);
1355     ret = AC3_encode_frame(&ctx, frame, samples);
1356     printf("ret=%d\n", ret);
1357 }
1358 #endif
1359
1360 AVCodec ac3_encoder = {
1361     "ac3",
1362     CODEC_TYPE_AUDIO,
1363     CODEC_ID_AC3,
1364     sizeof(AC3EncodeContext),
1365     AC3_encode_init,
1366     AC3_encode_frame,
1367     AC3_encode_close,
1368     NULL,
1369 };