王 膂 伍家松 Senhadji Lotfi 舒華忠

(1東南大學(xué)影像科學(xué)與技術(shù)實(shí)驗(yàn)室，南京 210096)

(2雷恩第一大學(xué)信號與圖像處理實(shí)驗(yàn)室，法國雷恩 35042)

(3東南大學(xué)中法生物醫(yī)學(xué)信息研究中心，南京 210096)

音頻壓縮中3種整數(shù)型MDCT變換的比較

王 膂1，3伍家松1，2，3Senhadji Lotfi2，3舒華忠1，3

(1東南大學(xué)影像科學(xué)與技術(shù)實(shí)驗(yàn)室，南京 210096)

(2雷恩第一大學(xué)信號與圖像處理實(shí)驗(yàn)室，法國雷恩 35042)

(3東南大學(xué)中法生物醫(yī)學(xué)信息研究中心，南京 210096)

為了快速計(jì)算整數(shù)型改進(jìn)的離散余弦變換(IntMDCT)，構(gòu)造了基于提升變換、模變換以及無窮范數(shù)旋轉(zhuǎn)變換的3種計(jì)算12點(diǎn)IntMDCT的算法.首先將12點(diǎn)MDCT轉(zhuǎn)化為6點(diǎn)Ⅳ型離散余弦變換(DCT－Ⅳ)，并將后者分解為7個(gè)Givens旋轉(zhuǎn)變換的乘積;然后分別利用提升變換算法、模變換算法和無窮范數(shù)旋轉(zhuǎn)變換算法實(shí)現(xiàn)Givens旋轉(zhuǎn)變換的整數(shù)型近似計(jì)算;最后，對這3種算法在語音信號無損和有損壓縮中的運(yùn)行速度和計(jì)算精確度進(jìn)行比較.實(shí)驗(yàn)結(jié)果表明，在這3種算法中，基于模變換的IntMDCT算法的運(yùn)行速度最快;基于無窮范數(shù)旋轉(zhuǎn)變換的IntMDCT算法的計(jì)算精度最高，并在有損音頻壓縮中獲得的信噪比最高.

提升算法;模變換;無窮范數(shù)旋轉(zhuǎn)變換;整數(shù)型MDCT;音頻壓縮

可逆的整數(shù)型變換可實(shí)現(xiàn)整數(shù)間的映射，因而被廣泛應(yīng)用于音頻與圖像的無損和有損壓縮中.提升算法(lifting scheme)作為一種成熟的整數(shù)變換算法，普遍應(yīng)用于整數(shù)型離散余弦變換(IntDCT)［1］以及整數(shù)型改進(jìn)離散余弦變換(integer modified discrete cosine transform，IntMDCT)［2］中.

改進(jìn)的離散余弦變換(MDCT)是基于第Ⅳ型離散余弦變換(DCT－Ⅳ)的重疊變換.MDCT的重疊性以及與DCT相似的能量壓縮性使其能夠有效地消除信號分塊壓縮時(shí)產(chǎn)生的邊界偽影，因而被廣泛應(yīng)用于音頻壓縮標(biāo)準(zhǔn)中.例如，MPEG－1音頻第3層編碼(MP3)［3］采用 12點(diǎn)和 36點(diǎn) MDCT;MPEG－4高級音頻編碼(AAC)［4］采用 256 點(diǎn)和2 048點(diǎn)MDCT.

