/* ***** BEGIN LICENSE BLOCK *****  
 * Source last modified: $Id: sbrqmf.c,v 1.2 2005/05/19 20:45:20 jrecker Exp $ 
 *   
 * Portions Copyright (c) 1995-2005 RealNetworks, Inc. All Rights Reserved.  
 *       
 * The contents of this file, and the files included with this file, 
 * are subject to the current version of the RealNetworks Public 
 * Source License (the "RPSL") available at 
 * http://www.helixcommunity.org/content/rpsl unless you have licensed 
 * the file under the current version of the RealNetworks Community 
 * Source License (the "RCSL") available at 
 * http://www.helixcommunity.org/content/rcsl, in which case the RCSL 
 * will apply. You may also obtain the license terms directly from 
 * RealNetworks.  You may not use this file except in compliance with 
 * the RPSL or, if you have a valid RCSL with RealNetworks applicable 
 * to this file, the RCSL.  Please see the applicable RPSL or RCSL for 
 * the rights, obligations and limitations governing use of the 
 * contents of the file. 
 *   
 * This file is part of the Helix DNA Technology. RealNetworks is the 
 * developer of the Original Code and owns the copyrights in the 
 * portions it created. 
 *   
 * This file, and the files included with this file, is distributed 
 * and made available on an 'AS IS' basis, WITHOUT WARRANTY OF ANY 
 * KIND, EITHER EXPRESS OR IMPLIED, AND REALNETWORKS HEREBY DISCLAIMS 
 * ALL SUCH WARRANTIES, INCLUDING WITHOUT LIMITATION, ANY WARRANTIES 
 * OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE, QUIET 
 * ENJOYMENT OR NON-INFRINGEMENT. 
 *  
 * Technology Compatibility Kit Test Suite(s) Location:  
 *    http://www.helixcommunity.org/content/tck  
 *  
 * Contributor(s):  
 *   
 * ***** END LICENSE BLOCK ***** */  

/**************************************************************************************
 * Fixed-point HE-AAC decoder
 * Jon Recker (jrecker@real.com)
 * February 2005
 *
 * sbrqmf.c - analysis and synthesis QMF filters for SBR
 **************************************************************************************/

#include "sbr.h"
#include "assembly.h"

/* PreMultiply64() table
 * format = Q30
 * reordered for sequential access
 *
 * for (i = 0; i < 64/4; i++) {
 *   angle = (i + 0.25) * M_PI / nmdct;
 *   x = (cos(angle) + sin(angle));
 *   x =  sin(angle);
 * 
 *   angle = (nmdct/2 - 1 - i + 0.25) * M_PI / nmdct;
 *   x = (cos(angle) + sin(angle));
 *   x =  sin(angle);
 * }
 */
static const int cos4sin4tab64[64] PROGMEM = {
	0x40c7d2bd, 0x00c90e90, 0x424ff28f, 0x3ff4e5e0, 0x43cdd89a, 0x03ecadcf, 0x454149fc, 0x3fc395f9,
	0x46aa0d6d, 0x070de172, 0x4807eb4b, 0x3f6af2e3, 0x495aada2, 0x0a2abb59, 0x4aa22036, 0x3eeb3347,
	0x4bde1089, 0x0d415013, 0x4d0e4de2, 0x3e44a5ef, 0x4e32a956, 0x104fb80e, 0x4f4af5d1, 0x3d77b192,
	0x50570819, 0x135410c3, 0x5156b6d9, 0x3c84d496, 0x5249daa2, 0x164c7ddd, 0x53304df6, 0x3b6ca4c4,
	0x5409ed4b, 0x19372a64, 0x54d69714, 0x3a2fcee8, 0x55962bc0, 0x1c1249d8, 0x56488dc5, 0x38cf1669,
	0x56eda1a0, 0x1edc1953, 0x57854ddd, 0x374b54ce, 0x580f7b19, 0x2192e09b, 0x588c1404, 0x35a5793c,
	0x58fb0568, 0x2434f332, 0x595c3e2a, 0x33de87de, 0x59afaf4c, 0x26c0b162, 0x59f54bee, 0x31f79948,
	0x5a2d0957, 0x29348937, 0x5a56deec, 0x2ff1d9c7, 0x5a72c63b, 0x2b8ef77d, 0x5a80baf6, 0x2dce88aa,
};

