dec_util.c

ffmpeg源码分析

/*
 *===================================================================
 *  3GPP AMR Wideband Floating-point Speech Codec
 *===================================================================
 */
#include <math.h>
#include <memory.h>
#include "typedef.h"
#include "dec_main.h"
#include "dec_lpc.h"

#define MAX_16       (Word16)0x7FFF
#define MIN_16       (Word16)0x8000
#define L_SUBFR      64       /* Subframe size                    */
#define L_SUBFR16k   80       /* Subframe size at 16kHz           */
#define M16k         20       /* Order of LP filter               */
#define PREEMPH_FAC  22282    /* preemphasis factor (0.68 in Q15) */
#define FAC4         4
#define FAC5         5
#define UP_FAC       20480    /* 5/4 in Q14                       */
#define INV_FAC5     6554     /* 1/5 in Q15                       */
#define NB_COEF_UP   12
#define L_FIR        31
#define MODE_7k      0
#define MODE_24k     8


extern const Word16 D_ROM_pow2[];
extern const Word16 D_ROM_isqrt[];
extern const Word16 D_ROM_log2[];
extern const Word16 D_ROM_fir_up[];
extern const Word16 D_ROM_fir_6k_7k[];
extern const Word16 D_ROM_fir_7k[];
extern const Word16 D_ROM_hp_gain[];

#ifdef WIN32
#pragma warning( disable : 4310)
#endif
/*
 * D_UTIL_random
 *
 * Parameters:
 *    seed        I/O: seed for random number
 *
 * Function:
 *    Signed 16 bits random generator.
 *
 * Returns:
 *    random number
 */
Word16 D_UTIL_random(Word16 *seed)
{
   /*static Word16 seed = 21845;*/
   *seed = (Word16)(*seed * 31821L + 13849L);
   return(*seed);
}


/*
 * D_UTIL_pow2
 *
 * Parameters:
 *    exponant    I: (Q0) Integer part.      (range: 0 <= val <= 30)
 *    fraction    I: (Q15) Fractionnal part. (range: 0.0 <= val < 1.0)
 *
 * Function:
 *    L_x = pow(2.0, exponant.fraction)         (exponant = interger part)
 *        = pow(2.0, 0.fraction) << exponant
 *
 *    Algorithm:
 *
 *    The function Pow2(L_x) is approximated by a table and linear
 *    interpolation.
 *
 *    1 - i = bit10 - b15 of fraction,   0 <= i <= 31
 *    2 - a = bit0 - b9   of fraction
 *    3 - L_x = table[i] << 16 - (table[i] - table[i + 1]) * a * 2
 *    4 - L_x = L_x >> (30-exponant)     (with rounding)
 *
 * Returns:
 *    range 0 <= val <= 0x7fffffff
 */
Word32 D_UTIL_pow2(Word16 exponant, Word16 fraction)
{
	Word32 L_x, tmp, i, exp;
	Word16 a;

	L_x = fraction * 32;          /* L_x = fraction<<6             */
	i = L_x >> 15;                /* Extract b10-b16 of fraction   */
	a = (Word16)(L_x);            /* Extract b0-b9   of fraction   */
	a = (Word16)(a & (Word16)0x7fff);
	L_x = D_ROM_pow2[i] << 16;    /* table[i] << 16                */
	tmp = D_ROM_pow2[i] - D_ROM_pow2[i + 1];  /* table[i] - table[i+1] */
	tmp = L_x - ((tmp * a) << 1); /* L_x -= tmp*a*2                */
	exp = 30 - exponant;
	if (exp <= 31)
	{
		L_x = tmp >> exp;

		if ((1 << (exp - 1)) & tmp)
		{
			L_x++;
		}
	}
	else
	{
		L_x = 0;
	}

	return(L_x);
}


/*
 * D_UTIL_norm_l
 *
 * Parameters:
 *    L_var1      I: 32 bit Word32 signed integer (Word32) whose value
 *                   falls in the range 0x8000 0000 <= var1 <= 0x7fff ffff.
 *
 * Function:
 *    Produces the number of left shifts needed to normalize the 32 bit
 *    variable L_var1 for positive values on the interval with minimum of
 *    1073741824 and maximum of 2147483647, and for negative values on
 *    the interval with minimum of -2147483648 and maximum of -1073741824;
 *    in order to normalize the result, the following operation must be done :
 *    norm_L_var1 = L_shl(L_var1,norm_l(L_var1)).
 *
 * Returns:
 *    16 bit Word16 signed integer (Word16) whose value falls in the range
 *    0x0000 0000 <= var_out <= 0x0000 001f.
 */
