sp_dec.c

FLOAT PINT

   for ( i = 0; i < M; i++ ) {
      /* ind = b8-b15 of lsf[i] */
      ind = lsf[i] >> 8;

      /* offset = b0-b7  of lsf[i] */
      offset = lsf[i] & 0x00ff;

      /* lsp[i] = table[ind]+ ((table[ind+1]-table[ind])*offset) / 256 */
      tmp = ( ( cos_table[ind+1]-cos_table[ind] )*offset ) << 1;
      lsp[i] = cos_table[ind] + ( tmp >> 9 );
   }
   return;
}


/*
 * D_plsf_3
 *
 *
 * Parameters:
 *    st->past_lsf_q    I: Past dequantized LFSs
 *    st->past_r_q      B: past quantized residual
 *    mode              I: AMR mode
 *    bfi               B: bad frame indicator
 *    indice            I: quantization indices of 3 submatrices, Q0
 *    lsp1_q            O: quantized 1st LSP vector
 *
 * Function:
 *    Decodes the LSP parameters using the received quantization indices.
 *    1st order MA prediction and split by 3 vector quantization (split-VQ)
 *
 * Returns:
 *    void
 */
static void D_plsf_3( D_plsfState *st, int mode, Word16 bfi, Word16 *
      indice, Word32 *lsp1_q )
{
   Word32 lsf1_r[M], lsf1_q[M];
   Word32 i, index, temp;
   const Word32 *p_cb1, *p_cb2, *p_cb3, *p_dico;


   /* if bad frame */
   if ( bfi != 0 ) {
      /* use the past LSFs slightly shifted towards their mean */
      for ( i = 0; i < M; i++ ) {
         /* lsfi_q[i] = ALPHA*past_lsf_q[i] + ONE_ALPHA*meanLsf[i]; */
         lsf1_q[i] = ( ( st->past_lsf_q[i] * ALPHA ) >> 15 ) + ( ( mean_lsf_3[i]
               * ONE_ALPHA ) >> 15 );
      }

      /* estimate past quantized residual to be used in next frame */
      if ( mode != MRDTX ) {
         for ( i = 0; i < M; i++ ) {
            /* temp  = meanLsf[i] +  pastR2_q[i] * pred_fac; */
            temp = mean_lsf_3[i] + ( ( st->past_r_q[i] * pred_fac[i] ) >> 15 );
            st->past_r_q[i] = lsf1_q[i] - temp;
         }
      }
      else {
         for ( i = 0; i < M; i++ ) {
            /* temp  = meanLsf[i] +  pastR2_q[i]; */
            temp = mean_lsf_3[i] + st->past_r_q[i];
            st->past_r_q[i] = lsf1_q[i] - temp;
         }
      }
   }

   /* if good LSFs received */
   else {
      if ( ( mode == MR475 ) | ( mode == MR515 ) ) {
         /* MR475, MR515 */
         p_cb1 = dico1_lsf_3;
         p_cb2 = dico2_lsf_3;
         p_cb3 = mr515_3_lsf;
      }
      else if ( mode == MR795 ) {
         /* MR795 */
         p_cb1 = mr795_1_lsf;
         p_cb2 = dico2_lsf_3;
         p_cb3 = dico3_lsf_3;
      }
      else {
         /* MR59, MR67, MR74, MR102, MRDTX */
         p_cb1 = dico1_lsf_3;
         p_cb2 = dico2_lsf_3;
         p_cb3 = dico3_lsf_3;
      }

      /* decode prediction residuals from 3 received indices */
      index = *indice++;
      p_dico = &p_cb1[index + index + index];
      index = *indice++;
      lsf1_r[0] = *p_dico++;
      lsf1_r[1] = *p_dico++;
      lsf1_r[2] = *p_dico++;

      if ( ( mode == MR475 ) | ( mode == MR515 ) ) {
         /* MR475, MR515 only using every second entry */
         index = index << 1;
      }
      p_dico = &p_cb2[index + index + index];
      index = *indice++;
      lsf1_r[3] = *p_dico++;
      lsf1_r[4] = *p_dico++;
      lsf1_r[5] = *p_dico++;
      p_dico = &p_cb3[index << 2];
      lsf1_r[6] = *p_dico++;
      lsf1_r[7] = *p_dico++;
      lsf1_r[8] = *p_dico++;
      lsf1_r[9] = *p_dico++;

      /* Compute quantized LSFs and update the past quantized residual */
      if ( mode != MRDTX ) {
         for ( i = 0; i < M; i++ ) {
            lsf1_q[i] = lsf1_r[i] + ( mean_lsf_3[i] + ( ( st->past_r_q[i] *
                  pred_fac[i] ) >> 15 ) );
         }
         memcpy( st->past_r_q, lsf1_r, M <<2 );
      }
      else {
         for ( i = 0; i < M; i++ ) {
            lsf1_q[i] = lsf1_r[i] + ( mean_lsf_3[i] + st->past_r_q[i] );
         }
         memcpy( st->past_r_q, lsf1_r, M <<2 );
      }
   }

   /* verification that LSFs has minimum distance of LSF_GAP Hz */
   temp = LSF_GAP;

   for ( i = 0; i < M; i++ ) {
      if ( lsf1_q[i] < temp ) {
         lsf1_q[i] = temp;
      }
      temp = lsf1_q[i] + LSF_GAP;
   }
   memcpy( st->past_lsf_q, lsf1_q, M <<2 );

   /*  convert LSFs to the cosine domain */
   Lsf_lsp( lsf1_q, lsp1_q );
   return;
}


