libavcodec/ra288.c

   1 /*
   2  * RealAudio 2.0 (28.8K)
   3  * Copyright (c) 2003 the ffmpeg project
   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 #include "avcodec.h"
  23 #define ALT_BITSTREAM_READER_LE
  24 #include "bitstream.h"
  25 #include "ra288.h"
  26
  27 typedef struct {
  28     float history[8];
  29     float output[40];
  30     float pr1[36];
  31     float pr2[10];
  32     int   phase, phasep;
  33
  34     float st1a[111], st1b[37], st1[37];
  35     float st2a[38], st2b[11], st2[11];
  36     float sb[41];
  37     float lhist[10];
  38 } Real288_internal;
  39
  40 /* Decode and produce output */
  41 static void decode(Real288_internal *glob, int amp_coef, int cb_coef)
  42 {
  43     unsigned int x, y;
  44     float f;
  45     double sum, sumsum;
  46     float *p1, *p2;
  47     float buffer[5];
  48     const float *table;
  49
  50     for (x=36; x--; glob->sb[x+5] = glob->sb[x]);
  51
  52     for (x=5; x--;) {
  53         p1 = glob->sb+x;
  54         p2 = glob->pr1;
  55         for (sum=0, y=36; y--; sum -= (*(++p1))*(*(p2++)));
  56
  57         glob->sb[x] = sum;
  58     }
  59
  60     f = amptable[amp_coef];
  61     table = codetable[cb_coef];
  62
  63     /* convert log and do rms */
  64     for (sum=32, x=10; x--; sum -= glob->pr2[x] * glob->lhist[x]);
  65
  66     if (sum < 0)
  67         sum = 0;
  68     else if (sum > 60)
  69         sum = 60;
  70
  71     sumsum = exp(sum * 0.1151292546497) * f;    /* pow(10.0,sum/20)*f */
  72
  73     for (sum=0, x=5; x--;) {
  74         buffer[x] = table[x] * sumsum;
  75         sum += buffer[x] * buffer[x];
  76     }
  77
  78     if ((sum /= 5) < 1)
  79         sum = 1;
  80
  81     /* shift and store */
  82     for (x=10; --x; glob->lhist[x] = glob->lhist[x-1]);
  83
  84     *glob->lhist = glob->history[glob->phase] = 10 * log10(sum) - 32;
  85
  86     for (x=1; x < 5; x++)
  87         for (y=x; y--; buffer[x] -= glob->pr1[x-y-1] * buffer[y]);
  88
  89     /* output */
  90     for (x=0; x < 5; x++) {
  91         f = glob->sb[4-x] + buffer[x];
  92
  93         if (f > 4095)
  94             f = 4095;
  95         else if (f < -4095)
  96             f = -4095;
  97
  98         glob->output[glob->phasep+x] = glob->sb[4-x] = f;
  99     }
 100 }
 101
 102 /* column multiply */
 103 static void colmult(float *tgt, float *m1, const float *m2, int n)
 104 {
 105     while (n--)
 106         *(tgt++) = (*(m1++)) * (*(m2++));
 107 }
 108
 109 static int pred(float *in, float *tgt, int n)
 110 {
 111     int x, y;
 112     float *p1, *p2;
 113     double f0, f1, f2;
 114     float temp;
 115
 116     if (in[n] == 0)
 117         return 0;
 118
 119     if ((f0 = *in) <= 0)
 120         return 0;
 121
 122     for (x=1 ; ; x++) {
 123         if (n < x)
 124             return 1;
 125
 126         p1 = in + x;
 127         p2 = tgt;
 128         f1 = *(p1--);
 129         for (y=x; --y; f1 += (*(p1--))*(*(p2++)));
 130
 131         p1 = tgt + x - 1;
 132         p2 = tgt;
 133         *(p1--) = f2 = -f1/f0;
 134         for (y=x >> 1; y--;) {
 135             temp = *p2 + *p1 * f2;
 136             *(p1--) += *p2 * f2;
 137             *(p2++) = temp;
 138         }
 139         if ((f0 += f1*f2) < 0)
 140             return 0;
 141     }
 142 }
 143
 144 /* product sum (lsf) */
 145 static void prodsum(float *tgt, float *src, int len, int n)
 146 {
 147     unsigned int x;
 148     float *p1, *p2;
 149     double sum;
 150
 151     while (n >= 0) {
 152         p1 = (p2 = src) - n;
 153         for (sum=0, x=len; x--; sum += (*p1++) * (*p2++));
 154         tgt[n--] = sum;
 155     }
 156 }
 157
 158 static void co(int n, int i, int j, float *in, float *out, float *st1,
 159                float *st2, const float *table)
 160 {
 161     int a, b, c;
 162     unsigned int x;
 163     float *fp;
 164     float buffer1[37];
 165     float buffer2[37];
 166     float work[111];
 167
 168     /* rotate and multiply */
 169     c = (b = (a = n + i) + j) - i;
 170     fp = st1 + i;
 171     for (x=0; x < b; x++) {
 172         if (x == c)
 173             fp=in;
 174         work[x] = *(table++) * (*(st1++) = *(fp++));
 175     }
 176
 177     prodsum(buffer1, work + n, i, n);
 178     prodsum(buffer2, work + a, j, n);
 179
 180     for (x=0;x<=n;x++) {
 181         *st2 = *st2 * (0.5625) + buffer1[x];
 182         out[x] = *(st2++) + buffer2[x];
 183     }
 184     *out *= 1.00390625; /* to prevent clipping */
 185 }
 186
 187 static void update(Real288_internal *glob)
 188 {
 189     int x,y;
 190     float buffer1[40], temp1[37];
 191     float buffer2[8], temp2[11];
 192
 193     for (x=0, y=glob->phasep+5; x < 40; buffer1[x++] = glob->output[(y++)%40]);
 194
 195     co(36, 40, 35, buffer1, temp1, glob->st1a, glob->st1b, table1);
 196
 197     if (pred(temp1, glob->st1, 36))
 198         colmult(glob->pr1, glob->st1, table1a, 36);
 199
 200     for (x=0, y=glob->phase + 1; x < 8; buffer2[x++] = glob->history[(y++) % 8]);
 201
 202     co(10, 8, 20, buffer2, temp2, glob->st2a, glob->st2b, table2);
 203
 204     if (pred(temp2, glob->st2, 10))
 205         colmult(glob->pr2, glob->st2, table2a, 10);
 206 }
 207
 208 /* Decode a block (celp) */
 209 static int ra288_decode_frame(AVCodecContext * avctx, void *data,
 210                               int *data_size, const uint8_t * buf,
 211                               int buf_size)
 212 {
 213     int16_t *out = data;
 214     int x, y;
 215     Real288_internal *glob = avctx->priv_data;
 216     GetBitContext gb;
 217
 218     if (buf_size < avctx->block_align) {
 219         av_log(avctx, AV_LOG_ERROR,
 220                "Error! Input buffer is too small [%d<%d]\n",
 221                buf_size, avctx->block_align);
 222         return 0;
 223     }
 224
 225     init_get_bits(&gb, buf, avctx->block_align * 8);
 226
 227     for (x=0; x < 32; x++) {
 228         int amp_coef = get_bits(&gb, 3);
 229         int cb_coef = get_bits(&gb, 6 + (x&1));
 230         glob->phasep = (glob->phase = x & 7) * 5;
 231         decode(glob, amp_coef, cb_coef);
 232
 233         for (y=0; y<5; *(out++) = 8 * glob->output[glob->phasep+(y++)]);
 234
 235         if (glob->phase == 3)
 236             update(glob);
 237     }
 238
 239     *data_size = (char *)out - (char *)data;
 240     return avctx->block_align;
 241 }
 242
 243 AVCodec ra_288_decoder =
 244 {
 245     "real_288",
 246     CODEC_TYPE_AUDIO,
 247     CODEC_ID_RA_288,
 248     sizeof(Real288_internal),
 249     NULL,
 250     NULL,
 251     NULL,
 252     ra288_decode_frame,
 253     .long_name = NULL_IF_CONFIG_SMALL("RealAudio 2.0 (28.8K)"),
 254 };