delay.c
上传用户:szhypcb168
上传日期:2007-01-06
资源大小:2187k
文件大小:5k
源码类别:

语音压缩

开发平台:

Unix_Linux

  1. #define NFRAC 5
  2. #define TRUE 1
  3. #define FALSE 0
  4. #define M1 -4
  5. #define M2  3
  6. /* five fractional delays calculated over an 8 point interpolation  */
  7. /* (-4 to 3) */
  8. static float frac[NFRAC] = {0.25, 0.33333333, 0.5, 0.66666667, 0.75};
  9. static int twelfths[NFRAC] = {3, 4, 6, 8, 9};
  10. /**************************************************************************
  11. *
  12. * NAME
  13. * sinc
  14. *
  15. * FUNCTION
  16. *
  17. *
  18. * SYNOPSIS
  19. * subroutine sinc(arg)
  20. *
  21. *   formal 
  22. *                       data    I/O
  23. *       name            type    type    function
  24. *       -------------------------------------------------------------------
  25. * arg float i
  26. * y(n) float func
  27. *
  28. ***************************************************************************
  29. *
  30. * DESCRIPTION
  31. *
  32. * This interpolating (sinc) function is Hamming windowed to bandlimit
  33. * the signal to reduce aliasing.
  34. *
  35. ***************************************************************************
  36. *
  37. * CALLED BY
  38. *
  39. * pitchvq  psearch
  40. *
  41. * CALLS
  42. *
  43. *
  44. ***************************************************************************/
  45. #include <math.h>
  46. #include <stdio.h>
  47. float
  48. sinc(arg)
  49. float arg;
  50. {
  51.   float pi;
  52.   pi = 4.0 * atan(1.0);
  53.   if (arg == 0.0)
  54.     return(1.0);
  55.   else
  56.     return(sin(pi * arg) / (pi * arg));
  57. /**************************************************************************
  58. *
  59. * NAME
  60. * qd
  61. *
  62. * FUNCTION
  63. *
  64. *
  65. * SYNOPSIS
  66. * subroutine qd(d)
  67. *
  68. *   formal 
  69. *                       data    I/O
  70. *       name            type    type    function
  71. *       -------------------------------------------------------------------
  72. * d float i
  73. * qd(d) int o quantize d function
  74. *
  75. ***************************************************************************
  76. *
  77. * DESCRIPTION
  78. *
  79. * Quantize d function
  80. *
  81. ***************************************************************************
  82. *
  83. * CALLED BY
  84. *
  85. * delay
  86. *
  87. * CALLS
  88. *
  89. *
  90. ***************************************************************************/
  91. int 
  92. qd(d)
  93. float d;
  94. {
  95.   int i, index, ok;
  96.   ok = FALSE;
  97.   for (i = 0; i < NFRAC; i++)
  98.   {
  99.     if (fabs(d - frac[i]) < 1.e-2)
  100.     {
  101.       index = i;
  102.       ok = TRUE;
  103.     }
  104.   }
  105.   if (!ok)
  106.   {
  107.     fprintf(stderr, "delay: Invalid pitch delay = %fn", d);
  108.     exit(1);
  109.   }
  110.   return(index);
  111. }
  112.   
  113. /**************************************************************************
  114. *
  115. * NAME
  116. * delay
  117. *
  118. * FUNCTION
  119. * Time delay a bandlimited signal
  120. * using point-by-point recursion
  121. *
  122. * SYNOPSIS
  123. * subroutine delay(x, start, n, d, m, y)
  124. *
  125. *   formal 
  126. *                       data    I/O
  127. *       name            type    type    function
  128. *       -------------------------------------------------------------------
  129. * x[n] float i signal input (output in last 60)
  130. * start int i beginning of output sequence
  131. * n int i length of input sequence
  132. * d float i fractional pitch
  133. * m int i integer pitch
  134. * y[n] float o delayed input signal
  135. *
  136. ***************************************************************************
  137. *
  138. * DESCRIPTION
  139. *
  140. * Subroutine to time delay a bandlimited signal by resampling the
  141. * reconstructed data (aka sinc interpolation).  The well known
  142. * reconstruction formula is:
  143. *
  144. *              |    M2      sin[pi(t-nT)/T]    M2
  145. *   y(n) = X(t)| = SUM x(n) --------------- = SUM x(n) sinc[(t-nT)/T]
  146. *              |   n=M1         pi(t-nT)/T    n=M1
  147. *              |t=n+d
  148. *
  149. * The interpolating (sinc) function is Hamming windowed to bandlimit
  150. * the signal to reduce aliasing.
  151. *
  152. * Multiple simultaneous time delays may be efficiently calculated
  153. * by polyphase filtering.  Polyphase filters decimate the unused
  154. * filter coefficients.  See Chapter 6 in C&R for details. 
  155. *
  156. ***************************************************************************
  157. *
  158. * CALLED BY
  159. *
  160. * pitchvq   psearch
  161. *
  162. * CALLS
  163. *
  164. * ham
  165. *
  166. ***************************************************************************
  167. *
  168. * REFERENCES
  169. *
  170. * Crochiere & Rabiner, Multirate Digital Signal Processing,
  171. * P-H 1983, Chapters 2 and 6.
  172. *
  173. *
  174. * Kroon & Atal, "Pitch Predictors with High Temporal Resolution,"
  175. * ICASSP '90, S12.6
  176. *
  177. **************************************************************************/
  178. #define SIZE (M2 - M1 + 1)
  179. delay(x, start, n, d, m, y)
  180. float x[], d, y[];
  181. int m, n, start;
  182. {
  183.   static float wsinc[SIZE][NFRAC], hwin[12*SIZE+1];
  184.   static int first = TRUE;
  185.   int i, j, k, index;
  186.   /* Generate Hamming windowed sinc interpolating function for each */
  187.   /* allowable fraction.  The interpolating functions are stored in   */
  188.   /* time reverse order (i.e., delay appears as advance) to align     */
  189.   /* with the adaptive code book v0 array.  The interp filters are:   */
  190.   /*  wsinc[.,0] frac = 1/4 (3/12)       */
  191.   /*  wsinc[.,1] frac = 1/3 (4/12)       */
  192.   /*  . .       */
  193.   /*  wsinc[.,4] frac = 3/4 (9/12)       */
  194.   if (first)
  195.   {
  196.     ham(hwin, 12*SIZE+1);
  197.     for (i = 0; i < NFRAC; i++)
  198.       for (k = M1, j = 0; j < SIZE; j++)
  199.       {
  200.         wsinc[j][i] = sinc(frac[i] + k) * hwin[12*j+twelfths[i]];
  201.         k++;
  202.       }
  203.     first = FALSE;
  204.   }
  205.   index = qd(d);
  206.   /* *Resample:  */
  207.   for (i = 0; i < n; i++)
  208.   {
  209.     x[start+i-1] = 0.0;
  210.     for (k = M1, j = 0; j < SIZE; j++)
  211.     {
  212.       x[start+i-1] += x[start-m+i+k-1] * wsinc[j][index];
  213.       k++;
  214.     }
  215.   }
  216.   /* *The v0 array in psearch.c/pgain.c must be zero above "start" */
  217.   /* *because of overlap and add convolution techniques used in pgain. */
  218.   for (i = 0; i < n; i++)
  219.   {
  220.     y[i] = x[start+i-1];
  221.     x[start+i-1] = 0.0;
  222.   }
  223. }