/*
 * pseudonoise
 *
 *
 * Parameters:
 *    shift_reg         B: Old CN generator shift register state
 *    no_bits           I: Number of bits
 *
 * Function:
 *    pseudonoise
 *
 * Returns:
 *    noise_bits
 */
static Word32 pseudonoise( Word32 *shift_reg, Word32 no_bits )
{
   Word32 noise_bits, Sn, i;
   Word32 s_reg;


   s_reg = *shift_reg;
   noise_bits = 0;

   for ( i = 0; i < no_bits; i++ ) {
      /* State n == 31 */
      Sn = s_reg & 0x00000001L;

      /* State n == 3 */
      if ( s_reg & 0x10000000L ) {
         Sn = Sn ^ 0x1L;
      }
      else {
         Sn = Sn ^ 0x0L;
      }
      noise_bits = ( noise_bits << 1 ) | ( s_reg & 1 );
      s_reg = s_reg >> 1;

      if ( Sn & 1 ) {
         s_reg = s_reg | 0x40000000L;
      }
   }
   *shift_reg = s_reg;
   return noise_bits;
}


/*
 * Lsp_lsf
 *
 *
 * Parameters:
 *    lsp               I: LSP vector (range: -1<=val<1)
 *    lsf               O: LSF vector Old CN generator shift register state
 *
 * Function:
 *    Transformation lsp to lsf, LPC order M
 *    lsf[i] = arccos(lsp[i])/(2*pi)
 *
 * Returns:
 *    void
 */
static void Lsp_lsf( Word32 lsp[], Word32 lsf[] )
{
   Word32 i, ind = 63;   /* begin at end of table -1 */


   for ( i = M - 1; i >= 0; i-- ) {
      /* find value in table that is just greater than lsp[i] */
      while ( cos_table[ind] < lsp[i] ) {
         ind--;
      }
      lsf[i] = ( ( ( ( lsp[i] - cos_table[ind] ) * acos_slope[ind] ) + 0x800 )
            >> 12 ) + ( ind << 8 );
   }
   return;
}


/*
 * Reorder_lsf
 *
 *
 * Parameters:
 *    lsf            B: vector of LSFs (range: 0<=val<=0.5)
 *    min_dist       I: minimum required distance
 *
 * Function:
 *    Make sure that the LSFs are properly ordered and to keep a certain minimum
 *    distance between adjacent LSFs. LPC order = M.
 *
 * Returns:
 *    void
 */
static void Reorder_lsf( Word32 *lsf, Word32 min_dist )
{
   Word32 lsf_min, i;


   lsf_min = min_dist;

   for ( i = 0; i < M; i++ ) {
      if ( lsf[i] < lsf_min ) {
         lsf[i] = lsf_min;
      }
      lsf_min = lsf[i] + min_dist;
   }
}