Word16 D_UTIL_norm_l(Word32 L_var1)
{
   Word16 var_out;

   if(L_var1 == 0)
   {
      var_out = 0;
   }
   else
   {
      if(L_var1 == (Word32)0xffffffffL)
      {
         var_out = 31;
      }
      else
      {
         if(L_var1 < 0)
         {
            L_var1 = ~L_var1;
         }

         for(var_out = 0; L_var1 < (Word32)0x40000000L; var_out++)
         {
            L_var1 <<= 1;
         }
      }
   }

   return(var_out);
}


/*
 * D_UTIL_norm_s
 *
 * Parameters:
 *    L_var1      I: 32 bit Word32 signed integer (Word32) whose value
 *                   falls in the range 0xffff 8000 <= var1 <= 0x0000 7fff.
 *
 * Function:
 *    Produces the number of left shift needed to normalize the 16 bit
 *    variable var1 for positive values on the interval with minimum
 *    of 16384 and maximum of 32767, and for negative values on
 *    the interval with minimum of -32768 and maximum of -16384.
 *
 * Returns:
 *    16 bit Word16 signed integer (Word16) whose value falls in the range
 *    0x0000 0000 <= var_out <= 0x0000 000f.
 */
Word16 D_UTIL_norm_s(Word16 var1)
{
   Word16 var_out;

   if(var1 == 0)
   {
      var_out = 0;
   }
   else
   {
      if(var1 == -1)
      {
         var_out = 15;
      }
      else
      {
         if(var1 < 0)
         {
            var1 = (Word16)~var1;
         }

         for(var_out = 0; var1 < 0x4000; var_out++)
         {
            var1 <<= 1;
         }
      }
   }
   return(var_out);
}


/*
 * D_UTIL_dot_product12
 *
 * Parameters:
 *    x        I: 12bit x vector
 *    y        I: 12bit y vector
 *    lg       I: vector length
 *    exp      O: exponent of result (0..+30)
 *
 * Function:
 *    Compute scalar product of <x[],y[]> using accumulator.
 *    The result is normalized (in Q31) with exponent (0..30).
 *
 * Returns:
 *    Q31 normalised result (1 < val <= -1)
 */
Word32 D_UTIL_dot_product12(Word16 x[], Word16 y[], Word16 lg, Word16 *exp)
{
   Word32 sum, i, sft;

   sum = 0L;

   for(i = 0; i < lg; i++)
   {
      sum += x[i] * y[i];
   }
   sum = (sum << 1) + 1;

   /* Normalize acc in Q31 */
   sft = D_UTIL_norm_l(sum);
   sum = sum << sft;
   *exp = (Word16)(30 - sft);   /* exponent = 0..30 */

   return(sum);
}


/*
 * D_UTIL_normalised_inverse_sqrt
 *
 * Parameters:
 *    frac     I/O: (Q31) normalized value (1.0 < frac <= 0.5)
 *    exp      I/O: exponent (value = frac x 2^exponent)
 *
 * Function:
 *    Compute 1/sqrt(value).
 *    If value is negative or zero, result is 1 (frac=7fffffff, exp=0).
 *
 *    The function 1/sqrt(value) is approximated by a table and linear
 *    interpolation.
 *    1. If exponant is odd then shift fraction right once.
 *    2. exponant = -((exponant - 1) >> 1)
 *    3. i = bit25 - b30 of fraction, 16 <= i <= 63 ->because of normalization.
 *    4. a = bit10 - b24
 *    5. i -= 16
 *    6. fraction = table[i]<<16 - (table[i] - table[i+1]) * a * 2
 *
 * Returns:
 *    void
 */
void D_UTIL_normalised_inverse_sqrt(Word32 *frac, Word16 *exp)
{
   Word32 i, tmp;
   Word16 a;

   if(*frac <= (Word32)0)
   {
      *exp = 0;
      *frac = 0x7fffffffL;
      return;
   }

   if((*exp & 0x1) == 1)   /* If exponant odd -> shift right */
   {
      *frac = *frac >> 1;
   }
   *exp = (Word16)(-((*exp - 1) >> 1));
   *frac = *frac >> 9;
   i = *frac >>16;      /* Extract b25-b31   */
   *frac = *frac >> 1;
   a = (Word16)(*frac); /* Extract b10-b24   */
   a = (Word16)(a & (Word16)0x7fff);
   i = i - 16;
   *frac = D_ROM_isqrt[i] << 16; /* table[i] << 16    */
   tmp = D_ROM_isqrt[i] - D_ROM_isqrt[i + 1];   /* table[i] - table[i+1]) */
   *frac = *frac - ((tmp * a) << 1);   /* frac -=  tmp*a*2  */

   return;
}