近年來，MDCT的快速算法得到了迅速發(fā)展［5－9］.作為 MPEG 無損音頻編碼標(biāo)準(zhǔn)的關(guān)鍵部分，IntMDCT 受到高度重視［2，10－11］.現(xiàn)有的一維IntMDCT一般是先將MDCT矩陣或相應(yīng)的DCT－Ⅳ矩陣分解為若干Givens旋轉(zhuǎn)變換的乘積，再通過提升算法進(jìn)行整數(shù)近似實(shí)現(xiàn).然而，提升算法存在一定的缺陷，如由級聯(lián)產(chǎn)生的延時(shí)和遞增的取整誤差等［12］.針對這些缺陷，Srinivasan［12］提出利用模變換(modulo transform)算法來取代提升算法.模變換算法利用勾股數(shù)三元組代替三角函數(shù)，對Givens旋轉(zhuǎn)變換進(jìn)行整數(shù)型近似，因而較提升算法具有更高的精度.Givens旋轉(zhuǎn)變換在二維平面內(nèi)的2范數(shù)不變性有利于擴(kuò)展輸出信號的動(dòng)態(tài)范圍.Yang等［13］提出了一種保持無窮范數(shù)不變性的無窮范數(shù)旋轉(zhuǎn)變換(infinity norm rotation transform)，該變換可以有效控制輸出信號的動(dòng)態(tài)范圍，且通過取整運(yùn)算來得到整數(shù)至整數(shù)的映射，從而實(shí)現(xiàn)Givens旋轉(zhuǎn)變換的整數(shù)型近似變換.

本文分別將提升算法、模變換和無窮范數(shù)旋轉(zhuǎn)變換用于計(jì)算12點(diǎn) IntMDCT，并對其在有損和無損音頻壓縮中的速度和精確度進(jìn)行比較.

1 改進(jìn)的離散余弦變換

對于長度為N的輸入序列x(n)，其一維MDCT/IMDCT 定義如下［14］:

式中，N為4的倍數(shù).

式(1)中MDCT的系數(shù)X(k)具有反對稱性，即

因此，只有N/2個(gè)線性獨(dú)立的系數(shù)需要計(jì)算.式(2)中得到的(n)是時(shí)域混疊的，并非直接等于原始信號x(n)［14］.

