g723源码详细分析(-)
?
完成了g723源代碼的分析,現(xiàn)作一些整理
1 信號高通濾波 Rem_Dc:
?
這個函數(shù)做高通濾波用的.將低頻噪聲濾除
濾波器的系統(tǒng)函數(shù)為
H(z)=(1-z^(-1)) / (1 - (127/128)*z^(-1))
將 單位圓上的值代入(即:cos(a)?+?sin(a)*i)
可以看出這個函數(shù)在π達到峰值,在 0 和 2π處是谷值,明顯是高通濾波
?
迭代公式如下
y[n] = x[n] - x[n-1] + (127/128) * y[n-1]
為了避免出現(xiàn)除法運算,所有的值都擴大了
?
x[n] - x[n - 1] 即為itu代碼中的fir部分
(127/128) * y[n-1] 即iir部分
?
我們來看看濾波代碼:
?? ? ? ? ? ?/* Do the Fir and scale by 2 */
?? ? ? ? ? ?Acc0 = L_mult( Dpnt[i], (Word16) 0x4000 ) ; //lsc 乘 128 * 128 * 2
?? ? ? ? ? ?Acc0 = L_mac ( Acc0, CodStat.HpfZdl, (Word16) 0xc000 ) ; //lsc 乘-128 * 128 * 2 + acc0 ? 0xc000是 -128 * 128的補碼
?? ? ? ? ? ?CodStat.HpfZdl = Dpnt[i] ; //lsc 完成濾波器零點部分運算 acc0 = 256 * 128 * x[n] - 256 * 128 * x[n-1]?
?
?? ? ? ? ? ?/* Do the Iir part */
?? ? ? ? ? ?Acc1 = L_mls( CodStat.HpfPdl, (Word16) 0x7f00 ) ; //lsc y[n-1] * (127 * 256 / 128 * 256) ?y[n - 1]在第一輪計算中已經(jīng)和x[n]處于一個數(shù)量級了
?? ? ? ? ? ?Acc0 = L_add( Acc0, Acc1 ) ; ?//lsc 完成零極點和,完成了整個濾波運算, 顯然,這里的運算結果擴大了 256 * 128倍
?? ? ? ? ? ?CodStat.HpfPdl = Acc0 ; ? //lsc 保存y[n] 至寄存器,形成下一步運算的y[n-1]
?? ? ? ? ? ?Dpnt[i] = round(Acc0) ; ? //lsc 取出運算結果的高16位
?
最終的結果是 輸入信號縮小了一半
?
從另一個分支也可以看出
?? ? ? Dpnt[i] = shr( Dpnt[i], (Word16) 1 ) ;
如果不做高通濾波,就直接把信號縮小一半
?
?
補充說明?
L_mult 計算兩個相乘，并將結擴大2倍
L_mac(var1, var2, var3) 令后兩個參數(shù)var2, var3相乘，并擴大兩倍,再與第一個參數(shù)var1相加,
為了使加法運算在同一個數(shù)量級 var1必須預先擴大2倍
L_mls(var1,var2) 完成var1 var2的乘法,運算是這么進行的,先計算var1低16與 var2的乘法,再計算var1高16與var2乘法結果
對低16的計算結果做一些舍棄,最終的結果是縮小至真實乘法值的2^(-15)
round 取出一個數(shù)的高16位
?
?
?
?
2 計算lpc系數(shù) Comp_Lpc
?
程序會保留之前幀的120個樣點,與當前幀的240個樣點,組合成一個新數(shù)組,lpc系統(tǒng)是針對這個數(shù)組進行計算的
lpc系數(shù)分四組運算,每組取180個樣點,組與組之間存在著120個樣點的交疊
?
ShAcf_sf[k] = Vec_Norm( Vect, (Word16) LpcFrame ) ;
這段代碼對信號進行歸一化處理,即找出絕對值最大的輸入信號,將其歸一化,即計算該輸入信號移到"最大"值所需要
左移次數(shù)n,然后將180個樣點進行左移n-3,完成了歸一化,該函數(shù)同時返回 n-3,(實際上是保留了"3空位",可能是為了防上后繼計算出現(xiàn)溢出)
?
可以看到itu使用了海明窗對信號進行截取,
窗函數(shù)的作用:截取信號,并盡量不影響信號的頻域特性,所以窗函數(shù)的頻域特性,應該是主瓣窄,并且旁瓣低,由于
同時滿足兩者有沖突,比如頻域上是個沖激(這是最理想的,不會改變信號的任何頻域特性),但在時域則為無窮個點,這種窗顯然沒有任何實際意義
于是數(shù)學家們設計了一些窗函數(shù),分別適用于某種特定應用場合,海明窗是其中的一個
?
itu中的海明窗是放大了2^15的,對應的乘法則縮小2^(-15),所以加窗過程沒有數(shù)值做任何放大
加窗代碼如下
?? ? ? ?/* Apply the Hamming window */ //對信號加一個海明窗
?? ? ? ?for ( i = 0 ; i < LpcFrame ; i ++ )
?? ? ? ? ? ?Vect[i] = mult_r(Vect[i], HammingWindowTable[i]) ;
?
然后計算零時自相關,即能量,根據(jù)柯西定理,也可以知道,零時自相關是所有自相關里絕對值最大的
代碼如下
?? ? ? ?/* Compute the zeroth-order coefficient (energy) */ //lsc 零時自相關,即能量
?? ? ? ?Acc1 = (Word32) 0 ;
?? ? ? ?for ( i = 0 ; i < LpcFrame ; i ++ ) {
?? ? ? ? ? ?Acc0 = L_mult( Vect[i], Vect[i] ) ;
?? ? ? ? ? ?Acc0 = L_shr( Acc0, (Word16) 1 ) ;
?? ? ? ? ? ?Acc1 = L_add( Acc1, Acc0 ) ;
?? ? ? ?}
可以看出,計算結果沒有做任何放大,能量還會加1/1024噪聲能量
?
然后將能量歸一化
?? ? ? ?/* Normalize the energy */ ?//acc1,能量歸一化,即左移至頂,需要多少次乘2
?? ? ? ?Exp = norm_l( Acc1 ) ;
?? ? ? ?Acc1 = L_shl( Acc1, Exp ) ;//能量歸一化
取出歸一化后能量的高16位
?? ? ? ?curAcf[0] = round( Acc1 ) ; //取出能量的高16位
如果能量的高16位為零(即能量為零),所有的自相關都認為是零
?
計算10個自相關系數(shù),用于萊文森-德賓遞推公式使用,代碼片段如下:
?? ? ? ? ? ?for ( i = 1 ; i <= LpcOrder ; i ++ ) {
?? ? ? ? ? ? ? ?Acc1 = (Word32) 0 ;
?? ? ? ? ? ? ? ?for ( j = i ; j < LpcFrame ; j ++ ) {
?? ? ? ? ? ? ? ? ? ?Acc0 = L_mult( Vect[j], Vect[j-i] ) ;
?? ? ? ? ? ? ? ? ? ?Acc0 = L_shr( Acc0, (Word16) 1 ) ;
?? ? ? ? ? ? ? ? ? ?Acc1 = L_add( Acc1, Acc0 ) ;//lsc 計算自相關,沒有任何放大
?? ? ? ? ? ? ? ?}
?? ? ? ? ? ? ? ?Acc0 = L_shl( Acc1, Exp ) ;/* lsc 已經(jīng)歸一化 */
?? ? ? ? ? ? ? ?Acc0 = L_mls( Acc0, BinomialWindowTable[i-1] ) ;/* 一個二項式窗函數(shù) 為何要加??? --- 加窗后的值應該沒有做改變的 */
?? ? ? ? ? ? ? ?curAcf[i] = round(Acc0) ;//lsc 取出高16位
?? ? ? ? ? ?}
?? ? ? ? ? ?/* Save Acf scaling factor */
?
?? ? ? ? ? ?ShAcf_sf[k] = add(Exp, shl(ShAcf_sf[k], 1));
?? ? ? ? ? ?//lsc歸一化的時候,本身移動了ShAcf_sf[k],而自相關為歸一化后的信號相乘,意味著
?? ? ? ? ? ?//得到的結果本身已經(jīng)擴了2^(ShAcf_sf[k]*2),歸一化又移動了Exp,自然要更新相應的ShAcf_sf[k] = ShAcf_sf[k] * 2 + Exp
?
?
?
Durbin --- 萊文森-德賓遞推算法,求出10個lpc系數(shù)
呼叫時傳入的參數(shù)
?? ? ? ?Durbin( &UnqLpc[k*LpcOrder], &curAcf[1], curAcf[0], &Pk2 );
第一個參數(shù)不解釋了,將記錄Durbin函數(shù)計算得到的10個lpc系數(shù)
&curAcf[1] 為自相關系數(shù)組,用于萊文森-德賓遞推,curAcf[0]為零時自相關,即能量,做為
遞推時分母的第一個因子,Durbin函數(shù)的返回值,為殘差信號的能量 &Pk2記錄部分相關系數(shù)k(該參數(shù)似乎是在生成舒適噪音時才使用的)
?
?
萊文森-德賓遞推公式證明:
?
首先從lpc系數(shù)的求解開始說吧
?
我們假定 s[n]是輸入的語音信號
s`[n]是10階預測信號
s`[n] = a10 * s[n - 10] + a9 * s[n - 9] + ... + a1 * s[n - 1]
?
我們取s`[n] 與 s[n] 方差最小值
?
(s[n] - s`[n])^2 ? ^2表示平方
然后對每個 a[i]求偏導，自然就得到了一個10元一次方程組,表示
?
10
Σa[i] * R(|k-i|) = R(k) ?k = 1,2,...,10 ?R(k)表示輸入信號的自相關 --- ?方程組1
i=1
?
即萊文森-德賓遞推公式就是一種適合于計算機實現(xiàn)的解這個10元一次方程組的算法
?
?
下面提到內積和正交的概念
A(z) B(z)分別表示前向預測與后向預測的逆濾波器的系統(tǒng)函數(shù)
則它們關于輸入s[n]的內積這么定義
?10 11 ? ? ? ? ??
?Σ Σ a[i]*b[k] R(|i-k|) ? R(n)同樣表示輸入信號s[n]的自相關函數(shù)
i=1 k=1
?
記作 <A(z),B(z)>
如果內積為零，則被稱之為正交
?
這里可以立即得出這么一個結論
?
<A(z), z^(-l)>一定為零,這個參見方程組1?
也就是說A(z)與z^(-l) l=1,2,...10正交時,A(z)是該階次下最優(yōu)估計的逆濾波器
?
開始遞推,首先從零開始,零階時,啥都沒有
最優(yōu)的逆濾波器自然就是
A0(z)=1
B0(z)=z^-1
?
構造出遞推公式
A[i](z) = A[i-1](z) + ki * B[i - 1](z)
B[i](z) = z^-1{ B[i - 1](z) + ki * A[i-1](z) } ? ? ? ? ?----- 等式2
?
其中 ki = -{ <A[i-1](z),B[i - 1](z)> } / { B[i - 1](z), B[i - 1](z) } ?---- 等式1
?
現(xiàn)在只需要證 A[i](z) 與 z^-i正交即可 (當然B[i](z)也要與z^-i,兩個證明的過程差不多)
將A[i](z) = A[i-1](z) + ki * B[i - 1](z)代入 <A[i](z), z^-i>
我們立該得到
<A[i-i](z), z^-i> + ki<B[i - 1](z), z^-i> = 0;
求出滿足這個條件的ki就是了
<A[i-i](z), z^-i> 實際上與 <A[i-i](z), B[i-1](z)>是相等的,為什么呢？
因為A[i-1](z)是最優(yōu)估計，那它一定與 z^-l (l=1,2,...,i-i)正交,而 B[i-i](z)的最高階系數(shù)是1...
同理<B[i - 1](z), z^-i> 與<B[i-i](z), B[i-1](z)>也是相等的...
自然ki就是等式1的那種形式,
?
萊文森-德賓遞推公式至此證明完畢
?
ki遞推公式化簡
?? ? ? ? ? ? ? ? ? ? ?i-1
分子可以化簡為 R[i] + Σ a(m-1)[n] * R(|l-i|) 這個可以用<A[i-i](z), z^-i> 直接推導出
?? ? ? ? ? ? ? ? ? ? ?l=1
?
分母可以化簡為 (1-k(l)^2)*(1-k(l-1)^2) ... (1-k(1)^2)*R(0)
這是由于 <z^-1 * F(z), z^-1 * G(z)> = <F(z), G(z)> = <F(1/z), G(1/z)> 由內積的定理可以直接證出
用 <B[i](z), B[i](z)> ?等式2進行代換,再進行相應的合并同類項處理,就可以得到化簡后的分子
?
遞推過程中,相應的更新每一輪的自相關系數(shù):
a(i+1)(j) = ai(j) + k ai(m-j)?
這個是從 A(i+1)(z) = Ai(z) + k Bi(z),推出的
我們可以證明 Bi(z) 的系數(shù) 與 Ai(z)的系數(shù)時"相反"的, 即把 Ai(z)的系數(shù)逆序排序其實就是Bi(z)的系數(shù)
?
?
Durbin函數(shù)代碼如下:
?
Word16 ?Durbin( Word16 *Lpc, Word16 *Corr, Word16 Err, Word16 *Pk2 )
{
?? ?int ? i,j ? ;
?
?? ?Word16 ? Temp[LpcOrder] ;
?? ?Word16 ? Pk ;
?
?? ?Word32 ? Acc0,Acc1,Acc2 ;
?
?/*
??* Initialize the LPC vector
??*/
?? ?for ( i = 0 ; i < LpcOrder ; i ++ )
?? ? ? ?Lpc[i] = (Word16) 0 ;
?
?/*
??* Do the recursion. ?At the ith step, the algorithm computes the
??* (i+1)th - order MMSE linear prediction filter.
??*/
?? ?for ( i = 0 ; i < LpcOrder ; i ++ ) {
?
/*
?* Compute the partial correlation (parcor) coefficient
?*/
?
?? ? ? ?/* Start parcor computation */
?? ? ? ?Acc0 = L_deposit_h( Corr[i] ) ;/* 將值擴至32位,對應遞推中的R(m)值 */
?? ? ? ?Acc0 = L_shr( Acc0, (Word16) 2 ) ;/* 右移2位,縮小4倍,因為L_msu擴了2 ??? */
?? ? ? ?for ( j = 0 ; j < i ; j ++ )
?? ? ? ? ? ?Acc0 = L_msu( Acc0, Lpc[j], Corr[i-j-1] ) ;/* 即ai * R(|i-m|),acc0就是k的分子 */
?? ? ? ?Acc0 = L_shl( Acc0, (Word16) 2 ) ;
?
?? ? ? ?/* Save sign */
?? ? ? ?Acc1 = Acc0 ;
?? ? ? ?Acc0 = L_abs( Acc0 ) ;
?
?? ? ? ?/* Finish parcor computation */
?? ? ? ?Acc2 = L_deposit_h( Err ) ;
?? ? ? ?if ( Acc0 >= Acc2 ) {
?? ? ? ? ? ?*Pk2 = 32767;
?? ? ? ? ? ?break ;
?? ? ? ?}
?
?? ? ? ?Pk = div_l( Acc0, Err ) ;/* lsc k的分子除分母 得出k pk就是k */
?
?? ? ? ?if ( Acc1 >= 0 )
?? ? ? ? ? ?Pk = negate(Pk) ;//lsc 負負得正,正負得負
?
?/*
??* Sine detector
??*/
?? ? ? ?if ( i == 1 ) *Pk2 = Pk;
?
?/*
??* Compute the ith LPC coefficient
??*/
?? ? ? ?Acc0 = L_deposit_h( negate(Pk) ) ;
?? ? ? ?Acc0 = L_shr( Acc0, (Word16) 2 ) ;
?? ? ? ?Lpc[i] = round( Acc0 ) ;/* a[m]=k */
?
?/*
??* Update the prediction error
??*/ //lsc 這段一比較繞,為了節(jié)省計算量itu使用了一些技巧,讀者需要注意
?? ? ? ?Acc1 = L_mls( Acc1, Pk ) ;//lsc (-分子原值) * k 得出的是 -k^2 * Err ? k * Err = 分子
?? ? ? ?Acc1 = L_add( Acc1, Acc2 ) ;//lsc Err(Acc2) + ???
?? ? ? ?Err = round( Acc1 ) ;/* lsc 更新了k分母 */
?
?
?/*
??* Compute the remaining LPC coefficients
??*/
?? ? ? ?for ( j = 0 ; j < i ; j ++ )
?? ? ? ? ? ?Temp[j] = Lpc[j] ;
?
?? ? ? ?for ( j = 0 ; j < i ; j ++ ) { /* lsc 更新lpc系數(shù) anew(i) = ai + ka(m-i), 這個是由A(i+1)(z) = Ai(z) + kBi(z) */
?? ? ? ? ? ?Acc0 = L_deposit_h( Lpc[j] ) ;
?? ? ? ? ? ?Acc0 = L_mac( Acc0, Pk, Temp[i-j-1] ) ;
?? ? ? ? ? ?Lpc[j] = round( Acc0 ) ;
?? ? ? ?}
?? ?}
?
?? ?return Err ;//lsc 返回的是最后一個分母,即殘差信號的能量,這個返回值被用于計算cng時的增益估計
}
?
?
ok 到這里已完成了10 lpc系數(shù)計算了,
由于一個lpc系數(shù)的誤差,會影響信號的每個頻段(這一點由A(z)的形式很容易想出)
畢竟 A(z) = P1 * Z^(-10) + P2 * z^(-9) + ....
這個系統(tǒng)函數(shù)不是以因式分解形式給出的,必然有些問題
我們必須把它轉化為因式分解的形式,然后對其求根,而每個根的誤差只會影響某頻域的能量,這樣
就適合進行矢量量化了
?
于是構造了lsf系數(shù),它實際是A(z)的變形,并且能把所有的根者限定的單位圓上,這樣求lsf系數(shù)就可以
簡單地沿著單位圓進行暴力搜索.(解一元10次方程是很復雜的,所以itu將求根過程簡化了)
?
至此,已經(jīng)分析完了g723的整個lpc求解過程,
筆者將在后繼章節(jié)中講述 lpc轉成lsf的過程,
?
完成了聲道系數(shù)求解,剩下的就是激勵編碼了,這些筆者均會在后繼章節(jié)中進行詳細分析
?
?
未完待繼...
?
版權木有,筆者不對本文轉載產生的任何后果負責(如耳麥炸裂,顯示器花屏等)
?
?? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? 林紹川
?? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? ? 2011-04-15 于杭州
?
?
?
?
總結
以上是生活随笔為你收集整理的g723源码详细分析(-)的全部內容，希望文章能夠幫你解決所遇到的問題。
- 上一篇： Play framework如何使用pl
- 下一篇： 使用redis做消息队列mq的总结