/* PostMultiply64() table
 * format = Q30
 * reordered for sequential access
 *
 * for (i = 0; i <= (32/2); i++) {
 *   angle = i * M_PI / 64;
 *   x = (cos(angle) + sin(angle));
 *   x = sin(angle);
 * }
 */
static const int cos1sin1tab64[34] PROGMEM = {
	0x40000000, 0x00000000, 0x43103085, 0x0323ecbe, 0x45f704f7, 0x0645e9af, 0x48b2b335, 0x09640837, 
	0x4b418bbe, 0x0c7c5c1e, 0x4da1fab5, 0x0f8cfcbe, 0x4fd288dc, 0x1294062f, 0x51d1dc80, 0x158f9a76, 
	0x539eba45, 0x187de2a7, 0x553805f2, 0x1b5d100a, 0x569cc31b, 0x1e2b5d38, 0x57cc15bc, 0x20e70f32, 
	0x58c542c5, 0x238e7673, 0x5987b08a, 0x261feffa, 0x5a12e720, 0x2899e64a, 0x5a6690ae, 0x2afad269, 
	0x5a82799a, 0x2d413ccd,
};

/**************************************************************************************
 * Function:    PreMultiply64
 *
 * Description: pre-twiddle stage of 64-point DCT-IV
 *
 * Inputs:      buffer of 64 samples
 *
 * Outputs:     processed samples in same buffer
 *
 * Return:      none
 *
 * Notes:       minimum 1 GB in, 2 GB out, gains 2 int bits
 *              gbOut = gbIn + 1
 *              output is limited to sqrt(2)/2 plus GB in full GB
 *              uses 3-mul, 3-add butterflies instead of 4-mul, 2-add
 **************************************************************************************/
static void PreMultiply64(int *zbuf1)
{
	int i, ar1, ai1, ar2, ai2, z1, z2;
	int t, cms2, cps2a, sin2a, cps2b, sin2b;
	int *zbuf2;
	const int *csptr;

	zbuf2 = zbuf1 + 64 - 1;
	csptr = cos4sin4tab64;

	/* whole thing should fit in registers - verify that compiler does this */
	for (i = 64 >> 2; i != 0; i--) {
		/* cps2 = (cos+sin), sin2 = sin, cms2 = (cos-sin) */
		cps2a = *csptr++;	
		sin2a = *csptr++;
		cps2b = *csptr++;	
		sin2b = *csptr++;

		ar1 = *(zbuf1 + 0);
		ai2 = *(zbuf1 + 1);
		ai1 = *(zbuf2 + 0);
		ar2 = *(zbuf2 - 1);

		/* gain 2 ints bit from MULSHIFT32 by Q30
		 * max per-sample gain (ignoring implicit scaling) = MAX(sin(angle)+cos(angle)) = 1.414
		 * i.e. gain 1 GB since worst case is sin(angle) = cos(angle) = 0.707 (Q30), gain 2 from
		 *   extra sign bits, and eat one in adding
		 */
		t  = MULSHIFT32(sin2a, ar1 + ai1);
		z2 = MULSHIFT32(cps2a, ai1) - t;
		cms2 = cps2a - 2*sin2a;
		z1 = MULSHIFT32(cms2, ar1) + t;
		*zbuf1++ = z1;	/* cos*ar1 + sin*ai1 */
		*zbuf1++ = z2;	/* cos*ai1 - sin*ar1 */

		t  = MULSHIFT32(sin2b, ar2 + ai2);
		z2 = MULSHIFT32(cps2b, ai2) - t;
		cms2 = cps2b - 2*sin2b;
		z1 = MULSHIFT32(cms2, ar2) + t;
		*zbuf2-- = z2;	/* cos*ai2 - sin*ar2 */
		*zbuf2-- = z1;	/* cos*ar2 + sin*ai2 */
	}
}