利用文獻(xiàn)［9－10］中的12點(diǎn)MDCT/IMDCT 將MDCT正交分解，信號流圖見圖1［10］.該算法先通過2級疊加將12點(diǎn)MDCT轉(zhuǎn)化為6點(diǎn)DCT－Ⅳ，再利用3次Givens旋轉(zhuǎn)變換將6點(diǎn)DCT－Ⅳ分解為一個(gè)3點(diǎn)第Ⅱ型離散余弦變換(DCT－Ⅱ)和一個(gè)3點(diǎn)第Ⅱ型離散正弦變換(DST－Ⅱ)，其中后者可分別由2次 Givens旋轉(zhuǎn)變換實(shí)現(xiàn).因此，12點(diǎn)MDCT可分解為7個(gè) Givens旋轉(zhuǎn)變換.逆變換IMDCT可通過信號流圖轉(zhuǎn)置得到.圖1中的旋轉(zhuǎn)系數(shù)(Ci，Si)均為無理數(shù)(見表1)，由計(jì)算機(jī)有限精度實(shí)現(xiàn)時(shí)必然產(chǎn)生截?cái)嗾`差，因此，由浮點(diǎn)型MDCT/IMDCT處理并重建的信號與原信號之間存在誤差，無法做到完全重建.

圖112點(diǎn)MDCT/IMDCT的信號流圖

文獻(xiàn)［9－10］中的12點(diǎn)MDCT/IMDCT算法適用于計(jì)算MP3編碼中的短窗MDCT.若直接將上述算法推廣至長窗36點(diǎn)MDCT，會(huì)引入較多的旋轉(zhuǎn)角度，使算法實(shí)現(xiàn)較為困難.為了解決這一問題，可以先采用文獻(xiàn)［6－7］中的算法，將 36 點(diǎn)MDCT轉(zhuǎn)化為3個(gè)長度為12點(diǎn)的MDCT，再對后者運(yùn)用前面介紹的12點(diǎn)算法，便可得到36點(diǎn)MDCT.

表112點(diǎn)MDCT/IMDCT的旋轉(zhuǎn)系數(shù)

2 12點(diǎn)MDCT的整數(shù)實(shí)現(xiàn)

2.1 提升算法

Givens旋轉(zhuǎn)矩陣可以分解為3個(gè)三角矩陣的乘積.將三角矩陣與輸入信號相乘稱為提升運(yùn)算，可用如下公式表示:

式(4)中，(cosθ－1)/sinθ和 sinθ即為提升系數(shù)(見圖2).在每次提升過程中，通過取整可實(shí)現(xiàn)Givens旋轉(zhuǎn)變換的整數(shù)型近似計(jì)算［15］.例如，第2次提升過程對輸入(x1，x2)進(jìn)行映射可表示為

圖2 Givens旋轉(zhuǎn)變換與提升結(jié)構(gòu)實(shí)現(xiàn)

式(5)可近似為

式中，round(·)表示取整.

更新后x1的值并未改變，因此仍能根據(jù)更新后的值計(jì)算出round(x1sinθ).式(6)的逆運(yùn)算為

可見提升結(jié)構(gòu)的整數(shù)近似是可逆的，且逆運(yùn)算無誤差.因此，通過三步提升方法實(shí)現(xiàn)的Givens旋轉(zhuǎn)變換的整數(shù)近似是完全可逆的.

圖1中的7個(gè)Givens旋轉(zhuǎn)變換均可由提升算法實(shí)現(xiàn)［10］，總計(jì)算量為21次整數(shù)加法、21次實(shí)數(shù)乘法和21次取整運(yùn)算.相應(yīng)的提升系數(shù)可事先由有限精度求得，計(jì)算時(shí)查表使用.

2.2 模變換

模變換是通過整數(shù)矩陣]和整數(shù)尺度對Givens旋轉(zhuǎn)矩陣]進(jìn)行整數(shù)型近似的［12］，其中勾股數(shù)三元組滿足.通過限定，可使模變換完全可逆，同時(shí)可減少計(jì)算量［12］.勾股數(shù)三元組(s，c，c+1)所對應(yīng)的旋轉(zhuǎn)角度θ=arctan(2s/(s2－1)).因此，旋轉(zhuǎn)角度為θ的旋轉(zhuǎn)矩陣可以進(jìn)行如下整數(shù)型近似:

模變換的信號流圖及相應(yīng)的逆變換如圖3所示.圖中，運(yùn)算符 sdiv 定義如下［12］:

式中，floor(x)表示取不大于x的最大整數(shù).

圖3 基于勾股數(shù)三元組(s，c，c+1)的模變換及逆變換

根據(jù)式(9)，模變換可表示為［12］

計(jì)算式(10)時(shí)需要4次整數(shù)加法、2次整數(shù)乘法、2次整數(shù)除法和2次取整運(yùn)算.

勾股數(shù)三元組能近似表示的角度見表2［12］.同一個(gè)旋轉(zhuǎn)角度θ可利用不同的勾股數(shù)組進(jìn)行不同精度的近似.此外，還可將θ分解為多個(gè)較小的角度之和，利用多個(gè)小角度模變換的級聯(lián)近似表示角度為θ的模變換.因此，表1中的旋轉(zhuǎn)系數(shù)可由表3中對應(yīng)的勾股數(shù)三元組進(jìn)行近似.圖1中7個(gè)Givens旋轉(zhuǎn)變換可用模變換進(jìn)行整數(shù)型近似，從而代替?zhèn)鹘y(tǒng)的三步提升方法.用模變換計(jì)算7個(gè)Givens旋轉(zhuǎn)變換的總運(yùn)算量為44次整數(shù)加法、22次整數(shù)乘法、22次整數(shù)除法和22次取整運(yùn)算.

表2 勾股數(shù)組及其所能近似表示的角度

表312點(diǎn)MDCT/IMDCT的旋轉(zhuǎn)系數(shù)對應(yīng)的勾股數(shù)三元組

2.3 無窮范數(shù)旋轉(zhuǎn)變換

提升算法和模變換都是對歐幾里得空間內(nèi)具有2范數(shù)不變性的Givens旋轉(zhuǎn)變換的整數(shù)型近似計(jì)算.具有無窮范數(shù)不變性的旋轉(zhuǎn)被稱為無窮范數(shù)旋轉(zhuǎn)變換［13］.在二維空間內(nèi)，2范數(shù)和無窮范數(shù)旋轉(zhuǎn)變換旋轉(zhuǎn)一周的軌跡形狀分別是一個(gè)圓形和一個(gè)正方形(見圖4).二維空間內(nèi)p范數(shù)旋轉(zhuǎn)變換的旋轉(zhuǎn)角θp可定義為游程與半徑的比值［13］.當(dāng)p=2時(shí)，θ2為Givens旋轉(zhuǎn)角，游程為旋轉(zhuǎn)劃過的弧長，即圖4中點(diǎn)A與點(diǎn)B'間的弧長，半徑為圓半徑;當(dāng)p=∞時(shí)，θ∞為無窮范數(shù)旋轉(zhuǎn)角，游程為旋轉(zhuǎn)在正方形邊上劃過的距離，即圖4中點(diǎn)A與點(diǎn)B間的線段長度，半徑為正方形邊長除以2.旋轉(zhuǎn)一周后，θ2=2π，θ∞=8.因此，θ2和 θ∞之間存在如下的對應(yīng)關(guān)系:

大于π/4的角度可以由多個(gè)小角度級聯(lián)構(gòu)成.

圖4 無窮范數(shù)旋轉(zhuǎn)變換的分段線性特性

為了計(jì)算無窮范數(shù)旋轉(zhuǎn)變換，采用4條相交虛線將xoy平面劃分為8個(gè)區(qū)域.每個(gè)區(qū)域內(nèi)的無窮范數(shù)旋轉(zhuǎn)變換均為基本線性變換，即錯(cuò)切變換、反射變換或錯(cuò)切變換(見圖4).這種無窮范數(shù)旋轉(zhuǎn)變換算法可表示為

式中，a為半徑，即L為游程;∧為且運(yùn)算.整數(shù)輸入x，y與實(shí)數(shù)θ∞相乘可得到實(shí)數(shù)輸出.為了實(shí)現(xiàn)整數(shù)到整數(shù)的映射，需要引入一次取整運(yùn)算.整數(shù)型無窮范數(shù)旋轉(zhuǎn)變換的單次計(jì)算量為1次整數(shù)加法、1次實(shí)數(shù)乘法和1次取整運(yùn)算，或2次整數(shù)加法、1次實(shí)數(shù)乘法、1次整數(shù)移位和1次取整運(yùn)算.與提升算法和模變換算法不同的是，這里的取整運(yùn)算會(huì)使該旋轉(zhuǎn)的逆運(yùn)算產(chǎn)生誤差.

由此可見，圖1中的7個(gè)Givens旋轉(zhuǎn)變換可以轉(zhuǎn)化為無窮范數(shù)旋轉(zhuǎn)變換進(jìn)行整數(shù)型近似.

3 算法效果比較

分別對基于提升算法、模變換和無窮范數(shù)旋轉(zhuǎn)變換的12點(diǎn)IntMDCT/IntIMDCT在語音信號無損和有損壓縮中的效果進(jìn)行比較.測試信號為一段時(shí)長為2.9 s、采樣頻率為8 kHz、8 bit整數(shù)編碼的真實(shí)語音序列，其波形圖如圖5(a)所示.效果評價(jià)指標(biāo)為均方誤差(mean square error，MSE)和信噪比(signal to noise ratio，SNR).測試環(huán)境為AMD 3.0 GHz CPU，Windows Server 2008系統(tǒng)和Matlab 7.9.基于提升變換［10］、模變換［12］、無窮范數(shù)旋轉(zhuǎn)變換［13］的IntMDCT分別簡記為算法1、算法2和算法3.

圖5 3種算法在有損語音壓縮中的效果比較

一維信號的均方誤差定義為

式中，em表示第m時(shí)刻時(shí)理想信號Ym與重建信號m之差，即em=Ym－m;M表示信號的長度.

一維信號的信噪比定義為

3.1 無損壓縮

分別對測試音頻信號運(yùn)行這3種算法，得到的變換系數(shù)不經(jīng)量化，直接進(jìn)行12點(diǎn) IntIMDCT并重建原始信號.對3種算法的運(yùn)行時(shí)間和所得重建信號的均方誤差進(jìn)行比較，結(jié)果見表4.由于模變換只涉及整數(shù)計(jì)算，因此算法2的速度最快:其正變換時(shí)間較算法1快92.9%，較算法3快66.0%;其反變換時(shí)間較算法1快87.3%，較算法3快56.6%.算法1和算法2在無損壓縮中不產(chǎn)生任何誤差，但算法3會(huì)產(chǎn)生少量誤差.

表4 3種算法在無損音頻壓縮中的運(yùn)行時(shí)間和均方誤差

3.2 有損壓縮

分別對測試信號運(yùn)行這3種算法，得到的變換系數(shù)進(jìn)行3～8 bit均勻量化.量化公式為

式中，q為量化比特?cái)?shù).對量化后的系數(shù)進(jìn)行逆量化和12點(diǎn)IntIMDCT，并重建信號.

對這3種算法的運(yùn)行時(shí)間和所得重建信號的信噪比進(jìn)行比較，結(jié)果見表5.由于有損壓縮的運(yùn)行時(shí)間與表4相同，此處不再一一列出.由表5可知，由于算法3能有效控制輸出的動(dòng)態(tài)范圍，因此在有損壓縮中得到了最高的信噪比.例如，在低比特(5 bit以下)均勻壓縮時(shí)，算法3實(shí)現(xiàn)的信噪比較算法1和算法2提高了1.5～2.2 dB.此外，算法2的取整運(yùn)算次數(shù)少于算法1，導(dǎo)致其取整誤差略小于算法1.因此，在大部分情況下，算法2在有損壓縮中的信噪比較算法1略高.圖5(b)～(d)給出了5 bit均勻量化時(shí)利用這3種算法得到的重建信號及誤差.圖中黑色部分為信號，灰色部分為誤差.由圖可知，算法3的誤差小于其他2種算法.

表5 3種算法在有損音頻壓縮中的信噪比 dB

4 結(jié)語

本文將提升算法、模變換和無窮范數(shù)旋轉(zhuǎn)應(yīng)用于12點(diǎn)整數(shù)型MDCT，并對基于這3種算法的IntMDCT在語音信號無損和有損壓縮中的速度和精確度進(jìn)行了比較.實(shí)驗(yàn)結(jié)果表明，基于模變換的IntMDCT運(yùn)算速度最快，基于無窮范數(shù)旋轉(zhuǎn)變換的IntMDCT的運(yùn)算速度次之，基于提升算法的IntMDCT算法的運(yùn)算速度最慢.在無損壓縮中，提升算法和模變換均不產(chǎn)生誤差，而無窮范數(shù)旋轉(zhuǎn)變換產(chǎn)生少量誤差.在有損壓縮中，無窮范數(shù)旋轉(zhuǎn)變換由于能夠有效控制輸出的動(dòng)態(tài)范圍，因而能實(shí)現(xiàn)最高的信噪比.

［1］ Zeng Y，Cheng L，Bi G.Integer DCTs and fast algorithms［J］.IEEE Transactions on Signal Processing，2001，49(19):2774－2782.

［2］ Huang H，Rahardja S，Yu R.Integer MDCT with enhanced approximation of the DCT－Ⅳ［J］.IEEE Transactions on Signal Processing，2006，54(11):1156－1159.

［3］林福宗.多媒體技術(shù)基礎(chǔ)［M］.北京:清華大學(xué)出版社，2000.

［4］ Pereira F，Ebrahimi T.The MPEG-4book［M］.Upper Saddle River，USA:Prentice Hall，2002.

［5］Britanak V.New universal rotation — based fast computational structures for an efficient implementation of the DCT－Ⅳ/DST－Ⅳ and analysis/synthesis MDCT/MDST filter banks［J］.Signal Processing，2009，89(11):2213－2232.

［6］ Shu H，Bao X，Toumoulin C.Radix－3 algorithm for the fast computation of forward and inverse MDCT［J］.IEEE Signal Processing Letters，2007，14(10):93－96.

［7］ Wu J，Shu H.Mixed－radix algorithm for the computation of forward and inverse MDCTs［J］.IEEE Transactions on Circuits and SystemsⅠ:Regular Papers，2009，56(4):784－794.

［8］ Britanak V.New fast computational structures for an efficient implementation of the forward/backward MDCT in MP3 audio coding standard ［J］.Signal Processing，2010，90(2):536－547.

［9］ Britanak V.A survey of efficient MDCT implementations in MP3 audio coding standard:retrospective and state－of－the－art［J］.Signal Processing，2011，91(4):624－672.

［10］ Krishnan T，Oraintara S.Fast and lossless implementation of the forward and inverse MDCT computation in MPEG audio coding［C］//IEEE International Symposium on Circuits and Systems.Phoenix，AZ，USA，2002:181－184.

［11］ Li J.Low noise reversible MDCT(RMDCT)and its application in progressive－to－lossless embedded audio coding ［J］.IEEE Transactions on Signal Processing，2005，53(5):1870－1880.

［12］ Srinivasan S.Modulo transforms— an alternative to lifting ［J］.IEEE Transactions on Signal Processing，2006，54(13):1864－1874.

［13］ Yang L，Hao P.Infinity－norm rotation transforms［J］.IEEE Transactions on Signal Processing，2009，57(7):2594－2603.

［14］ Princen J P，Johnson A W，Bradley A B.Subband/transform coding using filter bank designs based on time domain aliasing cancellation［C］//IEEE International Conference on Acoustics，Speech，and Signal Processing.Dallas，Texas，USA，1987:2161－2164.

［15］ Geiger R，Sporer T，Koller J.Audio coding based on integer transforms［C］//The111th Audio Engineering Society Convention.New York，USA，2001:5471－5479.

Comparison of three IntMDCT algorithms in audio compression

Wang Lü1，3Wu Jiasong1，2，3Senhadji Lotfi2，3Shu Huazhong1，3

(1Laboratory of Image Science and Technology，Southeast University，Nanjing 210096，China)
(2Laboratoire Traitement du Signal et de I'Image，Université de Rennes 1，Rennes 35042，F(xiàn)rance)
(3Centre de Recherche en Information Biomédicale Sino－Fran?ais，Nanjing 210096，China)

In order to improve the computation efficiency of the integer modified discrete cosine transform(IntMDCT)，three algorithms based on the lifting scheme，modulo transform and infinity norm rotation transform are formulated respectively for computing the 12－point IntMDCT.First，the 12－point IntMDCT is converted into the 6－point type－Ⅳ discrete cosine transform(DCT－Ⅳ)，which is then factorized into a product of 7 Givens rotation matrices.The integer type Givens rotation matrices are approximated by lifting scheme，modulo transform and infinity norm rotation transform，respectively.Finally，the speed and accuracy of these three IntMDCT algorithms are compared in both lossless and lossy audio compression.The experimental results show that in the three algorithms，the IntMDCT algorithm based on the modulo transform has the highest computation speed.The IntMDCT algorithm based on the infinity norm rotation transform has the highest accuracy，and can achieve the highest signal to noise ratio(SNR)in lossy audio compression.

lifting scheme;modulo transform;infinity norm rotation;integer modified discrete cosine transform;audio compression

TP391

A

1001－0505(2012)02－0259－06

10.3969/j.issn.1001－0505.2012.02.013

2011－09－08.

王膂(1986—)，男，碩士生;舒華忠(聯(lián)系人)，男，博士，教授，博士生導(dǎo)師，shu.list@seu.edu.cn.

國家自然科學(xué)基金資助項(xiàng)目(61073138，60873048)、國家自然科學(xué)基金國際合作與交流資助項(xiàng)目(60911130370).

王膂，伍家松，Senhadji Lotfi，等.音頻壓縮中3種整數(shù)型MDCT變換的比較［J］.東南大學(xué)學(xué)報(bào):自然科學(xué)版，2012，42(2):259－264.［doi:10.3969/j.issn.1001－0505.2012.02.013］