修改型離散餘弦轉換

**修改型離散餘弦轉換**（英語：**m**odified **d**iscrete **c**osine **t**ransform，下稱縮寫：**MDCT**）是一種衍生自傅立葉轉換的轉換。MDCT 以第四型離散餘弦轉換（DCT-IV）為基礎，但具有重疊轉換的性質：它在較大資料集合中，用於處理連續資料區塊，其中當前資料區塊的後半段與下一個資料區塊的前半段相重疊。這種重疊情形有助於避免資料區塊邊界產生多餘資料，除了具有離散餘弦轉換的能量壓縮特性外，也使 MDCT 尤其有利於資料壓縮方面應用。藉由以上優勢，MDCT 在音訊壓縮領域中成為應用最為廣泛的破壞性資料壓縮手段，應用於當前大多數音訊編碼格式，如 MP3、杜比數位、Ogg Vorbis、AAC 等等。

MDCT 是由 Princen、Johnson 和 Bradley 於 1987 年提出，其基本原理即**時域混疊消除法**（英語：time-domain aliasing cancellation，縮寫：**TDAC**）是由 Princen 和 Bradley 於 1986 年提出。其他類似的轉換還有如以離散正弦轉換為基礎的修改型離散正弦轉換（英語：modified discrete sine transform，縮寫：MDST）以及其他較少使用的轉換，例如以其他不同類型的 DCT 或 DCT/DST 的組合為基礎的 MDCT。

## 定義

作為一種重疊轉換，MDCT 與其他傅立葉相關的轉換比起有些不一樣，因為它的輸出數為輸入數的一半（而非相同）。特別的是，它是一個線性函式 _F_ : **R**2_N_ -> **R**N（其中 **R** 代表實數集合）。根據公式：

　
: $X_{k}=\sum _{n=0}^{2N-1}x_{n}\cos \left[{\frac {\pi }{N}}\left(n+{\frac {1}{2}}+{\frac {N}{2}}\right)\left(k+{\frac {1}{2}}\right)\right]$

將 2_N_ 個實數 _x_0, ..., _x_2_N_-1 轉換成 _N_ 個實數 _X_0, ..., _X_N-1。（在這個轉換前面的正規化係數，單位不為一定且根據不同的解釋方法而有所改變。下面只考慮 MDCT 跟 IMDCT 的正規化乘積。）

## 反轉換

MDCT 之反轉換為離散餘弦反轉換（英語：**i**nverse **m**odified **d**iscrete **c**osine **t**ransform，縮寫：**IMDCT**）。

MDCT 由於輸入數與輸出數並不相同，乍看之下應不可逆。然而，透過加入造成錯誤並且必須要回復原始資料的子帶重疊區塊中重疊的 IMDCT 部份，MDCT 可具備完全可逆性。這種技術稱為時域混疊消除法。

IMDCT 根據公式：

　
: $y_{n}={\frac {1}{N}}\sum _{k=0}^{N-1}X_{k}\cos \left[{\frac {\pi }{N}}\left(n+{\frac {1}{2}}+{\frac {N}{2}}\right)\left(k+{\frac {1}{2}}\right)\right]$

將 _N_ 個實數 _X_0, ..., _X_N-1 轉換成 2_N_ 個實數 _y_0, ..., _y_2_N_-1。（如同第四型離散餘弦轉換為一正交轉換，其逆轉換具有與正轉換一樣的形式）。具有一般窗式正規化（usual window normalization, 如下所示）的窗式 MDCT（windowed MDCT），其逆轉換的正規化係數必須乘 2（即變成 2/_N_）。

## 計算

雖然直接應用 MDCT 公式的計算複雜度為 O（_N_2），但是透過遞迴的方式將計算分解化，可以將計算複雜度降到 O（_N_ log _N_），並且得出一樣的結果，就如同快速傅立葉轉換一樣。另外也可以透過其他轉換計算 MDCT。通常一個離散傅立葉轉換（或快速傅立葉轉換）或是一個離散餘弦轉換包含的前置處理跟後置處理的複雜度為 O (_N_)。此外，如下所述，任何第四型離散餘弦轉換的演算法都能衍生出一種方法來計算偶數項的 MDCT 和 IMDCT。