/* VC5.0 Global optimization does not work with this function */
#if _MSC_VER == 1100
#pragma optimize( "g", off )
#endif
/*
 * Get_lsp_pol
 *
 *
 * Parameters:
 *    lsp               I: line spectral frequencies
 *    f                 O: polynomial F1(z) or F2(z)
 *
 * Function:
 *    Find the polynomial F1(z) or F2(z) from the LSPs.
 *
 *    F1(z) = product ( 1 - 2 lsp[i] z^-1 + z^-2 )
 *             i=0,2,4,6,8
 *    F2(z) = product   ( 1 - 2 lsp[i] z^-1 + z^-2 )
 *             i=1,3,5,7,9
 *
 *    where lsp[] is the LSP vector in the cosine domain.
 *
 *    The expansion is performed using the following recursion:
 *
 *    f[0] = 1
 *    b = -2.0 * lsp[0]
 *    f[1] = b
 *    for i=2 to 5 do
 *       b = -2.0 * lsp[2*i-2];
 *       f[i] = b*f[i-1] + 2.0*f[i-2];
 *       for j=i-1 down to 2 do
 *          f[j] = f[j] + b*f[j-1] + f[j-2];
 *       f[1] = f[1] + b;
 *
 * Returns:
 *    void
 */
static void Get_lsp_pol( Word32 *lsp, Word32 *f )
{
   volatile Word32 f0, f1, f2, f3, f4, f5;
   Word32 l1, l2, l3, l4;


   /* f[0] = 1.0; */
   f0 = 16777216L;

   /* f1 = *lsp * -1024; */
   f1 = -lsp[0] << 10;
   l1 = lsp[2];
   l2 = lsp[4];
   l3 = lsp[6];
   l4 = lsp[8];
   f2 = f0 << 1;
   f2 -= ( ( ( f1 >> 16 ) * l1 ) + ( ( ( f1 & 0xFFFE ) * l1 ) >> 16 ) ) << 2;
   f1 -= l1 << 10;
   f3 = f1 << 1;
   f3 -= ( ( ( f2 >> 16 ) * l2 ) + ( ( ( f2 & 0xFFFE ) * l2 ) >> 16 ) ) << 2;
   f2 += f0;
   f2 -= ( ( ( f1 >> 16 ) * l2 ) + ( ( ( f1 & 0xFFFE ) * l2 ) >> 16 ) ) << 2;
   f1 -= l2 << 10;
   f4 = f2 << 1;
   f4 -= ( ( ( f3 >> 16 ) * l3 ) + ( ( ( f3 & 0xFFFE ) * l3 ) >> 16 ) ) << 2;
   f3 += f1;
   f3 -= ( ( ( f2 >> 16 ) * l3 ) + ( ( ( f2 & 0xFFFE ) * l3 ) >> 16 ) ) << 2;
   f2 += f0;
   f2 -= ( ( ( f1 >> 16 ) * l3 ) + ( ( ( f1 & 0xFFFE ) * l3 ) >> 16 ) ) << 2;
   f1 -= l3 << 10;
   f5 = f3 << 1;
   f5 -= ( ( ( f4 >> 16 ) * l4 ) + ( ( ( f4 & 0xFFFE ) * l4 ) >> 16 ) ) << 2;
   f4 += f2;
   f4 -= ( ( ( f3 >> 16 ) * l4 ) + ( ( ( f3 & 0xFFFE ) * l4 ) >> 16 ) ) << 2;
   f3 += f1;
   f3 -= ( ( ( f2 >> 16 ) * l4 ) + ( ( ( f2 & 0xFFFE ) * l4 ) >> 16 ) ) << 2;
   f2 += f0;
   f2 -= ( ( ( f1 >> 16 ) * l4 ) + ( ( ( f1 & 0xFFFE ) * l4 ) >> 16 ) ) << 2;
   f1 -= l4 << 10;
   f[0] = f0;
   f[1] = f1;
   f[2] = f2;
   f[3] = f3;
   f[4] = f4;
   f[5] = f5;
   return;
}
#if _MSC_VER == 1100
#pragma optimize( "", on )
#endif


/*
 * Lsp_Az
 *
 *
 * Parameters:
 *    lsp                 I: Line spectral frequencies
 *    a                   O: Predictor coefficients
 *
 * Function:
 *    Converts from the line spectral pairs (LSP) to LP coefficients,
 *    for a 10th order filter.
 *
 *    Find the coefficients of F1(z) and F2(z)
 *    Multiply F1(z) by 1+z^{-1} and F2(z) by 1-z^{-1}
 *    A(z) = ( F1(z) + F2(z) ) / 2
 *
 * Returns:
 *    void
 */
