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 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: icdct.c,v 1.4 2000/10/13 14:29:02 ijr Exp $
25 ____________________________________________________________________________*/
26
27 /**** icdct.c ***************************************************
28
29
30 MPEG audio decoder, dct
31 portable C integer dct
32
33 mod 1/8/97 warnings
34
35 ******************************************************************/
36
37 #include <stdlib.h>
38 #include <stdio.h>
39 #include <float.h>
40 #include <math.h>
41 #include "itype.h"
42
43 /*-------------------------------------------------------------------*/
44 static DCTCOEF coef32[32]; /* 32 pt dct coefs */
45
46
47
48
49 #define forward_bf idx_forward_bf
50 /*--- #define forward_bf ptr_forward_bf ---*/
51 /*------------------------------------------------------------*/
i_dct_coef_addr()52 DCTCOEF *i_dct_coef_addr()
53 {
54 return coef32;
55 }
56 /*------------------------------------------------------------*/
idx_forward_bf(int m,int n,INT32 x[],INT32 f[],DCTCOEF coef[])57 static void idx_forward_bf(int m, int n, INT32 x[], INT32 f[], DCTCOEF coef[])
58 {
59 int i, j, n2;
60 int p, q, p0, k;
61
62 p0 = 0;
63 n2 = n >> 1;
64 for (i = 0; i < m; i++, p0 += n)
65 {
66 k = 0;
67 p = p0;
68 q = p + n - 1;
69 for (j = 0; j < n2; j++, p++, q--, k++)
70 {
71 f[p] = x[p] + x[q];
72 f[n2 + p] = ((x[p] - x[q]) * coef[k]) >> DCTBITS;
73 }
74 }
75 }
76 /*------------------------------------------------------------*/
77 /*--
78 static void ptr_forward_bf(int m, int n, INT32 x[], INT32 f[], DCTCOEF coef[])
79 {
80 int i, j, n2;
81 DCTCOEF *c;
82 INT32 *y;
83
84 n2 = n >> 1;
85 for(i=0; i<m; i++) {
86 c = coef;
87 y = x+n;
88 for(j=0; j<n2; j++) {
89 *f = *x + *--y;
90 *((f++)+n2) = ( (*x++ - *y) * (*c++) ) >> DCTBITS;
91 }
92 f+=n2;
93 x+=n2;
94 }
95 }
96 ---*/
97 /*------------------------------------------------------------*/
forward_bfm(int m,INT32 x[],INT32 f[])98 static void forward_bfm(int m, INT32 x[], INT32 f[])
99 {
100 int i;
101 int p;
102
103 /*--- special case last fwd stage ----*/
104 for (p = 0, i = 0; i < m; i++, p += 2)
105 {
106 f[p] = x[p] + x[p + 1];
107 f[p + 1] = ((x[p] - x[p + 1]) * coef32[30]) >> DCTBITS;
108 }
109 }
110 /*------------------------------------------------------------*/
back_bf(int m,int n,INT32 x[],INT32 f[])111 static void back_bf(int m, int n, INT32 x[], INT32 f[])
112 {
113 int i, j, n2, n21;
114 int p, q, p0;
115
116 p0 = 0;
117 n2 = n >> 1;
118 n21 = n2 - 1;
119 for (i = 0; i < m; i++, p0 += n)
120 {
121 p = p0;
122 q = p0;
123 for (j = 0; j < n2; j++, p += 2, q++)
124 f[p] = x[q];
125 p = p0 + 1;
126 for (j = 0; j < n21; j++, p += 2, q++)
127 f[p] = x[q] + x[q + 1];
128 f[p] = x[q];
129 }
130 }
131 /*------------------------------------------------------------*/
back_bf0(int n,INT32 x[],WININT f[])132 static void back_bf0(int n, INT32 x[], WININT f[])
133 {
134 int p, q;
135
136 n--;
137 #if DCTSATURATE
138 for (p = 0, q = 0; p < n; p += 2, q++)
139 {
140 tmp = x[q];
141 if (tmp > 32767)
142 tmp = 32767;
143 else if (tmp < -32768)
144 tmp = -32768;
145 f[p] = tmp;
146 }
147 for (p = 1; q < n; p += 2, q++)
148 {
149 tmp = x[q] + x[q + 1];
150 if (tmp > 32767)
151 tmp = 32767;
152 else if (tmp < -32768)
153 tmp = -32768;
154 f[p] = tmp;
155 }
156 tmp = x[q];
157 if (tmp > 32767)
158 tmp = 32767;
159 else if (tmp < -32768)
160 tmp = -32768;
161 f[p] = tmp;
162 #else
163 for (p = 0, q = 0; p < n; p += 2, q++)
164 f[p] = x[q];
165 for (p = 1; q < n; p += 2, q++)
166 f[p] = x[q] + x[q + 1];
167 f[p] = x[q];
168 #endif
169
170 }
171 /*------------------------------------------------------------*/
i_dct32(SAMPLEINT x[],WININT c[])172 void i_dct32(SAMPLEINT x[], WININT c[])
173 {
174 INT32 a[32]; /* ping pong buffers */
175 INT32 b[32];
176 int p, q;
177
178 /* special first stage */
179 for (p = 0, q = 31; p < 16; p++, q--)
180 {
181 a[p] = (INT32) x[p] + x[q];
182 a[16 + p] = (coef32[p] * ((INT32) x[p] - x[q])) >> DCTBITS;
183 }
184
185 forward_bf(2, 16, a, b, coef32 + 16);
186 forward_bf(4, 8, b, a, coef32 + 16 + 8);
187 forward_bf(8, 4, a, b, coef32 + 16 + 8 + 4);
188 forward_bfm(16, b, a);
189 back_bf(8, 4, a, b);
190 back_bf(4, 8, b, a);
191 back_bf(2, 16, a, b);
192 back_bf0(32, b, c);
193 }
194 /*------------------------------------------------------------*/
i_dct32_dual(SAMPLEINT x[],WININT c[])195 void i_dct32_dual(SAMPLEINT x[], WININT c[])
196 {
197 INT32 a[32]; /* ping pong buffers */
198 INT32 b[32];
199 int p, pp, qq;
200
201 /* special first stage for dual chan (interleaved x) */
202 pp = 0;
203 qq = 2 * 31;
204 for (p = 0; p < 16; p++, pp += 2, qq -= 2)
205 {
206 a[p] = (INT32) x[pp] + x[qq];
207 a[16 + p] = (coef32[p] * ((INT32) x[pp] - x[qq])) >> DCTBITS;
208 }
209 forward_bf(2, 16, a, b, coef32 + 16);
210 forward_bf(4, 8, b, a, coef32 + 16 + 8);
211 forward_bf(8, 4, a, b, coef32 + 16 + 8 + 4);
212 forward_bfm(16, b, a);
213 back_bf(8, 4, a, b);
214 back_bf(4, 8, b, a);
215 back_bf(2, 16, a, b);
216 back_bf0(32, b, c);
217 }
218 /*---------------convert dual to mono------------------------------*/
i_dct32_dual_mono(SAMPLEINT x[],WININT c[])219 void i_dct32_dual_mono(SAMPLEINT x[], WININT c[])
220 {
221 INT32 a[32]; /* ping pong buffers */
222 INT32 b[32];
223 INT32 t1, t2;
224 int p, pp, qq;
225
226 /* special first stage */
227 pp = 0;
228 qq = 2 * 31;
229 for (p = 0; p < 16; p++, pp += 2, qq -= 2)
230 {
231 t1 = ((INT32) x[pp] + x[pp + 1]);
232 t2 = ((INT32) x[qq] + x[qq + 1]);
233 a[p] = (t1 + t2) >> 1;
234 a[16 + p] = coef32[p] * (t1 - t2) >> (DCTBITS + 1);
235 }
236 forward_bf(2, 16, a, b, coef32 + 16);
237 forward_bf(4, 8, b, a, coef32 + 16 + 8);
238 forward_bf(8, 4, a, b, coef32 + 16 + 8 + 4);
239 forward_bfm(16, b, a);
240 back_bf(8, 4, a, b);
241 back_bf(4, 8, b, a);
242 back_bf(2, 16, a, b);
243 back_bf0(32, b, c);
244 }
245 /*------------------------------------------------------------*/
246 /*---------------- 16 pt dct -------------------------------*/
i_dct16(SAMPLEINT x[],WININT c[])247 void i_dct16(SAMPLEINT x[], WININT c[])
248 {
249 INT32 a[16]; /* ping pong buffers */
250 INT32 b[16];
251 int p, q;
252
253 /* special first stage (drop highest sb) */
254 a[0] = x[0];
255 a[8] = (a[0] * coef32[16]) >> DCTBITS;
256 for (p = 1, q = 14; p < 8; p++, q--)
257 {
258 a[p] = (INT32) x[p] + x[q];
259 a[8 + p] = (((INT32) x[p] - x[q]) * coef32[16 + p]) >> DCTBITS;
260 }
261 forward_bf(2, 8, a, b, coef32 + 16 + 8);
262 forward_bf(4, 4, b, a, coef32 + 16 + 8 + 4);
263 forward_bfm(8, a, b);
264 back_bf(4, 4, b, a);
265 back_bf(2, 8, a, b);
266 back_bf0(16, b, c);
267 }
268 /*------------------------------------------------------------*/
269 /*---------------- 16 pt dct dual chan---------------------*/
i_dct16_dual(SAMPLEINT x[],WININT c[])270 void i_dct16_dual(SAMPLEINT x[], WININT c[])
271 {
272 int p, pp, qq;
273 INT32 a[16]; /* ping pong buffers */
274 INT32 b[16];
275
276 /* special first stage for interleaved input */
277 a[0] = x[0];
278 a[8] = (coef32[16] * a[0]) >> DCTBITS;
279 pp = 2;
280 qq = 2 * 14;
281 for (p = 1; p < 8; p++, pp += 2, qq -= 2)
282 {
283 a[p] = (INT32) x[pp] + x[qq];
284 a[8 + p] = (coef32[16 + p] * ((INT32) x[pp] - x[qq])) >> DCTBITS;
285 }
286 forward_bf(2, 8, a, b, coef32 + 16 + 8);
287 forward_bf(4, 4, b, a, coef32 + 16 + 8 + 4);
288 forward_bfm(8, a, b);
289 back_bf(4, 4, b, a);
290 back_bf(2, 8, a, b);
291 back_bf0(16, b, c);
292 }
293 /*------------------------------------------------------------*/
294 /*---------------- 16 pt dct dual to mono-------------------*/
i_dct16_dual_mono(SAMPLEINT x[],WININT c[])295 void i_dct16_dual_mono(SAMPLEINT x[], WININT c[])
296 {
297 INT32 a[16]; /* ping pong buffers */
298 INT32 b[16];
299 INT32 t1, t2;
300 int p, pp, qq;
301
302 /* special first stage */
303 a[0] = ((INT32) x[0] + x[1]) >> 1;
304 a[8] = (coef32[16] * a[0]) >> DCTBITS;
305 pp = 2;
306 qq = 2 * 14;
307 for (p = 1; p < 8; p++, pp += 2, qq -= 2)
308 {
309 t1 = (INT32) x[pp] + x[pp + 1];
310 t2 = (INT32) x[qq] + x[qq + 1];
311 a[p] = (t1 + t2) >> 1;
312 a[8 + p] = (coef32[16 + p] * (t1 - t2)) >> (DCTBITS + 1);
313 }
314 forward_bf(2, 8, a, b, coef32 + 16 + 8);
315 forward_bf(4, 4, b, a, coef32 + 16 + 8 + 4);
316 forward_bfm(8, a, b);
317 back_bf(4, 4, b, a);
318 back_bf(2, 8, a, b);
319 back_bf0(16, b, c);
320 }
321 /*------------------------------------------------------------*/
322 /*---------------- 8 pt dct -------------------------------*/
i_dct8(SAMPLEINT x[],WININT c[])323 void i_dct8(SAMPLEINT x[], WININT c[])
324 {
325 int p, q;
326 INT32 a[8]; /* ping pong buffers */
327 INT32 b[8];
328
329 /* special first stage */
330
331 for (p = 0, q = 7; p < 4; p++, q--)
332 {
333 b[p] = (INT32) x[p] + x[q];
334 b[4 + p] = (coef32[16 + 8 + p] * ((INT32) x[p] - x[q])) >> DCTBITS;
335 }
336
337 forward_bf(2, 4, b, a, coef32 + 16 + 8 + 4);
338 forward_bfm(4, a, b);
339 back_bf(2, 4, b, a);
340 back_bf0(8, a, c);
341 }
342 /*------------------------------------------------------------*/
343 /*---------------- 8 pt dct dual chan---------------------*/
i_dct8_dual(SAMPLEINT x[],WININT c[])344 void i_dct8_dual(SAMPLEINT x[], WININT c[])
345 {
346 int p, pp, qq;
347 INT32 a[8]; /* ping pong buffers */
348 INT32 b[8];
349
350 /* special first stage for interleaved input */
351 for (p = 0, pp = 0, qq = 14; p < 4; p++, pp += 2, qq -= 2)
352 {
353 b[p] = (INT32) x[pp] + x[qq];
354 b[4 + p] = (coef32[16 + 8 + p] * ((INT32) x[pp] - x[qq])) >> DCTBITS;
355 }
356 forward_bf(2, 4, b, a, coef32 + 16 + 8 + 4);
357 forward_bfm(4, a, b);
358 back_bf(2, 4, b, a);
359 back_bf0(8, a, c);
360 }
361 /*------------------------------------------------------------*/
362 /*---------------- 8 pt dct dual to mono---------------------*/
i_dct8_dual_mono(SAMPLEINT x[],WININT c[])363 void i_dct8_dual_mono(SAMPLEINT x[], WININT c[])
364 {
365 int p, pp, qq;
366 INT32 a[8]; /* ping pong buffers */
367 INT32 b[8];
368 INT32 t1, t2;
369
370 /* special first stage */
371 for (p = 0, pp = 0, qq = 14; p < 4; p++, pp += 2, qq -= 2)
372 {
373 t1 = (INT32) x[pp] + x[pp + 1];
374 t2 = (INT32) x[qq] + x[qq + 1];
375 b[p] = (t1 + t2) >> 1;
376 b[4 + p] = (coef32[16 + 8 + p] * (t1 - t2)) >> (DCTBITS + 1);
377 }
378 forward_bf(2, 4, b, a, coef32 + 16 + 8 + 4);
379 forward_bfm(4, a, b);
380 back_bf(2, 4, b, a);
381 back_bf0(8, a, c);
382 }
383 /*------------------------------------------------------------*/
384