## 窗函式

在一般的訊號壓縮應用中，轉換特性可透過利用將窗函式 _w_n（_n_ = 0, ..., 2_N_-1）乘以 MDCT 中的 _x_n 及 IMDCT 中的 _y_n 作進一步的改良，為了避免在 _n_ = 0 and 2_N_ 的邊界點產生不連續並且使函式於這些點可以變得平滑而趨近零（即 window the data _before_ the MDCT and _after_ the IMDCT）。原則上，_x_ 和 _y_ 可以有不同的窗函式，而窗函式再資料區塊轉換時可以改變（尤其是當兩個不同尺寸的資料區塊結合在一起時）。但為簡單起見，我們考慮給相同尺寸的資料區塊應用的窗函式。 
轉換仍然是可逆的（即 TDAC 依然可以正常作用），考慮一個對稱的窗 _w_n = _w_2_N_-1-_n_，只要 _w_ 滿足 Princen-Bradley condition：

　
: $w_{n}^{2}+w_{n+N}^{2}=1$

各種不同應用的窗函式都是一樣的，例如：

　
: $w_{n}=\sin \left[{\frac {\pi }{2N}}\left(n+{\frac {1}{2}}\right)\right]$

為應用於 MP3 和 MPEG - 2 的 AAC 的窗函式，而

　
: $w_{n}=\sin \left({\frac {\pi }{2}}\sin ^{2}\left[{\frac {\pi }{2N}}\left(n+{\frac {1}{2}}\right)\right]\right)$

為 Vorbis 所應用。此外 AC-3 與 MPEG - 4 AAC 均採用了 Kaiser-Bessel derived (KBD) window。值得注意的是，應用在 MDCT 的窗函式與應用於其他類型的訊號分析的窗函式是不同的。因為它們必須滿足 Princen – Bradley condition。原因之一，MDCT 應用了兩次窗函式，包括 MDCT（分析）和 IMDCT（合成）。

## 第四型離散餘弦轉換與原始 TDAC 的關係

由定義可以看出，一個偶數點 _N_ 的 MDCT 基本上就相當於一個輸入位移了 _N_/2 且兩個 N 點資料區塊同時轉換的第四型離散餘弦轉換。而再更進一步的檢視這個等效關係，則重要的特性例如 TDAC 可以很容易地推導出來。為了更明確的定義與第四型離散餘弦轉換的關係，我們必須了解第四型離散餘弦轉換分別對應到偶數點跟奇數點的邊界條件：偶數點在左邊邊界，（約在 _n_=–1/2），奇數點在右邊邊界（約在 _n_=_N_–1/2）諸如此類（而不是 DFT 那樣的週期性邊界）。這些遵從以下所述：

　
: $\cos \left[{\frac {\pi }{N}}\left(-n-1+{\frac {1}{2}}\right)\left(k+{\frac {1}{2}}\right)\right]=\cos \left[{\frac {\pi }{N}}\left(n+{\frac {1}{2}}\right)\left(k+{\frac {1}{2}}\right)\right]$

和

　
: $\cos \left[{\frac {\pi }{N}}\left(2N-n-1+{\frac {1}{2}}\right)\left(k+{\frac {1}{2}}\right)\right]=-\cos \left[{\frac {\pi }{N}}\left(n+{\frac {1}{2}}\right)\left(k+{\frac {1}{2}}\right)\right]$

