pitch_a.c
上传用户:zhouyunkk
上传日期:2013-01-10
资源大小:59k
文件大小:13k
源码类别:

语音压缩

开发平台:

C/C++

  1. /*
  2.    ITU-T G.729 Annex C - Reference C code for floating point
  3.                          implementation of G.729 Annex A
  4.                          Version 1.01 of 15.September.98
  5. */
  6. /*
  7. ----------------------------------------------------------------------
  8.                     COPYRIGHT NOTICE
  9. ----------------------------------------------------------------------
  10.    ITU-T G.729 Annex C ANSI C source code
  11.    Copyright (C) 1998, AT&T, France Telecom, NTT, University of
  12.    Sherbrooke.  All rights reserved.
  13. ----------------------------------------------------------------------
  14. */
  15. /***********************************************************************/
  16. /*           Long Term (Pitch) Prediction Functions                    */
  17. /***********************************************************************/
  18. #include <math.h>
  19. #include "typedef.h"
  20. #include "ld8a.h"
  21. /* prototypes for local functions */
  22. static FLOAT dot_product(FLOAT x[], FLOAT y[], int lg);
  23. /*----------------------------------------------------------------------*
  24.  * pitch_ol_fast -> compute the open loop pitch lag -> fast version     *
  25.  *----------------------------------------------------------------------*/
  26. int pitch_ol_fast(         /* output: open-loop pitch lag       */
  27.     FLOAT signal[],        /* input : signal to compute pitch   */
  28.                            /*         s[-PIT_MAX : l_frame-1]   */
  29.     int l_frame            /* input : error minimization window */
  30. )
  31. {
  32.     int  i, j;
  33.     int  T1, T2, T3;
  34.     FLOAT  max1, max2, max3;
  35.     FLOAT  *p, *p1, sum;
  36.    /*--------------------------------------------------------------------*
  37.     *  The pitch lag search is divided in three sections.                *
  38.     *  Each section cannot have a pitch multiple.                        *
  39.     *  A maximum is find for each section.                               *
  40.     *  The final lag is selected by taking into account the multiple.    *
  41.     *                                                                    *
  42.     *  First section:  lag delay = 20 to 39                              *
  43.     *  Second section: lag delay = 40 to 79                              *
  44.     *  Third section:  lag delay = 80 to 143                             *
  45.     *--------------------------------------------------------------------*/
  46.     /* First section */
  47.     max1 = FLT_MIN_G729;
  48.     for (i = 20; i < 40; i++) {
  49.         p  = signal;
  50.         p1 = &signal[-i];
  51.         sum = (F)0.0;
  52.         /* Dot product with decimation by 2 */
  53.         for (j=0; j<l_frame; j+=2, p+=2, p1+=2)
  54.             sum += *p * *p1;
  55.         if (sum > max1) { max1 = sum; T1 = i;}
  56.     }
  57.     /* compute energy of maximum1 */
  58.     sum = (F)0.01;                   /* to avoid division by zero */
  59.     p = &signal[-T1];
  60.     for(i=0; i<l_frame; i+=2, p+=2)
  61.         sum += *p * *p;
  62.     sum = (F)1.0/(FLOAT)sqrt(sum);          /* 1/sqrt(energy)    */
  63.     max1 *= sum;                  /* max/sqrt(energy)  */
  64.     /* Second section */
  65.     max2 = FLT_MIN_G729;
  66.     for (i = 40; i < 80; i++) {
  67.         p  = signal;
  68.         p1 = &signal[-i];
  69.         sum = (F)0.0;
  70.         /* Dot product with decimation by 2 */
  71.         for (j=0; j<l_frame; j+=2, p+=2, p1+=2)
  72.             sum += *p * *p1;
  73.         if (sum > max2) { max2 = sum; T2 = i;}
  74.     }
  75.     /* compute energy of maximum2 */
  76.     sum = (F)0.01;                   /* to avoid division by zero */
  77.     p = &signal[-T2];
  78.     for(i=0; i<l_frame; i+=2, p+=2)
  79.         sum += *p * *p;
  80.     sum  = (F)1.0/(FLOAT)sqrt(sum);         /* 1/sqrt(energy)    */
  81.     max2 *= sum;                  /* max/sqrt(energy)  */
  82.     /* Third section */
  83.     max3 = FLT_MIN_G729;
  84.     /* decimation by 2 for the possible delay */
  85.     for (i = 80; i < 143; i+=2) {
  86.         p  = signal;
  87.         p1 = &signal[-i];
  88.         sum = (F)0.0;
  89.         /* Dot product with decimation by 2 */
  90.         for (j=0; j<l_frame; j+=2, p+=2, p1+=2)
  91.             sum += *p * *p1;
  92.         if (sum > max3) { max3 = sum; T3 = i;}
  93.     }
  94.      /* Test around max3 */
  95.      i = T3;
  96.      p  = signal;
  97.      p1 = &signal[-(i+1)];
  98.      sum = (F)0.0;
  99.      for (j=0; j<l_frame; j+=2, p+=2, p1+=2)
  100.          sum += *p * *p1;
  101.      if (sum > max3) { max3 = sum; T3 = i+1;}
  102.      p  = signal;
  103.      p1 = &signal[-(i-1)];
  104.      sum = (F)0.0;
  105.      for (j=0; j<l_frame; j+=2, p+=2, p1+=2)
  106.          sum += *p * *p1;
  107.      if (sum > max3) { max3 = sum; T3 = i-1;}
  108.     /* compute energy of maximum3 */
  109.     sum = (F)0.01;                   /* to avoid division by zero */
  110.     p = &signal[-T3];
  111.     for(i=0; i<l_frame; i+=2, p+=2)
  112.         sum += *p * *p;
  113.     sum  = (F)1.0/(FLOAT)sqrt(sum);         /* 1/sqrt(energy)    */
  114.     max3 *= sum;                  /* max/sqrt(energy)  */
  115.    /*-----------------------*
  116.     * Test for multiple.    *
  117.     *-----------------------*/
  118.     if( abs(T2*2 - T3) < 5)
  119.        max2 += max3 * (F)0.25;
  120.     if( abs(T2*3 - T3) < 7)
  121.        max2 += max3 * (F)0.25;
  122.     if( abs(T1*2 - T2) < 5)
  123.        max1 += max2 * (F)0.20;
  124.     if( abs(T1*3 - T2) < 7)
  125.        max1 += max2 * (F)0.20;
  126.    /*--------------------------------------------------------------------*
  127.     * Compare the 3 sections maxima.                                     *
  128.     *--------------------------------------------------------------------*/
  129.     if( max1 < max2 ) {max1 = max2; T1 = T2;}
  130.     if( max1 < max3 )  T1 = T3;
  131.     return (T1);
  132. }
  133. /*------------------------------------------------------------------*
  134.  * dot_product()                                                    *
  135.  *  dot product between vector x[] and y[] of lenght lg             *
  136.  *------------------------------------------------------------------*/
  137. static FLOAT dot_product(   /* Return the dot product between x[] an y[] */
  138.   FLOAT x[],                /* First vector.                             */
  139.   FLOAT y[],                /* Second vector.                            */
  140.   int lg                    /* Lenght of the product.                    */
  141. )
  142. {
  143.   int i;
  144.   FLOAT sum;
  145.   sum = (F)0.1;
  146.   for(i=0; i<lg; i++)
  147.      sum +=  x[i] * y[i];
  148.   return sum;
  149. }
  150. /*-------------------------------------------------------------------------*
  151.  * pitch_fr3_fast()                                                        *
  152.  *  Find the pitch period in close loop with 1/3 subsample resolution      *
  153.  *  Fast version                                                           *
  154.  *-------------------------------------------------------------------------*/
  155. int pitch_fr3_fast(     /* output: integer part of pitch period        */
  156.   FLOAT exc[],          /* input : excitation buffer                   */
  157.   FLOAT xn[],           /* input : target vector                       */
  158.   FLOAT h[],            /* input : impulse response.                   */
  159.   int l_subfr,          /* input : Length of subframe                  */
  160.   int t0_min,           /* input : minimum value in the searched range */
  161.   int t0_max,           /* input : maximum value in the searched range */
  162.   int i_subfr,          /* input : indicator for first subframe        */
  163.   int *pit_frac         /* output: chosen fraction                     */
  164. )
  165. {
  166.   int  t, t0;
  167.   FLOAT  dn[L_SUBFR];
  168.   FLOAT  exc_tmp[L_SUBFR];
  169.   FLOAT  corr, max;
  170.   /* Compute  correlations of input response h[] with the target vector xn[].*/
  171.   cor_h_x (h, xn, dn);
  172.   /* Find maximum integer delay */
  173.   max = FLT_MIN_G729;
  174.   for(t=t0_min; t<=t0_max; t++)
  175.   {
  176.     corr = dot_product(dn, &exc[-t], l_subfr);
  177.     if(corr > max) {max = corr; t0 = t;}
  178.   }
  179.   /* Test fractions */
  180.   /* Fraction 0 */
  181.   pred_lt_3(exc, t0, 0, l_subfr);
  182.   max = dot_product(dn, exc, l_subfr);
  183.   *pit_frac = 0;
  184.   /* If first subframe and lag > 84 do not search fractional pitch */
  185.   if( (i_subfr == 0) && (t0 > 84) )
  186.     return t0;
  187.   copy(exc, exc_tmp, l_subfr);
  188.   /* Fraction -1/3 */
  189.   pred_lt_3(exc, t0, -1, l_subfr);
  190.   corr = dot_product(dn, exc, l_subfr);
  191.   if(corr > max){
  192.     max = corr;
  193.     *pit_frac = -1;
  194.     copy(exc, exc_tmp, l_subfr);
  195.   }
  196.   /* Fraction +1/3 */
  197.   pred_lt_3(exc, t0, 1, l_subfr);
  198.   corr = dot_product(dn, exc, l_subfr);
  199.   if(corr > max){
  200.     max = corr;
  201.     *pit_frac = 1;
  202.   }
  203.   else
  204.     copy(exc_tmp, exc, l_subfr);
  205.   return t0;
  206. }
  207. /*---------------------------------------------------------------------------*
  208.  * g_pitch  - compute adaptive codebook gain and compute <y1,y1> , -2<xn,y1> *
  209.  *---------------------------------------------------------------------------*/
  210. FLOAT g_pitch(          /* output: pitch gain                        */
  211.   FLOAT xn[],           /* input : target vector                     */
  212.   FLOAT y1[],           /* input : filtered adaptive codebook vector */
  213.   FLOAT g_coeff[],      /* output: <y1,y1> and -2<xn,y1>             */
  214.   int l_subfr           /* input : vector dimension                  */
  215. )
  216. {
  217.     FLOAT xy, yy, gain;
  218.     int   i;
  219.     yy = (F)0.01;
  220.     for (i = 0; i < l_subfr; i++) {
  221.         yy += y1[i] * y1[i];          /* energy of filtered excitation */
  222.     }
  223.     xy = (F)0.0;
  224.     for (i = 0; i < l_subfr; i++) {
  225.         xy += xn[i] * y1[i];
  226.     }
  227.     g_coeff[0] = yy;
  228.     g_coeff[1] = (F)-2.0*xy +(F)0.01;
  229.     /* find pitch gain and bound it by [0,1.2] */
  230.     gain = xy/yy;
  231.     if (gain<(F)0.0)  gain = (F)0.0;
  232.     if (gain>GAIN_PIT_MAX) gain = GAIN_PIT_MAX;
  233.     return gain;
  234. }
  235. /*----------------------------------------------------------------------*
  236.  *    Functions enc_lag3()                                              *
  237.  *              ~~~~~~~~~~                                              *
  238.  *   Encoding of fractional pitch lag with 1/3 resolution.              *
  239.  *----------------------------------------------------------------------*
  240.  * The pitch range for the first subframe is divided as follows:        *
  241.  *   19 1/3  to   84 2/3   resolution 1/3                               *
  242.  *   85      to   143      resolution 1                                 *
  243.  *                                                                      *
  244.  * The period in the first subframe is encoded with 8 bits.             *
  245.  * For the range with fractions:                                        *
  246.  *   index = (T-19)*3 + frac - 1;   where T=[19..85] and frac=[-1,0,1]  *
  247.  * and for the integer only range                                       *
  248.  *   index = (T - 85) + 197;        where T=[86..143]                   *
  249.  *----------------------------------------------------------------------*
  250.  * For the second subframe a resolution of 1/3 is always used, and the  *
  251.  * search range is relative to the lag in the first subframe.           *
  252.  * If t0 is the lag in the first subframe then                          *
  253.  *  t_min=t0-5   and  t_max=t0+4   and  the range is given by           *
  254.  *       t_min - 2/3   to  t_max + 2/3                                  *
  255.  *                                                                      *
  256.  * The period in the 2nd subframe is encoded with 5 bits:               *
  257.  *   index = (T-(t_min-1))*3 + frac - 1;    where T[t_min-1 .. t_max+1] *
  258.  *----------------------------------------------------------------------*/
  259. int  enc_lag3(     /* output: Return index of encoding */
  260.   int  T0,         /* input : Pitch delay              */
  261.   int  T0_frac,    /* input : Fractional pitch delay   */
  262.   int  *T0_min,    /* in/out: Minimum search delay     */
  263.   int  *T0_max,    /* in/out: Maximum search delay     */
  264.   int pit_min,     /* input : Minimum pitch delay      */
  265.   int pit_max,     /* input : Maximum pitch delay      */
  266.   int  i_subfr     /* input : Flag for 1st subframe    */
  267. )
  268. {
  269.   int index;
  270.   if (i_subfr == 0)   /* if 1st subframe */
  271.   {
  272.      /* encode pitch delay (with fraction) */
  273.      if (T0 <= 85)
  274.        index = T0*3 - 58 + T0_frac;
  275.      else
  276.        index = T0 + 112;
  277.      /* find T0_min and T0_max for second subframe */
  278.      *T0_min = T0 - 5;
  279.      if (*T0_min < pit_min) *T0_min = pit_min;
  280.      *T0_max = *T0_min + 9;
  281.      if (*T0_max > pit_max)
  282.      {
  283.          *T0_max = pit_max;
  284.          *T0_min = *T0_max - 9;
  285.      }
  286.   }
  287.   else                    /* second subframe */
  288.   {
  289.      index = T0 - *T0_min;
  290.      index = index*3 + 2 + T0_frac;
  291.   }
  292.   return index;
  293. }