Commit unrolled version of WMA's exponent decode taken from latest ffmpeg. Gives...
[kugel-rb.git] / apps / codecs / libwma / wmadeci.c
blob87af30b5180e4d92121d2801566267a6b3270738
1 /*
2 * WMA compatible decoder
3 * Copyright (c) 2002 The FFmpeg Project.
5 * This library is free software; you can redistribute it and/or
6 * modify it under the terms of the GNU Lesser General Public
7 * License as published by the Free Software Foundation; either
8 * version 2 of the License, or (at your option) any later version.
10 * This library is distributed in the hope that it will be useful,
11 * but WITHOUT ANY WARRANTY; without even the implied warranty of
12 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 * Lesser General Public License for more details.
15 * You should have received a copy of the GNU Lesser General Public
16 * License along with this library; if not, write to the Free Software
17 * Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
20 /**
21 * @file wmadec.c
22 * WMA compatible decoder.
25 #include <codecs.h>
26 #include <codecs/lib/codeclib.h>
27 #include "asf.h"
28 #include "wmadec.h"
29 #include "wmafixed.h"
30 #include "wmadata.h"
32 static const uint8_t ff_log2_tab[256]={
33 0,0,1,1,2,2,2,2,3,3,3,3,3,3,3,3,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,4,
34 5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,
35 6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
36 6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
37 7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
38 7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
39 7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
40 7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7
43 static inline int av_log2(unsigned int v)
45 int n;
47 n = 0;
48 if (v & 0xffff0000) {
49 v >>= 16;
50 n += 16;
52 if (v & 0xff00) {
53 v >>= 8;
54 n += 8;
56 n += ff_log2_tab[v];
58 return n;
61 static void wma_lsp_to_curve_init(WMADecodeContext *s, int frame_len);
62 inline void vector_fmul_add_add(fixed32 *dst, const fixed32 *data,
63 const fixed32 *window, int n);
64 inline void vector_fmul_reverse(fixed32 *dst, const fixed32 *src0,
65 const fixed32 *src1, int len);
67 /*declarations of statically allocated variables used to remove malloc calls*/
69 fixed32 coefsarray[MAX_CHANNELS][BLOCK_MAX_SIZE] IBSS_ATTR;
70 /*decode and window into IRAM on targets with at least 80KB of codec IRAM*/
71 fixed32 frame_out_buf[MAX_CHANNELS][BLOCK_MAX_SIZE * 2] IBSS_ATTR_WMA_LARGE_IRAM;
73 /*MDCT reconstruction windows*/
74 fixed32 stat0[2048], stat1[1024], stat2[512], stat3[256], stat4[128];
76 /*VLC lookup tables*/
77 uint16_t *runtabarray[2], *levtabarray[2];
79 /*these could be made smaller since only one can be 1336*/
80 uint16_t runtab0[1336], runtab1[1336], levtab0[1336], levtab1[1336];
82 #define VLCBUF1SIZE 4598
83 #define VLCBUF2SIZE 3574
84 #define VLCBUF3SIZE 360
85 #define VLCBUF4SIZE 540
87 /*putting these in IRAM actually makes PP slower*/
89 VLC_TYPE vlcbuf1[VLCBUF1SIZE][2];
90 VLC_TYPE vlcbuf2[VLCBUF2SIZE][2];
91 VLC_TYPE vlcbuf3[VLCBUF3SIZE][2];
92 VLC_TYPE vlcbuf4[VLCBUF4SIZE][2];
97 /**
98 * Apply MDCT window and add into output.
100 * We ensure that when the windows overlap their squared sum
101 * is always 1 (MDCT reconstruction rule).
103 * The Vorbis I spec has a great diagram explaining this process.
104 * See section 1.3.2.3 of http://xiph.org/vorbis/doc/Vorbis_I_spec.html
106 static void wma_window(WMADecodeContext *s, fixed32 *in, fixed32 *out)
108 //float *in = s->output;
109 int block_len, bsize, n;
111 /* left part */
113 /* previous block was larger, so we'll use the size of the current
114 * block to set the window size*/
115 if (s->block_len_bits <= s->prev_block_len_bits) {
116 block_len = s->block_len;
117 bsize = s->frame_len_bits - s->block_len_bits;
119 vector_fmul_add_add(out, in, s->windows[bsize], block_len);
121 } else {
122 /*previous block was smaller or the same size, so use it's size to set the window length*/
123 block_len = 1 << s->prev_block_len_bits;
124 /*find the middle of the two overlapped blocks, this will be the first overlapped sample*/
125 n = (s->block_len - block_len) / 2;
126 bsize = s->frame_len_bits - s->prev_block_len_bits;
128 vector_fmul_add_add(out+n, in+n, s->windows[bsize], block_len);
130 memcpy(out+n+block_len, in+n+block_len, n*sizeof(fixed32));
132 /* Advance to the end of the current block and prepare to window it for the next block.
133 * Since the window function needs to be reversed, we do it backwards starting with the
134 * last sample and moving towards the first
136 out += s->block_len;
137 in += s->block_len;
139 /* right part */
140 if (s->block_len_bits <= s->next_block_len_bits) {
141 block_len = s->block_len;
142 bsize = s->frame_len_bits - s->block_len_bits;
144 vector_fmul_reverse(out, in, s->windows[bsize], block_len);
146 } else {
147 block_len = 1 << s->next_block_len_bits;
148 n = (s->block_len - block_len) / 2;
149 bsize = s->frame_len_bits - s->next_block_len_bits;
151 memcpy(out, in, n*sizeof(fixed32));
153 vector_fmul_reverse(out+n, in+n, s->windows[bsize], block_len);
155 memset(out+n+block_len, 0, n*sizeof(fixed32));
162 /* XXX: use same run/length optimization as mpeg decoders */
163 static void init_coef_vlc(VLC *vlc,
164 uint16_t **prun_table, uint16_t **plevel_table,
165 const CoefVLCTable *vlc_table, int tab)
167 int n = vlc_table->n;
168 const uint8_t *table_bits = vlc_table->huffbits;
169 const uint32_t *table_codes = vlc_table->huffcodes;
170 const uint16_t *levels_table = vlc_table->levels;
171 uint16_t *run_table, *level_table;
172 const uint16_t *p;
173 int i, l, j, level;
176 init_vlc(vlc, VLCBITS, n, table_bits, 1, 1, table_codes, 4, 4, 0);
178 run_table = runtabarray[tab];
179 level_table= levtabarray[tab];
181 p = levels_table;
182 i = 2;
183 level = 1;
184 while (i < n)
186 l = *p++;
187 for(j=0;j<l;++j)
189 run_table[i] = j;
190 level_table[i] = level;
191 ++i;
193 ++level;
195 *prun_table = run_table;
196 *plevel_table = level_table;
199 int wma_decode_init(WMADecodeContext* s, asf_waveformatex_t *wfx)
202 int i, flags1, flags2;
203 fixed32 *window;
204 uint8_t *extradata;
205 fixed64 bps1;
206 fixed32 high_freq;
207 fixed64 bps;
208 int sample_rate1;
209 int coef_vlc_table;
210 // int filehandle;
211 #ifdef CPU_COLDFIRE
212 coldfire_set_macsr(EMAC_FRACTIONAL | EMAC_SATURATE);
213 #endif
215 /*clear stereo setting to avoid glitches when switching stereo->mono*/
216 s->channel_coded[0]=0;
217 s->channel_coded[1]=0;
218 s->ms_stereo=0;
220 s->sample_rate = wfx->rate;
221 s->nb_channels = wfx->channels;
222 s->bit_rate = wfx->bitrate;
223 s->block_align = wfx->blockalign;
225 s->coefs = &coefsarray;
226 s->frame_out = &frame_out_buf;
228 if (wfx->codec_id == ASF_CODEC_ID_WMAV1) {
229 s->version = 1;
230 } else if (wfx->codec_id == ASF_CODEC_ID_WMAV2 ) {
231 s->version = 2;
232 } else {
233 /*one of those other wma flavors that don't have GPLed decoders */
234 return -1;
237 /* extract flag infos */
238 flags1 = 0;
239 flags2 = 0;
240 extradata = wfx->data;
241 if (s->version == 1 && wfx->datalen >= 4) {
242 flags1 = extradata[0] | (extradata[1] << 8);
243 flags2 = extradata[2] | (extradata[3] << 8);
244 }else if (s->version == 2 && wfx->datalen >= 6){
245 flags1 = extradata[0] | (extradata[1] << 8) |
246 (extradata[2] << 16) | (extradata[3] << 24);
247 flags2 = extradata[4] | (extradata[5] << 8);
249 s->use_exp_vlc = flags2 & 0x0001;
250 s->use_bit_reservoir = flags2 & 0x0002;
251 s->use_variable_block_len = flags2 & 0x0004;
253 /* compute MDCT block size */
254 if (s->sample_rate <= 16000){
255 s->frame_len_bits = 9;
256 }else if (s->sample_rate <= 22050 ||
257 (s->sample_rate <= 32000 && s->version == 1)){
258 s->frame_len_bits = 10;
259 }else{
260 s->frame_len_bits = 11;
262 s->frame_len = 1 << s->frame_len_bits;
263 if (s-> use_variable_block_len)
265 int nb_max, nb;
266 nb = ((flags2 >> 3) & 3) + 1;
267 if ((s->bit_rate / s->nb_channels) >= 32000)
269 nb += 2;
271 nb_max = s->frame_len_bits - BLOCK_MIN_BITS; //max is 11-7
272 if (nb > nb_max)
273 nb = nb_max;
274 s->nb_block_sizes = nb + 1;
276 else
278 s->nb_block_sizes = 1;
281 /* init rate dependant parameters */
282 s->use_noise_coding = 1;
283 high_freq = itofix64(s->sample_rate) >> 1;
286 /* if version 2, then the rates are normalized */
287 sample_rate1 = s->sample_rate;
288 if (s->version == 2)
290 if (sample_rate1 >= 44100)
291 sample_rate1 = 44100;
292 else if (sample_rate1 >= 22050)
293 sample_rate1 = 22050;
294 else if (sample_rate1 >= 16000)
295 sample_rate1 = 16000;
296 else if (sample_rate1 >= 11025)
297 sample_rate1 = 11025;
298 else if (sample_rate1 >= 8000)
299 sample_rate1 = 8000;
302 fixed64 tmp = itofix64(s->bit_rate);
303 fixed64 tmp2 = itofix64(s->nb_channels * s->sample_rate);
304 bps = fixdiv64(tmp, tmp2);
305 fixed64 tim = bps * s->frame_len;
306 fixed64 tmpi = fixdiv64(tim,itofix64(8));
307 s->byte_offset_bits = av_log2(fixtoi64(tmpi+0x8000)) + 2;
309 /* compute high frequency value and choose if noise coding should
310 be activated */
311 bps1 = bps;
312 if (s->nb_channels == 2)
313 bps1 = fixmul32(bps,0x1999a);
314 if (sample_rate1 == 44100)
316 if (bps1 >= 0x9c29)
317 s->use_noise_coding = 0;
318 else
319 high_freq = fixmul32(high_freq,0x6666);
321 else if (sample_rate1 == 22050)
323 if (bps1 >= 0x128f6)
324 s->use_noise_coding = 0;
325 else if (bps1 >= 0xb852)
326 high_freq = fixmul32(high_freq,0xb333);
327 else
328 high_freq = fixmul32(high_freq,0x999a);
330 else if (sample_rate1 == 16000)
332 if (bps > 0x8000)
333 high_freq = fixmul32(high_freq,0x8000);
334 else
335 high_freq = fixmul32(high_freq,0x4ccd);
337 else if (sample_rate1 == 11025)
339 high_freq = fixmul32(high_freq,0xb333);
341 else if (sample_rate1 == 8000)
343 if (bps <= 0xa000)
345 high_freq = fixmul32(high_freq,0x8000);
347 else if (bps > 0xc000)
349 s->use_noise_coding = 0;
351 else
353 high_freq = fixmul32(high_freq,0xa666);
356 else
358 if (bps >= 0xcccd)
360 high_freq = fixmul32(high_freq,0xc000);
362 else if (bps >= 0x999a)
364 high_freq = fixmul32(high_freq,0x999a);
366 else
368 high_freq = fixmul32(high_freq,0x8000);
372 /* compute the scale factor band sizes for each MDCT block size */
374 int a, b, pos, lpos, k, block_len, i, j, n;
375 const uint8_t *table;
377 if (s->version == 1)
379 s->coefs_start = 3;
381 else
383 s->coefs_start = 0;
385 for(k = 0; k < s->nb_block_sizes; ++k)
387 block_len = s->frame_len >> k;
389 if (s->version == 1)
391 lpos = 0;
392 for(i=0;i<25;++i)
394 a = wma_critical_freqs[i];
395 b = s->sample_rate;
396 pos = ((block_len * 2 * a) + (b >> 1)) / b;
397 if (pos > block_len)
398 pos = block_len;
399 s->exponent_bands[0][i] = pos - lpos;
400 if (pos >= block_len)
402 ++i;
403 break;
405 lpos = pos;
407 s->exponent_sizes[0] = i;
409 else
411 /* hardcoded tables */
412 table = NULL;
413 a = s->frame_len_bits - BLOCK_MIN_BITS - k;
414 if (a < 3)
416 if (s->sample_rate >= 44100)
417 table = exponent_band_44100[a];
418 else if (s->sample_rate >= 32000)
419 table = exponent_band_32000[a];
420 else if (s->sample_rate >= 22050)
421 table = exponent_band_22050[a];
423 if (table)
425 n = *table++;
426 for(i=0;i<n;++i)
427 s->exponent_bands[k][i] = table[i];
428 s->exponent_sizes[k] = n;
430 else
432 j = 0;
433 lpos = 0;
434 for(i=0;i<25;++i)
436 a = wma_critical_freqs[i];
437 b = s->sample_rate;
438 pos = ((block_len * 2 * a) + (b << 1)) / (4 * b);
439 pos <<= 2;
440 if (pos > block_len)
441 pos = block_len;
442 if (pos > lpos)
443 s->exponent_bands[k][j++] = pos - lpos;
444 if (pos >= block_len)
445 break;
446 lpos = pos;
448 s->exponent_sizes[k] = j;
452 /* max number of coefs */
453 s->coefs_end[k] = (s->frame_len - ((s->frame_len * 9) / 100)) >> k;
454 /* high freq computation */
456 fixed32 tmp1 = high_freq*2; /* high_freq is a fixed32!*/
457 fixed32 tmp2=itofix32(s->sample_rate>>1);
458 s->high_band_start[k] = fixtoi32( fixdiv32(tmp1, tmp2) * (block_len>>1) +0x8000);
461 s->high_band_start[k] = (int)((block_len * 2 * high_freq) /
462 s->sample_rate + 0.5);*/
464 n = s->exponent_sizes[k];
465 j = 0;
466 pos = 0;
467 for(i=0;i<n;++i)
469 int start, end;
470 start = pos;
471 pos += s->exponent_bands[k][i];
472 end = pos;
473 if (start < s->high_band_start[k])
474 start = s->high_band_start[k];
475 if (end > s->coefs_end[k])
476 end = s->coefs_end[k];
477 if (end > start)
478 s->exponent_high_bands[k][j++] = end - start;
480 s->exponent_high_sizes[k] = j;
484 /*Not using the ffmpeg IMDCT anymore*/
486 /* mdct_init_global();
488 for(i = 0; i < s->nb_block_sizes; ++i)
490 ff_mdct_init(&s->mdct_ctx[i], s->frame_len_bits - i + 1, 1);
495 /* ffmpeg uses malloc to only allocate as many window sizes as needed.
496 * However, we're really only interested in the worst case memory usage.
497 * In the worst case you can have 5 window sizes, 128 doubling up 2048
498 * Smaller windows are handled differently.
499 * Since we don't have malloc, just statically allocate this
501 fixed32 *temp[5];
502 temp[0] = stat0;
503 temp[1] = stat1;
504 temp[2] = stat2;
505 temp[3] = stat3;
506 temp[4] = stat4;
508 /* init MDCT windows : simple sinus window */
509 for(i = 0; i < s->nb_block_sizes; i++)
511 int n, j;
512 fixed32 alpha;
513 n = 1 << (s->frame_len_bits - i);
514 window = temp[i];
516 /* this calculates 0.5/(2*n) */
517 alpha = (1<<15)>>(s->frame_len_bits - i+1);
518 for(j=0;j<n;++j)
520 fixed32 j2 = itofix32(j) + 0x8000;
521 /*alpha between 0 and pi/2*/
522 window[j] = fsincos(fixmul32(j2,alpha)<<16, 0);
524 s->windows[i] = window;
528 s->reset_block_lengths = 1;
530 if (s->use_noise_coding)
532 /* init the noise generator */
533 if (s->use_exp_vlc)
535 s->noise_mult = 0x51f;
536 s->noise_table = noisetable_exp;
538 else
540 s->noise_mult = 0xa3d;
541 /* LSP values are simply 2x the EXP values */
542 for (i=0;i<NOISE_TAB_SIZE;++i)
543 noisetable_exp[i] = noisetable_exp[i]<< 1;
544 s->noise_table = noisetable_exp;
546 #if 0
547 /* We use a lookup table computered in advance, so no need to do this*/
549 unsigned int seed;
550 fixed32 norm;
551 seed = 1;
552 norm = 0; // PJJ: near as makes any diff to 0!
553 for (i=0;i<NOISE_TAB_SIZE;++i)
555 seed = seed * 314159 + 1;
556 s->noise_table[i] = itofix32((int)seed) * norm;
559 #endif
561 s->hgain_vlc.table = vlcbuf4;
562 s->hgain_vlc.table_allocated = VLCBUF4SIZE;
563 init_vlc(&s->hgain_vlc, HGAINVLCBITS, sizeof(hgain_huffbits),
564 hgain_huffbits, 1, 1,
565 hgain_huffcodes, 2, 2, 0);
568 if (s->use_exp_vlc)
571 s->exp_vlc.table = vlcbuf3;
572 s->exp_vlc.table_allocated = VLCBUF3SIZE;
574 init_vlc(&s->exp_vlc, EXPVLCBITS, sizeof(scale_huffbits),
575 scale_huffbits, 1, 1,
576 scale_huffcodes, 4, 4, 0);
578 else
580 wma_lsp_to_curve_init(s, s->frame_len);
583 /* choose the VLC tables for the coefficients */
584 coef_vlc_table = 2;
585 if (s->sample_rate >= 32000)
587 if (bps1 < 0xb852)
588 coef_vlc_table = 0;
589 else if (bps1 < 0x128f6)
590 coef_vlc_table = 1;
593 runtabarray[0] = runtab0; runtabarray[1] = runtab1;
594 levtabarray[0] = levtab0; levtabarray[1] = levtab1;
596 s->coef_vlc[0].table = vlcbuf1;
597 s->coef_vlc[0].table_allocated = VLCBUF1SIZE;
598 s->coef_vlc[1].table = vlcbuf2;
599 s->coef_vlc[1].table_allocated = VLCBUF2SIZE;
602 init_coef_vlc(&s->coef_vlc[0], &s->run_table[0], &s->level_table[0],
603 &coef_vlcs[coef_vlc_table * 2], 0);
604 init_coef_vlc(&s->coef_vlc[1], &s->run_table[1], &s->level_table[1],
605 &coef_vlcs[coef_vlc_table * 2 + 1], 1);
607 s->last_superframe_len = 0;
608 s->last_bitoffset = 0;
610 return 0;
614 /* compute x^-0.25 with an exponent and mantissa table. We use linear
615 interpolation to reduce the mantissa table size at a small speed
616 expense (linear interpolation approximately doubles the number of
617 bits of precision). */
618 static inline fixed32 pow_m1_4(WMADecodeContext *s, fixed32 x)
620 union {
621 float f;
622 unsigned int v;
623 } u, t;
624 unsigned int e, m;
625 fixed32 a, b;
627 u.f = fixtof64(x);
628 e = u.v >> 23;
629 m = (u.v >> (23 - LSP_POW_BITS)) & ((1 << LSP_POW_BITS) - 1);
630 /* build interpolation scale: 1 <= t < 2. */
631 t.v = ((u.v << LSP_POW_BITS) & ((1 << 23) - 1)) | (127 << 23);
632 a = s->lsp_pow_m_table1[m];
633 b = s->lsp_pow_m_table2[m];
635 /* lsp_pow_e_table contains 32.32 format */
636 /* TODO: Since we're unlikely have value that cover the whole
637 * IEEE754 range, we probably don't need to have all possible exponents */
639 return (lsp_pow_e_table[e] * (a + fixmul32(b, ftofix32(t.f))) >>32);
642 static void wma_lsp_to_curve_init(WMADecodeContext *s, int frame_len)
644 fixed32 wdel, a, b, temp, temp2;
645 int i, m;
647 wdel = fixdiv32(M_PI_F, itofix32(frame_len));
648 temp = fixdiv32(itofix32(1), itofix32(frame_len));
649 for (i=0; i<frame_len; ++i)
651 /* TODO: can probably reuse the trig_init values here */
652 fsincos((temp*i)<<15, &temp2);
653 /* get 3 bits headroom + 1 bit from not doubleing the values */
654 s->lsp_cos_table[i] = temp2>>3;
657 /* NOTE: these two tables are needed to avoid two operations in
658 pow_m1_4 */
659 b = itofix32(1);
660 int ix = 0;
662 /*double check this later*/
663 for(i=(1 << LSP_POW_BITS) - 1;i>=0;i--)
665 m = (1 << LSP_POW_BITS) + i;
666 a = pow_a_table[ix++]<<4;
667 s->lsp_pow_m_table1[i] = 2 * a - b;
668 s->lsp_pow_m_table2[i] = b - a;
669 b = a;
674 /* NOTE: We use the same code as Vorbis here */
675 /* XXX: optimize it further with SSE/3Dnow */
676 static void wma_lsp_to_curve(WMADecodeContext *s,
677 fixed32 *out,
678 fixed32 *val_max_ptr,
679 int n,
680 fixed32 *lsp)
682 int i, j;
683 fixed32 p, q, w, v, val_max, temp, temp2;
685 val_max = 0;
686 for(i=0;i<n;++i)
688 /* shift by 2 now to reduce rounding error,
689 * we can renormalize right before pow_m1_4
692 p = 0x8000<<5;
693 q = 0x8000<<5;
694 w = s->lsp_cos_table[i];
696 for (j=1;j<NB_LSP_COEFS;j+=2)
698 /* w is 5.27 format, lsp is in 16.16, temp2 becomes 5.27 format */
699 temp2 = ((w - (lsp[j - 1]<<11)));
700 temp = q;
702 /* q is 16.16 format, temp2 is 5.27, q becomes 16.16 */
703 q = fixmul32b(q, temp2 )<<4;
704 p = fixmul32b(p, (w - (lsp[j]<<11)))<<4;
707 /* 2 in 5.27 format is 0x10000000 */
708 p = fixmul32(p, fixmul32b(p, (0x10000000 - w)))<<3;
709 q = fixmul32(q, fixmul32b(q, (0x10000000 + w)))<<3;
711 v = (p + q) >>9; /* p/q end up as 16.16 */
712 v = pow_m1_4(s, v);
713 if (v > val_max)
714 val_max = v;
715 out[i] = v;
718 *val_max_ptr = val_max;
721 /* decode exponents coded with LSP coefficients (same idea as Vorbis)
722 * only used for low bitrate (< 16kbps) files
724 static void decode_exp_lsp(WMADecodeContext *s, int ch)
726 fixed32 lsp_coefs[NB_LSP_COEFS];
727 int val, i;
729 for (i = 0; i < NB_LSP_COEFS; ++i)
731 if (i == 0 || i >= 8)
732 val = get_bits(&s->gb, 3);
733 else
734 val = get_bits(&s->gb, 4);
735 lsp_coefs[i] = lsp_codebook[i][val];
738 wma_lsp_to_curve(s,
739 s->exponents[ch],
740 &s->max_exponent[ch],
741 s->block_len,
742 lsp_coefs);
745 /* decode exponents coded with VLC codes - used for bitrate >= 32kbps*/
746 static int decode_exp_vlc(WMADecodeContext *s, int ch)
748 int last_exp, n, code;
749 const uint16_t *ptr, *band_ptr;
750 fixed32 v, max_scale;
751 fixed32 *q,*q_end;
753 /*accommodate the 60 negative indices */
754 const fixed32 *pow_10_to_yover16_ptr = &pow_10_to_yover16[61];
756 band_ptr = s->exponent_bands[s->frame_len_bits - s->block_len_bits];
757 ptr = band_ptr;
758 q = s->exponents[ch];
759 q_end = q + s->block_len;
760 max_scale = 0;
763 if (s->version == 1) //wmav1 only
765 last_exp = get_bits(&s->gb, 5) + 10;
767 v = pow_10_to_yover16_ptr[last_exp];
768 max_scale = v;
769 n = *ptr++;
770 switch (n & 3) do {
771 case 0: *q++ = v;
772 case 3: *q++ = v;
773 case 2: *q++ = v;
774 case 1: *q++ = v;
775 } while ((n -= 4) > 0);
776 } else {
777 last_exp = 36;
780 while (q < q_end)
782 code = get_vlc2(&s->gb, s->exp_vlc.table, EXPVLCBITS, EXPMAX);
783 if (code < 0)
785 return -1;
787 /* NOTE: this offset is the same as MPEG4 AAC ! */
788 last_exp += code - 60;
790 v = pow_10_to_yover16_ptr[last_exp];
791 if (v > max_scale)
793 max_scale = v;
795 n = *ptr++;
796 switch (n & 3) do {
797 case 0: *q++ = v;
798 case 3: *q++ = v;
799 case 2: *q++ = v;
800 case 1: *q++ = v;
801 } while ((n -= 4) > 0);
804 s->max_exponent[ch] = max_scale;
805 return 0;
808 /* return 0 if OK. return 1 if last block of frame. return -1 if
809 unrecorrable error. */
810 static int wma_decode_block(WMADecodeContext *s, int32_t *scratch_buffer)
812 int n, v, a, ch, code, bsize;
813 int coef_nb_bits, total_gain;
814 int nb_coefs[MAX_CHANNELS];
815 fixed32 mdct_norm;
817 /*DEBUGF("***decode_block: %d (%d samples of %d in frame)\n", s->block_num, s->block_len, s->frame_len);*/
819 /* compute current block length */
820 if (s->use_variable_block_len)
822 n = av_log2(s->nb_block_sizes - 1) + 1;
824 if (s->reset_block_lengths)
826 s->reset_block_lengths = 0;
827 v = get_bits(&s->gb, n);
828 if (v >= s->nb_block_sizes)
830 return -2;
832 s->prev_block_len_bits = s->frame_len_bits - v;
833 v = get_bits(&s->gb, n);
834 if (v >= s->nb_block_sizes)
836 return -3;
838 s->block_len_bits = s->frame_len_bits - v;
840 else
842 /* update block lengths */
843 s->prev_block_len_bits = s->block_len_bits;
844 s->block_len_bits = s->next_block_len_bits;
846 v = get_bits(&s->gb, n);
848 if (v >= s->nb_block_sizes)
850 // rb->splash(HZ*4, "v was %d", v); //5, 7
851 return -4; //this is it
853 else{
854 //rb->splash(HZ, "passed v block (%d)!", v);
856 s->next_block_len_bits = s->frame_len_bits - v;
858 else
860 /* fixed block len */
861 s->next_block_len_bits = s->frame_len_bits;
862 s->prev_block_len_bits = s->frame_len_bits;
863 s->block_len_bits = s->frame_len_bits;
865 /* now check if the block length is coherent with the frame length */
866 s->block_len = 1 << s->block_len_bits;
868 if ((s->block_pos + s->block_len) > s->frame_len)
870 return -5; //oddly 32k sample from tracker fails here
873 if (s->nb_channels == 2)
875 s->ms_stereo = get_bits1(&s->gb);
877 v = 0;
878 for (ch = 0; ch < s->nb_channels; ++ch)
880 a = get_bits1(&s->gb);
881 s->channel_coded[ch] = a;
882 v |= a;
884 /* if no channel coded, no need to go further */
885 /* XXX: fix potential framing problems */
886 if (!v)
888 goto next;
891 bsize = s->frame_len_bits - s->block_len_bits;
893 /* read total gain and extract corresponding number of bits for
894 coef escape coding */
895 total_gain = 1;
896 for(;;)
898 a = get_bits(&s->gb, 7);
899 total_gain += a;
900 if (a != 127)
902 break;
906 if (total_gain < 15)
907 coef_nb_bits = 13;
908 else if (total_gain < 32)
909 coef_nb_bits = 12;
910 else if (total_gain < 40)
911 coef_nb_bits = 11;
912 else if (total_gain < 45)
913 coef_nb_bits = 10;
914 else
915 coef_nb_bits = 9;
917 /* compute number of coefficients */
918 n = s->coefs_end[bsize] - s->coefs_start;
920 for(ch = 0; ch < s->nb_channels; ++ch)
922 nb_coefs[ch] = n;
924 /* complex coding */
925 if (s->use_noise_coding)
928 for(ch = 0; ch < s->nb_channels; ++ch)
930 if (s->channel_coded[ch])
932 int i, n, a;
933 n = s->exponent_high_sizes[bsize];
934 for(i=0;i<n;++i)
936 a = get_bits1(&s->gb);
937 s->high_band_coded[ch][i] = a;
938 /* if noise coding, the coefficients are not transmitted */
939 if (a)
940 nb_coefs[ch] -= s->exponent_high_bands[bsize][i];
944 for(ch = 0; ch < s->nb_channels; ++ch)
946 if (s->channel_coded[ch])
948 int i, n, val, code;
950 n = s->exponent_high_sizes[bsize];
951 val = (int)0x80000000;
952 for(i=0;i<n;++i)
954 if (s->high_band_coded[ch][i])
956 if (val == (int)0x80000000)
958 val = get_bits(&s->gb, 7) - 19;
960 else
962 //code = get_vlc(&s->gb, &s->hgain_vlc);
963 code = get_vlc2(&s->gb, s->hgain_vlc.table, HGAINVLCBITS, HGAINMAX);
964 if (code < 0)
966 return -6;
968 val += code - 18;
970 s->high_band_values[ch][i] = val;
977 /* exponents can be reused in short blocks. */
978 if ((s->block_len_bits == s->frame_len_bits) || get_bits1(&s->gb))
980 for(ch = 0; ch < s->nb_channels; ++ch)
982 if (s->channel_coded[ch])
984 if (s->use_exp_vlc)
986 if (decode_exp_vlc(s, ch) < 0)
988 return -7;
991 else
993 decode_exp_lsp(s, ch);
995 s->exponents_bsize[ch] = bsize;
1000 /* parse spectral coefficients : just RLE encoding */
1001 for(ch = 0; ch < s->nb_channels; ++ch)
1003 if (s->channel_coded[ch])
1005 VLC *coef_vlc;
1006 int level, run, sign, tindex;
1007 int16_t *ptr, *eptr;
1008 const int16_t *level_table, *run_table;
1010 /* special VLC tables are used for ms stereo because
1011 there is potentially less energy there */
1012 tindex = (ch == 1 && s->ms_stereo);
1013 coef_vlc = &s->coef_vlc[tindex];
1014 run_table = s->run_table[tindex];
1015 level_table = s->level_table[tindex];
1016 /* XXX: optimize */
1017 ptr = &s->coefs1[ch][0];
1018 eptr = ptr + nb_coefs[ch];
1019 memset(ptr, 0, s->block_len * sizeof(int16_t));
1021 for(;;)
1023 code = get_vlc2(&s->gb, coef_vlc->table, VLCBITS, VLCMAX);
1025 if (code < 0)
1027 return -8;
1029 if (code == 1)
1031 /* EOB */
1032 break;
1034 else if (code == 0)
1036 /* escape */
1037 level = get_bits(&s->gb, coef_nb_bits);
1038 /* NOTE: this is rather suboptimal. reading
1039 block_len_bits would be better */
1040 run = get_bits(&s->gb, s->frame_len_bits);
1042 else
1044 /* normal code */
1045 run = run_table[code];
1046 level = level_table[code];
1048 sign = get_bits1(&s->gb);
1049 if (!sign)
1050 level = -level;
1051 ptr += run;
1052 if (ptr >= eptr)
1054 break;
1056 *ptr++ = level;
1059 /* NOTE: EOB can be omitted */
1060 if (ptr >= eptr)
1061 break;
1064 if (s->version == 1 && s->nb_channels >= 2)
1066 align_get_bits(&s->gb);
1071 int n4 = s->block_len >> 1;
1074 mdct_norm = 0x10000>>(s->block_len_bits-1);
1076 if (s->version == 1)
1078 mdct_norm *= fixtoi32(fixsqrt32(itofix32(n4)));
1083 /* finally compute the MDCT coefficients */
1084 for(ch = 0; ch < s->nb_channels; ++ch)
1086 if (s->channel_coded[ch])
1088 int16_t *coefs1;
1089 fixed32 *exponents;
1090 fixed32 *coefs, atemp;
1091 fixed64 mult;
1092 fixed64 mult1;
1093 fixed32 noise, temp1, temp2, mult2;
1094 int i, j, n, n1, last_high_band, esize;
1095 fixed32 exp_power[HIGH_BAND_MAX_SIZE];
1097 //total_gain, coefs1, mdctnorm are lossless
1099 coefs1 = s->coefs1[ch];
1100 exponents = s->exponents[ch];
1101 esize = s->exponents_bsize[ch];
1102 coefs = (*(s->coefs))[ch];
1103 n=0;
1106 * The calculation of coefs has a shift right by 2 built in. This
1107 * prepares samples for the Tremor IMDCT which uses a slightly
1108 * different fixed format then the ffmpeg one. If the old ffmpeg
1109 * imdct is used, each shift storing into coefs should be reduced
1110 * by 1.
1111 * See SVN logs for details.
1115 if (s->use_noise_coding)
1117 /*This case is only used for low bitrates (typically less then 32kbps)*/
1119 /*TODO: mult should be converted to 32 bit to speed up noise coding*/
1121 mult = fixdiv64(pow_table[total_gain+20],Fixed32To64(s->max_exponent[ch]));
1122 mult = mult* mdct_norm;
1123 mult1 = mult;
1125 /* very low freqs : noise */
1126 for(i = 0;i < s->coefs_start; ++i)
1128 *coefs++ = fixmul32( (fixmul32(s->noise_table[s->noise_index],
1129 exponents[i<<bsize>>esize])>>4),Fixed32From64(mult1)) >>2;
1130 s->noise_index = (s->noise_index + 1) & (NOISE_TAB_SIZE - 1);
1133 n1 = s->exponent_high_sizes[bsize];
1135 /* compute power of high bands */
1136 exponents = s->exponents[ch] +(s->high_band_start[bsize]<<bsize);
1137 last_high_band = 0; /* avoid warning */
1138 for (j=0;j<n1;++j)
1140 n = s->exponent_high_bands[s->frame_len_bits -
1141 s->block_len_bits][j];
1142 if (s->high_band_coded[ch][j])
1144 fixed32 e2, v;
1145 e2 = 0;
1146 for(i = 0;i < n; ++i)
1148 /*v is normalized later on so its fixed format is irrelevant*/
1149 v = exponents[i<<bsize>>esize]>>4;
1150 e2 += fixmul32(v, v)>>3;
1152 exp_power[j] = e2/n; /*n is an int...*/
1153 last_high_band = j;
1155 exponents += n<<bsize;
1158 /* main freqs and high freqs */
1159 exponents = s->exponents[ch] + (s->coefs_start<<bsize);
1160 for(j=-1;j<n1;++j)
1162 if (j < 0)
1164 n = s->high_band_start[bsize] -
1165 s->coefs_start;
1167 else
1169 n = s->exponent_high_bands[s->frame_len_bits -
1170 s->block_len_bits][j];
1172 if (j >= 0 && s->high_band_coded[ch][j])
1174 /* use noise with specified power */
1175 fixed32 tmp = fixdiv32(exp_power[j],exp_power[last_high_band]);
1177 /*mult1 is 48.16, pow_table is 48.16*/
1178 mult1 = fixmul32(fixsqrt32(tmp),
1179 pow_table[s->high_band_values[ch][j]+20]) >> 16;
1181 /*this step has a fairly high degree of error for some reason*/
1182 mult1 = fixdiv64(mult1,fixmul32(s->max_exponent[ch],s->noise_mult));
1183 mult1 = mult1*mdct_norm>>PRECISION;
1184 for(i = 0;i < n; ++i)
1186 noise = s->noise_table[s->noise_index];
1187 s->noise_index = (s->noise_index + 1) & (NOISE_TAB_SIZE - 1);
1188 *coefs++ = fixmul32((fixmul32(exponents[i<<bsize>>esize],noise)>>4),
1189 Fixed32From64(mult1)) >>2;
1192 exponents += n<<bsize;
1194 else
1196 /* coded values + small noise */
1197 for(i = 0;i < n; ++i)
1199 noise = s->noise_table[s->noise_index];
1200 s->noise_index = (s->noise_index + 1) & (NOISE_TAB_SIZE - 1);
1202 /*don't forget to renormalize the noise*/
1203 temp1 = (((int32_t)*coefs1++)<<16) + (noise>>4);
1204 temp2 = fixmul32(exponents[i<<bsize>>esize], mult>>18);
1205 *coefs++ = fixmul32(temp1, temp2);
1207 exponents += n<<bsize;
1211 /* very high freqs : noise */
1212 n = s->block_len - s->coefs_end[bsize];
1213 mult2 = fixmul32(mult>>16,exponents[((-1<<bsize))>>esize]) ;
1214 for (i = 0; i < n; ++i)
1216 /*renormalize the noise product and then reduce to 14.18 precison*/
1217 *coefs++ = fixmul32(s->noise_table[s->noise_index],mult2) >>6;
1219 s->noise_index = (s->noise_index + 1) & (NOISE_TAB_SIZE - 1);
1222 else
1224 /*Noise coding not used, simply convert from exp to fixed representation*/
1226 fixed32 mult3 = (fixed32)(fixdiv64(pow_table[total_gain+20],
1227 Fixed32To64(s->max_exponent[ch])));
1228 mult3 = fixmul32(mult3, mdct_norm);
1230 /*zero the first 3 coefficients for WMA V1, does nothing otherwise*/
1231 for(i=0; i<s->coefs_start; i++)
1232 *coefs++=0;
1234 n = nb_coefs[ch];
1236 /* XXX: optimize more, unrolling this loop in asm
1237 might be a good idea */
1239 for(i = 0;i < n; ++i)
1241 /*ffmpeg imdct needs 15.17, while tremor 14.18*/
1242 atemp = (coefs1[i] * mult3)>>2;
1243 *coefs++=fixmul32(atemp,exponents[i<<bsize>>esize]);
1245 n = s->block_len - s->coefs_end[bsize];
1246 memset(coefs, 0, n*sizeof(fixed32));
1253 if (s->ms_stereo && s->channel_coded[1])
1255 fixed32 a, b;
1256 int i;
1257 fixed32 (*coefs)[MAX_CHANNELS][BLOCK_MAX_SIZE] = (s->coefs);
1259 /* nominal case for ms stereo: we do it before mdct */
1260 /* no need to optimize this case because it should almost
1261 never happen */
1262 if (!s->channel_coded[0])
1264 memset((*(s->coefs))[0], 0, sizeof(fixed32) * s->block_len);
1265 s->channel_coded[0] = 1;
1268 for(i = 0; i < s->block_len; ++i)
1270 a = (*coefs)[0][i];
1271 b = (*coefs)[1][i];
1272 (*coefs)[0][i] = a + b;
1273 (*coefs)[1][i] = a - b;
1277 for(ch = 0; ch < s->nb_channels; ++ch)
1279 if (s->channel_coded[ch])
1281 int n4, index;
1283 n4 = s->block_len >>1;
1285 /*faster IMDCT from Vorbis*/
1286 mdct_backward( (1 << (s->block_len_bits+1)), (int32_t*)(*(s->coefs))[ch], (int32_t*)scratch_buffer);
1288 /*slower but more easily understood IMDCT from FFMPEG*/
1289 //ff_imdct_calc(&s->mdct_ctx[bsize],
1290 // output,
1291 // (*(s->coefs))[ch]);
1294 /* add in the frame */
1295 index = (s->frame_len / 2) + s->block_pos - n4;
1296 wma_window(s, scratch_buffer, &((*s->frame_out)[ch][index]));
1300 /* specific fast case for ms-stereo : add to second
1301 channel if it is not coded */
1302 if (s->ms_stereo && !s->channel_coded[1])
1304 wma_window(s, scratch_buffer, &((*s->frame_out)[1][index]));
1308 next:
1309 /* update block number */
1310 ++s->block_num;
1311 s->block_pos += s->block_len;
1312 if (s->block_pos >= s->frame_len)
1314 return 1;
1316 else
1318 return 0;
1322 /* decode a frame of frame_len samples */
1323 static int wma_decode_frame(WMADecodeContext *s, int32_t *samples)
1325 int ret, i, n, ch, incr;
1326 int32_t *ptr;
1327 fixed32 *iptr;
1329 /* read each block */
1330 s->block_num = 0;
1331 s->block_pos = 0;
1334 for(;;)
1336 ret = wma_decode_block(s, samples);
1337 if (ret < 0)
1340 DEBUGF("wma_decode_block failed with code %d\n", ret);
1341 return -1;
1343 if (ret)
1345 break;
1349 /* return frame with full 30-bit precision */
1350 n = s->frame_len;
1351 incr = s->nb_channels;
1352 for(ch = 0; ch < s->nb_channels; ++ch)
1354 ptr = samples + ch;
1355 iptr = &((*s->frame_out)[ch][0]);
1357 for (i=0;i<n;++i)
1359 *ptr = (*iptr++);
1360 ptr += incr;
1363 memmove(&((*s->frame_out)[ch][0]), &((*s->frame_out)[ch][s->frame_len]),
1364 s->frame_len * sizeof(fixed32));
1367 return 0;
1370 /* Initialise the superframe decoding */
1372 int wma_decode_superframe_init(WMADecodeContext* s,
1373 const uint8_t *buf, /*input*/
1374 int buf_size)
1376 if (buf_size==0)
1378 s->last_superframe_len = 0;
1379 return 0;
1382 s->current_frame = 0;
1384 init_get_bits(&s->gb, buf, buf_size*8);
1386 if (s->use_bit_reservoir)
1388 /* read super frame header */
1389 skip_bits(&s->gb, 4); /* super frame index */
1390 s->nb_frames = get_bits(&s->gb, 4);
1392 if (s->last_superframe_len == 0)
1393 s->nb_frames --;
1394 else if (s->nb_frames == 0)
1395 s->nb_frames++;
1397 s->bit_offset = get_bits(&s->gb, s->byte_offset_bits + 3);
1398 } else {
1399 s->nb_frames = 1;
1402 return 1;
1406 /* Decode a single frame in the current superframe - return -1 if
1407 there was a decoding error, or the number of samples decoded.
1410 int wma_decode_superframe_frame(WMADecodeContext* s,
1411 int32_t* samples, /*output*/
1412 const uint8_t *buf, /*input*/
1413 int buf_size)
1415 int pos, len;
1416 uint8_t *q;
1417 int done = 0;
1418 if ((s->use_bit_reservoir) && (s->current_frame == 0))
1420 if (s->last_superframe_len > 0)
1422 /* add s->bit_offset bits to last frame */
1423 if ((s->last_superframe_len + ((s->bit_offset + 7) >> 3)) >
1424 MAX_CODED_SUPERFRAME_SIZE)
1426 DEBUGF("superframe size too large error\n");
1427 goto fail;
1429 q = s->last_superframe + s->last_superframe_len;
1430 len = s->bit_offset;
1431 while (len > 7)
1433 *q++ = (get_bits)(&s->gb, 8);
1434 len -= 8;
1436 if (len > 0)
1438 *q++ = (get_bits)(&s->gb, len) << (8 - len);
1441 /* XXX: s->bit_offset bits into last frame */
1442 init_get_bits(&s->gb, s->last_superframe, MAX_CODED_SUPERFRAME_SIZE*8);
1443 /* skip unused bits */
1444 if (s->last_bitoffset > 0)
1445 skip_bits(&s->gb, s->last_bitoffset);
1447 /* this frame is stored in the last superframe and in the
1448 current one */
1449 if (wma_decode_frame(s, samples) < 0)
1451 goto fail;
1453 done = 1;
1456 /* read each frame starting from s->bit_offset */
1457 pos = s->bit_offset + 4 + 4 + s->byte_offset_bits + 3;
1458 init_get_bits(&s->gb, buf + (pos >> 3), (MAX_CODED_SUPERFRAME_SIZE - (pos >> 3))*8);
1459 len = pos & 7;
1460 if (len > 0)
1461 skip_bits(&s->gb, len);
1463 s->reset_block_lengths = 1;
1466 /* If we haven't decoded a frame yet, do it now */
1467 if (!done)
1469 if (wma_decode_frame(s, samples) < 0)
1471 goto fail;
1475 s->current_frame++;
1477 if ((s->use_bit_reservoir) && (s->current_frame == s->nb_frames))
1479 /* we copy the end of the frame in the last frame buffer */
1480 pos = get_bits_count(&s->gb) + ((s->bit_offset + 4 + 4 + s->byte_offset_bits + 3) & ~7);
1481 s->last_bitoffset = pos & 7;
1482 pos >>= 3;
1483 len = buf_size - pos;
1484 if (len > MAX_CODED_SUPERFRAME_SIZE || len < 0)
1486 DEBUGF("superframe size too large error after decoding\n");
1487 goto fail;
1489 s->last_superframe_len = len;
1490 memcpy(s->last_superframe, buf + pos, len);
1493 return s->frame_len;
1495 fail:
1496 /* when error, we reset the bit reservoir */
1498 s->last_superframe_len = 0;
1499 return -1;