因此，如果輸入為一組長度為 _N_ 的 _x_ 陣列，我們可以想見這個陣列會擴增至（_x_, –_x_R, –_x_, _x_R, ...），其中 _x_R 代表 _x_ 的反向序列。考慮一個有 2_N_ 個輸入與 _N_ 個輸出的 MDCT。我們將輸入分成 4 組區塊（_a_, _b_, _c_, _d_），每組大小為 _N_/2。如果我們把這些資料位移 _N_/2（即 MDCT 的定義中的 +_N_/2 項），然後（_b_, _c_, _d_）擴增到 _N_ 點第四型離散餘弦轉換的輸入結尾。因此，我們必須將這些資料根據上述邊界條件「疊」（「fold」）回來。是故，一個具有 2_N_ 點輸入（_a_, _b_, _c_, _d_）的 MDCT，就確實的相當於一個具有 _N_ 點輸入（–_c_R–_d_, _a_–_b_R）的第四型離散餘弦轉換，其中 _R_ 的意義與前述相同。（以這種方式，任何用來計算第四型離散餘弦轉換的演算法都可以很直觀的套用在 MDCT 上）。同樣地，上述 IMDCT 公式恰好正是第四型離散餘弦轉換的 1/2，其中輸出位移了 _N_/2 且擴增長度為 2_N_（透過邊界條件）。而第四型離散餘弦逆轉換可以簡單的給定如上述之輸入（–_c_R–_d_, _a_–_b_R）。當透過邊界條件完成位移跟擴增的動作後，我們可以得到：

　
: IMDCT(MDCT(_a_, _b_, _c_, _d_)) = (_a_–_b_R, _b_–_a_R, _c_+_d_R, _c_R+_d_) / 2

（因此有一半的 IMDCT 輸出是多餘的。）

我們現在明白了 TDAC 的運作原理。考慮計算一個 MDCT 的子帶，有 50％的重疊，2_N_ 個資料區塊（_c_, _d_, _e_, _f_）。則得到的 IMDCT 會跟上述類似：(_c_–_d_R, _d_–_c_R, _e_+_f_R, _e_R+_f_) / 2。當加入這些到先前由重疊到一半所產生的 IMDCT 時，反向的項被消除了，並獲得一簡單的已回復的原始資料（_c_, _d_）。

## TDAC 的起源

」時間域混疊消除法」的起源現在很清楚了。邏輯上第四型離散餘弦轉換中根據邊界條件所擴增的輸入資料造成資料發生混疊的現象，正如同頻譜中奈奎斯特頻率會與較低的頻率發生混疊一樣。因此，當 _c_–_d_R 及其他被加進來會有準確的符號（right sign）提供作消除之用。對於**奇數點**_N_（事實上很少使用），_N_/2 不是整數，所以 MDCT 並不是由第四型離散餘弦轉換透過一個簡單的位移就可以互相對應。在這種情況下，額外的位移一半的樣本，代表著 MDCT/IMDCT 等同於第三 / 二型離散餘弦轉換，其分析類似上述。

## 窗式 MDCT 的 TDAC

如上所述，TDAC 性質已被證明適用於普通 MDCT，顯示加入子帶資料區塊的 IMDCT 於重疊的一半部份可回復原始的資料。推導這項窗式 MDCT 逆性質較上述稍微複雜一些，敘述如下。根據前述當 $(a,b,c,d)$ 和 $(c,d,e,f)$ 經過 MDCT 與 IMDCT，並且加上其重疊一半時，我們得到 $(c+d_{R},c_{R}+d)/2+(c-d_{R},d-c_{R})/2=(c,d)$ ，即原始的資料。 
現在我們假設我們將 MDCT 的輸入和 IMDCT 的輸出兩者均乘上一個長度為 2_N_ 的窗函式。如前所述，我們假設一個對稱的窗函式其形式為 $(w,z,z_{R},w_{R})$ 其中 _w_ 跟 _z_ 是長度 _N_/2 的向量，而 _R_ 代表反向。如此則 Princen – Bradley condition 可以寫成： $w^{2}+z_{R}^{2}=(1,1,\ldots )$ ，其中乘法和加法為點對點的運算，或等效成 $w_{R}^{2}+z^{2}=(1,1,\ldots )$ （將 _w_ 跟 _z_ 反向）。 
因此，我們對 $(wa,zb,z_{R}c,w_{R}d)$ 作 MDCT 取代對 $(a,b,c,d)$ 作 MDCT（所有的乘法跟加法均為點對點）。當完成 IMDCT 再乘以窗函式（點對點相乘），最後 _N_ 點的一半，成為：

　
: $(z_{R},w_{R})\cdot (z_{R}c+wd_{R},zc_{R}+w_{R}d)=(z_{R}^{2}c+wz_{R}d_{R},w_{R}zc_{R}+w_{R}^{2}d)$