/*
 * D_UTIL_inverse_sqrt
 *
 * Parameters:
 *    L_x     I/O: (Q0) input value (range: 0<=val<=7fffffff)
 *
 * Function:
 *    Compute 1/sqrt(L_x).
 *    If value is negative or zero, result is 1 (7fffffff).
 *
 *    The function 1/sqrt(value) is approximated by a table and linear
 *    interpolation.
 *    1. Normalization of L_x
 *    2. call Normalised_Inverse_sqrt(L_x, exponant)
 *    3. L_y = L_x << exponant
 *
 * Returns:
 *    (Q31) output value (range: 0 <= val < 1)
 */
Word32 D_UTIL_inverse_sqrt(Word32 L_x)
{
   Word32 L_y;
   Word16 exp;

   exp = D_UTIL_norm_l(L_x);
   L_x = (L_x << exp);   /* L_x is normalized */
   exp = (Word16)(31 - exp);
   D_UTIL_normalised_inverse_sqrt(&L_x, &exp);

   if(exp < 0)
   {
      L_y = (L_x >> -exp);   /* denormalization   */
   }
   else
   {
      L_y = (L_x << exp);   /* denormalization   */
   }

   return(L_y);
}


/*
 * D_UTIL_normalised_log2
 *
 * Parameters:
 *    L_x      I: input value (normalized)
 *    exp      I: norm_l (L_x)
 *    exponent O: Integer part of Log2.   (range: 0<=val<=30)
 *    fraction O: Fractional part of Log2. (range: 0<=val<1)
 *
 * Function:
 *    Computes log2(L_x, exp),  where   L_x is positive and
 *    normalized, and exp is the normalisation exponent
 *    If L_x is negative or zero, the result is 0.
 *
 *    The function Log2(L_x) is approximated by a table and linear
 *    interpolation. The following steps are used to compute Log2(L_x)
 *
 *    1. exponent = 30 - norm_exponent
 *    2. i = bit25 - b31 of L_x;  32 <= i <= 63  (because of normalization).
 *    3. a = bit10 - b24
 *    4. i -= 32
 *    5. fraction = table[i] << 16 - (table[i] - table[i + 1]) * a * 2
 *
 *
 * Returns:
 *    void
 */
static void D_UTIL_normalised_log2(Word32 L_x, Word16 exp, Word16 *exponent,
                                   Word16 *fraction)
{
   Word32 i, a, tmp;
   Word32 L_y;

   if (L_x <= 0)
   {
      *exponent = 0;
      *fraction = 0;
      return;
   }

   *exponent = (Word16)(30 - exp);

   L_x = L_x >> 10;
   i = L_x >> 15;         /* Extract b25-b31               */
   a = L_x;               /* Extract b10-b24 of fraction   */
   a = a & 0x00007fff;
   i = i - 32;
   L_y = D_ROM_log2[i] << 16;               /* table[i] << 16        */
   tmp = D_ROM_log2[i] - D_ROM_log2[i + 1]; /* table[i] - table[i+1] */
   L_y = L_y - ((tmp * a) << 1);            /* L_y -= tmp*a*2        */
   *fraction = (Word16)(L_y >> 16);

   return;
}


/*
 * D_UTIL_log2
 *
 * Parameters:
 *    L_x      I: input value
 *    exponent O: Integer part of Log2.   (range: 0<=val<=30)
 *    fraction O: Fractional part of Log2. (range: 0<=val<1)
 *
 * Function:
 *    Computes log2(L_x),  where   L_x is positive.
 *    If L_x is negative or zero, the result is 0.
 *
 * Returns:
 *    void
 */
void D_UTIL_log2(Word32 L_x, Word16 *exponent, Word16 *fraction)
{
   Word16 exp;

   exp = D_UTIL_norm_l(L_x);
   D_UTIL_normalised_log2((L_x <<exp), exp, exponent, fraction);
}


/*
 * D_UTIL_l_extract
 *
 * Parameters:
 *    L_32        I: 32 bit integer.
 *    hi          O: b16 to b31 of L_32
 *    lo          O: (L_32 - hi<<16)>>1
 *
 * Function:
 *    Extract from a 32 bit integer two 16 bit DPF.
 *
 * Returns:
 *    void
 */
void D_UTIL_l_extract(Word32 L_32, Word16 *hi, Word16 *lo)
{
   *hi = (Word16)(L_32 >> 16);
   *lo = (Word16)((L_32 >> 1) - (*hi * 32768));

   return;
}


/*
 * D_UTIL_mpy_32_16
 *
 * Parameters:
 *    hi          I: hi part of 32 bit number
 *    lo          I: lo part of 32 bit number
 *    n           I: 16 bit number
 *
 * Function:
 *    Multiply a 16 bit integer by a 32 bit (DPF). The result is divided
 *    by 2^15.
 *
 *    L_32 = (hi1*lo2)<<1 + ((lo1*lo2)>>15)<<1