git.sesse.net Git - vlc/blob - plugins/imdct/ac3_imdct_c.c

   1 /*****************************************************************************
   2  * ac3_imdct_c.c: ac3 DCT in C
   3  *****************************************************************************
   4  * Copyright (C) 1999-2001 VideoLAN
   5  * $Id: ac3_imdct_c.c,v 1.5 2001/12/30 07:09:55 sam Exp $
   6  *
   7  * Authors: Renaud Dartus <reno@videolan.org>
   8  *          Aaron Holtzman <aholtzma@engr.uvic.ca>
   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., 59 Temple Place - Suite 330, Boston, MA  02111, USA.
  23  *****************************************************************************/
  24
  25 /*****************************************************************************
  26  * Preamble
  27  *****************************************************************************/
  28 #include <string.h>                                              /* memcpy() */
  29
  30 #include <math.h>
  31 #include <stdio.h>
  32
  33 #include <videolan/vlc.h>
  34
  35 #include "ac3_imdct.h"
  36 #include "ac3_imdct_common.h"
  37 #include "ac3_retables.h"
  38
  39 #ifndef M_PI
  40 #   define M_PI 3.14159265358979323846
  41 #endif
  42
  43 void _M( fft_64p )  ( complex_t *x );
  44 void _M( fft_128p ) ( complex_t *x );
  45
  46 void _M( imdct_init ) (imdct_t * p_imdct)
  47 {
  48     int i;
  49     float scale = 181.019;
  50
  51     /* Twiddle factors to turn IFFT into IMDCT */
  52     for (i=0; i < 128; i++) {
  53         p_imdct->xcos1[i] = cos(2.0f * M_PI * (8*i+1)/(8*N)) * scale;
  54         p_imdct->xsin1[i] = sin(2.0f * M_PI * (8*i+1)/(8*N)) * scale;
  55     }
  56 }
  57
  58 void _M( imdct_do_512 ) (imdct_t * p_imdct, float data[], float delay[])
  59 {
  60     int i, j;
  61     float tmp_a_r, tmp_a_i;
  62     float *data_ptr;
  63     float *delay_ptr;
  64     float *window_ptr;
  65
  66     /* 512 IMDCT with source and dest data in 'data'
  67      * Pre IFFT complex multiply plus IFFT complex conjugate */
  68
  69     for( i=0; i < 128; i++) {
  70         j = pm128[i];
  71         /* a = (data[256-2*j-1] - data[2*j]) * (xcos1[j] + xsin1[j]);
  72          * c = data[2*j] * xcos1[j];
  73          * b = data[256-2*j-1] * xsin1[j];
  74          * buf1[i].real = a - b + c;
  75          * buf1[i].imag = b + c; */
  76         p_imdct->buf[i].real = (data[256-2*j-1] * p_imdct->xcos1[j]) - (data[2*j] * p_imdct->xsin1[j]);
  77         p_imdct->buf[i].imag = -1.0 * (data[2*j] * p_imdct->xcos1[j] + data[256-2*j-1] * p_imdct->xsin1[j]);
  78     }
  79
  80     _M( fft_128p ) ( &p_imdct->buf[0] );
  81
  82     /* Post IFFT complex multiply  plus IFFT complex conjugate */
  83     for (i=0; i < 128; i++) {
  84         tmp_a_r = p_imdct->buf[i].real;
  85         tmp_a_i = p_imdct->buf[i].imag;
  86         /* a = (tmp_a_r - tmp_a_i) * (xcos1[j] + xsin1[j]);
  87          * b = tmp_a_r * xsin1[j];
  88          * c = tmp_a_i * xcos1[j];
  89          * buf[j].real = a - b + c;
  90          * buf[j].imag = b + c; */
  91         p_imdct->buf[i].real =(tmp_a_r * p_imdct->xcos1[i])  +  (tmp_a_i  * p_imdct->xsin1[i]);
  92         p_imdct->buf[i].imag =(tmp_a_r * p_imdct->xsin1[i])  -  (tmp_a_i  * p_imdct->xcos1[i]);
  93     }
  94
  95     data_ptr = data;
  96     delay_ptr = delay;
  97     window_ptr = window;
  98
  99     /* Window and convert to real valued signal */
 100     for (i=0; i< 64; i++) {
 101         *data_ptr++ = -p_imdct->buf[64+i].imag  * *window_ptr++ + *delay_ptr++;
 102         *data_ptr++ = p_imdct->buf[64-i-1].real * *window_ptr++ + *delay_ptr++;
 103     }
 104
 105     for(i=0; i< 64; i++) {
 106         *data_ptr++ = -p_imdct->buf[i].real      * *window_ptr++ + *delay_ptr++;
 107         *data_ptr++ = p_imdct->buf[128-i-1].imag * *window_ptr++ + *delay_ptr++;
 108     }
 109
 110     /* The trailing edge of the window goes into the delay line */
 111     delay_ptr = delay;
 112
 113     for(i=0; i< 64; i++) {
 114         *delay_ptr++ = -p_imdct->buf[64+i].real   * *--window_ptr;
 115         *delay_ptr++ =  p_imdct->buf[64-i-1].imag * *--window_ptr;
 116     }
 117
 118     for(i=0; i<64; i++) {
 119         *delay_ptr++ =  p_imdct->buf[i].imag       * *--window_ptr;
 120         *delay_ptr++ = -p_imdct->buf[128-i-1].real * *--window_ptr;
 121     }
 122 }
 123
 124
 125 void _M( imdct_do_512_nol ) (imdct_t * p_imdct, float data[], float delay[])
 126 {
 127     int i, j;
 128
 129     float tmp_a_i;
 130     float tmp_a_r;
 131
 132     float *data_ptr;
 133     float *delay_ptr;
 134     float *window_ptr;
 135
 136     /* 512 IMDCT with source and dest data in 'data'
 137      * Pre IFFT complex multiply plus IFFT cmplx conjugate */
 138
 139     for( i=0; i < 128; i++) {
 140         /* z[i] = (X[256-2*i-1] + j * X[2*i]) * (xcos1[i] + j * xsin1[i]) */
 141         j = pm128[i];
 142         /* a = (data[256-2*j-1] - data[2*j]) * (xcos1[j] + xsin1[j]);
 143          * c = data[2*j] * xcos1[j];
 144          * b = data[256-2*j-1] * xsin1[j];
 145          * buf1[i].real = a - b + c;
 146          * buf1[i].imag = b + c; */
 147         p_imdct->buf[i].real = (data[256-2*j-1] * p_imdct->xcos1[j]) - (data[2*j] * p_imdct->xsin1[j]);
 148         p_imdct->buf[i].imag = -1.0 * (data[2*j] * p_imdct->xcos1[j] + data[256-2*j-1] * p_imdct->xsin1[j]);
 149     }
 150
 151     _M( fft_128p ) ( &p_imdct->buf[0] );
 152
 153     /* Post IFFT complex multiply  plus IFFT complex conjugate*/
 154     for (i=0; i < 128; i++) {
 155         /* y[n] = z[n] * (xcos1[n] + j * xsin1[n]) ;
 156          * int j1 = i; */
 157         tmp_a_r = p_imdct->buf[i].real;
 158         tmp_a_i = p_imdct->buf[i].imag;
 159         /* a = (tmp_a_r - tmp_a_i) * (xcos1[j] + xsin1[j]);
 160          * b = tmp_a_r * xsin1[j];
 161          * c = tmp_a_i * xcos1[j];
 162          * buf[j].real = a - b + c;
 163          * buf[j].imag = b + c; */
 164         p_imdct->buf[i].real =(tmp_a_r * p_imdct->xcos1[i]) + (tmp_a_i  * p_imdct->xsin1[i]);
 165         p_imdct->buf[i].imag =(tmp_a_r * p_imdct->xsin1[i]) - (tmp_a_i  * p_imdct->xcos1[i]);
 166     }
 167
 168     data_ptr = data;
 169     delay_ptr = delay;
 170     window_ptr = window;
 171
 172     /* Window and convert to real valued signal, no overlap here*/
 173     for (i=0; i< 64; i++) {
 174         *data_ptr++ = -p_imdct->buf[64+i].imag  * *window_ptr++;
 175         *data_ptr++ = p_imdct->buf[64-i-1].real * *window_ptr++;
 176     }
 177
 178     for(i=0; i< 64; i++) {
 179         *data_ptr++ = -p_imdct->buf[i].real      * *window_ptr++;
 180         *data_ptr++ = p_imdct->buf[128-i-1].imag * *window_ptr++;
 181     }
 182
 183     /* The trailing edge of the window goes into the delay line */
 184     delay_ptr = delay;
 185
 186     for(i=0; i< 64; i++) {
 187         *delay_ptr++ = -p_imdct->buf[64+i].real   * *--window_ptr;
 188         *delay_ptr++ =  p_imdct->buf[64-i-1].imag * *--window_ptr;
 189     }
 190
 191     for(i=0; i<64; i++) {
 192         *delay_ptr++ =  p_imdct->buf[i].imag       * *--window_ptr;
 193         *delay_ptr++ = -p_imdct->buf[128-i-1].real * *--window_ptr;
 194     }
 195 }
 196