（請注意，我們不再乘 1/2，因為 IMDCT 的正規化與窗式的例子中的因數 2 有所不同）。同樣的，對於產生 $(c,d,e,f)$ 的窗式 MDCT 和窗式 IMDCT，是在其前 _N_ 點的一半：

　
: $(w,z)\cdot (wc-z_{R}d_{R},zd-w_{R}c_{R})=(w^{2}c-wz_{R}d_{R},z^{2}d-w_{R}zc_{R})$

當我們這兩個各半的結果加總起來，我們得到：

　
: $(z_{R}^{2}c+wz_{R}d_{R},w_{R}zc_{R}+w_{R}^{2}d)+(w^{2}c-wz_{R}d_{R},z^{2}d-w_{R}zc_{R})$

　
: $=\left([z_{R}^{2}+w^{2}]c+[wz_{R}-wz_{R}]d_{R},[w_{R}^{2}+z^{2}]d+[w_{R}z-w_{R}z]c_{R}\right)=(c,d)$

即回復的原始資料。

## 音訊編碼應用

就 MP3 格式而言，MDCT 並不直接用於處理音訊訊號，而是用於處理 32 波段多相正交濾波器（英語：polyphase quadrature filter，縮寫：PQF）陣列輸出端訊號。MDCT 輸出結果會再經由一條混疊削減公式作後處理，以減弱 PQF 當中典型的混疊。這種轉換與濾波器陣列組合，被稱作混合濾波器陣列或子帶 MDCT。

AAC 則相反，通常純粹採用 MDCT；其中唯一例外是索尼的罕見格式 MPEG – 4 AAC - SSR—— 既採用 MDCT，又在其後用到四波段 PQF 陣列。自適應聽覺轉換編碼（英語：adaptive transform acoustic coding，縮寫：ATRAC）利用運用 MDCT 的堆疊型正交鏡像濾波器（英語：quadrature mirror filter，縮寫：QMF）。

## 參考

* Henrique S. Malvar, _Signal Processing with Lapped Transforms_ (Artech House: Norwood MA, 1992).
* John P. Princen and Alan B. Bradley, "Analysis/synthesis filter bank design based on time domain aliasing cancellation," _IEEE Trans. Acoust. Speech Sig. Proc._ **ASSP-34** (5), 1153-1161 (1986).（Described a precursor to the MDCT using a combination of discrete cosine and sine transforms.）
* J. P. Princen and A. W. Johnson and A. B. Bradley, "Subband/transform coding using filter bank designs based on time domain aliasing cancellation," _IEEE Proc. Intl. Conf. on Acoustics, Speech, and Signal Processing（ICASSP）_ **12**, 2161-2164 (1987).（Initial description of what is now called the MDCT.）
* A. W. Johnson and A. B. Bradley, "Adaptive transform coding incorporating time domain aliasing cancellation," _Speech Comm._ **6**, 299-308 (1987).
* For algorithms, see e.g.:
  * Chi-Min Liu and Wen-Chieh Lee, "A unified fast algorithm for cosine modulated filterbanks in current audio standards", _J. Audio Engineering_ **47** (12), 1061-1075 (1999).
  * V. Britanak and K. R. Rao, "A new fast algorithm for the unified forward and inverse MDCT/MDST computation," _Signal Processing_ **82**, 433-459（2002）
  * Vladimir Nikolajevic and Gerhard Fettweis, "Computation of forward and inverse MDCT using Clenshaw's recurrence formula," _IEEE Trans. Sig. Proc._ **51** (5), 1439-1444（2003）
  * Che-Hong Chen, Bin-Da Liu, and Jar-Ferr Yang, "Recursive architectures for realizing modified discrete cosine transform and its inverse," _IEEE Trans. Circuits Syst. II: Analog Dig. Sig. Proc._ **50** (1), 38-45（2003）
  * J.S. Wu, H.Z. Shu, L. Senhadji, and L.M. Luo, "Mixed-radix algorithm for the computation of forward and inverse MDCTs," _IEEE Trans. Circuits Syst. I: Reg. Papers._ **56** (4), 784-794（2009）
  * V. Britanak, "A survey of efficient MDCT implementations in MP3 audio coding standard: retrospective and state-of-the-art," Signal. Process. _**91** (4), 624-672(2011)_
  * ...and references thereof.