codecs/mp3/src/icdct.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 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 /****  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 /*------------------------------------------------------------*/
  52 DCTCOEF *i_dct_coef_addr(void)
  53 {
  54    return coef32;
  55 }
  56 /*------------------------------------------------------------*/
  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 /*------------------------------------------------------------*/
  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 /*------------------------------------------------------------*/
 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 /*------------------------------------------------------------*/
 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 /*------------------------------------------------------------*/
 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 /*------------------------------------------------------------*/
 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------------------------------*/
 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 -------------------------------*/
 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---------------------*/
 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-------------------*/
 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 -------------------------------*/
 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---------------------*/
 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---------------------*/
 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 /*------------------------------------------------------------*/