Version 0.1.1 from FTP
[asterisk/asterisk.git] / codecs / mp3 / src / cdct.c
1 /*____________________________________________________________________________
2         
3         FreeAmp - The Free MP3 Player
4
5         MP3 Decoder originally Copyright (C) 1995-1997 Xing Technology
6         Corp.  http://www.xingtech.com
7
8         Portions Copyright (C) 1998-1999 EMusic.com
9
10         This program is free software; you can redistribute it and/or modify
11         it under the terms of the GNU General Public License as published by
12         the Free Software Foundation; either version 2 of the License, or
13         (at your option) any later version.
14
15         This program is distributed in the hope that it will be useful,
16         but WITHOUT ANY WARRANTY; without even the implied warranty of
17         MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
18         GNU General Public License for more details.
19
20         You should have received a copy of the GNU General Public License
21         along with this program; if not, write to the Free Software
22         Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
23         
24         $Id$
25 ____________________________________________________________________________*/
26
27 /****  cdct.c  ***************************************************
28
29 mod 5/16/95 first stage in 8 pt dct does not drop last sb mono
30
31
32 MPEG audio decoder, dct
33 portable C
34
35 ******************************************************************/
36
37 #include <stdlib.h>
38 #include <stdio.h>
39 #include <float.h>
40 #include <math.h>
41 #include "L3.h"
42 #include "mhead.h"
43
44 #ifdef ASM_X86
45 extern void fdct32_asm(float*a, float*b);
46 extern void fdct32_dual_asm(float*a, float*b);
47 #endif /* ASM_X86 */
48
49 /*------------------------------------------------------------*/
50 float *dct_coef_addr(MPEG *m)
51 {
52    return m->cdct.coef32;
53 }
54 /*------------------------------------------------------------*/
55 static void forward_bf(int m, int n, float x[], float f[], float coef[])
56 {
57    int i, j, n2;
58    int p, q, p0, k;
59
60    p0 = 0;
61    n2 = n >> 1;
62    for (i = 0; i < m; i++, p0 += n)
63    {
64       k = 0;
65       p = p0;
66       q = p + n - 1;
67       for (j = 0; j < n2; j++, p++, q--, k++)
68       {
69          f[p] = x[p] + x[q];
70          f[n2 + p] = coef[k] * (x[p] - x[q]);
71       }
72    }
73 }
74 /*------------------------------------------------------------*/
75 static void back_bf(int m, int n, float x[], float f[])
76 {
77    int i, j, n2, n21;
78    int p, q, p0;
79
80    p0 = 0;
81    n2 = n >> 1;
82    n21 = n2 - 1;
83    for (i = 0; i < m; i++, p0 += n)
84    {
85       p = p0;
86       q = p0;
87       for (j = 0; j < n2; j++, p += 2, q++)
88          f[p] = x[q];
89       p = p0 + 1;
90       for (j = 0; j < n21; j++, p += 2, q++)
91          f[p] = x[q] + x[q + 1];
92       f[p] = x[q];
93    }
94 }
95 /*------------------------------------------------------------*/
96
97 #ifdef  _EQUALIZER_ENABLE_
98 extern float equalizer[32];
99 extern int enableEQ;
100 #endif
101
102 void fdct32(MPEG *m, float x[], float c[])
103 {
104 #if (!defined(ASM_X86) && !defined(ASM_X86_OLD) || defined(_EQUALIZER_ENABLE_))
105    float a[32];                 /* ping pong buffers */
106    float b[32];
107    int p, q;
108 #endif
109
110    float *src = x;
111
112 #ifdef  _EQUALIZER_ENABLE_
113    int i;
114    float b[32];
115    if (enableEQ) {
116        for(i=0; i<32; i++)
117            b[i] = x[i] * equalizer[i];
118        src = b;
119    }
120 #endif  /* _EQUALIZER_ENABLE_ */
121 #undef  _EQUALIZER_ENABLE_
122
123 #ifdef ASM_X86
124    fdct32_asm(src, c);
125 #elif defined(ASM_X86_OLD)
126    asm_fdct32(src, c);
127 #else
128 /* special first stage */
129    for (p = 0, q = 31; p < 16; p++, q--)
130    {
131       a[p] = src[p] + src[q];
132       a[16 + p] = m->cdct.coef32[p] * (src[p] - src[q]);
133    }
134    forward_bf(2, 16, a, b, m->cdct.coef32 + 16);
135    forward_bf(4, 8, b, a, m->cdct.coef32 + 16 + 8);
136    forward_bf(8, 4, a, b, m->cdct.coef32 + 16 + 8 + 4);
137    forward_bf(16, 2, b, a, m->cdct.coef32 + 16 + 8 + 4 + 2);
138    back_bf(8, 4, a, b);
139    back_bf(4, 8, b, a);
140    back_bf(2, 16, a, b);
141    back_bf(1, 32, b, c);
142 #endif
143 }
144 /*------------------------------------------------------------*/
145 void fdct32_dual(MPEG *m, float x[], float c[])
146 {
147 #ifdef ASM_X86
148    fdct32_dual_asm(x, c);
149 #else
150    float a[32];                 /* ping pong buffers */
151    float b[32];
152    int p, pp, qq;
153
154 /* special first stage for dual chan (interleaved x) */
155    pp = 0;
156    qq = 2 * 31;
157    for (p = 0; p < 16; p++, pp += 2, qq -= 2)
158    {
159       a[p] = x[pp] + x[qq];
160       a[16 + p] = m->cdct.coef32[p] * (x[pp] - x[qq]);
161    }
162    forward_bf(2, 16, a, b, m->cdct.coef32 + 16);
163    forward_bf(4, 8, b, a, m->cdct.coef32 + 16 + 8);
164    forward_bf(8, 4, a, b, m->cdct.coef32 + 16 + 8 + 4);
165    forward_bf(16, 2, b, a, m->cdct.coef32 + 16 + 8 + 4 + 2);
166    back_bf(8, 4, a, b);
167    back_bf(4, 8, b, a);
168    back_bf(2, 16, a, b);
169    back_bf(1, 32, b, c);
170 #endif
171 }
172 /*---------------convert dual to mono------------------------------*/
173 void fdct32_dual_mono(MPEG *m, float x[], float c[])
174 {
175    float a[32];                 /* ping pong buffers */
176    float b[32];
177    float t1, t2;
178    int p, pp, qq;
179
180 /* special first stage  */
181    pp = 0;
182    qq = 2 * 31;
183    for (p = 0; p < 16; p++, pp += 2, qq -= 2)
184    {
185       t1 = 0.5F * (x[pp] + x[pp + 1]);
186       t2 = 0.5F * (x[qq] + x[qq + 1]);
187       a[p] = t1 + t2;
188       a[16 + p] = m->cdct.coef32[p] * (t1 - t2);
189    }
190    forward_bf(2, 16, a, b, m->cdct.coef32 + 16);
191    forward_bf(4, 8, b, a, m->cdct.coef32 + 16 + 8);
192    forward_bf(8, 4, a, b, m->cdct.coef32 + 16 + 8 + 4);
193    forward_bf(16, 2, b, a, m->cdct.coef32 + 16 + 8 + 4 + 2);
194    back_bf(8, 4, a, b);
195    back_bf(4, 8, b, a);
196    back_bf(2, 16, a, b);
197    back_bf(1, 32, b, c);
198 }
199 /*------------------------------------------------------------*/
200 /*---------------- 16 pt fdct -------------------------------*/
201 void fdct16(MPEG *m, float x[], float c[])
202 {
203    float a[16];                 /* ping pong buffers */
204    float b[16];
205    int p, q;
206
207 /* special first stage (drop highest sb) */
208    a[0] = x[0];
209    a[8] = m->cdct.coef32[16] * x[0];
210    for (p = 1, q = 14; p < 8; p++, q--)
211    {
212       a[p] = x[p] + x[q];
213       a[8 + p] = m->cdct.coef32[16 + p] * (x[p] - x[q]);
214    }
215    forward_bf(2, 8, a, b, m->cdct.coef32 + 16 + 8);
216    forward_bf(4, 4, b, a, m->cdct.coef32 + 16 + 8 + 4);
217    forward_bf(8, 2, a, b, m->cdct.coef32 + 16 + 8 + 4 + 2);
218    back_bf(4, 4, b, a);
219    back_bf(2, 8, a, b);
220    back_bf(1, 16, b, c);
221 }
222 /*------------------------------------------------------------*/
223 /*---------------- 16 pt fdct dual chan---------------------*/
224 void fdct16_dual(MPEG *m, float x[], float c[])
225 {
226    float a[16];                 /* ping pong buffers */
227    float b[16];
228    int p, pp, qq;
229
230 /* special first stage for interleaved input */
231    a[0] = x[0];
232    a[8] = m->cdct.coef32[16] * x[0];
233    pp = 2;
234    qq = 2 * 14;
235    for (p = 1; p < 8; p++, pp += 2, qq -= 2)
236    {
237       a[p] = x[pp] + x[qq];
238       a[8 + p] = m->cdct.coef32[16 + p] * (x[pp] - x[qq]);
239    }
240    forward_bf(2, 8, a, b, m->cdct.coef32 + 16 + 8);
241    forward_bf(4, 4, b, a, m->cdct.coef32 + 16 + 8 + 4);
242    forward_bf(8, 2, a, b, m->cdct.coef32 + 16 + 8 + 4 + 2);
243    back_bf(4, 4, b, a);
244    back_bf(2, 8, a, b);
245    back_bf(1, 16, b, c);
246 }
247 /*------------------------------------------------------------*/
248 /*---------------- 16 pt fdct dual to mono-------------------*/
249 void fdct16_dual_mono(MPEG *m, float x[], float c[])
250 {
251    float a[16];                 /* ping pong buffers */
252    float b[16];
253    float t1, t2;
254    int p, pp, qq;
255
256 /* special first stage  */
257    a[0] = 0.5F * (x[0] + x[1]);
258    a[8] = m->cdct.coef32[16] * a[0];
259    pp = 2;
260    qq = 2 * 14;
261    for (p = 1; p < 8; p++, pp += 2, qq -= 2)
262    {
263       t1 = 0.5F * (x[pp] + x[pp + 1]);
264       t2 = 0.5F * (x[qq] + x[qq + 1]);
265       a[p] = t1 + t2;
266       a[8 + p] = m->cdct.coef32[16 + p] * (t1 - t2);
267    }
268    forward_bf(2, 8, a, b, m->cdct.coef32 + 16 + 8);
269    forward_bf(4, 4, b, a, m->cdct.coef32 + 16 + 8 + 4);
270    forward_bf(8, 2, a, b, m->cdct.coef32 + 16 + 8 + 4 + 2);
271    back_bf(4, 4, b, a);
272    back_bf(2, 8, a, b);
273    back_bf(1, 16, b, c);
274 }
275 /*------------------------------------------------------------*/
276 /*---------------- 8 pt fdct -------------------------------*/
277 void fdct8(MPEG *m, float x[], float c[])
278 {
279    float a[8];                  /* ping pong buffers */
280    float b[8];
281    int p, q;
282
283 /* special first stage  */
284
285    b[0] = x[0] + x[7];
286    b[4] = m->cdct.coef32[16 + 8] * (x[0] - x[7]);
287    for (p = 1, q = 6; p < 4; p++, q--)
288    {
289       b[p] = x[p] + x[q];
290       b[4 + p] = m->cdct.coef32[16 + 8 + p] * (x[p] - x[q]);
291    }
292
293    forward_bf(2, 4, b, a, m->cdct.coef32 + 16 + 8 + 4);
294    forward_bf(4, 2, a, b, m->cdct.coef32 + 16 + 8 + 4 + 2);
295    back_bf(2, 4, b, a);
296    back_bf(1, 8, a, c);
297 }
298 /*------------------------------------------------------------*/
299 /*---------------- 8 pt fdct dual chan---------------------*/
300 void fdct8_dual(MPEG *m, float x[], float c[])
301 {
302    float a[8];                  /* ping pong buffers */
303    float b[8];
304    int p, pp, qq;
305
306 /* special first stage for interleaved input */
307    b[0] = x[0] + x[14];
308    b[4] = m->cdct.coef32[16 + 8] * (x[0] - x[14]);
309    pp = 2;
310    qq = 2 * 6;
311    for (p = 1; p < 4; p++, pp += 2, qq -= 2)
312    {
313       b[p] = x[pp] + x[qq];
314       b[4 + p] = m->cdct.coef32[16 + 8 + p] * (x[pp] - x[qq]);
315    }
316    forward_bf(2, 4, b, a, m->cdct.coef32 + 16 + 8 + 4);
317    forward_bf(4, 2, a, b, m->cdct.coef32 + 16 + 8 + 4 + 2);
318    back_bf(2, 4, b, a);
319    back_bf(1, 8, a, c);
320 }
321 /*------------------------------------------------------------*/
322 /*---------------- 8 pt fdct dual to mono---------------------*/
323 void fdct8_dual_mono(MPEG *m, float x[], float c[])
324 {
325    float a[8];                  /* ping pong buffers */
326    float b[8];
327    float t1, t2;
328    int p, pp, qq;
329
330 /* special first stage  */
331    t1 = 0.5F * (x[0] + x[1]);
332    t2 = 0.5F * (x[14] + x[15]);
333    b[0] = t1 + t2;
334    b[4] = m->cdct.coef32[16 + 8] * (t1 - t2);
335    pp = 2;
336    qq = 2 * 6;
337    for (p = 1; p < 4; p++, pp += 2, qq -= 2)
338    {
339       t1 = 0.5F * (x[pp] + x[pp + 1]);
340       t2 = 0.5F * (x[qq] + x[qq + 1]);
341       b[p] = t1 + t2;
342       b[4 + p] = m->cdct.coef32[16 + 8 + p] * (t1 - t2);
343    }
344    forward_bf(2, 4, b, a, m->cdct.coef32 + 16 + 8 + 4);
345    forward_bf(4, 2, a, b, m->cdct.coef32 + 16 + 8 + 4 + 2);
346    back_bf(2, 4, b, a);
347    back_bf(1, 8, a, c);
348 }
349 /*------------------------------------------------------------*/