线性预测编码LPC

第1页LPC系数预测实验目的语音线性预测的基本思想是：一个语音信号的抽样值可以用过去若干个取样值的线性组合来逼近。通过使实际语音抽样值与线性预测抽样值的均方误差达到最小，可以确定唯一的一组线性预测系数。本实验要求掌握LPC原理，利用自相关法，将语音序列加窗，然后对加窗语音进行LP分析，编写程序求12阶线性预测系数。实验原理1、线性预测编码LPC算法由于语音样点之间存在相关性，所以可以用过去的样点值来预测现在或未来的样点值，从而可以通过使实际语音和线性预测结果之间的误差在某个准则下达到最小值来决定唯一的一组预测系数。而这组系数就能反映语音信号的特性，可以作为语音信号特征参数来用于语音编码、语音合成和语音识别等应用中去。假设)(ny是一实数据列，n，我们可以用过去时刻的N个数据来预测当前时刻的数据)(ny，即：n-,)()()(ˆ1NkNknykany（1)这里)(kaN即为预测系数。定义预测误差)(ne为)(ˆ)()(nynyne(2)我们将采用最小均方误差准则来选择)(kaN的值，使得式(3)总误差NE最小。nNkNnNknykanyneE212)()()()((3)这种优化参数)(kaN的方法导致了求解如下的正则方程组NllRlkRkayyNkyyN,,2,,)()()(1(4)这里的)(kRyy是序列)(ny的自相关系数。式（4）可写成如下的矩阵形式：NNNraR(5)其中TNNNNNaaa)(,),2(),1(a(6)TyyyyyyNNRRR)(,),2(),1(r(7)注意到)(kRyy具有)()(jiRjiRyyyy的性质，式(5)中的NR可写成如下形式第2页)0()2()1()2()0()1()1()1()0(yyyyyyyyyyyyyyyyyyNRNRNRNRRRNRRRR(8)2、Levinson-Durbin算法Levinson-Durbin算法是求解正则方程组中的预测系数Na的有效算法。这种算法利用了自相关矩阵中特殊的对称性。注意到)(),(jiRjiRyyN，即对角线上的元素都相等，所以这个自相关矩阵是Toeplitz矩阵。Levinson-Durbin算法利用了Toeplitz矩阵的特点来进行迭代计算。即首先由一阶预测器(1N)开始，计算预测系数)1(1a。然后增加阶数，利用低阶的结果得到下一个高阶的计算结果。根据(4)式求解得到的一阶预测器的预测系数)1(1a是：)0()1()1(1yyyyRRa(9)其最小均方误差是：212111)1(1)0()0()1()1()1(2)0(aRRaRaREyyyyyyyy(10)这里11)1(a是格形滤波器的第一反射系数。下一步是求解二阶预测器的系数)1(2a和)2(2a，并将结果用)1(1a表示，根据式(5)得到的两个方程是:)1()1()2()0()1(22yyyyyyRRaRa)2()0()2()1()1(22yyyyyyRRaRa（11）通过用(9)的解来消去)1(yyR，我们得到解：112112)1()1()2()1(1)0()1()1()2()2(ERaRaRRaRayyyyyyyyyy)1()2()1()1(1212aaaa（12）这样我们得到了二阶预测器的预测系数，我们再次注意到22)2(a是格形滤波器中的第二反射系数。据此类推，我们可以用)1(m阶预测器的预测系数来表示m阶预测器的系数。这样，我们可将m阶预测系数矢量ma写成两矢量的和，也就是mmmmmmmmaaa110)()2()1(daa(13)这里1ma矢量是第)1(m阶预测器的预测系数，)1(m维的矢量1md和标量m是待定的。我们将mm自相关矩阵mR分区如下：第3页)0(111yybtmbmmmRrrRR(14)这里TbmyyyyyybtmRmRmR)()1(,),2(),1(11rr，bm1r的上标b表示TyyyyyymmRRR)1(,),2(),1(1r元素的倒序排列。根据式(13)和(14)，式(5)可以写成如下形式)(0)0(111111mRRyymmmmyybtmbmmrdarrR(15)这是Levinson-Durbin算法中的关键一步，从式(15)中我们得到两个方程111111mbmmmmmmrrdRaR(16))()0(1111mRRyyyymmbtmmbtmdrar(17)由于111mmmraR，由式(16)得到bmmmm1111rRd(18)又由于bm1r仅是1mr倒序排列，且1mR是Toeplitz，因此可得bmbmmraR1(19)即：bmmbmrRa11(20)因此式(18)可写成)1()2()1(11111mmmmbmmmamamaad(21)现在可用式(17)这个方程来求解m，如果我们用式(21)来消去式(17)中的1md，可得)()0(1111mRRyymbtmbmbtmyymarar（22）注意到1111mtmbmbtmarar，由式(22)可得1111)0()(mtmyymbtmyymRmRarar(23)因此，通过用式(21)和(23)的结果，替换式(13)中的值，我们得到了求解预测器系数的Levinson-Durbin迭代算法。实验内容有一段语音信号，采样率为8kHz；加长度为120样点的汉宁窗。编写程序，求12阶LPC系数。实验方法及结果第4页本实验开发工具为MicrosoftVisualStudio2010，采用语音信号为es01_8k_mono.snd,其信号的采样率为8kHz，对其信号数据进行预处理，采用120个采样点长度的汉宁窗进行窗处理，并取整个语音文件中的其中一帧作为实验对象，用C语言编写程序，实现Lenvinson-Durbin算法，从而求解语音序列的12阶LPC系数。VS2010程序运行结构如下：实验代码#includestdio.h#includeiostream#includemath.h#includesample.hintmain(){float*sample,sample_H[N];floatR0,*R;floatE0,EP[M];floatK[M];floatA[M][M];floatCoefficient[M+1];floatsum;第5页inti,j;FILE*fp;fp=fopen(C:\\Users\\Sun\\Desktop\\新建文件夹(2)\\es01_8k_mono.snd,rb);sample=Read_Data(fp);for(i=0;iN;i++){sample_H[i]=sample[i]*(float)(0.5-0.5*cos(2*PI*i/(N-1)));}R0=*correlation(sample_H);R=(correlation(sample_H)+1);E0=R0;K[0]=R[0]/R0;A[0][0]=K[0];EP[0]=(1-K[0]*K[0])*E0;for(i=1;iM;i++){sum=0;for(j=0;j=i-1;j++){sum+=A[j][i-1]*R[i-j-1];}K[i]=(R[i]-sum)/EP[i-1];A[i][i]=K[i];EP[i]=(1-K[i]*K[i])*EP[i-1];for(j=0;j=i-1;j++){A[j][i]=A[j][i-1]-K[i]*A[i-j-1][i-1];}}Coefficient[0]=1.0;printf(LPC系数:\n);printf(a0=%f\t,Coefficient[0]);for(j=1;j=M;j++){Coefficient[j]=-A[j-1][M-1];printf(a%d=%f\t,j,Coefficient[j]);if((j+1)%3==0)printf(\n\n);}printf(\n\n%s\n\n,__FILE__);system(pause);return0;}第6页Sample.h#definePI3.1415926#defineM12#defineN120float*correlation(floata[]){floatstaticrelate[M+1];inti,k;floatsum;for(k=0;k=M;k++){sum=0;for(i=k;iN;i++){sum+=a[i]*a[i-k];}relate[k]=sum;}returnrelate;}float*Read_Data(FILE*fp_speech){shortdata;floatstatica[N];for(inti=0;i=N;i++){fread(&data,sizeof(short),1,fp_speech);a[i]=(float)data/(float)(1024*32);}returna+1;}