/**************************************************************************************
 * Function:    PostMultiply64
 *
 * Description: post-twiddle stage of 64-point type-IV DCT
 *
 * Inputs:      buffer of 64 samples
 *              number of output samples to calculate
 *
 * Outputs:     processed samples in same buffer
 *
 * Return:      none
 *
 * Notes:       minimum 1 GB in, 2 GB out, gains 2 int bits
 *              gbOut = gbIn + 1
 *              output is limited to sqrt(2)/2 plus GB in full GB
 *              nSampsOut is rounded up to next multiple of 4, since we calculate
 *                4 samples per loop
 **************************************************************************************/
static void PostMultiply64(int *fft1, int nSampsOut)
{
	int i, ar1, ai1, ar2, ai2;
	int t, cms2, cps2, sin2;
	int *fft2;
	const int *csptr;

	csptr = cos1sin1tab64;
	fft2 = fft1 + 64 - 1;

	/* load coeffs for first pass
	 * cps2 = (cos+sin)/2, sin2 = sin/2, cms2 = (cos-sin)/2
	 */
	cps2 = *csptr++;
	sin2 = *csptr++;
	cms2 = cps2 - 2*sin2;

	for (i = (nSampsOut + 3) >> 2; i != 0; i--) {
		ar1 = *(fft1 + 0);
		ai1 = *(fft1 + 1);
		ar2 = *(fft2 - 1);
		ai2 = *(fft2 + 0);

		/* gain 2 int bits (multiplying by Q30), max gain = sqrt(2) */
		t = MULSHIFT32(sin2, ar1 + ai1);
		*fft2-- = t - MULSHIFT32(cps2, ai1);
		*fft1++ = t + MULSHIFT32(cms2, ar1);

		cps2 = *csptr++;
		sin2 = *csptr++;

		ai2 = -ai2;
		t = MULSHIFT32(sin2, ar2 + ai2);
		*fft2-- = t - MULSHIFT32(cps2, ai2);
		cms2 = cps2 - 2*sin2;
		*fft1++ = t + MULSHIFT32(cms2, ar2);
	}
}

/**************************************************************************************
 * Function:    QMFAnalysisConv
 *
 * Description: convolution kernel for analysis QMF 
 *
 * Inputs:      pointer to coefficient table, reordered for sequential access
 *              delay buffer of size 32*10 = 320 real-valued PCM samples
 *              index for delay ring buffer (range = [0, 9])
 *
 * Outputs:     64 consecutive 32-bit samples
 *
 * Return:      none
 *
 * Notes:       this is carefully written to be efficient on ARM
 *              use the assembly code version in sbrqmfak.s when building for ARM!
 **************************************************************************************/
