`
oraclestudy
  • 浏览: 497608 次
文章分类
社区版块
存档分类

快速离散余弦变换代码实现(FDCT)

 
阅读更多

网上常见的快速离散余弦变换的代码如下:

#include <iostream>
using namespace std;

#define DCTSIZE 8


void FDCT(double* lpBuff)
{
double tmp0, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7;
double tmp10, tmp11, tmp12, tmp13;
double z1, z2, z3, z4, z5, z11, z13;
double* dataptr;
int ctr;

/* 第一部分,对行进行计算 */
dataptr = lpBuff;
for (ctr = DCTSIZE-1; ctr >= 0; ctr--)
{
tmp0 = dataptr[0] + dataptr[7];
tmp7 = dataptr[0] - dataptr[7];
tmp1 = dataptr[1] + dataptr[6];
tmp6 = dataptr[1] - dataptr[6];
tmp2 = dataptr[2] + dataptr[5];
tmp5 = dataptr[2] - dataptr[5];
tmp3 = dataptr[3] + dataptr[4];
tmp4 = dataptr[3] - dataptr[4];

/* 对偶数项进行运算 */
tmp10 = tmp0 + tmp3; /* phase 2 */
tmp13 = tmp0 - tmp3;
tmp11 = tmp1 + tmp2;
tmp12 = tmp1 - tmp2;

dataptr[0] = tmp10 + tmp11; /* phase 3 */
dataptr[4] = tmp10 - tmp11;

z1 = (tmp12 + tmp13) * (0.707106781); /* c4 */
dataptr[2] = tmp13 + z1; /* phase 5 */
dataptr[6] = tmp13 - z1;

/* 对奇数项进行计算 */
tmp10 = tmp4 + tmp5; /* phase 2 */
tmp11 = tmp5 + tmp6;
tmp12 = tmp6 + tmp7;

z5 = (tmp10 - tmp12) * ( 0.382683433); /* c6 */
z2 = (0.541196100) * tmp10 + z5; /* c2-c6 */
z4 = (1.306562965) * tmp12 + z5; /* c2+c6 */
z3 = tmp11 * (0.707106781); /* c4 */

z11 = tmp7 + z3; /* phase 5 */
z13 = tmp7 - z3;

dataptr[5] = z13 + z2; /* phase 6 */
dataptr[3] = z13 - z2;
dataptr[1] = z11 + z4;
dataptr[7] = z11 - z4;

dataptr += DCTSIZE; /* 将指针指向下一行 */
}

/* 第二部分,对列进行计算 */
dataptr = lpBuff;
for (ctr = DCTSIZE-1; ctr >= 0; ctr--)
{
tmp0 = dataptr[DCTSIZE*0] + dataptr[DCTSIZE*7];
tmp7 = dataptr[DCTSIZE*0] - dataptr[DCTSIZE*7];
tmp1 = dataptr[DCTSIZE*1] + dataptr[DCTSIZE*6];
tmp6 = dataptr[DCTSIZE*1] - dataptr[DCTSIZE*6];
tmp2 = dataptr[DCTSIZE*2] + dataptr[DCTSIZE*5];
tmp5 = dataptr[DCTSIZE*2] - dataptr[DCTSIZE*5];
tmp3 = dataptr[DCTSIZE*3] + dataptr[DCTSIZE*4];
tmp4 = dataptr[DCTSIZE*3] - dataptr[DCTSIZE*4];

/* 对偶数项进行运算 */
tmp10 = tmp0 + tmp3; /* phase 2 */
tmp13 = tmp0 - tmp3;
tmp11 = tmp1 + tmp2;
tmp12 = tmp1 - tmp2;

dataptr[DCTSIZE*0] = tmp10 + tmp11; /* phase 3 */
dataptr[DCTSIZE*4] = tmp10 - tmp11;

z1 = (tmp12 + tmp13) * (0.707106781); /* c4 */
dataptr[DCTSIZE*2] = tmp13 + z1; /* phase 5 */
dataptr[DCTSIZE*6] = tmp13 - z1;

/* 对奇数项进行计算 */
tmp10 = tmp4 + tmp5; /* phase 2 */
tmp11 = tmp5 + tmp6;
tmp12 = tmp6 + tmp7;

z5 = (tmp10 - tmp12) * (0.382683433); /* c6 */
z2 = (0.541196100) * tmp10 + z5; /* c2-c6 */
z4 = (1.306562965) * tmp12 + z5; /* c2+c6 */
z3 = tmp11 * (0.707106781); /* c4 */

z11 = tmp7 + z3; /* phase 5 */
z13 = tmp7 - z3;

dataptr[DCTSIZE*5] = z13 + z2; /* phase 6 */
dataptr[DCTSIZE*3] = z13 - z2;
dataptr[DCTSIZE*1] = z11 + z4;
dataptr[DCTSIZE*7] = z11 - z4;

++dataptr; /* 将指针指向下一列 */
}
}

int main(void)
{
const short NUM = 8;
short i = 0;
short j = 0;
short u = 0;
short v = 0;

double input[NUM][NUM] =
{
{89.00, 101.00, 114.00, 125.00, 126.00, 115.00, 105.00, 96.00},
{97.00, 115.00, 131.00, 147.00, 149.00, 135.00, 123.00, 113.00},
{114.00, 134.00, 159.00, 178.00, 175.00, 164.00, 149.00, 137.00},
{121.00, 143.00, 177.00, 196.00, 201.00, 189.00, 165.00, 150.00},
{119.00, 141.00, 175.00, 201.00, 207.00, 186.00, 162.00, 144.00},
{107.00, 130.00, 165.00, 189.00, 192.00, 171.00, 144.00, 125.00},
{97.00, 119.00, 149.00, 171.00, 172.00, 145.00, 117.00, 96.00},
{88.00, 107.00, 136.00, 156.00, 155.00, 129.00, 97.00, 75.00}
};

double *data = new double[NUM * NUM];
for(i = 0; i < NUM; i++)
{
for(j = 0; j < NUM; j++)
{
data[i * NUM + j] = input[i][j];
}
}
FDCT(data);
// print the result of FDCT
for(u = 0; u < NUM; u++)
{
for(v = 0; v < NUM; v++)
{
printf("%7.2f ", data[u*NUM + v]);
}
cout << endl;
}
return 0;
}

运算得到的结果:
9000.00 -353.39 -1933.67 -61.84 16.00 -4.74 -10.33 3.96
-241.34 -243.40 650.60 -42.79 -24.64 -1.12 -9.34 -5.43
-1725.67 464.18 230.84 25.81 13.49 -6.50 -15.02 -0.25
-66.65 -56.94 2.18 17.28 15.93 -5.71 -5.37 4.21
-18.00 2.53 -5.07 28.89 -2.00 -1.57 9.07 2.16
20.23 8.66 11.55 -5.36 -11.20 2.24 5.68 0.07
-0.33 -1.85 9.02 -6.61 -3.49 6.67 3.16 -1.46
-0.24 8.50 3.08 -2.01 3.92 1.92 -2.39 -0.12

这个结果和离散余弦正逆变换中最原始的或者说是真正的离散余弦变换的结果是不同的。原因请看下面程序中的红色字体部分。

=================================================================

// FDCT.h

void fdct(double *block);

// FDCT.cpp

#define USE_ACCURATE_ROUNDING

#define RIGHT_SHIFT(x, shft) ((x) >> (shft))

#ifdef USE_ACCURATE_ROUNDING
#define ONE ((int) 1)
#define DESCALE(x, n) RIGHT_SHIFT((x) + (ONE << ((n) - 1)), n)
#else
#define DESCALE(x, n) RIGHT_SHIFT(x, n)
#endif

#define CONST_BITS 13
#define PASS1_BITS 2

#define FIX_0_298631336 ((int) 2446)/* FIX(0.298631336) */
#define FIX_0_390180644 ((int) 3196)/* FIX(0.390180644) */
#define FIX_0_541196100 ((int) 4433)/* FIX(0.541196100) */
#define FIX_0_765366865 ((int) 6270)/* FIX(0.765366865) */
#define FIX_0_899976223 ((int) 7373)/* FIX(0.899976223) */
#define FIX_1_175875602 ((int) 9633)/* FIX(1.175875602) */
#define FIX_1_501321110 ((int) 12299)/* FIX(1.501321110) */
#define FIX_1_847759065 ((int) 15137)/* FIX(1.847759065) */
#define FIX_1_961570560 ((int) 16069)/* FIX(1.961570560) */
#define FIX_2_053119869 ((int) 16819)/* FIX(2.053119869) */
#define FIX_2_562915447 ((int) 20995)/* FIX(2.562915447) */
#define FIX_3_072711026 ((int) 25172)/* FIX(3.072711026) */

void fdct(double *block)
{
int tmp0, tmp1, tmp2, tmp3, tmp4, tmp5, tmp6, tmp7;
int tmp10, tmp11, tmp12, tmp13;
int z1, z2, z3, z4, z5;
double *blkptr;
int *dataptr;
int data[64];
int i;

/* Pass 1: process rows. */
/* Note results are scaled up by sqrt(8) compared to a true DCT; */
/* furthermore, we scale the results by 2**PASS1_BITS. */

dataptr = data;
blkptr = block;
for (i = 0; i < 8; i++) {
tmp0 = (int)(blkptr[0] + blkptr[7]);
tmp7 = (int)(blkptr[0] - blkptr[7]);
tmp1 = (int)(blkptr[1] + blkptr[6]);
tmp6 = (int)(blkptr[1] - blkptr[6]);
tmp2 = (int)(blkptr[2] + blkptr[5]);
tmp5 = (int)(blkptr[2] - blkptr[5]);
tmp3 = (int)(blkptr[3] + blkptr[4]);
tmp4 = (int)(blkptr[3] - blkptr[4]);

tmp10 = tmp0 + tmp3;
tmp13 = tmp0 - tmp3;
tmp11 = tmp1 + tmp2;
tmp12 = tmp1 - tmp2;

dataptr[0] = (tmp10 + tmp11) << PASS1_BITS;
dataptr[4] = (tmp10 - tmp11) << PASS1_BITS;

z1 = (tmp12 + tmp13) * FIX_0_541196100;
dataptr[2] = DESCALE(z1 + tmp13 * FIX_0_765366865, CONST_BITS - PASS1_BITS);
dataptr[6] =DESCALE(z1 + tmp12 * (-FIX_1_847759065), CONST_BITS - PASS1_BITS);

z1 = tmp4 + tmp7;
z2 = tmp5 + tmp6;
z3 = tmp4 + tmp6;
z4 = tmp5 + tmp7;
z5 = (z3 + z4) * FIX_1_175875602;/* sqrt(2) * c3 */

tmp4 *= FIX_0_298631336;/* sqrt(2) * (-c1+c3+c5-c7) */
tmp5 *= FIX_2_053119869;/* sqrt(2) * ( c1+c3-c5+c7) */
tmp6 *= FIX_3_072711026;/* sqrt(2) * ( c1+c3+c5-c7) */
tmp7 *= FIX_1_501321110;/* sqrt(2) * ( c1+c3-c5-c7) */
z1 *= -FIX_0_899976223;/* sqrt(2) * (c7-c3) */
z2 *= -FIX_2_562915447;/* sqrt(2) * (-c1-c3) */
z3 *= -FIX_1_961570560;/* sqrt(2) * (-c3-c5) */
z4 *= -FIX_0_390180644;/* sqrt(2) * (c5-c3) */

z3 += z5;
z4 += z5;

dataptr[7] = DESCALE(tmp4 + z1 + z3, CONST_BITS - PASS1_BITS);
dataptr[5] = DESCALE(tmp5 + z2 + z4, CONST_BITS - PASS1_BITS);
dataptr[3] = DESCALE(tmp6 + z2 + z3, CONST_BITS - PASS1_BITS);
dataptr[1] = DESCALE(tmp7 + z1 + z4, CONST_BITS - PASS1_BITS);

dataptr += 8;/* advance pointer to next row */
blkptr += 8;
}

/* Pass 2: process columns.
* We remove the PASS1_BITS scaling, but leave the results scaled up
* by an overall factor of 8.
*/

dataptr = data;
for (i = 0; i < 8; i++) {
tmp0 = dataptr[0] + dataptr[56];
tmp7 = dataptr[0] - dataptr[56];
tmp1 = dataptr[8] + dataptr[48];
tmp6 = dataptr[8] - dataptr[48];
tmp2 = dataptr[16] + dataptr[40];
tmp5 = dataptr[16] - dataptr[40];
tmp3 = dataptr[24] + dataptr[32];
tmp4 = dataptr[24] - dataptr[32];

tmp10 = tmp0 + tmp3;
tmp13 = tmp0 - tmp3;
tmp11 = tmp1 + tmp2;
tmp12 = tmp1 - tmp2;

dataptr[0] = DESCALE(tmp10 + tmp11, PASS1_BITS);
dataptr[32] = DESCALE(tmp10 - tmp11, PASS1_BITS);

z1 = (tmp12 + tmp13) * FIX_0_541196100;
dataptr[16] = DESCALE(z1 + tmp13 * FIX_0_765366865, CONST_BITS + PASS1_BITS);
dataptr[48] = DESCALE(z1 + tmp12 * (-FIX_1_847759065), CONST_BITS + PASS1_BITS);

z1 = tmp4 + tmp7;
z2 = tmp5 + tmp6;
z3 = tmp4 + tmp6;
z4 = tmp5 + tmp7;
z5 = (z3 + z4) * FIX_1_175875602;/* sqrt(2) * c3 */

tmp4 *= FIX_0_298631336;/* sqrt(2) * (-c1+c3+c5-c7) */
tmp5 *= FIX_2_053119869;/* sqrt(2) * ( c1+c3-c5+c7) */
tmp6 *= FIX_3_072711026;/* sqrt(2) * ( c1+c3+c5-c7) */
tmp7 *= FIX_1_501321110;/* sqrt(2) * ( c1+c3-c5-c7) */
z1 *= -FIX_0_899976223;/* sqrt(2) * (c7-c3) */
z2 *= -FIX_2_562915447;/* sqrt(2) * (-c1-c3) */
z3 *= -FIX_1_961570560;/* sqrt(2) * (-c3-c5) */
z4 *= -FIX_0_390180644;/* sqrt(2) * (c5-c3) */

z3 += z5;
z4 += z5;

dataptr[56] = DESCALE(tmp4 + z1 + z3, CONST_BITS + PASS1_BITS);
dataptr[40] = DESCALE(tmp5 + z2 + z4, CONST_BITS + PASS1_BITS);
dataptr[24] = DESCALE(tmp6 + z2 + z3, CONST_BITS + PASS1_BITS);
dataptr[8] = DESCALE(tmp7 + z1 + z4, CONST_BITS + PASS1_BITS);

dataptr++;/* advance pointer to next column */
}
/* descale */
for (i = 0; i < 64; i++)
block[i] = (double)DESCALE(data[i], 3);
}

// 测试程序:calculate.cpp

#include "FDCT.h"
#include <iostream>
using namespace std;

int main(void)
{
const short NUM = 8;
short i = 0;
short j = 0;
short u = 0;
short v = 0;

double input[NUM][NUM] =
{
{89.00, 101.00, 114.00, 125.00, 126.00, 115.00, 105.00, 96.00},
{97.00, 115.00, 131.00, 147.00, 149.00, 135.00, 123.00, 113.00},
{114.00, 134.00, 159.00, 178.00, 175.00, 164.00, 149.00, 137.00},
{121.00, 143.00, 177.00, 196.00, 201.00, 189.00, 165.00, 150.00},
{119.00, 141.00, 175.00, 201.00, 207.00, 186.00, 162.00, 144.00},
{107.00, 130.00, 165.00, 189.00, 192.00, 171.00, 144.00, 125.00},
{97.00, 119.00, 149.00, 171.00, 172.00, 145.00, 117.00, 96.00},
{88.00, 107.00, 136.00, 156.00, 155.00, 129.00, 97.00, 75.00}
};

double *data = new double[NUM * NUM];
for(i = 0; i < NUM; i++)
{
for(j = 0; j < NUM; j++)
{
data[i * NUM + j] = input[i][j];
}
}
fdct(data);
// print the result of FDCT
for(u = 0; u < NUM; u++)
{
for(v = 0; v < NUM; v++)
{
printf("%7.2f ", data[u*NUM + v]);
//cout << input[u][v] << '/t';
}
cout << endl;
}
return 0;
}

运算结果:

1125.00 -32.00 -185.00 -7.00 2.00 -1.00 -2.00 2.00
-22.00 -16.00 45.00 -3.00 -2.00 0.00 -1.00 -2.00
-165.00 32.00 17.00 2.00 1.00 -1.00 -3.00 0.00
-7.00 -4.00 0.00 2.00 2.00 -1.00 -1.00 2.00
-2.00 0.00 0.00 3.00 0.00 0.00 2.00 1.00
3.00 1.00 1.00 -1.00 -2.00 1.00 2.00 0.00
0.00 0.00 2.00 -1.00 -1.00 2.00 1.00 -1.00
0.00 3.00 1.00 -1.00 2.00 1.00 -2.00 0.00

上述运算结果和离散余弦正逆变换一文中的结果保持完全一致。

尽管两种快速算法的结果有很大的不同,但都反应了信号的频谱特性,因此在实际应用中都可以认为是对的。

分享到:
评论

相关推荐

    FDCT_快速离散余弦变换_

    快速离散余弦变换(Fast Discrete Cosine Transform, FDCT)是数字信号处理和图像压缩领域中的一个重要算法。它是离散余弦变换(DCT)的一种优化实现,旨在提高计算效率,降低计算复杂度。在诸如JPEG图像压缩标准中...

    离散余弦变换代码

    在文件"fdct"中,可能包含了一个离散余弦变换的实现,这可能是用某种编程语言编写的,例如C、C++或Python。这个实现可能涵盖了DCT的计算过程,以及可能的优化技巧,比如使用快速算法来提高效率。此外,还可能包括...

    DCT.rar_余弦变换_离散余弦变换

    本压缩包包含两个C++源代码文件:FDCT.cpp和DCT.cpp,分别可能对应前向离散余弦变换(Forward DCT)和逆向离散余弦变换(Inverse DCT)的实现。 离散余弦变换的基本概念是将一个信号或图像的像素值表示为一组频率...

    DCT(离散余弦变换)C#源代码

    ### DCT(离散余弦变换)C#源代码解析及知识点详解 #### 一、基础知识简介 **离散余弦变换(Discrete Cosine Transform, DCT)**是一种广泛应用于数字信号处理领域的数学工具,尤其在图像压缩、音频编码等多媒体技术...

    FDCT_Vc_

    2. **快速离散余弦变换(FDCT)**: 是DCT的一种优化算法,通过减少计算量和提高效率来实现快速计算。 3. **逆离散余弦变换(IDCT)**: 与DCT相反,它将频率域的系数转换回空间域,恢复原始图像。 4. **VC++编程**:...

    图像正交变换

    在实际应用中,为了提高计算效率,我们通常采用快速傅里叶变换(FFT)算法来执行DFT,以及快速离散余弦变换(FDCT)来执行DCT。这些快速算法大大减少了计算量,使得大规模图像处理成为可能。 通过学习和理解傅里叶...

    基于VC++实现的DCT变换

    FDCT是快速离散余弦变换的缩写,是快速算法实现DCT的方式。在图像编码标准如JPEG中,FDCT被广泛采用。在VC++环境下,FDCT的实现通常会涉及以下步骤: - 分块处理:将大图像划分为小块,分别进行DCT。 - 二维DCT:对...

    dct_sse2.rar_short

    标题 "dct_sse2.rar_short" 暗示了这是一个与数字信号处理相关的压缩文件,特别是针对视频编码中的快速离散余弦变换(FDCT)优化的实现,使用了SSE2指令集。SSE2是Intel在 Pentium 4 和 Xeon 处理器中引入的一种 ...

    一个C编写的快速DCT程序

    【快速离散余弦变换(DCT)与C语言实现】 快速离散余弦变换(Fast Discrete Cosine Transform,简称FDCT)是数字信号处理领域中的一个重要算法,广泛应用于图像压缩,如JPEG等标准中。它将图像数据转换到频率域,...

    DCT算法研究与总结

    正向离散余弦变换(FDCT)是数字信号处理中的重要算法,主要用于图像和音频的压缩编码。FDCT-1处理一维序列,而FDCT-2则应用于二维矩阵,通常用于图像数据。在二维DCT中,图像被划分为8x8的子块进行处理,以减少计算...

    VC++DCT变换源码

    VC++实现的DCT(离散余弦变换)和IDCT(逆离散余弦变换)源码是一个非常实用的学习资源,对于深入理解图像处理和数字信号处理领域的核心算法具有重要意义。DCT是一种在频域分析中广泛使用的线性变换方法,尤其在图像...

    VC++实现图像的DCT,IDCT变换

    在图像处理领域,离散余弦变换(Discrete Cosine Transform, DCT)和反离散余弦变换(Inverted Discrete Cosine Transform, IDCT)是两种非常重要的技术,常用于图像压缩,如JPEG格式的图像编码。本文将深入探讨如何...

    JPEG压缩编码系统源代码.rar

    - **DCT函数**:实现离散余弦变换,通常包括正向DCT(FDCT)和反向DCT(IDCT)。 - **量化函数**:对DCT系数应用量化表。 - **熵编码模块**:包括哈夫曼编码器或算术编码器,用于将量化后的系数编码为位流。 - **...

    VC++开发的快速DCT变换算法

    快速离散余弦变换(Fast Discrete Cosine Transform,简称FDCT)是图像处理和数字信号处理中的一个重要算法,主要用于图像和音频数据的压缩,如JPEG图像压缩标准和MP3音频压缩格式。在VC++环境下,开发快速DCT算法...

    DCT变换的快速算法

    离散余弦变换(Discrete Cosine Transform, DCT)是一种广泛应用在图像处理、音频编码、信号处理等领域的数学变换技术。它将数据从时域或空间域转换到频域,便于对信号进行分析和压缩。传统的DCT计算方法虽然理论...

    C++ 实现DCT算法

    在IT领域,离散余弦变换(Discrete Cosine Transform, DCT)是一种广泛应用的信号处理技术,尤其在图像和音频压缩中,如JPEG和MP3等格式就依赖于DCT来实现数据的高效编码。C++作为一种强大的编程语言,可以用来实现...

    Fbst.rar_dct

    标题"Fbst.rar_dct"可能指的是一个包含快速离散余弦变换(Fast Discrete Cosine Transform,简称FDCT)算法的压缩文件。描述中的“DCT快速算法Fast DCT algorithm good”进一步证实了这一点,表明这个压缩包内含有...

    dct_mmx.rar_short

    标题 "dct_mmx.rar_short" 提供的信息表明,这是一个与数字信号处理相关的压缩文件,特别是针对VP8视频编码标准的快速离散余弦变换(FDCT)算法。VP8是一种广泛使用的网络视频编解码器,常用于WebM开源视频格式。MMX...

    模式识别DCT快速变换源码

    而离散余弦变换(Discrete Cosine Transform,简称DCT)作为模式识别中的关键算法,被广泛应用于图像压缩和信号处理。本文将详细探讨DCT快速变换的原理,并结合提供的源码,解析其在VC++环境下的实现细节。 DCT是一...

    2-D DCT快速算法

    二维离散余弦变换(2-D Discrete Cosine Transform, DCT)是图像处理和信号处理领域中的一个重要工具,尤其在图像压缩技术如JPEG中扮演着核心角色。DCT能够将图像数据转换到频域,使得大部分能量集中在低频部分,...

Global site tag (gtag.js) - Google Analytics