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