MDCT/MDST, IMDCT与DCT-IV的关系以及快速FFT实现

本文主要是MDCT/MDST/IMDCT与DCT-IV的关系，然后是采用N/4点快速FFT傅里叶变换来实现N点的MDCT/IMDCT的公式推导以及重叠相加的基本原理，并根据MDCT/IMDCT的特点推导满足完全重建的窗函数的特性。

MDCT/MDST/IMDCT与DCT-IV的关系

如下是DCT-IV，MDCT/MDST/IMDCT/IMDST的公式

MDCT/MDST变换的基本特性

相同长度的输入序列

拼接成一个序列

,然后经过前向和后向MDCT变换，则时域重建的信号输出为：


。

这里的

表示序列X1的逆序排列。那么如果输入的时域序列是有50%重叠的，那么



50%重叠相加序列

，那么重叠部分的序列

就得到完全重建

.

从IMDCT(MDCT)的输出来看，MDCT具有前半部分奇对称，后半部分偶对称的特性，而MDST刚好相反，前半部分偶对称，后半部分奇对称。

MDST与MDCT的关系

通过对序列

的MDCT变换，然后对结果做一个重排序，就得到了MDST的变换了。

使用DCT-IV来实现MDCT/MDST, IMDCT

令N=2M



令

,那么根据三角函数的周期特性，有


那么可以得到MDCT的如下结果


那么可以得到

, 即N点的MDCT可以先创建一个重排序的序列，然后通过一个M=N/2点的DCT-IV来实现。

对于IMDCT, 很容易由DCT-IV实现


然后根据MDCT的前半部分的奇对称，后半部分的偶对称特性，有如下结果：



即N点的IMDCT可以转化为N/2点的DCT-IV，然后对计算结果做一个重排序。

对于 MDST，其可以转化为MDCT，然后也就能转化为DCT-IV。



记得到

，然后是一个后排序。



i.e. N点的MDST 可以转化为一个重排序，然后M=N/2点的DCT-IV，最后是一个后排序。

MDCT的N/4点FFT实现：

MDCT的表达式：





利用三角函数关系


得到



将n划分为奇数和偶数部分，然后将奇数部分做一个变量代换写成如下的形式：



其中

，进一步化简



得到：



表达成指数形式：



把指数项拆开：


把n分为

和

两部分：



进一步的令


对

做

点FFT得到

，若令

，继续进行下去，可以得到

，进一步求得

用N/4点FFT实现IMDCT变换

计算N点IMDCT
需要做N×N/2次乘加运算，当N较大时运算量非常庞大，为此采用快速实现算法。利用IMDCT的对称性，将输入分为奇、偶序列，然后用N/4点FFT实现。实验结果表明，经过快速变换以后，IMDCT算法速度提高200~300倍，大大降低了运算时间。

为了表达式的清晰，下面的公式中忽略了块标号i。

首先把输入序列组合成N/4点的复数序列。


（1）

由于

，可以得到


（2）

将k划分为奇数和偶数部分，然后将奇数部分写成如下的形式：


…………………………………………………（3）

而


所以上面（3）式进一步的化简，得到


………………………………（4）

由

，所以使用旋转因子做进一步的代换，可以进一步的化简如下：


（5）

上面的推导中

。对（5）两端取共轭得到



……………………………………（6）

然后令

（7）

令

，同理

，

，对组合然后预旋转的复序列

做

点的IFFT，得到


（8）

由傅立叶变换性质得到


，于是得到


（9）

将模运算去掉，对m分情况求解，得到


（10）

所以若令

（11）

则有


（12）

则由

，可以得到下面的表达式：


（13）

这样就得到了全部的偶数部分的表达式

然后由x的对称关系：


（14）

得到的

奇数系数的求解表达式如下：


（15）

进而有


（15）


奇数系数的求解表达式如下



……………………（16）

重叠相加原理