static void Lsp_Az( Word32 lsp[], Word32 a[] )
{
   Word32 f1[6], f2[6];
   Word32 T0, i, j;


   Get_lsp_pol( &lsp[0], f1 );
   Get_lsp_pol( &lsp[1], f2 );

   for ( i = 5; i > 0; i-- ) {
      f1[i] += f1[i - 1];
      f2[i] -= f2[i - 1];
   }
   a[0] = 4096;

   for ( i = 1, j = 10; i <= 5; i++, j-- ) {
      T0 = f1[i] + f2[i];
      a[i] = (Word16)(T0 >> 13);  /* emulate fixed point bug */
      if ( ( T0 & 4096 ) != 0 ) {
         a[i]++;
      }
      T0 = f1[i] - f2[i];
      a[j] = (Word16)(T0 >> 13);   /* emulate fixed point bug */

      if ( ( T0 & 4096 ) != 0 ) {
         a[j]++;
      }
   }
   return;
}


/*
 * A_Refl
 *
 *
 * Parameters:
 *    a                 I: Directform coefficients
 *    refl              O: Reflection coefficients
 *
 * Function:
 *    Converts from the directform coefficients to reflection coefficients
 *
 * Returns:
 *    void
 */
static void A_Refl( Word32 a[], Word32 refl[] )
{
   /* local variables */
   int normShift;
   Word32 aState[M], bState[M];
   Word32 normProd, acc, temp, mult, scale, i, j;


   /* initialize states */
   memcpy( aState, a, M <<2 );

   /* backward Levinson recursion */
   for ( i = M - 1; i >= 0; i-- ) {
      if ( labs( aState[i] ) >= 4096 ) {
         goto ExitRefl;
      }
      refl[i] = aState[i] << 3;
      temp = ( refl[i] * refl[i] ) << 1;
      acc = ( MAX_32 - temp );
      normShift=0;
      if (acc != 0){
         temp = acc;
         while (!(temp & 0x40000000))
         {
            normShift++;
            temp = temp << 1;
         }
      }
      else{
         normShift = 0;
      }
      scale = 15 - normShift;
      acc = ( acc << normShift );
      temp = ( acc + ( Word32 )0x00008000L );

      if ( temp > 0 ) {
         normProd = temp >> 16;
         mult = 0x20000000L / normProd;
      }
      else
         mult = 16384;

      for ( j = 0; j < i; j++ ) {
         acc = aState[j] << 16;
         acc -= ( refl[i] * aState[i - j - 1] ) << 1;
         temp = ( acc + ( Word32 )0x00008000L ) >> 16;
         temp = ( mult * temp ) << 1;

         if ( scale > 0 ) {
            if ( ( temp & ( ( Word32 )1 << ( scale - 1 ) ) ) != 0 ) {
               temp = ( temp >> scale ) + 1;
            }
            else
               temp = ( temp >> scale );
         }
         else
            temp = ( temp >> scale );

         if ( labs( temp ) > 32767 ) {
            goto ExitRefl;
         }
         bState[j] = temp;
      }
      memcpy( aState, bState, i <<2 );
   }
   return;
ExitRefl:
   memset( refl, 0, M <<2 );
}


/*
 * Log2_norm
 *
 *
 * Parameters:
 *    x                 I: input value
 *    exp               I: exponent
 *    exponent          O: Integer part of Log2. (range: 0<=val<=30)
 *    fraction          O: Fractional part of Log2. (range: 0<=val<1)
 *
 * Function:
 *    Computes log2
 *
 *    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)
 *
 *    exponent = 30-normExponent
 *    i = bit25-b31 of L_x;  32<=i<=63  (because of normalization).
 *    a = bit10-b24
 *    i -=32
 *    fraction = table[i]<<16 - (table[i] - table[i+1]) * a * 2
 *
 * Returns:
 *    void
 */
static void Log2_norm( Word32 x, Word32 exp, Word32 *exponent, Word32 *
      fraction )
{
   Word32 y, i, a;


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

   /* Extract b25-b31 */
   i = x >> 25;
   i = i - 32;

   /* Extract b10-b24 of fraction */
   a = x >> 9;
   a = a & 0xFFFE;   /* 2a */

   /* fraction */
   y = ( log2_table[i] << 16 ) - a * ( log2_table[i] - log2_table[i + 1] );
   *fraction = y >> 16;
   *exponent = 30 - exp;
   return;
}


/*
 * Log2
 *
 *
 * Parameters:
 *    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)
 *    If x is negative or zero, the result is 0.
 *
 * Returns:
 *    void
 */
static void Log2( Word32 x, Word32 *exponent, Word32 *fraction )
{
   int tmp, exp=0;

   if (x != 0){
         tmp = x;
         while (!((tmp & 0x80000000) ^ ((tmp & 0x40000000) << 1)))
         {
            exp++;
            tmp = tmp << 1;