/****************************************************************************//* *	mpegimdc.c -- MPEG Audio IMDCT hybrid filter * *	Author  :   St閜hane TAVENARD * *	(C) Copyright 1997-1999 St閜hane TAVENARD *	    All Rights Reserved * *	This program is free software; you can redistribute it and/or modify *	it under the terms of the GNU General Public License as published by *	the Free Software Foundation; either version 2 of the License, or *	(at your option) any later version. *	 *	This program is distributed in the hope that it will be useful, *	but WITHOUT ANY WARRANTY; without even the implied warranty of *	MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the *	GNU General Public License for more details. *	 *	You should have received a copy of the GNU General Public License *	along with this program; if not, write to the Free Software *	Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA. *//****************************************************************************/#include "defs.h"#include "mpegaud.h"#include "mpegimdc.h"#include <math.h>#ifdef MPEGAUD_INT#define IMDCT_INT#ifdef USE_IMDCT_TABLE#include "imdct.h"#endif#else// This one is faster for float (but not accurate for int math)#define NEW_IMDCT	#endif#ifdef IMDCT_TABLE_GEN#include <stdlib.h>#endif/****************************************************************************/#ifdef NEW_IMDCT#define SBLIMIT	MPA_SBLIMIT#define SSLIMIT	MPA_SSLIMIT#ifndef M_PI#define M_PI	3.14159265358979#endif#ifdef IMDCT_INT#define WIN_BITS 14#define TFC_BITS 14#define COS_BITS 14#define WIN_MULT(val, windex)	(((val) * w[windex]) >> WIN_BITS)#define T36_MULT(val, tindex)	(((val) * tfcos36[tindex]) >> TFC_BITS)#define T12_MULT(val, tindex)	(((val) * tfcos12[tindex]) >> TFC_BITS)#define COS_MULT(val, cosval)	(((val) * (cosval)) >> COS_BITS)#define CM(val, cindex)		((val) * c[cindex]))#define CM1(val)		((val) << COS_BITS)#define CF(val)			((val) >> (COS_BITS+SH_BITS))#define IMDCT_TYPE	INT32#define WIN_TYPE	INT32#define COS_TYPE	INT16#else#define WIN_MULT(val, windex)	((val) * w[windex])#define T36_MULT(val, tindex)	((val) * tfcos36[tindex ])#define T12_MULT(val, tindex)	((val) * tfcos12[tindex ])#define COS_MULT(val, cosval)	((val) * (cosval))#define CM(val, cindex)		((val) * c[cindex])#define CM1(val)		(val)#define CF(val)			(val)#define IMDCT_TYPE REAL#define WIN_TYPE REAL#define COS_TYPE REAL#endif/****************************************************************************/static COS_TYPE COS1[12][6];static WIN_TYPE win[4][36];static WIN_TYPE win1[4][36];static COS_TYPE COS9[9];static COS_TYPE COS6_1, COS6_2;static IMDCT_TYPE tfcos36[9];static IMDCT_TYPE tfcos12[3];/****************************************************************************/static void dct36(MPEGAUD_FRACT_TYPE *inbuf, MPEGIMDCT_BLOCK_TYPE *prev,	WIN_TYPE *wintab, MPEGAUD_FRACT_TYPE *tsbuf){	MPEGAUD_FRACT_TYPE *in = inbuf;#define I0(i)	in[2*(i)+0]#define I1(i)	in[2*(i)+1]	in[17] += in[16];	in[16] += in[15];	in[15] += in[14];	in[14] += in[13];	in[13] += in[12];	in[12] += in[11];	in[11] += in[10];	in[10] += in[9];	in[9] += in[8];	in[8] += in[7];	in[7] += in[6];	in[6] += in[5];	in[5] += in[4];	in[4] += in[3];	in[3] += in[2];	in[2] += in[1];	in[1] += in[0];	in[17] += in[15];	in[15] += in[13];	in[13] += in[11];	in[11] += in[9];	in[9] += in[7];	in[7] += in[5];	in[5] += in[3];	in[3] += in[1];	{#define MACRO0(v) { \	IMDCT_TYPE tmp = sum0 + sum1; \	sum0 -= sum1; \	ts[SBLIMIT*(8-(v))] = prev[8-(v)] + WIN_MULT(sum0, 8-(v)); \	prev[8-(v)] = WIN_MULT(tmp, 26-(v)); \	ts[SBLIMIT*(9+(v))] = prev[9+(v)] + WIN_MULT(sum0, 9+(v)); \	prev[9+(v)] = WIN_MULT(tmp, 27+(v));}#define MACRO1(v) { \	IMDCT_TYPE sum0, sum1; \	sum0 = tmp1a + tmp2a; \	sum1 = T36_MULT(tmp1b + tmp2b, v); \	MACRO0(v);}#define MACRO2(v) { \	IMDCT_TYPE sum0, sum1; \	sum0 = tmp2a - tmp1a; \	sum1 = T36_MULT(tmp2b - tmp1b, v); \	MACRO0(v);}		const COS_TYPE *c = COS9;		WIN_TYPE *w = wintab;		MPEGAUD_FRACT_TYPE *ts = tsbuf;		IMDCT_TYPE ta33, ta66, tb33, tb66;		ta33 = CM(I0(3), 3);		ta66 = CM(I0(6), 6);		tb33 = CM(I1(3), 3);		tb66 = CM(I1(6), 6);		{			IMDCT_TYPE tmp1a, tmp2a, tmp1b, tmp2b;			tmp1a = CF(CM(I0(1), 1) + ta33 + CM(I0(5), 5) +				   CM(I0(7), 7));			tmp1b = CF(CM(I1(1), 1) + tb33 + CM(I1(5), 5) +				   CM(I1(7), 7));			tmp2a = CF(CM1(I0(0)) + CM(I0(2), 2) +				   CM(I0(4), 4) + ta66 + CM(I0(8), 8));			tmp2b = CF(CM1(I1(0)) + CM(I1(2), 2) +				   CM(I1(4), 4) + tb66 + CM(I1(8), 8));			MACRO1(0);			MACRO2(8);		}		{			IMDCT_TYPE tmp1a, tmp2a, tmp1b, tmp2b;			tmp1a = CF(CM(I0(1) - I0(5) - I0(7), 3));			tmp1b = CF(CM(I1(1) - I1(5) - I1(7), 3));			tmp2a = CF(CM(I0(2) - I0(4) - I0(8), 6) +				   CM1(I0(0) - I0(6)));			tmp2b = CF(CM(I1(2) - I1(4) - I1(8), 6) +				   CM1(I1(0) - I1(6)));			MACRO1(1);			MACRO2(7);		}		{			IMDCT_TYPE tmp1a, tmp2a, tmp1b, tmp2b;			tmp1a = CF(CM(I0(1), 5) - ta33 - CM(I0(5), 7) +				   CM(I0(7), 1));			tmp1b = CF(CM(I1(1), 5) - tb33 - CM(I1(5), 7) +				   CM(I1(7), 1));			tmp2a = CF(CM1(I0(0)) - CM(I0(2), 8) -				   CM(I0(4), 2) + ta66 + CM(I0(8), 4));			tmp2b = CF(CM1(I1(0)) - CM(I1(2), 8) -				   CM(I1(4), 2) + tb66 + CM(I1(8), 4));			MACRO1(2);			MACRO2(6);		}		{			IMDCT_TYPE tmp1a, tmp2a, tmp1b, tmp2b;			tmp1a = CF(CM(I0(1), 7) - ta33 + CM(I0(5), 1) -				   CM(I0(7), 5));			tmp1b = CF(CM(I1(1), 7) - tb33 + CM(I1(5), 1) -				   CM(I1(7), 5));			tmp2a = CF(CM1(I0(0)) - CM(I0(2), 4) +				   CM(I0(4), 8) + ta66 - CM(I0(8), 2));			tmp2b = CF(CM1(I1(0)) - CM(I1(2), 4) +				   CM(I1(4), 8) + tb66 - CM(I1(8), 2));			MACRO1(3);			MACRO2(5);		}		{			IMDCT_TYPE sum0, sum1;			sum0 = CF(CM1 (I0(0) - I0(2) + I0(4) - I0(6) + I0(8)));			sum1 = T36_MULT(CF(CM1(I1(0) - I1(2) + I1(4) - I1(6) +				I1(8))), 4);			MACRO0(4);		}	}}/****************************************************************************/static void dct12(MPEGAUD_FRACT_TYPE *in, MPEGIMDCT_BLOCK_TYPE *prev,	WIN_TYPE *w, MPEGAUD_FRACT_TYPE *ts){#define DCT12_PART1 \	in5 = in[5*3];  \	in5 += (in4 = in[4*3]); \	in4 += (in3 = in[3*3]); \	in3 += (in2 = in[2*3]); \	in2 += (in1 = in[1*3]); \	in1 += (in0 = in[0*3]); \	\	in5 += in3; in3 += in1; \	\	in2 = COS_MULT(in2, COS6_1); \	in3 = COS_MULT(in3, COS6_1);#define DCT12_PART2 \	in0 += COS_MULT(in4, COS6_2); \	\	in4 = in0 + in2; \	in0 -= in2; \	\	in1 += COS_MULT(in5, COS6_2); \	\	in5 = T12_MULT((in1 + in3), 0); \	in1 = T12_MULT((in1 - in3), 2); \	\	in3 = in4 + in5; \	in4 -= in5; \	\	in2 = in0 + in1; \	in0 -= in1;	{		IMDCT_TYPE in0, in1, in2, in3, in4, in5;		ts[SBLIMIT * 0] = prev[0];		ts[SBLIMIT * 1] = prev[1];		ts[SBLIMIT * 2] = prev[2];		ts[SBLIMIT * 3] = prev[3];		ts[SBLIMIT * 4] = prev[4];		ts[SBLIMIT * 5] = prev[5];		DCT12_PART1		{			IMDCT_TYPE tmp0, tmp1 = (in0 - in4);			{				IMDCT_TYPE tmp2 = T12_MULT((in1 - in5), 1);				tmp0 = tmp1 + tmp2;				tmp1 -= tmp2;			}			ts[(17 - 1) * SBLIMIT] =				prev[17 - 1] + WIN_MULT(tmp0, 11 - 1);			ts[(12 + 1) * SBLIMIT] =				prev[12 + 1] + WIN_MULT(tmp0, 6 + 1);			ts[(6 + 1) * SBLIMIT] =				prev[6 + 1] + WIN_MULT(tmp1, 1);			ts[(11 - 1) * SBLIMIT] =				prev[11 - 1] + WIN_MULT(tmp1, 5 - 1);		}		DCT12_PART2		ts[(17 - 0) * SBLIMIT] = prev[17 - 0] + WIN_MULT(in2, 11 - 0);		ts[(12 + 0) * SBLIMIT] = prev[12 + 0] + WIN_MULT(in2, 6 + 0);		ts[(12 + 2) * SBLIMIT] = prev[12 + 2] + WIN_MULT(in3, 6 + 2);		ts[(17 - 2) * SBLIMIT] = prev[17 - 2] + WIN_MULT(in3, 11 - 2);		ts[(6 + 0) * SBLIMIT] = prev[6 + 0] + WIN_MULT(in0, 0);		ts[(11 - 0) * SBLIMIT] = prev[11 - 0] + WIN_MULT(in0, 5 - 0);		ts[(6 + 2) * SBLIMIT] = prev[6 + 2] + WIN_MULT(in4, 2);		ts[(11 - 2) * SBLIMIT] = prev[11 - 2] + WIN_MULT(in4, 5 - 2);	}	in++;	{		IMDCT_TYPE in0, in1, in2, in3, in4, in5;		DCT12_PART1		{			IMDCT_TYPE tmp0, tmp1 = (in0 - in4);			{				IMDCT_TYPE tmp2 = T12_MULT((in1 - in5), 1);				tmp0 = tmp1 + tmp2;				tmp1 -= tmp2;			}			prev[5 - 1] = WIN_MULT(tmp0, 11 - 1);			prev[0 + 1] = WIN_MULT(tmp0, 6 + 1);			ts[(12 + 1) * SBLIMIT] += WIN_MULT(tmp1, 1);			ts[(17 - 1) * SBLIMIT] += WIN_MULT(tmp1, 5 - 1);		}		DCT12_PART2		prev[5 - 0] = WIN_MULT(in2, 11 - 0);		prev[0 + 0] = WIN_MULT(in2, 6 + 0);		prev[0 + 2] = WIN_MULT(in3, 6 + 2);		prev[5 - 2] = WIN_MULT(in3, 11 - 2);		ts[(12 + 0) * SBLIMIT] += WIN_MULT(in0, 0);		ts[(17 - 0) * SBLIMIT] += WIN_MULT(in0, 5 - 0);		ts[(12 + 2) * SBLIMIT] += WIN_MULT(in4, 2);		ts[(17 - 2) * SBLIMIT] += WIN_MULT(in4, 5 - 2);	}	in++;	{		IMDCT_TYPE in0, in1, in2, in3, in4, in5;		prev[12] = prev[13] = prev[14] = prev[15] = prev[16] =			prev[17] = (MPEGIMDCT_BLOCK_TYPE) 0;		DCT12_PART1		{			IMDCT_TYPE tmp0, tmp1 = (in0 - in4);			{				IMDCT_TYPE tmp2 = T12_MULT((in1 - in5), 1);				tmp0 = tmp1 + tmp2;				tmp1 -= tmp2;			}			prev[11 - 1] = WIN_MULT(tmp0, 11 - 1);			prev[6 + 1] = WIN_MULT(tmp0, 6 + 1);			prev[0 + 1] += WIN_MULT(tmp1, 1);			prev[5 - 1] += WIN_MULT(tmp1, 5 - 1);		}		DCT12_PART2		prev[11 - 0] = WIN_MULT(in2, 11 - 0);		prev[6 + 0] = WIN_MULT(in2, 6 + 0);		prev[6 + 2] = WIN_MULT(in3, 6 + 2);		prev[11 - 2] = WIN_MULT(in3, 11 - 2);		prev[0 + 0] += WIN_MULT(in0, 0);		prev[5 - 0] += WIN_MULT(in0, 5 - 0);		prev[0 + 2] += WIN_MULT(in4, 2);		prev[5 - 2] += WIN_MULT(in4, 5 - 2);	}}/****************************************************************************//* *	Apply the hybrid imdct to a granule *	Return 0 if Ok*/int MPEGIMDCT_hybrid(MPEGIMDCT *mpegimdct, MPEGAUD_FRACT_TYPE *in,	MPEGAUD_FRACT_TYPE *out, INT16 block_type,	BOOL mixed, INT16 ch, INT16 sb_max){	static BOOL init = FALSE;	MPEGIMDCT_BLOCK_TYPE *prev;	INT16 bt1, bt2;	INT16 sb;#define O(i)	out[i*32]	if (!init) {		int i, j;#ifdef IMDCT_INT#define KW	((double)(1<<(WIN_BITS-1))+0)#define KC	(double)(1<<COS_BITS)#define KT	((double)(1<<(TFC_BITS-1))+0)#else#define KW	0.5#define KC	1#define KT	0.5#endif		for (i = 0; i < 18; i++) {			win[0][i] = win[1][i] = KW * sin(M_PI / 72.0 *				(double) (2 * (i + 0) + 1)) / cos(M_PI *				(double) (2 * (i + 0) + 19) / 72.0);			win[0][i + 18] = win[3][i + 18] = KW * sin(M_PI / 72.0 *				(double) (2 * (i + 18) + 1)) / cos(M_PI *				(double) (2 * (i + 18) + 19) / 72.0);		}		for (i = 0; i < 6; i++) {			win[1][i + 18] = KW / cos(M_PI *				(double) (2 * (i + 18) + 19) / 72.0);			win[3][i + 12] = KW / cos(M_PI *				(double) (2 * (i + 12) + 19) / 72.0);			win[1][i + 24] = KW * sin(M_PI / 24.0 *				(double) (2 * i + 13)) / cos(M_PI *				(double) (2 * (i + 24) + 19) / 72.0);			win[1][i + 30] = win[3][i] = 0.0;			win[3][i + 6] = KW * sin(M_PI / 24.0 *				(double) (2 * i + 1)) / cos(M_PI *				(double) (2 * (i + 6) + 19) / 72.0);		}		for (i = 0; i < 9; i++)			COS9[i] = KC * cos(M_PI / 18.0 * (double) i);		for (i = 0; i < 9; i++)			tfcos36[i] =				KT / cos(M_PI * (double) (i * 2 + 1) / 36.0);		for (i = 0; i < 3; i++)			tfcos12[i] =				KT / cos(M_PI * (double) (i * 2 + 1) / 12.0);		COS6_1 = KC * cos(M_PI / 6.0 * (double) 1);		COS6_2 = KC * cos(M_PI / 6.0 * (double) 2);		for (i = 0; i < 12; i++) {			win[2][i] = KW * sin(M_PI / 24.0 *				(double) (2 * i + 1)) / cos(M_PI *				(double) (2 * i + 7) / 24.0);			for (j = 0; j < 6; j++)				COS1[i][j] = KC * cos(M_PI / 24.0 *					(double) ((2 * i + 7) * (2 * j + 1)));		}		for (j = 0; j < 4; j++) {			static const int len[4] = { 36, 36, 12, 36 };			for (i = 0; i < len[j]; i += 2)				win1[j][i] = +win[j][i];			for (i = 1; i < len[j]; i += 2)				win1[j][i] = -win[j][i];		}		init = TRUE;	}	prev = mpegimdct->prevblk[ch];	bt1 = (mixed) ? 0 : block_type;	bt2 = block_type;	if (bt2 == 2) {		if (!bt1) {			dct36(in, prev, win[0], out);			in += SSLIMIT;			dct36(in, prev + 18, win1[0], out + 1);			in += SSLIMIT;		} else {			dct12(in, prev, win[2], out);			in += SSLIMIT;			dct12(in, prev + 18, win1[2], out + 1);			in += SSLIMIT;		}		prev += 36;		out += 2;		for (sb = 2; sb < sb_max; sb += 2, out += 2, prev += 36) {			dct12(in, prev, win[2], out);			in += SSLIMIT;			dct12(in, prev + 18, win1[2], out + 1);			in += SSLIMIT;		}	} else {		dct36(in, prev, win[bt1], out);		in += SSLIMIT;		dct36(in, prev + 18, win1[bt1], out + 1);		in += SSLIMIT;		prev += 36;		out += 2;		for (sb = 2; sb < sb_max; sb += 2, out += 2, prev += 36) {			dct36(in, prev, win[bt2], out);			in += SSLIMIT;			dct36(in, prev + 18, win1[bt2], out + 1);			in += SSLIMIT;		}	}	for (; sb < SBLIMIT; sb++) {		O(0) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(1) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(2) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(3) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(4) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(5) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(6) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(7) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(8) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(9) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(10) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(11) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(12) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(13) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(14) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(15) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(16) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		O(17) = *prev;		*prev++ = (MPEGIMDCT_BLOCK_TYPE) 0;		out++;	}	return(0);}/****************************************************************************/#else /* !NEW_IMDCT *//****************************************************************************/#ifndef PI#define PI	3.14159265358979323846#endif#ifdef IMDCT_INT#define IMDCT_IO_TYPE	MPEGAUD_FRACT_TYPE#define IMDCT_TYPE	INT32#define IMDCT_BITS	14#define WIN_TYPE	INT16#define WIN_BITS	14#define WIN_MULT(t, w)	(((t) * (w))>>WIN_BITS)#else#define IMDCT_IO_TYPE	REAL#define IMDCT_TYPE	REAL#define IMDCT_BITS	0#define WIN_TYPE	REAL#define WIN_BITS	0#define WIN_MULT(t, w)	((t) * (w))#endif/****************************************************************************/static int imdct_l(IMDCT_IO_TYPE *x, IMDCT_IO_TYPE *out,	MPEGIMDCT_BLOCK_TYPE *prev, WIN_TYPE *win)