线性预测中的自相关系数

线性预测中的自相关系数1．原理线性预测是语音编码中的基本算法，其基本原理如下：设语音信号的样值序列{}()12,,,k k k X X x x x ==，第k 时刻的取样值x k 可以用之前的P 个样值的线性组合来预测。

1ˆPk i k i i xa x -==∑ 实际样值与预测值之间的误差为：1ˆPk k k k i k i i e x xx a x -==-=-∑ 因此预测系统的传递函数为：()()()()1111P ii i X Z H Z E Z A Z a z -====-∑ 其中H(Z)是一个全极点滤波器，称为综合滤波器。

A(Z)是H(Z)的逆滤波器，称为分析滤波器。

在语音线性预测编码中，A(Z)的系数反映了声道特性。

为了使预测误差最小，采用最小均方误差准则，即使误差的均方值[]{}2221211ˆPkk k k k P P k i k i k i E x x E e E x a x σ=-===-=⎛⎫=- ⎪⎝⎭∑∑∑最小。

在预测阶数P 给定后，2k σ就是所有预测系数{}i a 的函数，因此：[]2ˆˆ20k k k k i i e x E x x a a ⎧⎫∂∂=--=⎨⎬∂∂⎩⎭ []{}ˆ0k k k i E x xx -⇒-= 可见，要使k e 的预测误差最小，则k e 必须与所有数据k i x -正交，称为正交性原理。

将上式展开，可得：{}{}1Pk k i j k j k i j E x x a E x x ---==∑其中{}(),k j k i E x x R k j k i --=--，即信号的自相关系数。

对平稳信号（语音信号一般不是平稳信号，但对单独处理的每帧来说，可以近似认为是短时平稳信号）来说，()()()(),R k i k j R i j R i R i --=--=因此，可以得到：()()()()()()()()()()()()1210112122120P R R R R P a R R R R P a R P R P R P R a -⎡⎤⎡⎤⎡⎤⎢⎥⎢⎥⎢⎥-⎢⎥⎢⎥⎢⎥=⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥⎢⎥--⎢⎥⎢⎥⎣⎦⎣⎦⎣⎦ 解此线性方程组，即可得到各预测系数。

一般采用Levison-Durbin 算法递推求解。

在已知输入信号的情况下，必须先求各阶自相关系数。

2．G .729中的线性预测ITU-T G .729语音编码标准采用CS-ACELP 混合编码方法，输出码率为8kbps ，每帧语音长度10ms ，在8kHz 采样条件下，具有80个样点。

在线性预测时一般采用连续3帧进行加窗处理后进行，因此计算自相关系数的数组具有240个样点值。

下面是G.729参考代码中的相关程序部分，并附加说明。

