/**********************************
語音線性預測分析子程序
輸入：*original 待語音樣點
輸出： a        線性預測系數

注：Framelen   --200 幀長（對應采樣率為8Hz,幀長為
25ms)
    Fshift     --100 幀移
    alpha      --0.93預加重系數
*************************************/

void lpc(int *original,int a)
{
/*預加重*/
float  *data;
int j;
for(j=0;j<Framelen;j++)
   data[j]=original[j]-alpha*original[j-1];


/*乘窗長*/
MultiHamming(data);



/*自相關分析*/
AutoCorrel(Framelen,data,10,a,1);


void MultiHamming(int *data)
{register int i;
for(i=0;i<Framelen;i++)
   data[i]=0.54-0.46*cos(2*PI*i/(Framelen-1));
}



/*******************************************
FUNCTION  自相關

輸入：n  NO.OF POINTS  IN EVERY FRAME 
      m  ORDER OF FILTER (M<21, DEE NOTE *l)
      sl[n] INPUT POINTS SEWENCE IN EVERY FRAME
輸出：a[n]  FILTER  COFFICIECE
注：  r     TIME --DOMAIN AUTOCORRELATION COFFICENCE
      e     RESIDUAL "ENERGY"
      PROGRAME LIMITED M<21 BECAUSE OF DIMENSIONS OF R[.]
**********************************************/
int AutoCorrelate(int n,int *sl,int m,float *a,int mark)
{
   register int i,j;
   float r[21],e[21],k[21],aa[22][21],sum;
for(i=0;i<21;i++)
{
  r[i]=0;  e[i]=0; k[i]=0;
}
for(i=0;i<21;i++)
    for(j=0;j<21;j++)
        aa[i][j]=0;
for(j=0;j<=m;j++)
{
  r[j]=0;
  for(i=0;i<n-j;i++)
    r[j]+=(float)sl[i]*sl[i+j];
}

if(mark==1)
e[0]=r[0];
aa[1][1]=r[1]/(r[0]+0.00001);
k[1]=r[1]/(e[0]+0.00001)  ;
e[1]=(1-k[1]*k[1]*e[0];

for(i=2;i<=m;i++)
{
  sum=0.0；
  for(j=1;j<=i-1;j++)
        sum+=aa[i-1][j]*r[i-1];

  if(e[i-1]==0.0)  return(-1);
  k[i]=(r[i]-sum)/e[i-1];
  aa[i][i]=k[i];
  for(j=1;j<=i-1;j++)
     aa[i][j]=aa[i-1][j]-k[i]*aa[i-1][i-j];
  e[i]=(1.0-k[i]*k[i]*e[i-1];
}
for(j=1;j<=m;j++)
  a[j]=aa[m][j];
return(0);
}