// #if (defined (__arm) && defined (__ARMCC_VERSION)) || (defined (_WIN32) && defined (_WIN32_WCE) && defined (ARM)) || (defined(__GNUC__) && defined(__arm__)) || defined(ARDUINO)
// #ifdef __cplusplus
// extern "C"
// #endif
// void QMFAnalysisConv(int *cTab, int *delay, int dIdx, int *uBuf);
// #else
void QMFAnalysisConv(int *cTab, int *delay, int dIdx, int *uBuf)
{
	int k, dOff;
	int *cPtr0, *cPtr1;
	U64 u64lo, u64hi;

	dOff = dIdx*32 + 31;
	cPtr0 = cTab;
	cPtr1 = cTab + 33*5 - 1;

	/* special first pass since we need to flip sign to create cTab[384], cTab[512] */
	u64lo.w64 = 0;
	u64hi.w64 = 0;
	u64lo.w64 = MADD64(u64lo.w64,  *cPtr0++,   delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
	u64hi.w64 = MADD64(u64hi.w64,  *cPtr0++,   delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
	u64lo.w64 = MADD64(u64lo.w64,  *cPtr0++,   delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
	u64hi.w64 = MADD64(u64hi.w64,  *cPtr0++,   delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
	u64lo.w64 = MADD64(u64lo.w64,  *cPtr0++,   delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
	u64hi.w64 = MADD64(u64hi.w64,  *cPtr1--,   delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
	u64lo.w64 = MADD64(u64lo.w64, -(*cPtr1--), delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
	u64hi.w64 = MADD64(u64hi.w64,  *cPtr1--,   delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
	u64lo.w64 = MADD64(u64lo.w64, -(*cPtr1--), delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
	u64hi.w64 = MADD64(u64hi.w64,  *cPtr1--,   delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}

	uBuf[0]  = u64lo.r.hi32;
	uBuf[32] = u64hi.r.hi32;
	uBuf++;
	dOff--;

	/* max gain for any sample in uBuf, after scaling by cTab, ~= 0.99 
	 * so we can just sum the uBuf values with no overflow problems
	 */
	for (k = 1; k <= 31; k++) {
		u64lo.w64 = 0;
		u64hi.w64 = 0;
		u64lo.w64 = MADD64(u64lo.w64, *cPtr0++, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
		u64hi.w64 = MADD64(u64hi.w64, *cPtr0++, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
		u64lo.w64 = MADD64(u64lo.w64, *cPtr0++, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
		u64hi.w64 = MADD64(u64hi.w64, *cPtr0++, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
		u64lo.w64 = MADD64(u64lo.w64, *cPtr0++, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
		u64hi.w64 = MADD64(u64hi.w64, *cPtr1--, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
		u64lo.w64 = MADD64(u64lo.w64, *cPtr1--, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
		u64hi.w64 = MADD64(u64hi.w64, *cPtr1--, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
		u64lo.w64 = MADD64(u64lo.w64, *cPtr1--, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}
		u64hi.w64 = MADD64(u64hi.w64, *cPtr1--, delay[dOff]);	dOff -= 32; if (dOff < 0) {dOff += 320;}

		uBuf[0]  = u64lo.r.hi32;
		uBuf[32] = u64hi.r.hi32;
		uBuf++;
		dOff--;
	}
}
//#endif

/**************************************************************************************
 * Function:    QMFAnalysis
 *
 * Description: 32-subband analysis QMF (4.6.18.4.1)
 *
 * Inputs:      32 consecutive samples of decoded 32-bit PCM, format = Q(fBitsIn)
 *              delay buffer of size 32*10 = 320 PCM samples
 *              number of fraction bits in input PCM
 *              index for delay ring buffer (range = [0, 9])
 *              number of subbands to calculate (range = [0, 32])
 *
 * Outputs:     qmfaBands complex subband samples, format = Q(FBITS_OUT_QMFA)
 *              updated delay buffer
 *              updated delay index
 *
 * Return:      guard bit mask
 *
 * Notes:       output stored as RE{X0}, IM{X0}, RE{X1}, IM{X1}, ... RE{X31}, IM{X31}
 *              output stored in int buffer of size 64*2 = 128 
 *                (zero-filled from XBuf[2*qmfaBands] to XBuf[127])
 **************************************************************************************/
int QMFAnalysis(int *inbuf, int *delay, int *XBuf, int fBitsIn, int *delayIdx, int qmfaBands)
{
	int n, y, shift, gbMask;
	int *delayPtr, *uBuf, *tBuf;

	/* use XBuf[128] as temp buffer for reordering */
	uBuf = XBuf;		/* first 64 samples */
	tBuf = XBuf + 64;	/* second 64 samples */

	/* overwrite oldest PCM with new PCM
	 * delay[n] has 1 GB after shifting (either << or >>)
	 */
	delayPtr = delay + (*delayIdx * 32);
	if (fBitsIn > FBITS_IN_QMFA) {
		shift = MIN(fBitsIn - FBITS_IN_QMFA, 31);
		for (n = 32; n != 0; n--) {
			y = (*inbuf) >> shift;
			inbuf++;
			*delayPtr++ = y;
		}
	} else {
		shift = MIN(FBITS_IN_QMFA - fBitsIn, 30);
		for (n = 32; n != 0; n--) {
			y = *inbuf++;
			CLIP_2N_SHIFT30(y, shift);
			*delayPtr++ = y;
		}
	}
	
	QMFAnalysisConv((int *)cTabA, delay, *delayIdx, uBuf);
	
	/* uBuf has at least 2 GB right now (1 from clipping to Q(FBITS_IN_QMFA), one from
	 *   the scaling by cTab (MULSHIFT32(*delayPtr--, *cPtr++), with net gain of < 1.0)
	 * TODO - fuse with QMFAnalysisConv to avoid separate reordering
	 */
    tBuf[2*0 + 0] = uBuf[0];
    tBuf[2*0 + 1] = uBuf[1];
    for (n = 1; n < 31; n++) {
        tBuf[2*n + 0] = -uBuf[64-n];
        tBuf[2*n + 1] =  uBuf[n+1];
    }
    tBuf[2*31 + 1] =  uBuf[32];
    tBuf[2*31 + 0] = -uBuf[33];
	
	/* fast in-place DCT-IV - only need 2*qmfaBands output samples */
	PreMultiply64(tBuf);	/* 2 GB in, 3 GB out */
	FFT32C(tBuf);			/* 3 GB in, 1 GB out */
	PostMultiply64(tBuf, qmfaBands*2);	/* 1 GB in, 2 GB out */

	/* TODO - roll into PostMultiply (if enough registers) */
	gbMask = 0;
	for (n = 0; n < qmfaBands; n++) {
		XBuf[2*n+0] =  tBuf[ n + 0];	/* implicit scaling of 2 in our output Q format */
		gbMask |= FASTABS(XBuf[2*n+0]);
		XBuf[2*n+1] = -tBuf[63 - n];
		gbMask |= FASTABS(XBuf[2*n+1]);
	}

	/* fill top section with zeros for HF generation */
	for (    ; n < 64; n++) {
		XBuf[2*n+0] = 0;
		XBuf[2*n+1] = 0;
	}

	*delayIdx = (*delayIdx == NUM_QMF_DELAY_BUFS - 1 ? 0 : *delayIdx + 1);

	/* minimum of 2 GB in output */
	return gbMask;
}

/* lose FBITS_LOST_DCT4_64 in DCT4, gain 6 for implicit scaling by 1/64, lose 1 for cTab multiply (Q31) */
#define FBITS_OUT_QMFS	(FBITS_IN_QMFS - FBITS_LOST_DCT4_64 + 6 - 1)
#define RND_VAL			(1 << (FBITS_OUT_QMFS-1))

/**************************************************************************************
 * Function:    QMFSynthesisConv
 *
 * Description: final convolution kernel for synthesis QMF 
 *
 * Inputs:      pointer to coefficient table, reordered for sequential access
 *              delay buffer of size 64*10 = 640 complex samples (1280 ints)
 *              index for delay ring buffer (range = [0, 9])
 *              number of QMF subbands to process (range = [0, 64])
 *              number of channels
 *
 * Outputs:     64 consecutive 16-bit PCM samples, interleaved by factor of nChans
 *
 * Return:      none
 *
 * Notes:       this is carefully written to be efficient on ARM
 *              use the assembly code version in sbrqmfsk.s when building for ARM!
 **************************************************************************************/
// #if (defined (__arm) && defined (__ARMCC_VERSION)) || (defined (_WIN32) && defined (_WIN32_WCE) && defined (ARM)) || (defined(__GNUC__) && defined(__arm__)) || (defined(ARDUINO))
// #ifdef __cplusplus
// extern "C"
// #endif
// void QMFSynthesisConv(int *cPtr, int *delay, int dIdx, short *outbuf, int nChans);
// #else
void QMFSynthesisConv(int *cPtr, int *delay, int dIdx, short *outbuf, int nChans)
{
	int k, dOff0, dOff1;
	U64 sum64;

	dOff0 = (dIdx)*128;
	dOff1 = dOff0 - 1;
	if (dOff1 < 0)
		dOff1 += 1280;

	/* scaling note: total gain of coefs (cPtr[0]-cPtr[9] for any k) is < 2.0, so 1 GB in delay values is adequate */
	for (k = 0; k <= 63; k++) {
		sum64.w64 = 0;
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff0]);	dOff0 -= 256; if (dOff0 < 0) {dOff0 += 1280;}
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff1]);	dOff1 -= 256; if (dOff1 < 0) {dOff1 += 1280;}
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff0]);	dOff0 -= 256; if (dOff0 < 0) {dOff0 += 1280;}
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff1]);	dOff1 -= 256; if (dOff1 < 0) {dOff1 += 1280;}
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff0]);	dOff0 -= 256; if (dOff0 < 0) {dOff0 += 1280;}
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff1]);	dOff1 -= 256; if (dOff1 < 0) {dOff1 += 1280;}
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff0]);	dOff0 -= 256; if (dOff0 < 0) {dOff0 += 1280;}
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff1]);	dOff1 -= 256; if (dOff1 < 0) {dOff1 += 1280;}
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff0]);	dOff0 -= 256; if (dOff0 < 0) {dOff0 += 1280;}
		sum64.w64 = MADD64(sum64.w64, *cPtr++, delay[dOff1]);	dOff1 -= 256; if (dOff1 < 0) {dOff1 += 1280;}

		dOff0++;
		dOff1--;
		*outbuf = CLIPTOSHORT((sum64.r.hi32 + RND_VAL) >> FBITS_OUT_QMFS);
		outbuf += nChans;
	}
}
//#endif

