sbrfreq.c

nds上大名鼎鼎的DSOrganize 2.8最新源代码。很值得研究。

/* ***** BEGIN LICENSE BLOCK *****  
 * Source last modified: $Id: sbrfreq.c,v 1.1.2.1.4.1 2005/05/20 23:40:55 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 ([e-mail])
 * February 2005
 *
 * sbrfreq.c - frequency band table calculation for SBR
 **************************************************************************************/

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

/**************************************************************************************
 * Function:    BubbleSort
 *
 * Description: in-place sort of unsigned chars
 *
 * Inputs:      buffer of elements to sort
 *              number of elements to sort
 *
 * Outputs:     sorted buffer
 *
 * Return:      none
 **************************************************************************************/
static void BubbleSort(unsigned char *v, int nItems)
{
	int i;
	unsigned char t;

	while (nItems >= 2) {
		for (i = 0; i < nItems-1; i++) {
			if (v[i+1] < v[i]) {
				t = v[i+1];	
				v[i+1] = v[i];	
				v[i] = t;
			}
		}
		nItems--;
	}
}

/**************************************************************************************
 * Function:    VMin
 *
 * Description: find smallest element in a buffer of unsigned chars
 *
 * Inputs:      buffer of elements to search
 *              number of elements to search
 *
 * Outputs:     none
 *
 * Return:      smallest element in buffer
 **************************************************************************************/
static unsigned char VMin(unsigned char *v, int nItems)
{
	int i;
	unsigned char vMin;

	vMin = v[0];
	for (i = 1; i < nItems; i++) {
		if (v[i] < vMin)
			vMin = v[i];
	}
	return vMin;
}

/**************************************************************************************
 * Function:    VMax
 *
 * Description: find largest element in a buffer of unsigned chars
 *
 * Inputs:      buffer of elements to search
 *              number of elements to search
 *
 * Outputs:     none
 *
 * Return:      largest element in buffer
 **************************************************************************************/
static unsigned char VMax(unsigned char *v, int nItems)
{
	int i;
	unsigned char vMax;

	vMax = v[0];
	for (i = 1; i < nItems; i++) {
		if (v[i] > vMax)
			vMax = v[i];
	}
	return vMax;
}

/**************************************************************************************
 * Function:    CalcFreqMasterScaleZero
 *
 * Description: calculate master frequency table when freqScale == 0 
 *                (4.6.18.3.2.1, figure 4.39)
 *
 * Inputs:      alterScale flag
 *              index of first QMF subband in master freq table (k0)
 *              index of last QMF subband (k2)
 *
 * Outputs:     master frequency table
 *
 * Return:      number of bands in master frequency table
 *
 * Notes:       assumes k2 - k0 <= 48 and k2 >= k0 (4.6.18.3.6)
 **************************************************************************************/
static int CalcFreqMasterScaleZero(unsigned char *freqMaster, int alterScale, int k0, int k2)
{
	int nMaster, k, nBands, k2Achieved, dk, vDk[64], k2Diff;

	if (alterScale) {
		dk = 2;
		nBands = 2 * ((k2 - k0 + 2) >> 2);
	} else {
		dk = 1;
		nBands = 2 * ((k2 - k0) >> 1);
	}

	if (nBands <= 0)
		return 0;

	k2Achieved = k0 + nBands * dk;
	k2Diff = k2 - k2Achieved;
	for (k = 0; k < nBands; k++)
		vDk[k] = dk;

	if (k2Diff > 0) {
		k = nBands - 1;
		while (k2Diff) {
			vDk[k]++;
			k--;
			k2Diff--;
		}
	} else if (k2Diff < 0) {
		k = 0;
		while (k2Diff) {
			vDk[k]--;
			k++;
			k2Diff++;
		}
	}

	nMaster = nBands;
	freqMaster[0] = k0;
	for (k = 1; k <= nBands; k++)
		freqMaster[k] = freqMaster[k-1] + vDk[k-1];

	return nMaster;
}

/* mBandTab[i] = temp1[i] / 2 */
static const int mBandTab[3] = {6, 5, 4};

/* invWarpTab[i] = 1.0 / temp2[i], Q30 (see 4.6.18.3.2.1) */
static const int invWarpTab[2] = {0x40000000, 0x313b13b1};

/**************************************************************************************
 * Function:    CalcFreqMasterScale
 *
 * Description: calculate master frequency table when freqScale > 0 
 *                (4.6.18.3.2.1, figure 4.39)
 *
 * Inputs:      alterScale flag
 *              freqScale flag
 *              index of first QMF subband in master freq table (k0)
 *              index of last QMF subband (k2)
 *
 * Outputs:     master frequency table
 *
 * Return:      number of bands in master frequency table
 *
 * Notes:       assumes k2 - k0 <= 48 and k2 >= k0 (4.6.18.3.6)
 **************************************************************************************/
