LPC解算的burg算法

声音的线性预测系数LPC的快速算法有很多种，其中Burg算法的流图如下：

C++语言实现的源代码如下：

 1 double LPCBurg(const double x[]/*in wave data*/, size_t N, double a[]/*out LPC*/, size_t M)
 2 {
 3     std::vector<double> b1 (N, 0);
 4     std::vector<double> b2 (N, 0);
 5     std::vector<double> aa (M, 0);
 6 
 7     // (3)
 8 
 9     double gain = 0.0;
10     for(size_t j = 0; j < N; j++)
11     {
12         gain += x[j] * x[j];
13     }
14 
15     gain /= N;
16     if(gain <= 0)
17     {
18         return 0.0;    // warning empty
19     }
20 
21     // (9)
22 
23     b1[0] = x[0];
24     b2[N - 2] = x[N - 1];
25     for(size_t j = 1; j < N - 1; j++)
26     {
27         b1[j] = b2[j - 1] = x[j];
28     }
29 
30     for(size_t i = 0; ; ++i)
31     {
32         // (7)
33 
34         double num = 0.0, denum = 0.0;
35         for(size_t j = 0; j < N - i - 1; j++)
36         {
37             num += b1[j] * b2[j];
38             denum += b1[j] * b1[j] + b2[j] * b2[j];
39         }
40 
41         if(denum <= 0)
42         {
43             return 0.0;    // warning ill-conditioned
44         }
45 
46         a[i] = 2.0 * num / denum;
47 
48         // (10)
49 
50         gain *= 1.0 - a[i] * a[i];
51 
52         // (5)
53 
54         for(size_t j = 0; j+1 <= i; j++)
55         {
56             a[j] = aa[j] - a[i] * aa[i - j - 1];
57         }
58 
59         if(i == M-1) break;
60 
61         // (8)  Watch out: i -> i+1
62 
63         for(size_t j = 0; j <= i; j++)
64         {
65             aa[j] = a[j];
66         }
67 
68         for(size_t j = 0; j <= N - i - 2; j++)
69         {
70             b1[j] -= aa[i] * b2[j];
71             b2[j] = b2[j + 1] - aa[i] * b1[j + 1];
72         }
73     }
74 
75     return gain;
76 }

其中x为输入的语音帧数据，N为语音帧的大小，a为输出的LPC系数，N为需要解算的LPC系数阶数，返回增益gain。该函数返回的gain为均方差平方的均值，真实的预测增益G等于sqrt(gain*N)

posted @ 2012-12-27 22:31  Goncely  阅读(2071)  评论(0)    收藏  举报