/**************************************************************************************
 * Function:    QMFSynthesis
 *
 * Description: 64-subband synthesis QMF (4.6.18.4.2)
 *
 * Inputs:      64 consecutive complex subband QMF samples, format = Q(FBITS_IN_QMFS)
 *              delay buffer of size 64*10 = 640 complex samples (1280 ints)
 *              index for delay ring buffer (range = [0, 9])
 *              number of QMF subbands to process (range = [0, 64])
 *              number of channels
 *
 * Outputs:     64 consecutive 16-bit PCM samples, interleaved by factor of nChans
 *              updated delay buffer
 *              updated delay index
 *
 * Return:      none
 *
 * Notes:       assumes MIN_GBITS_IN_QMFS guard bits in input, either from
 *                QMFAnalysis (if upsampling only) or from MapHF (if SBR on)
 **************************************************************************************/
void QMFSynthesis(int *inbuf, int *delay, int *delayIdx, int qmfsBands, short *outbuf, int nChans)
{
	int n, a0, a1, b0, b1, dOff0, dOff1, dIdx;
	int *tBufLo, *tBufHi;

	dIdx = *delayIdx;
	tBufLo = delay + dIdx*128 + 0;
	tBufHi = delay + dIdx*128 + 127;

	/* reorder inputs to DCT-IV, only use first qmfsBands (complex) samples 
	 * TODO - fuse with PreMultiply64 to avoid separate reordering steps
	 */
    for (n = 0; n < qmfsBands >> 1; n++) {
		a0 = *inbuf++;
		b0 = *inbuf++;
		a1 = *inbuf++;
		b1 = *inbuf++;
		*tBufLo++ = a0;
        *tBufLo++ = a1;
        *tBufHi-- = b0;
        *tBufHi-- = b1;
    }
	if (qmfsBands & 0x01) {
		a0 = *inbuf++;
		b0 = *inbuf++;
		*tBufLo++ = a0;
        *tBufHi-- = b0;
        *tBufLo++ = 0;
		*tBufHi-- = 0;
		n++;
	}
    for (     ; n < 32; n++) {
		*tBufLo++ = 0;
        *tBufHi-- = 0;
        *tBufLo++ = 0;
        *tBufHi-- = 0;
	}

	tBufLo = delay + dIdx*128 + 0;
	tBufHi = delay + dIdx*128 + 64;

	/* 2 GB in, 3 GB out */
	PreMultiply64(tBufLo);
	PreMultiply64(tBufHi);

	/* 3 GB in, 1 GB out */
	FFT32C(tBufLo);
	FFT32C(tBufHi);

	/* 1 GB in, 2 GB out */
	PostMultiply64(tBufLo, 64);
	PostMultiply64(tBufHi, 64);

	/* could fuse with PostMultiply64 to avoid separate pass */
	dOff0 = dIdx*128;
	dOff1 = dIdx*128 + 64;
	for (n = 32; n != 0; n--) {
		a0 =  (*tBufLo++);
		a1 =  (*tBufLo++);
		b0 =  (*tBufHi++);
		b1 = -(*tBufHi++);

		delay[dOff0++] = (b0 - a0);
		delay[dOff0++] = (b1 - a1);
		delay[dOff1++] = (b0 + a0);
		delay[dOff1++] = (b1 + a1);
	}

	QMFSynthesisConv((int *)cTabS, delay, dIdx, outbuf, nChans);

	*delayIdx = (*delayIdx == NUM_QMF_DELAY_BUFS - 1 ? 0 : *delayIdx + 1);
}
