g7231_lsp.c

G723.1语音压缩解压在tms320c54系列上的实现代码,本人已在CCS上仿真通过. 包含全部源代码,主函数请自已写(本人的就不奉送了:

#include "typedef.h"
#include "G7231_basop.h"
#include "G7231_cst_lbc.h"
#include "G7231_tab_lbc.h"
#include "G7231_lsp.h"
#include "intrindefs.h"

void G7231AtoLsp( Word16 *LspVect, Word16 *Lpc, Word16 *PrevLsp )
{

    int   i,j,k ;

    Word32   Lpq[G7231LpcOrder+2] ;
    Word16   Spq[G7231LpcOrder+2] ;

    Word16   Exp   ;
    Word16   LspCnt ;

    Word32   PrevVal,CurrVal   ;
    Word32   Acc0,Acc1   ;


 /*
  * Perform a bandwidth expansion on the LPC coefficients.  This
  * scales the poles of the LPC synthesis filter by a factor of
  * 0.994.
  */
    for ( i = 0 ; i < G7231LpcOrder ; i ++ )
        LspVect[i] = G7231mult_r( Lpc[i], G7231BandExpTable[i] ) ;


 /*
  * Compute the sum and difference polynomials with the roots at z =
  * -1 (sum) or z = +1 (difference) removed.  Let these polynomials
  * be P(z) and Q(z) respectively, and let their coefficients be
  * {p_i} amd {q_i}.  The coefficients are stored in the array Lpq[]
  * as follows: p_0, q_0, p_1, q_1, ..., p_5, q_5.  There is no need
  * to store the other coefficients because of symmetry.
  */


 /*
  * Set p_0 = q_0 = 1.  The LPC coefficients are already scaled by
  *  1/4.  P(z) and Q(z) are scaled by an additional scaling factor of
  *  1/16, for an overall factor of 1/64 = 0x02000000L.
  */

    Lpq[0] = Lpq[1] = (Word32) 0x02000000L ;

 /*
  * This loop computes the coefficients of P(z) and Q(z).  The long
  * division (to remove the real zeros) is done recursively.
  */
    for ( i = 0 ; i < G7231LpcOrder/2 ; i ++ ) {

        /* P(z) */
        Acc0 = G7231L_negate( Lpq[2*i+0] ) ;
        Acc1 = G7231L_deposit_h( LspVect[i] ) ;
        Acc1 = G7231L_shr( Acc1, (Word16) 4 ) ;
        Acc0 = G7231L_sub( Acc0, Acc1 ) ;
        Acc1 = G7231L_deposit_h( LspVect[G7231LpcOrder-1-i] ) ;
        Acc1 = G7231L_shr( Acc1, (Word16) 4 ) ;
        Acc0 = G7231L_sub( Acc0, Acc1 ) ;
        Lpq[2*i+2] = Acc0 ;

        /* Q(z) */
        Acc0 = Lpq[2*i+1] ;
        Acc1 = G7231L_deposit_h( LspVect[i] ) ;
        Acc1 = G7231L_shr( Acc1, (Word16) 4 ) ;

        Acc0 = G7231L_sub( Acc0, Acc1 ) ;
        Acc1 = G7231L_deposit_h( LspVect[G7231LpcOrder-1-i] ) ;
        Acc1 = G7231L_shr( Acc1, (Word16) 4 ) ;
        Acc0 = G7231L_add( Acc0, Acc1 ) ;
        Lpq[2*i+3] = Acc0 ;
    }

 /*
  * Divide p_5 and q_5 by 2 for proper weighting during polynomial
  * evaluation.
  */
    Lpq[G7231LpcOrder+0] = G7231L_shr( Lpq[G7231LpcOrder+0], (Word16) 1 ) ;
    Lpq[G7231LpcOrder+1] = G7231L_shr( Lpq[G7231LpcOrder+1], (Word16) 1 ) ;

 /*
  * Normalize the polynomial coefficients and convert to shorts
  */

    /* Find the maximum */
    Acc1 = G7231L_abs( Lpq[0] ) ;
    for ( i = 1 ; i < G7231LpcOrder+2 ; i ++ ) {
        Acc0 = G7231L_abs( Lpq[i] ) ;
        if ( Acc0 > Acc1 )
            Acc1 = Acc0 ;
    }

    /* Compute the normalization factor */
    Exp = G7231norm_l( Acc1 ) ;


    /* Normalize and convert to shorts */
    for ( i = 0 ; i < G7231LpcOrder+2 ; i ++ ) {
        Acc0 = G7231L_shl( Lpq[i], Exp ) ;
        Spq[i] = G7231round( Acc0 ) ;
    }

 /*
  * Initialize the search loop
  */

 /*
  * The variable k is a flag that indicates which polynomial (sum or
  * difference) the algorithm is currently evaluating.  Start with
  * the sum.
  */
    k = 0 ;

    /* Evaluate the sum polynomial at frequency zero */
    PrevVal = (Word32) 0 ;
    for ( j = 0 ; j <= G7231LpcOrder/2 ; j ++ )
        PrevVal = G7231L_mac( PrevVal, Spq[2*j], G7231CosineTable[0] ) ;


 /*
  * Search loop.  Evaluate P(z) and Q(z) at uniform intervals of
  * pi/256 along the unit circle.  Check for zero crossings.  The
  * zeros of P(w) and Q(w) alternate, so only one of them need by
  * evaluated at any given step.
  */
    LspCnt = (Word16) 0 ;
    for ( i = 1 ; i < G7231CosineTableSize/2 ; i ++ ) {

        /* Evaluate the selected polynomial */
        CurrVal = (Word32) 0 ;
        for ( j = 0 ; j <= G7231LpcOrder/2 ; j ++ )
            CurrVal = G7231L_mac( CurrVal, Spq[G7231LpcOrder-2*j+k],
                                    G7231CosineTable[i*j%G7231CosineTableSize] ) ;

        /* Check for a sign change, indicating a zero crossing */
        if ( (CurrVal ^ PrevVal) < (Word32) 0 ) {

 /*
  * Interpolate to find the bottom 7 bits of the
  * zero-crossing frequency
  */
            Acc0 = G7231L_abs( CurrVal ) ;
            Acc1 = G7231L_abs( PrevVal ) ;
            Acc0 = G7231L_add( Acc0, Acc1 ) ;

            /* Normalize the sum */
            Exp = G7231norm_l( Acc0 ) ;
            Acc0 = G7231L_shl( Acc0, Exp ) ;
            Acc1 = G7231L_shl( Acc1, Exp ) ;

            Acc1 = G7231L_shr( Acc1, (Word16) 8 ) ;

            LspVect[LspCnt] = G7231div_l( Acc1, G7231extract_h( Acc0 ) ) ;

 /*
  * Add the upper part of the zero-crossing frequency,
  * i.e. bits 7-15
  */
            Exp = G7231shl( (Word16) (i-1), (Word16) 7 ) ;
            LspVect[LspCnt] = G7231add( LspVect[LspCnt], Exp ) ;
            LspCnt ++ ;

            /* Check if all zeros have been found */
            if ( LspCnt == (Word16) G7231LpcOrder )
                break ;

 /*
  * Switch the pointer between sum and difference polynomials
  */
            k ^= 1 ;

 /*
  * Evaluate the new polynomial at the current frequency
  */
            CurrVal = (Word32) 0 ;
            for ( j = 0 ; j <= G7231LpcOrder/2 ; j ++ )
                CurrVal = G7231L_mac( CurrVal, Spq[G7231LpcOrder-2*j+k],
                                    G7231CosineTable[i*j%G7231CosineTableSize] ) ;
        }

        /* Update the previous value */
        PrevVal = CurrVal ;
    }


 /*
  * Check if all 10 zeros were found.  If not, ignore the results of
  * the search and use the previous frame's LSP frequencies instead.
  */
    if ( LspCnt != (Word16) G7231LpcOrder ) {
        for ( j = 0 ; j < G7231LpcOrder ; j ++ )
            LspVect[j] = PrevLsp[j] ;
    }

    return ;
}

Word32   G7231Lsp_Qnt( Word16 *CurrLsp, Word16 *PrevLsp )
{
    int   i ;

    Word16   Wvect[G7231LpcOrder] ;

    Word16   Tmp0,Tmp1   ;
    Word16   Exp   ;


 /*
  * Compute the VQ weighting vector.  The weights assign greater
  * precision to those frequencies that are closer together.
  */

    /* Compute the end differences */
    Wvect[0] = G7231sub( CurrLsp[1], CurrLsp[0] ) ;
    Wvect[G7231LpcOrder-1] = G7231sub( CurrLsp[G7231LpcOrder-1], CurrLsp[G7231LpcOrder-2] ) ;

    /* Compute the rest of the differences */
    for ( i = 1 ; i < G7231LpcOrder-1 ; i ++ ) 
    {
        Tmp0 = G7231sub( CurrLsp[i+1], CurrLsp[i] ) ;
        Tmp1 = G7231sub( CurrLsp[i], CurrLsp[i-1] ) ;
        if ( Tmp0 > Tmp1 )
            Wvect[i] = Tmp1 ;
        else
            Wvect[i] = Tmp0 ;
    }

    /* Invert the differences */
    Tmp0 = (Word16) 0x0020 ;
    for ( i = 0 ; i < G7231LpcOrder ; i ++ ) 
    {

        if ( Wvect[i] > Tmp0 )
        {
            Wvect[i] = G7231div_s( Tmp0, Wvect[i] ) ;
        }
        else
            Wvect[i] = G7231MAX_16 ;
    }

    /* Normalize the weight vector */
    Tmp0 = (Word16) 0 ;
    for ( i = 0 ; i < G7231LpcOrder ; i ++ )
        if ( Wvect[i] > Tmp0 )
            Tmp0 = Wvect[i] ;

    Exp = G7231norm_s( Tmp0 ) ;
    for ( i = 0 ; i < G7231LpcOrder ; i ++ )
        Wvect[i] = G7231shl( Wvect[i], Exp ) ;


 /*
  * Compute the VQ target vector.  This is the residual that remains
  * after subtracting both the DC and predicted
  * components.
  */

 /*
  * Subtract the DC component from both the current and previous LSP
  * vectors.
  */
    for ( i = 0 ; i < G7231LpcOrder ; i ++ ) {
        CurrLsp[i] = G7231sub( CurrLsp[i], G7231LspDcTable[i] ) ;
        PrevLsp[i] = G7231sub( PrevLsp[i], G7231LspDcTable[i] ) ;
    }

 /*
  * Generate the prediction vector and subtract it.  Use a constant
  * first-order predictor based on the previous (DC-free) LSP
  * vector.
  */
    for ( i = 0 ; i < G7231LpcOrder ; i ++ ) {
        Tmp0 = G7231mult_r( PrevLsp[i], (Word16) G7231LspPrd0 ) ;
        CurrLsp[i] = G7231sub( CurrLsp[i], Tmp0 ) ;
    }

 /*
  * Add the DC component back to the previous LSP vector.  This
  * vector is needed in later routines.
  */
    for ( i = 0 ; i < G7231LpcOrder ; i ++ )
        PrevLsp[i] = G7231add( PrevLsp[i], G7231LspDcTable[i] ) ;

 /*
  * Do the vector quantization for all three bands
  */
    return G7231Lsp_Svq( CurrLsp, Wvect ) ;
}

Word32   G7231Lsp_Svq( Word16 *Tv, Word16 *Wvect )
{
    int   i,j,k ;

    Word32   Rez,Indx    ;
    Word32   Acc0,Acc1   ;

    Word16   Tmp[G7231LpcOrder] ;
    Word16  *LspQntPnt  ;


 /*
  * Initialize the return value
  */
    Rez = (Word32) 0 ;

 /*
  * Quantize each band separately
  */
    for ( k = 0 ; k < G7231LspQntBands ; k ++ ) {

 /*
  * Search over the entire VQ table to find the index that
  * minimizes the error.
  */

        /* Initialize the search */
        Acc1 = (Word32) -1 ;
        Indx = (Word32) 0 ;
        LspQntPnt = G7231BandQntTable[k] ;

        for ( i = 0 ; i < G7231LspCbSize ; i ++ ) {

 /*
  * Generate the metric, which is the negative error with the