void Autocorr(Word16 x[], /* (i) : Input signal */Word16 m, /* (i) : LPC order */ Word16 r_h[], /* (o) : Autocorrelations (msb) */Word16 r_l[] /* (o) : Autocorrelations (lsb) */){Word16 i, j, norm;Word16 y[L_WINDOW];Word32 sum;extern Flag Overflow;/* Windowing of signal */for(i=0; i<L_WINDOW; i++){y[i] = mult_r(x[i], hamwindow[i]);}/* Compute r[0] and test for overflow */do {Overflow = 0;sum = 1; /* Avoid case of all zeros */for(i=0; i<L_WINDOW; i++)sum = L_mac(sum, y[i], y[i]);/* If overflow divide y[] by 4 */if(Overflow != 0){for(i=0; i<L_WINDOW; i++){y[i] = shr(y[i], 2);}}}while (Overflow != 0);/* Normalization of r[0] */norm = norm_l(sum);sum = L_shl(sum, norm);L_Extract(sum, &r_h[0], &r_l[0]); /* Put in DPF format (see oper_32b) *//* r[1] to r[m] */for (i = 1; i <= m; i++){sum = 0;for(j=0; j<L_WINDOW-i; j++)sum = L_mac(sum, y[j], y[j+i]);sum = L_shl(sum, norm);L_Extract(sum, &r_h[i], &r_l[i]);}return;}/*___________________________________________________________________________ | | | Function Name : mult_r || | | Purpose : | | | | Same as mult with rounding, i.e.: || mult_r(var1,var2) = shr(((var1*var2) + 16384),15) and || mult_r(-32768,-32768) = 32767. || | | Complexity weight : 2 || | | Inputs : | | | | var1 | | 16 bit short signed integer (Word16) whose value falls in the || range : 0xffff 8000 <= var1 <= 0x0000 7fff. || | | var2 | | 16 bit short signed integer (Word16) whose value falls in the || range : 0xffff 8000 <= var1 <= 0x0000 7fff. || | | Outputs : | | | | none | | | | Return Value : || | | var_out || 16 bit short signed integer (Word16) whose value falls in the || range : 0xffff 8000 <= var_out <= 0x0000 7fff. ||___________________________________________________________________________| */Word16 mult_r(Word16 var1, Word16 var2){Word16 var_out;Word32 L_produit_arr;L_produit_arr = (Word32)var1 * (Word32)var2; /* product */L_produit_arr += (Word32) 0x00004000; /* round */L_produit_arr &= (Word32) 0xffff8000L;L_produit_arr >>= 15; /* shift */if (L_produit_arr & (Word32) 0x00010000L) /* sign extend when necessary */{L_produit_arr |= (Word32) 0xffff0000L;}var_out = sature(L_produit_arr);return(var_out);}/*___________________________________________________________________________ | | | Function Name : norm_l | | | | Purpose : | | | | Produces the number of left shift needed to normalize the 32 bit varia- || ble l_var1 for positive values on the interval with minimum of || 1073741824 and maximum of 2147483647, and for negative values on the in-|| terval with minimum of -2147483648 and maximum of -1073741824; in order || to normalize the result, the following operation must be done : || norm_L_var1 = L_shl(L_var1,norm_l(L_var1)). | | | | Complexity weight : 30 | | | | Inputs : | | | | L_var1 | | 32 bit long signed integer (Word32) whose value falls in the || range : 0x8000 0000 <= var1 <= 0x7fff ffff. || || | | none | | | | Return Value : || | | var_out | | 16 bit short signed integer (Word16) whose value falls in the || range : 0x0000 0000 <= var_out <= 0x0000 001f. ||___________________________________________________________________________| */Word16 norm_l(Word32 L_var1){Word16 var_out;if (L_var1 == 0){var_out = 0;}else{if (L_var1 == (Word32)0xffffffffL){var_out = 31;}else{if (L_var1 < 0){L_var1 = ~L_var1;}for(var_out = 0;L_var1 < (Word32)0x40000000L;var_out++){L_var1 <<= 1;}}}return(var_out);}/*___________________________________________________________________________| Function Name : L_mac | | | | Purpose : | | | | Multiply var1 by var2 and shift the result left by 1. Add the 32 bit || result to L_var3 with saturation, return a 32 bit result: || L_mac(L_var3,var1,var2) = L_add(L_var3,(L_mult(var1,var2)). || | | Complexity weight : 1 || | | Inputs : | | | | L_var3 32 bit long signed integer (Word32) whose value falls in the || range : 0x8000 0000 <= L_var3 <= 0x7fff ffff. || | | var1 | | 16 bit short signed integer (Word16) whose value falls in the || range : 0xffff 8000 <= var1 <= 0x0000 7fff. || | | var2 | | 16 bit short signed integer (Word16) whose value falls in the || range : 0xffff 8000 <= var1 <= 0x0000 7fff. || | | Outputs : | | | | none | | | | Return Value : || | | L_var_out | | 32 bit long signed integer (Word32) whose value falls in the || range : 0x8000 0000 <= L_var_out <= 0x7fff ffff. ||___________________________________________________________________________| */Word32 L_mac(Word32 L_var3, Word16 var1, Word16 var2){Word32 L_var_out;Word32 L_produit;L_produit = L_mult(var1,var2);L_var_out = L_add(L_var3,L_produit);return(L_var_out);}/* Hamming_cos window for LPC analysis. *//* Create with function ham_cos(window,200,40) */Word16 hamwindow[L_WINDOW] = {2621, 2623, 2629, 2638, 2651, 2668, 2689, 2713, 2741, 2772,2808, 2847, 2890, 2936, 2986, 3040, 3097, 3158, 3223, 3291,3363, 3438, 3517, 3599, 3685, 3774, 3867, 3963, 4063, 4166,4272, 4382, 4495, 4611, 4731, 4853, 4979, 5108, 5240, 5376,5514, 5655, 5800, 5947, 6097, 6250, 6406, 6565, 6726, 6890,7057, 7227, 7399, 7573, 7750, 7930, 8112, 8296, 8483, 8672,8863, 9057, 9252, 9450, 9650, 9852, 10055, 10261, 10468, 10677,10888, 11101, 11315, 11531, 11748, 11967, 12187, 12409, 12632, 12856,13082, 13308, 13536, 13764, 13994, 14225, 14456, 14688, 14921, 15155,15389, 15624, 15859, 16095, 16331, 16568, 16805, 17042, 17279, 17516,17754, 17991, 18228, 18465, 18702, 18939, 19175, 19411, 19647, 19882,20117, 20350, 20584, 20816, 21048, 21279, 21509, 21738, 21967, 22194,22420, 22644, 22868, 23090, 23311, 23531, 23749, 23965, 24181, 24394,24606, 24816, 25024, 25231, 25435, 25638, 25839, 26037, 26234, 26428,26621, 26811, 26999, 27184, 27368, 27548, 27727, 27903, 28076, 28247,28415, 28581, 28743, 28903, 29061, 29215, 29367, 29515, 29661, 29804,29944, 30081, 30214, 30345, 30472, 30597, 30718, 30836, 30950, 31062,31170, 31274, 31376, 31474, 31568, 31659, 31747, 31831, 31911, 31988,32062, 32132, 32198, 32261, 32320, 32376, 32428, 32476, 32521, 32561,32599, 32632, 32662, 32688, 32711, 32729, 32744, 32755, 32763, 32767,32767, 32741, 32665, 32537, 32359, 32129, 31850, 31521, 31143, 30716,30242, 29720, 29151, 28538, 27879, 27177, 26433, 25647, 24821, 23957,23055, 22117, 21145, 20139, 19102, 18036, 16941, 15820, 14674, 13505,12315, 11106, 9879, 8637, 7381, 6114, 4838, 3554, 2264, 971};3．设计要求利用TMS320C54x汇编实现自相关函数的计算。