语音压缩

开发平台：
Unix_Linux

find.f：源码内容
							c==========================================================================
c
c ROUTINE
c               find
c
c FUNCTION
c                
c               computes filter coefficients, cepstral coefficients, and
c		filter coefficient autocorrelations
c
c SYNOPSIS
c               subroutine find(m,nf,r,cep,ra,alpha,a,rc)
c
c   formal 
c
c                       data    I/O
c       name            type    type    function
c       -------------------------------------------------------------------
c       m		int	i	filter order
c	nf		int	i	number of terms to be found for
c					the cepstrum
c	r		real	i	auto correlation sequence
c	cep		real	o	cepstral coefficients
c	ra		real	o	filter autocorrelation sequence
c	a		real	o	filter coefficients
c	rc		real	o	reflection coefficients
c
c==========================================================================
c
c DESCRIPTION
c
c	See references.  For use with dist.f
c
c==========================================================================
c 
c REFERENCES
c
c	"Distance Measures for Speech Processing", A.H. Gray  
c	and J.D. Markel,IEEE Trans. on ASSP, Vol. ASSP-24, 
c	no. 5, Oct. 1976
c
c	"Quantization and Bit Allocation in Speech Processing",
c	A.H. Gray and J.D. Markel,IEEE Trans. on ASSP, Vol. ASSP-24
c	no. 6, Dec. 1976
c
c	"A Note on Quantization and Bit Allocation in Speech Processing",
c	A.H. Gray and J.D. Markel,IEEE Trans. on ASSP, Vol. ASSP-25
c	no. 3, June 1977
c		
c**************************************************************************
c
	subroutine find(m,nf,r,cep,ra,alpha,a,rc)
	implicit undefined(a-z)
	include 'ccsub.h'
convex	#include "ccsub.h"
	integer j, jl, jm, k, kb, l, lb, m, mh, mp, nf
	real alpha, at, q
	real r(2*maxno+1),cep(maxl*6),ra(2*maxno+1)
	real a(2*maxno+1),rc(2*maxno+1)
	mp=m+1
	a(1)=1.
	if (r(1) .eq. 0) then
	   print *,' find:  r(1)=0, resetting to 1e-6'
	   r(1)=1e-6
	end if
	a(2)=-r(2)/r(1)
	rc(1)=a(2)
	alpha=r(1)*(1.-a(2)*a(2))
	do 450 j=2,m
	   mh=j/2
	   jm=j-1
	   q=r(j+1)
	   do 420 l=1,jm
	      lb=j+1-l
	      q=q+a(l+1)*r(lb)
420	   continue
	   q=-q/alpha
	   rc(j)=q
	   do 430 k=1,mh
	      kb=j-k+1
	      at=a(k+1)+q*a(kb)
	      a(kb)=a(kb)+q*a(k+1)
	      a(k+1)=at
430	   continue
	   a(j+1)=q
	   alpha=alpha*(1.-q*q)
c	   *kill job if unstable filter
	   if (alpha .le. 0.) then
	      stop ' find:  unstable filter'
	   end if
450	continue
c
c	*** evaluation of cepstrum
c
	cep(1)=a(2)
	do 455 j=2,m
	   cep(j)=float(j)*a(j+1)
	   jm=j-1
	   do 455 k=1,jm
	      kb=j-k+1
	      cep(j)=cep(j)-cep(k)*a(j-k+1)
455	continue
	if(nf .le. m) goto 480
	do 460 j=mp,nf
	   cep(j)=0.
	   do 460 k=1,m
	      cep(j)=cep(j)-cep(j-k)*a(k+1)
460	continue
	do 470 j=1,nf
	   cep(j)=-cep(j)/float(j)
470	continue
c
c	*** evaluation of polynomial autocorrelation
c
480	do 500 l=1,mp
c			*bug fix (see last reference)
c	   k=mp+l-1
	   k=mp-l+1
	   ra(l)=0.
	   do 500 j=1,k
	      jl=l+j-1
	      ra(l)=ra(l)+a(j)*a(jl)
500	continue
c	write(*,555)(rc(j),j=1,m)
c555	format(4f10.6)
	return
	end

						