一維2N浮點FDCT變換
1、 fdct.h檔案
#ifndef __JENC__ //判斷宏定義是否被定義
#define __JENC__
#define PI 3.1415926
#include <math.h>
#include <string>
double *C=NULL;
class Fdct
{
public:
void Invoke()
{
double dataf[16] = {137,136, 133, 136, 138, 134, 134, 132, 132, 138, 129, 131, 132, 127, 129, 128};
fdct_1D(dataf, GetTwoIndex(16));
for(int i = 0; i < 16; i++){
printf("%f\n",dataf[i]);
}
}
private:
//*******************************************************************
//方法名稱:GetTwoIndex
//說明:算出deg的值
//日期:2012.2.29
//參數說明:
//num: 資料的數量
//*******************************************************************
int GetTwoIndex(int num)
{
int temp = 0;
if(num<=0) return -1;
while(1<<temp<num)
temp++;
return temp;
}
//*******************************************************************
//方法名稱:swap
//說明:交換temp1與temp2的值
//日期:2012.2.29
//參數說明:
//a:所要交換的資料
//b:所要交換的資料
//*******************************************************************
void swap(double &a,double &b)
{
double temp;
temp=a;
a=b;
b=temp;
}
//*******************************************************************
//方法名稱:bitrev
//說明:位反序函數
//日期:2012.2.29
//參數說明:
//bi:不清楚
//deg:為DCT變換資料的長度
//*******************************************************************
int bitrev(int bi,int deg)
{
int j=1,temp=0,degnum,halfnum;
degnum=deg;
//if(deg<0) return 0;
if(deg==0) return bi;
halfnum=1<<(deg-1);
while(halfnum)
{
if(halfnum&bi)
temp+=j;
j<<=1;
halfnum>>=1;
}
return temp;
}
//*******************************************************************
//方法名稱:fbitrev
//說明: 實作位反序排列(注意和上面函數的差別)
//日期:2012.2.29
//參數說明:
//f:不清楚
//deg:為DCT變換資料的長度
//*******************************************************************
void fbitrev(double *f,int deg)
{
if(deg==1) return;
int len=(1<<deg)-1;
int i=1;
int ii;
double temp;
while(i<len)
{ ii=bitrev(i,deg);
if(ii>i)
{ temp=f[ii];
f[ii]=f[i];
f[i]=temp;
}
i++;
}
}
//*******************************************************************
//方法名稱:initDCTParam
//說明:考慮到DCT變換中的系數要重複計算,可使用查找表來提高運作的效率,隻要開始DCT變換之前計算一次,DCT變換中就可以隻查找而無需計算系數。
//日期:2012.2.29
//參數說明:
//deg:為DCT變換資料的長度的幂
//*******************************************************************
void initDCTParam(int deg)
{
int total,halftotal,i,group,endstart,factor;
total=1<<deg; //total==16
if(C!=NULL) delete []C;
C=(double *)new double[total];
halftotal=total>>1; //右移一位相當于給該數除以2 halftotal == 8
for(i=0;i<halftotal;i++)
C[total-i-1]=(double)(2*i+1);
for(group=0;group<deg-1;group++)
{
endstart=1<<(deg-1-group); //endstart == 8
int len=endstart>>1; //len==4
factor=1<<(group+1);
for(int j=0;j<len;j++)
C[endstart-j-1]=factor*C[total-j-1];
}
for(i=1;i<total;i++)
C[i]=2.0*cos(C[i]*PI/(total<<1)); //C[0]空着,沒使用
}
//*******************************************************************
//方法名稱:dct_forward
//說明:DCT變換過程可分為兩部分:前向追底和後向回根,這個是實作前向追底。
//日期:2012.2.29
//參數說明:
//f:存儲DCT資料
//deg:為DCT變換資料的長度
//*******************************************************************
void dct_forward(double *f,int deg)
{
int i_deg,i_halfwing,total,wing,wings,winglen,halfwing;
double temp1,temp2;
total=1<<deg;
for(i_deg=0;i_deg<deg;i_deg++)
{
wings=1<<i_deg;
winglen=total>>i_deg;
halfwing=winglen>>1;
for(wing=0;wing<wings;wing++)
{
for(i_halfwing=0;i_halfwing<halfwing;i_halfwing++)
{
temp1=f[wing*winglen+i_halfwing];
temp2=f[(wing+1)*winglen-1-i_halfwing];
if(wing%2)
swap(temp1,temp2); // 交換temp1與temp2的值
f[wing*winglen+i_halfwing]=temp1+temp2;
f[(wing+1)*winglen-1-i_halfwing]=(temp1-temp2)*C[winglen-1-i_halfwing];
}
}
}
}
//*******************************************************************
//方法名稱:dct_backward
//說明:DCT變換過程可分為兩部分:前向追底和後向回根,這個是實作後向回根。
//日期:2012.2.29
//參數說明:
//f:存儲DCT資料
//deg:為DCT變換資料的長度
//*******************************************************************
void dct_backward(double *f,int deg)
{
int total,i_deg,wing,wings,halfwing,winglen,i_halfwing,temp1,temp2;
total=1<<deg;
for(i_deg=deg-1;i_deg>=0;i_deg--)
{
wings=1<<i_deg;
winglen=1<<(deg-i_deg);
halfwing=winglen>>1;
for(wing=0;wing<wings;wing++)
{
for(i_halfwing=0;i_halfwing<halfwing;i_halfwing++)
{ //f[i_halfwing+wing*winglen]=f[i_halfwing+wing*winglen];
if(i_halfwing==0)
f[halfwing+wing*winglen+i_halfwing]=0.5*f[halfwing+wing*winglen+i_halfwing];
else
{
temp1=bitrev(i_halfwing,deg-i_deg-1); // bitrev為位反序
temp2=bitrev(i_halfwing-1,deg-i_deg-1); // 第一參數為要變換的數
// 第二參數為二進制長度
f[halfwing+wing*winglen+temp1]=f[halfwing+wing*winglen+temp1]-f[halfwing+wing*winglen+temp2];
}
}
}
}
}
//*******************************************************************
//方法名稱:fdct_1D_no_param
//說明:完整實作一維DCT變換
//日期:2012.2.29
//參數說明:
//f:存儲DCT資料
//deg:為DCT變換資料的長度
//*******************************************************************
void fdct_1D_no_param(double *f,int deg) //deg == 4
{
initDCTParam(deg);
dct_forward(f,deg);
dct_backward(f,deg);
fbitrev(f,deg); // 實作位反序排列
f[0]=1/(sqrt(2.0))*f[0];
}
void fdct_1D(double *f,int deg)
{
fdct_1D_no_param(f,deg);
int total=1<<deg;
double param=sqrt(2.0/total);
for(int i=0;i<total;i++){
f[i]=param*f[i];
}
}
};
#endif // __JENC__
2、FDCT.c檔案
// NewFDCT.cpp : 定義控制台應用程式的入口點。
#include "stdafx.h"
#include "fdct.h" //引用頭檔案
//實作FDCT一維快速變換
int _tmain(int argc, _TCHAR* argv[])
{
Fdct fdct;
fdct.Invoke();
return 0;
}
一維2N浮點IFDCT變換
1、IFDCT.h
#ifndef __JENC__ //判斷宏定義是否被定義
#define __JENC__
#define PI 3.1415926
#include <math.h>
#include <string>
double *C=NULL;
class Fidct
{
public:
void Invoke()
{
double dataf[16] = {531.500000, 10.491843, -2.369268, 1.307582, -0.718939, -0.948647, 2.013496, 3.625035, 1.500000, -1.973283, 2.122850, -0.347127, -5.655363, 2.111544, 0.908794, -1.302960};
fidct_1D(dataf, 4);
for(int i = 0; i < 16; i++){
printf("%d\n",int(dataf[i]));
}
}
private:
//*******************************************************************
//方法名稱:initIDCTParam
//說明:考慮到IDCT變換中的系數要重複計算,可使用查找表來提高運作的效率,隻要開始IDCT變換之前計算一次,IDCT變換中就可以隻查找而無需計算系數。
//日期:2012.2.29
//參數說明:
//deg:為IDCT變換資料長度的幂
//*******************************************************************
void initIDCTParam(int deg)
{
int total,halftotal,i,group,endstart,factor;
total=1<<deg;
if(C!=NULL) delete []C;
C=(double *)new double[total];
halftotal=total>>1;
for(i=0;i<halftotal;i++)
C[total-i-1]=(double)(2*i+1);
for(group=0;group<deg-1;group++)
{
endstart=1<<(deg-1-group);
int len=endstart>>1;
factor=1<<(group+1);
for(int j=0;j<len;j++)
C[endstart-j-1]=factor*C[total-j-1];
}
for(i=1;i<total;i++)
C[i]=1.0/(2.0*cos(C[i]*PI/(total<<1))); //C[0]空着,沒使用
}
//*******************************************************************
//方法名稱:bitrev
//說明:位反序函數
//日期:2012.2.29
//參數說明:
//bi:不清楚
//deg:為IDCT變換資料長度的幂
//*******************************************************************
int bitrev(int bi,int deg)
{
int j=1,temp=0,degnum,halfnum;
degnum=deg;
//if(deg<0) return 0;
if(deg==0) return bi;
halfnum=1<<(deg-1);
while(halfnum)
{
if(halfnum&bi)
temp+=j;
j<<=1;
halfnum>>=1;
}
return temp;
}
//*******************************************************************
//方法名稱:fbitrev
//說明: 實作位反序排列(注意和上面函數的差別)
//日期:2012.2.29
//參數說明:
//f:不清楚
//deg:為DCT變換資料長度的幂
//*******************************************************************
void fbitrev(double *f,int deg)
{
if(deg==1) return;
int len=(1<<deg)-1;
int i=1;
int ii;
double temp;
while(i<len)
{ ii=bitrev(i,deg);
if(ii>i)
{ temp=f[ii];
f[ii]=f[i];
f[i]=temp;
}
i++;
}
}
//*******************************************************************
//方法名稱:idct_forward
//說明:IDCT變換過程可分為兩部分:前向追底和後向回根,這個是實作前向追底。
//日期:2012.2.29
//參數說明:
//f:存儲IDCT資料
//deg:為IDCT變換資料長度的幂
//*******************************************************************
void idct_forward(double *F,int deg)
{
int total,i_deg,wing,wings,halfwing,winglen,i_halfwing,temp1,temp2;
total=1<<deg;
for(i_deg=0;i_deg<deg;i_deg++)
{
wings=1<<i_deg;
winglen=1<<(deg-i_deg);
halfwing=winglen>>1;
for(wing=0;wing<wings;wing++)
{
for(i_halfwing=halfwing-1;i_halfwing>=0;i_halfwing--)
{
if(i_halfwing==0)
F[halfwing+wing*winglen+i_halfwing]=2.0*F[halfwing+wing*winglen+i_halfwing];
else
{
temp1=bitrev(i_halfwing,deg-i_deg-1); // bitrev為位反序
temp2=bitrev(i_halfwing-1,deg-i_deg-1);// 第一參數為要變換的數
// 第二參數為二進制長度
F[halfwing+wing*winglen+temp1]=F[halfwing+wing*winglen+temp1]+F[halfwing+wing*winglen+temp2];
}
}
}
}
}
//*******************************************************************
//方法名稱:idct_backward
//說明:IDCT變換過程可分為兩部分:前向追底和後向回根,這個是實作後向回根。
//日期:2012.2.29
//參數說明:
//F:存儲DCT資料
//deg:為DCT變換資料長度的幂
//*******************************************************************
void idct_backward(double *F,int deg)
{
int i_deg,i_halfwing,total,wing,wings,winglen,halfwing;
double temp1,temp2;
total=1<<deg;
for(i_deg=deg-1;i_deg>=0;i_deg--)
{
wings=1<<i_deg;
winglen=total>>i_deg;
halfwing=winglen>>1;
for(wing=0;wing<wings;wing++)
{
for(i_halfwing=0;i_halfwing<halfwing;i_halfwing++)
{
temp1=F[wing*winglen+i_halfwing];
temp2=F[(wing+1)*winglen-1-i_halfwing]*C[winglen-1-i_halfwing];
if(wing%2)
{
F[wing*winglen+i_halfwing]=(temp1-temp2)*0.5;
F[(wing+1)*winglen-1-i_halfwing]=(temp1+temp2)*0.5;
}
else
{
F[wing*winglen+i_halfwing]=(temp1+temp2)*0.5;
F[(wing+1)*winglen-1-i_halfwing]=(temp1-temp2)*0.5;
}
}
}
}
}
//*******************************************************************
//方法名稱:fidct_1D_no_param
//說明:完整實作一維DCT變換
//日期:2012.2.29
//參數說明:
//F:存儲DCT資料
//deg:為DCT變換資料長度的幂
//*******************************************************************
void fidct_1D_no_param(double *F,int deg)
{
initIDCTParam(deg);//111
F[0]=F[0]*sqrt(2.0);
fbitrev(F,deg); // 實作位反序排列
idct_forward(F,deg);
idct_backward(F,deg);
}
void fidct_1D(double *F,int deg)
{
int total=1<<deg;
double param=sqrt((double)total/2.0);
for(int i=0;i<total;i++)
{
F[i]=param*F[i];
}
fidct_1D_no_param(F,deg);
}
};
#endif // __JENC__
2、IFDCT.c
// NewIFDCT.cpp : 定義控制台應用程式的入口點。
//
#include "stdafx.h"
#include "fidct.h" //引用頭檔案
//實作一維反DCT快速變換
int _tmain(int argc, _TCHAR* argv[])
{
Fidct fidct;
fidct.Invoke();
return 0;
}