DCT(离散余弦变换)算法原理和源码(python)

时间:2022-04-29 12:37:30

  原理:

  离散余弦变换(DCT for Discrete Cosine Transform)是与傅里叶变换相关的一种变换,它类似于离散傅里叶变换(DFT for Discrete Fourier Transform),但是只使用实数。离散余弦变换相当于一个长度大概是它两倍的离散傅里叶变换,这个离散傅里叶变换是对一个实偶函数进行的(因为一个实偶函数的傅里叶变换仍然是一个实偶函数),在有些变形里面需要将输入或者输出的位置移动半个单位(DCT有8种标准类型,其中4种是常见的)。

  使用场景:

  离散余弦变换,尤其是它的第二种类型,经常被信号处理和图像处理使用,用于对信号和图像(包括静止图像和运动图像)进行有损数据压缩。一个类似的变换, 改进的离散余弦变换被用在高级音频编码(AAC for Advanced Audio Coding),Vorbis 和 MP3 音频压缩当中。

  python源码实现:  

# import numpy
import numpy as np # import dct
from scipy.fftpack import dct # numpy array
x = np.array([1.0, 2.0, 1.0, 2.0, -1.0])
print("x : ",x) # apply dct function on array
y = dct(x)
print("dct(x) : ",y)

  C源码实现:

#include <stdio.h>
#include <math.h>
#include <stdlib.h> void dct(float **DCTMatrix, float **Matrix, int N, int M);
void write_mat(FILE *fp, float **testRes, int N, int M);
void idct(float **Matrix, float **DCTMatrix, int N, int M);
float **calloc_mat(int dimX, int dimY);
void free_mat(float **p); float **calloc_mat(int dimX, int dimY){
float **m = calloc(dimX, sizeof(float*));
float *p = calloc(dimX*dimY, sizeof(float));
int i;
for(i=; i <dimX;i++){
m[i] = &p[i*dimY]; }
return m;
} void free_mat(float **m){
free(m[]);
free(m);
} void write_mat(FILE *fp, float **m, int N, int M){ int i, j;
for(i =; i< N; i++){
fprintf(fp, "%f", m[i][]);
for(j = ; j < M; j++){
fprintf(fp, "\t%f", m[i][j]);
}
fprintf(fp, "\n");
}
fprintf(fp, "\n");
} void dct(float **DCTMatrix, float **Matrix, int N, int M){ int i, j, u, v;
for (u = ; u < N; ++u) {
for (v = ; v < M; ++v) {
DCTMatrix[u][v] = ;
for (i = ; i < N; i++) {
for (j = ; j < M; j++) {
DCTMatrix[u][v] += Matrix[i][j] * cos(M_PI/((float)N)*(i+./.)*u)*cos(M_PI/((float)M)*(j+./.)*v);
}
}
}
}
} void idct(float **Matrix, float **DCTMatrix, int N, int M){
int i, j, u, v; for (u = ; u < N; ++u) {
for (v = ; v < M; ++v) {
Matrix[u][v] = /.*DCTMatrix[][];
for(i = ; i < N; i++){
Matrix[u][v] += /.*DCTMatrix[i][];
}
for(j = ; j < M; j++){
Matrix[u][v] += /.*DCTMatrix[][j];
} for (i = ; i < N; i++) {
for (j = ; j < M; j++) {
Matrix[u][v] += DCTMatrix[i][j] * cos(M_PI/((float)N)*(u+./.)*i)*cos(M_PI/((float)M)*(v+./.)*j);
}
}
Matrix[u][v] *= ./((float)N)*./((float)M);
}
}
} int main() { float
testBlockA[][] = { {, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , } }, testBlockB[][] = {{, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , },
{, , , , , , , } }; FILE * fp = fopen("mydata.csv", "w");
int dimX = , dimY = ;
int i, j; float **testBlock = calloc_mat(dimX, dimY);
float **testDCT = calloc_mat(dimX, dimY);
float **testiDCT = calloc_mat(dimX, dimY); for(i = ; i<dimX; i++){
for(j = ; j<dimY; j++){
testBlock[i][j] = testBlockB[i][j];
}
} dct(testDCT, testBlock, dimX, dimY);
write_mat(fp, testDCT, dimX, dimY); idct(testiDCT, testDCT, dimX, dimY);
write_mat(fp, testiDCT, dimX, dimY); fclose(fp);
free_mat(testBlock);
free_mat(testDCT);
free_mat(testiDCT); return ;
}

  代码路径: https://github.com/DyLanCao/DCT.git

参考文档:

https://*.com/questions/8310749/discrete-cosine-transform-dct-implementation-c