离散余弦变换(DCT)是一种广泛使用的变换方法,尤其在数据压缩和处理领域。DCT因其系数的实数特性而比傅里叶变换更受欢迎。尽管DCT的应用非常广泛,但目前只有一种广为人知的快速算法,即二维DCT算法。对于一维DCT,虽然存在Byeong Lee的DCT算法,但其递归性质使得在某些设备上(例如GPGPU)的使用变得不那么实用。为此,本文提出了一种迭代算法,它不依赖于递归调用,并且能够处理几乎无限大小的输入数据(仅需为2的幂)。
本文选择Byeong Lee的DCT算法作为基础,该算法广泛用于硬件实现,并且基于蝴蝶分布。算法的实现主要分为两个路径:输入蝴蝶和输出蝴蝶。迭代过程中,一个数据块用于保持数据,另一个数据块用于存储计算得到的余弦系数。每次迭代将前缓冲区分为偶数和奇数两部分,然后将它们复制到后缓冲区,并交换前后缓冲区。
以下是C++实现的非递归DCT算法的代码片段:
double *Front, *Back, *temp;
// 指向迭代过程中的前、后和临时数据块的指针。
double *Array = new double[N << 1];
// 存储迭代过程中数据的数组,大小为N的两倍。
double *Curent_Cos = new double[N >> 1];
// 用于预先计算余弦系数的数组。
Front = Array;
// 将Front指向数组的开始。
step = 0;
Back = Front + N;
// 将Back指向数组的中间。
for (int i = 0; i < N; i++) {
Front[i] = block[i];
// 加载输入数据
}
while (j > 1) { // 输入蝴蝶循环
half_N = half_N >> 1;
current_PI_half_By_N = PI_half / (double)prev_half_N;
current_PI_By_N = 0;
step_Phase = current_PI_half_By_N * 2.0;
for (int i = 0; i < half_N; i++) {
// 预先计算余弦系数
Curent_Cos[i] = 0.5 / cos(current_PI_By_N + current_PI_half_By_N);
current_PI_By_N += step_Phase;
}
shif_Block = 0;
for (int x = 0; x < step; x++) {
for (int i = 0; i < half_N; i++) {
shift = shif_Block + prev_half_N - i - 1;
Back[shif_Block + i] = Front[shif_Block + i] + Front[shift];
Back[shif_Block + half_N + i] = (Front[shif_Block + i] - Front[shift]) * Curent_Cos[i];
}
shif_Block += prev_half_N;
}
temp = Front;
Front = Back;
Back = temp;
j = j >> 1;
step = step << 1;
prev_half_N = half_N;
}
while (j < N) { // 输出蝴蝶循环
shif_Block = 0;
for (int x = 0; x < step; x++) {
for (int i = 0; i < half_N - 1; i++) {
Back[shif_Block + (i << 1)] = Front[shif_Block + i];
Back[shif_Block + (i << 1) + 1] = Front[shif_Block + half_N + i] + Front[shif_Block + half_N + i + 1];
}
Back[shif_Block + ((half_N - 1) << 1)] = Front[shif_Block + half_N - 1];
Back[shif_Block + (half_N << 1) - 1] = Front[shif_Block + (half_N << 1) - 1];
shif_Block += prev_half_N;
}
temp = Front;
Front = Back;
Back = temp;
j = j << 1;
step = step >> 1;
half_N = prev_half_N;
prev_half_N = prev_half_N >> 1;
}
for (int i = 0; i < N; i++) {
block[i] = Front[i];
// 卸载DCT系数
}
block[0] = block[0] / dN;
// 计算直流分量。
一个非常重要的事实是,在第一次迭代之后,只有一组奇数系数,但在第二次迭代之后,有两组奇数系数。对于这两组奇数系数,需要相同的余弦系数。预先计算的存储允许只为一组系数计算一次,另一组则使用存储的数据。这对于计算1024点及以上的数据非常有用。