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

语音压缩

开发平台:

Unix_Linux

  1. C==========================================================================
  2. C
  3. C ROUTINE
  4. C ldelay
  5. C
  6. C FUNCTION
  7. C Time delay a bandlimited signal
  8. C using point-by-point recursion
  9. C with a long interpolation filter (see delay.f)
  10. C
  11. C SYNOPSIS
  12. C subroutine ldelay(x, bufr, start, n, d, m, y)
  13. C
  14. C   formal 
  15. C                       data    I/O
  16. C       name            type    type    function
  17. C       -------------------------------------------------------------------
  18. C x(n) real i signal input (output in last 60)
  19. C bufr int i length of input sequence x(n)
  20. C start int i beginning of output sequence
  21. C n int i length of output sequence
  22. C d real i fractional pitch
  23. C m int i integer pitch
  24. C y(n) real o delayed input signal
  25. C==========================================================================
  26. C
  27. C DESCRIPTION
  28. C
  29. C Subroutine to time delay a bandlimited signal by resampling the
  30. C reconstructed data (aka sinc interpolation).  The well known
  31. C reconstruction formula is:
  32. C
  33. C              |    M2      sin[pi(t-nT)/T]    M2
  34. C   y(n) = X(t)| = SUM x(n) --------------- = SUM x(n) sinc[(t-nT)/T]
  35. C              |   n=M1         pi(t-nT)/T    n=M1
  36. C              |t=n+d
  37. C
  38. C The interpolating (sinc) function is Hamming windowed to bandlimit
  39. C the signal to reduce aliasing.
  40. C
  41. C Multiple simultaneous time delays may be efficiently calculated
  42. C by polyphase filtering.  Polyphase filters decimate the unused
  43. C filter coefficients.  See Chapter 6 in C&R for details. 
  44. C==========================================================================
  45. C
  46. C REFERENCES
  47. C
  48. C Crochiere & Rabiner, Multirate Digital Signal Processing,
  49. C P-H 1983, Chapters 2 and 6.
  50. C
  51. C Kroon & Atal, "Pitch Predictors with High Temporal Resolution,"
  52. C ICASSP '90, S12.6
  53. C
  54. C**************************************************************************
  55. C
  56. subroutine ldelay(x, bufr, start, n, d, m, y)
  57. implicit undefined(a-z)
  58. integer n, bufr, start
  59. real x(bufr), y(0:n-1), d
  60. integer nfrac, i, j, k, m, qd, M1, M2
  61. c * 5 fractional delays calculated over a 40 point interpolation (-20 to 19)
  62. parameter (nfrac = 5, M1 = -20, M2 = 19)
  63. integer twelfths(nfrac)
  64. real sinc, wsinc(M1:M2, nfrac), hwin(12*M1:12*(M2+1)), frac(nfrac)
  65. logical first
  66. c
  67. save first, hwin, wsinc
  68. c
  69. data first /.true./
  70. data twelfths/3,    4,          6,   8,          9/
  71. data frac/    0.25, 0.33333333, 0.5, 0.66666667, 0.75/
  72. c
  73. c Generate Hamming windowed sinc interpolating function
  74. c for each allowable fraction.  The interpolating functions
  75. c are stored in time reverse order (i.e., delay appears as
  76. c advance) to allign with the adaptive code book v0 array.
  77. c The interp filters are:
  78. c wsinc(.,1) frac = 1/4 (3/12)
  79. c wsinc(.,2) frac = 1/3 (4/12)
  80. c . .
  81. c wsinc(.,5) frac = 3/4 (9/12)
  82. c
  83. if (first) then
  84.    call ham(hwin, 12*(M2-M1+1)+1)
  85.    do 20 i = 1, nfrac
  86.       do 10 j = M1, M2
  87.          wsinc(j,i) = sinc(frac(i)+j)*hwin(12*j+twelfths(i))
  88. 10       continue
  89. 20    continue
  90.    first = .false.
  91. end if
  92. c
  93. k=qd(d)
  94. c
  95. c Resample:
  96. do 80 i = 0, n-1
  97.    x(start+i) = 0.0
  98.    do 70 j = M1, M2
  99.       x(start+i) = x(start+i) + x(start-m+i+j)*wsinc(j,k)
  100. 70    continue
  101. 80 continue
  102. c
  103. c The v0 array in psearch.FF/pgain.f must be zero above "start"
  104. c because of overlap and add convolution techniques used in pgain.
  105. c
  106. do 100 i = 0, n-1
  107.    y(i) = x(start + i)
  108.    x(start + i) = 0
  109. 100 continue
  110. return
  111. end