static int CalcFreqMaster(unsigned char *freqMaster, int freqScale, int alterScale, int k0, int k2)
{
	int bands, twoRegions, k, k1, t, vLast, vCurr, pCurr;
	int invWarp, nBands0, nBands1, change;
	unsigned char vDk1Min, vDk0Max;
	unsigned char *vDelta;

	if (freqScale < 1 || freqScale > 3)
		return -1;

	bands = mBandTab[freqScale - 1];
	invWarp = invWarpTab[alterScale];

	/* tested for all k0 = [5, 64], k2 = [k0, 64] */
	if (k2*10000 > 22449*k0) {
		twoRegions = 1;
		k1 = 2*k0;
	} else {
		twoRegions = 0;
		k1 = k2;
	}
	
	/* tested for all k0 = [5, 64], k1 = [k0, 64], freqScale = [1,3] */
	t = (log2Tab[k1] - log2Tab[k0]) >> 3;				/* log2(k1/k0), Q28 to Q25 */
	nBands0 = 2 * (((bands * t) + (1 << 24)) >> 25);	/* multiply by bands/2, round to nearest int (mBandTab has factor of 1/2 rolled in) */

	/* tested for all valid combinations of k0, k1, nBands (from sampRate, freqScale, alterScale) 
	 * roundoff error can be a problem with fixpt (e.g. pCurr = 12.499999 instead of 12.50003)
	 *   because successive multiplication always undershoots a little bit, but this
	 *   doesn't occur in any of the ratios we encounter from the valid k0/k1 bands in the spec
	 */
	t = RatioPowInv(k1, k0, nBands0);
	pCurr = k0 << 24;
	vLast = k0;
	vDelta = freqMaster + 1;	/* operate in-place */
	for (k = 0; k < nBands0; k++) {
		pCurr = MULSHIFT32(pCurr, t) << 8;	/* keep in Q24 */
		vCurr = (pCurr + (1 << 23)) >> 24;
		vDelta[k] = (vCurr - vLast);
		vLast = vCurr;
	}

	/* sort the deltas and find max delta for first region */
	BubbleSort(vDelta, nBands0);
	vDk0Max = VMax(vDelta, nBands0);

	/* fill master frequency table with bands from first region */
	freqMaster[0] = k0;
	for (k = 1; k <= nBands0; k++)
		freqMaster[k] += freqMaster[k-1];

	/* if only one region, then the table is complete */
	if (!twoRegions)
		return nBands0;

	/* tested for all k1 = [10, 64], k2 = [k0, 64], freqScale = [1,3] */
	t = (log2Tab[k2] - log2Tab[k1]) >> 3;		/* log2(k1/k0), Q28 to Q25 */
	t = MULSHIFT32(bands * t, invWarp) << 2;	/* multiply by bands/2, divide by warp factor, keep Q25 */
	nBands1 = 2 * ((t + (1 << 24)) >> 25);		/* round to nearest int */
				
	/* see comments above for calculations in first region */
	t = RatioPowInv(k2, k1, nBands1);
	pCurr = k1 << 24;
	vLast = k1;
	vDelta = freqMaster + nBands0 + 1;	/* operate in-place */
	for (k = 0; k < nBands1; k++) {
		pCurr = MULSHIFT32(pCurr, t) << 8;	/* keep in Q24 */
		vCurr = (pCurr + (1 << 23)) >> 24;
		vDelta[k] = (vCurr - vLast);
		vLast = vCurr;
	}

	/* sort the deltas, adjusting first and last if the second region has smaller deltas than the first */
	vDk1Min = VMin(vDelta, nBands1);
	if (vDk1Min < vDk0Max) {
		BubbleSort(vDelta, nBands1);
		change = vDk0Max - vDelta[0];
		if (change > ((vDelta[nBands1 - 1] - vDelta[0]) >> 1))
			 change = ((vDelta[nBands1 - 1] - vDelta[0]) >> 1);
		vDelta[0] += change;
		vDelta[nBands1-1] -= change;
	}
	BubbleSort(vDelta, nBands1);

	/* fill master frequency table with bands from second region 
	 * Note: freqMaster[nBands0] = k1
	 */
	for (k = 1; k <= nBands1; k++)
		freqMaster[k + nBands0] += freqMaster[k + nBands0 - 1];

	return (nBands0 + nBands1);
}

/**************************************************************************************
 * Function:    CalcFreqHigh
 *
 * Description: calculate high resolution frequency table (4.6.18.3.2.2)
 *
 * Inputs:      master frequency table
 *              number of bands in master frequency table 
 *              crossover band from header
 *
 * Outputs:     high resolution frequency table
 *
 * Return:      number of bands in high resolution frequency table
 **************************************************************************************/
static int CalcFreqHigh(unsigned char *freqHigh, unsigned char *freqMaster, int nMaster, int crossOverBand)
{
	int k, nHigh;

	nHigh = nMaster - crossOverBand;