无论是AAC 还是在其之前的音频编码方法(MPEG-1 Layer I 、II 、III 等)，编码器的对信号进行的各种处理绝大多数是在频域进行的，负责把信号从时域转换到频域的工具就是编码器中的滤波器组，AAC 之前的MPEG 音频编码器都是采用的子带编码的技术，就是按照心理声学的临界频带理论使用一组多相正交滤波器组把声音信号对应人耳的听觉特性分成不同的子带，在这些子带的基础上进行处理，这种做法在理论上无可非议，但是由于子带滤波器在分割频带的时候会引入噪声，该噪声会影响后级的编码组件的性能。AAC 编码标准摒弃了以前的多相正交滤波器组，
改用1024 点重叠的MDCT 滤波器组来对信号进行时-频转换，这种方法更接近人耳的工作原理， 从而进一步降低了在时-频转换时引入的噪声。

　在解码器端，量化后的频域数据被送入滤波器组，滤波器组输出时域的音频数据图2 是滤波器组的构成。

在滤波器组中，时域数据是经过MDCT 变换为频域数据的，关于MDCT 的特性在参考文献

[2 ]中有详细的描述，每一次进行MDCT 变换的数据块(其长度为2N 或N/4 ，N = 1024) 都和上一次进行MDCT 变换的数据块有50 %的重叠，所以滤波器组中有一个重叠缓冲区(overlap buffer) ，当前输入的数据帧(长度为N) 一方面要和重叠缓冲区里的数据拼接在一起进行加窗和MDCT 变换另一方面要存入重叠缓冲区，留待与下一帧数据拼接。AAC 中采用的窗函数有两种， 正弦窗(SIN E) 和KBD ( Kaiser 2 Bessel derived) 窗，分别用于语音成分丰富的情况和需要收敛信号频谱的情况，AAC
标准还规定了窗函数的长度可以是长窗(ONL Y LON G 2048 点) 或短窗(EIGHT SHORT 256 点)，分别用于信号的较为平稳的情况和信号高频分量较为丰富的情况，另外为了在长窗和短窗之间AAC 还采用了LONG START 窗和LONG STOP 窗作为过渡窗，对于当前输入的数据和重叠缓冲区里的数据拼接后采用什么形状以及什么长度的窗函数是在心理声学模型编码组件里确定的，这些信息作为控制信息传送给滤波器组。需要说明的是窗函数的作用是因为IMDCT这类块变换会造成因为量化噪声导致的在块边界的块效应，即块边界的不连续。因而窗函数的一个选择标准就是频率选择性和阻带的衰减。



图1. Mpeg-2 AAC中使用的长窗函数

滤波器组的原理

参考文献[4 ]指出，若将等长的2 个序列

顺序拼接后得到

进行MDCT 变换，然后再对变换结果进行IMDCT 变换， 那么输出的结果是


。

其中

表示序列x 的逆序排列，如果两次进行MDCT 变换的数据有50 %的重叠，那么就可以通过IMDCT 把重叠的部分恢复出来。设


把

重叠50 %相加，就得到序列

，可见重复的

被恢复出来了，然而在AAC
的滤波器组中，数据并不是简单的进行叠加。可以验证AAC 中的定义的正弦窗和KBD 窗满足以下两个特性(设窗长度是N) :



同时AAC 规定，重叠部分的数据必须用同一种形状的窗函数。其工作过程描述如下。在50 %重叠的情况下，编码器端重叠的数据长度是

，第i－1数据块中的后半部分是第i数据块的前半部分，设第i－1数据块的前半部分数据是

，i数据块的前半部分数据是

，则在编码时的MDCT
输入的第i－1数据块中对输入时域信号加窗后的结果是



在第二段中x[ k] 加窗后的结果是



利用前文所述的MDCT 的性质，在解码器滤波器组的IMDCT 输出端，第i－1块的输出数据中后半段的表达式为


第i块输出数据中前半段的表达式为



然后对得到的IMDCT变换的结果进行交叠加窗，首先计算这两段数据分别加窗后的结果

第i－1块：


第i块：


把两段结果叠加，根据窗函数的性质1、2 ， 相加的结果就是

。

Reference：

http://houh-1984.blog.163.com/

http://en.wikipedia.org/wiki/Modified_discrete_cosine_transform

http://www.zh.wikipedia.org/wiki/窗函数

本文主要是MDCT/MDST/IMDCT与DCT-IV的关系，然后是采用N/4点快速FFT傅里叶变换来实现N点的MDCT/IMDCT的公式推导以及重叠相加的基本原理，并根据MDCT/IMDCT的特点推导满足完全重建的窗函数的特性。