dct_64.cpp
来自「DCT快速算法64位的C语言实现」· C++ 代码 · 共 1,770 行 · 第 1/4 页
CPP
1,770 行
t5=(double *)malloc(32); //为t5 分配4 个double的空间
t5[0]=t4[0]+t4[1];
t5[1]=(t4[0]-t4[1])*c_pi_1_4;
}
//03_02_3_2_2_2_3 t4[2]--t4[3]乘以B2
// 4次加法 3次乘法
{
t4[2]*=c_pi_1_8; t4[3]*=c_pi_3_8;
t5[2]=t4[2]+t4[3]; t5[3]=(t4[2]-t4[3])*c_pi_1_4;
t5[3]=t5[3]+t5[3];
t5[3]=t5[3]-t5[2];
}// 此时 t5 记录的是P4之前的列向量
//03_02_3_2_2_2_4 t5乘以P4 用t4记录
{
t4[0]=t5[0];
t4[1]=t5[2];
t4[2]=t5[1];
t4[3]=t5[3];
}
}//此时 t4[0]--t4[3]记录的是乘以A4后的t3[0]--t3[3]元素 t5 已没用
//03_02_3_2_2_3 t3[4]--t3[7]乘以B4
//
{
//03_02_3_2_2_3_1 t3[4]--t3[7]乘以B4分解式最右侧的以 cos(PI*(2*i+1)/2/8)为对角元的对角阵
// 0次加法 4次乘法
{
t3[4]*=c_pi_1_16; t3[5]*=c_pi_3_16;
t3[6]*=c_pi_5_16; t3[7]*=c_pi_7_16;
}
//03_02_3_2_2_3_2 接下来 t3[4]--t3[7]乘以A4
{
//03_02_3_2_2_3_2_1 t3[4]--t3[7]乘以A4分解式右侧的矩阵 t5记录
// 4次加法 0次乘法
{
t5[0]=t3[4]+t3[7]; t5[2]=t3[4]-t3[7];
t5[1]=t3[5]+t3[6]; t5[3]=t3[5]-t3[6];
}
//03_02_3_2_2_3_2_2 t5[0]--t5[1]乘以乘以A2 这时可以直接乘了！！！！
// 2 次加法 1 次乘法
{
t6=(double *)malloc(32); //为t6 分配4 个double的空间
t6[0]=t5[0]+t5[1];
t6[1]=(t5[0]-t5[1])*c_pi_1_4;
}
//03_02_3_2_2_3_2_3 t5[2]--t5[3]乘以B2
// 4次加法 3次乘法
{
t5[2]*=c_pi_1_8; t5[3]*=c_pi_3_8;
t6[2]=t5[2]+t5[3]; t6[3]=(t5[2]-t5[3])*c_pi_1_4;
t6[3]=t6[3]+t6[3];
t6[3]=t6[3]-t6[2];
}// 此时 t6 记录的是P4之前的列向量
//03_02_3_2_2_3_2_4 t6乘以P4 用t5记录
//
{
t5[0]=t6[0];
t5[1]=t6[2];
t5[2]=t6[1];
t5[3]=t6[3];
}//此时 t5[0]--t5[3]记录的是乘以A4后的t3[4]--t3[7]元素 t6 已没用
}
//03_02_3_2_2_3_3 t5[0]--t5[3]乘以 一二阵 和 正负一阵
// 9次加法 0次乘法
{
t5[1]=t5[1]+t5[1]; t5[2]=t5[2]+t5[2]; t5[3]=t5[3]+t5[3];
t6[0]=t5[0];
t6[1]=t5[1]-t5[0];
t6[2]=t5[0]-t5[1]+t5[2];
t6[3]=t5[1]-t5[0]+t5[3]-t5[2];
}
}// 此时t6[0]--t6[3]记录的是t3[4]--t3[7]乘以B4后的元素
//03_02_3_2_2_4 P8之前的t3[0]--t3[7]乘以P8, 仍记录在t2[0]--t2[7]里
//
{
t2[0]=t4[0];
t2[1]=t6[0];
t2[2]=t4[1];
t2[3]=t6[1];
t2[4]=t4[2];
t2[5]=t6[2];
t2[6]=t4[3];
t2[7]=t6[3];
}
}//此时 t2[0]--t2[7]记录的是乘完A8的t2[0]--t2[7]
//03_02_3_2_3 t2[8]--t2[15]乘以B8 仍保存到t2[8]--t2[15]里
{
//03_02_3_2_3_1 t2[8]--t2[15]乘以B8分解式最右侧的以 cos(PI*(2*i+1)/2/16)为对角元的对角阵
// 0次加法 8次乘法
{
t2[8]*=c_pi_1_32; t2[9]*=c_pi_3_32;
t2[10]*=c_pi_5_32; t2[11]*=c_pi_7_32;
t2[12]*=c_pi_9_32; t2[13]*=c_pi_11_32;
t2[14]*=c_pi_13_32; t2[15]*=c_pi_15_32;
}
//03_02_3_2_3_2 接下来 t2[8]--t2[15]乘以A8
{
//03_02_3_2_3_2_1 tt2[8]--t2[15]乘以A8后分解式右侧的矩阵
// 8次加法 0次乘法
{
free(t4);
t4=(double *)malloc(64); //重新为t4 分配8 个double的空间
t4[0]=t2[8]+t2[15]; t4[4]=t2[8]-t2[15];
t4[1]=t2[9]+t2[14]; t4[5]=t2[9]-t2[14];
t4[2]=t2[10]+t2[13]; t4[6]=t2[10]-t2[13];
t4[3]=t2[11]+t2[12]; t4[7]=t2[11]-t2[12];
}
//03_02_3_2_3_2_2 t4[0]--t4[3]乘以A4
{
//03_02_3_2_3_2_2_1 t4[0]--t4[3]乘以A4分解式右侧的矩阵
// 4次加法 0次乘法
{
t5[0]=t4[0]+t4[3]; t5[2]=t4[0]-t4[3];
t5[1]=t4[1]+t4[2]; t5[3]=t4[1]-t4[2];
}
//03_02_3_2_3_2_2_2 t5[0]--t5[1]乘以A2 这时可以直接乘了！！！！
// 2 次加法 1 次乘法
{
t6[0]=t5[0]+t5[1];
t6[1]=(t5[0]-t5[1])*c_pi_1_4;
}
//03_02_3_2_3_2_2_3 t5[2]--t5[3]乘以B2
// 4次加法 3次乘法
{
t5[2]*=c_pi_1_8; t5[3]*=c_pi_3_8;
t6[2]=t5[2]+t5[3]; t6[3]=(t5[2]-t5[3])*c_pi_1_4;
t6[3]=t6[3]+t6[3];
t6[3]=t6[3]-t6[2];
}// 此时 t6 记录的是P4之前的列向量
//03_02_3_2_3_2_2_4 t6乘以P4 用t5记录
{
t5[0]=t6[0];
t5[1]=t6[2];
t5[2]=t6[1];
t5[3]=t6[3];
}
}//此时 t5[0]--t5[3]记录的是乘以A4后的t4[0]--t4[3]元素 t6 已没用
//03_02_3_2_3_2_3 t4[4]--t4[7]乘以B4
//
{
//03_02_3_2_3_2_3_1 t4[4]--t4[7]乘以B4分解式最右侧的以 cos(PI*(2*i+1)/2/8)为对角元的对角阵
// 0次加法 4次乘法
{
t4[4]*=c_pi_1_16; t4[5]*=c_pi_3_16;
t4[6]*=c_pi_5_16; t4[7]*=c_pi_7_16;
}
//03_02_3_2_3_2_3_2 接下来 t4[4]--t4[7]乘以A4
{
//03_02_3_2_3_2_3_2_1 t4[4]--t4[7]乘以A4分解式右侧的矩阵 t6记录
// 4次加法 0次乘法
{
t6[0]=t4[4]+t4[7]; t6[2]=t4[4]-t4[7];
t6[1]=t4[5]+t4[6]; t6[3]=t4[5]-t4[6];
}
//03_02_3_2_3_2_3_2_2 t6[0]--t6[1]乘以乘以A2 这时可以直接乘了！！！！
// 2 次加法 1 次乘法
{
t7=(double *)malloc(32); //为t7 分配4 个double的空间
t7[0]=t6[0]+t6[1];
t7[1]=(t6[0]-t6[1])*c_pi_1_4;
}
//03_02_3_2_3_2_3_2_3 t6[2]--t6[3]乘以B2
// 4次加法 3次乘法
{
t6[2]*=c_pi_1_8; t6[3]*=c_pi_3_8;
t7[2]=t6[2]+t6[3]; t7[3]=(t6[2]-t6[3])*c_pi_1_4;
t7[3]=t7[3]+t7[3];
t7[3]=t7[3]-t7[2];
}// 此时 t7 记录的是P4之前的列向量
//03_02_3_2_3_2_3_2_4 t7乘以P4 用t6记录
//
{
t6[0]=t7[0];
t6[1]=t7[2];
t6[2]=t7[1];
t6[3]=t7[3];
}//此时 t6[0]--t6[3]记录的是乘以A4后的t4[4]--t4[7]元素 t7 已没用
}
//03_02_3_2_3_2_3_3 t6[0]--t6[3]乘以 一二阵 和 正负一阵
// 9次加法 0次乘法
{
t6[1]=t6[1]+t6[1]; t6[2]=t6[2]+t6[2]; t6[3]=t6[3]+t6[3];
t7[0]=t6[0];
t7[1]=t6[1]-t6[0];
t7[2]=t6[0]-t6[1]+t6[2];
t7[3]=t6[1]-t6[0]+t6[3]-t6[2];
}// 此时t7[0]--t7[3]记录的是t4[4]--t4[7]乘以B4后的元素
}
//03_02_3_2_3_2_4 P8之前的43[0]--t4[7]乘以P8,记录在t4[0]--t4[7]里
{
t4[0]=t5[0];
t4[1]=t7[0];
t4[2]=t5[1];
t4[3]=t7[1];
t4[4]=t5[2];
t4[5]=t7[2];
t4[6]=t5[3];
t4[7]=t7[3];
}
}//此时t4[0]--t4[7]记录的是乘完A8后的t2[8]--t2[15]
//03_02_3_2_3_3 t4[0]--t4[7]乘以 一二阵 和 正负一阵，得到 ,记录在t2[8]--t2[15]里
{
t4[1]=t4[1]+t4[1];
t4[2]=t4[2]+t4[2];
t4[3]=t4[3]+t4[3];
t4[4]=t4[4]+t4[4];
t4[5]=t4[5]+t4[5];
t4[6]=t4[6]+t4[6];
t4[7]=t4[7]+t4[7];
t2[8]=t4[0];
t2[9]=t4[1]-t2[8];
t2[10]=t4[2]-t2[9];
t2[11]=t4[3]-t2[10];
t2[12]=t4[4]-t2[11];
t2[13]=t4[5]-t2[12];
t2[14]=t4[6]-t2[13];
t2[15]=t4[7]-t2[14];
}
}// 此时 t2[8]--t2[15]记录的是乘完B8的t2[8]--t2[15]
//03_02_3_2_4 乘完A8 B8的t2[0]--t2[15] 再乘以P16
{
t1[0]=t2[0]; t1[1]=t2[8];
t1[2]=t2[1]; t1[3]=t2[9];
t1[4]=t2[2]; t1[5]=t2[10];
t1[6]=t2[3]; t1[7]=t2[11];
t1[8]=t2[4]; t1[9]=t2[12];
t1[10]=t2[5]; t1[11]=t2[13];
t1[12]=t2[6]; t1[13]=t2[14];
t1[14]=t2[7]; t1[15]=t2[15];
}
}
//03_02_3_3 乘完A16的t1[16]--t1[31]乘以 一二阵 和 正负一阵
{
t1[1]=t1[1]+t1[1];
t1[2]=t1[2]+t1[2];
t1[3]=t1[3]+t1[3];
t1[4]=t1[4]+t1[4];
t1[5]=t1[5]+t1[5];
t1[6]=t1[6]+t1[6];
t1[7]=t1[7]+t1[7];
t1[8]=t1[8]+t1[8];
t1[9]=t1[9]+t1[9];
t1[10]=t1[10]+t1[10];
t1[11]=t1[11]+t1[11];
t1[12]=t1[12]+t1[12];
t1[13]=t1[13]+t1[13];
t1[14]=t1[14]+t1[14];
t1[15]=t1[15]+t1[15];
f[16]=t1[0];
f[17]=t1[1]-f[16];
f[18]=t1[2]-f[17];
f[19]=t1[3]-f[18];
f[20]=t1[4]-f[19];
f[21]=t1[5]-f[20];
f[22]=t1[6]-f[21];
f[23]=t1[7]-f[22];
f[24]=t1[8]-f[23];
f[25]=t1[9]-f[24];
f[26]=t1[10]-f[25];
f[27]=t1[11]-f[26];
f[28]=t1[12]-f[27];
f[29]=t1[13]-f[28];
f[30]=t1[14]-f[29];
f[31]=t1[15]-f[30];
}
}
//03_02_4 f[0]--f[31]乘以P32
{
t1[0]=f[0]; t1[1]=f[16];
t1[2]=f[1]; t1[3]=f[17];
t1[4]=f[2]; t1[5]=f[18];
t1[6]=f[3]; t1[7]=f[19];
t1[8]=f[4]; t1[9]=f[20];
t1[10]=f[5]; t1[11]=f[21];
t1[12]=f[6]; t1[13]=f[22];
t1[14]=f[7]; t1[15]=f[23];
t1[16]=f[8]; t1[17]=f[24];
t1[18]=f[9]; t1[19]=f[25];
t1[20]=f[10]; t1[21]=f[26];
t1[22]=f[11]; t1[23]=f[27];
t1[24]=f[12]; t1[25]=f[28];
t1[26]=f[13]; t1[27]=f[29];
t1[28]=f[14]; t1[29]=f[30];
t1[30]=f[15]; t1[31]=f[31];
}
}
//03_03 乘以一二阵和正负一阵
{
//03_03_01 乘以一二阵
{
t1[1]+=t1[1];
t1[2]+=t1[2]; t1[3]+=t1[3];
t1[4]+=t1[4]; t1[5]+=t1[5];
t1[6]+=t1[6]; t1[7]+=t1[7];
t1[8]+=t1[8]; t1[9]+=t1[9];
t1[10]+=t1[10]; t1[11]+=t1[11];
t1[12]+=t1[12]; t1[13]+=t1[13];
t1[14]+=t1[14]; t1[15]+=t1[15];
t1[16]+=t1[16]; t1[17]+=t1[17];
t1[18]+=t1[18]; t1[19]+=t1[19];
t1[20]+=t1[20]; t1[21]+=t1[21];
t1[22]+=t1[22]; t1[23]+=t1[23];
t1[24]+=t1[24]; t1[25]+=t1[25];
t1[26]+=t1[26]; t1[27]+=t1[27];
t1[28]+=t1[28]; t1[29]+=t1[29];
t1[30]+=t1[30]; t1[31]+=t1[31];
}
//03_03_02 乘以正负一阵
{
t0[32]=t1[0];
t0[33]=t1[1]-t0[32];
t0[34]=t1[2]-t0[33];
t0[35]=t1[3]-t0[34];
t0[36]=t1[4]-t0[35];
t0[37]=t1[5]-t0[36];
t0[38]=t1[6]-t0[37];
t0[39]=t1[7]-t0[38];
t0[40]=t1[8]-t0[39];
t0[41]=t1[9]-t0[40];
t0[42]=t1[10]-t0[41];
t0[43]=t1[11]-t0[42];
t0[44]=t1[12]-t0[43];
t0[45]=t1[13]-t0[44];
t0[46]=t1[14]-t0[45];
t0[47]=t1[15]-t0[46];
t0[48]=t1[16]-t0[47];
t0[49]=t1[17]-t0[48];
t0[50]=t1[18]-t0[49];
t0[51]=t1[19]-t0[50];
t0[52]=t1[20]-t0[51];
t0[53]=t1[21]-t0[52];
t0[54]=t1[22]-t0[53];
t0[55]=t1[23]-t0[54];
t0[56]=t1[24]-t0[55];
t0[57]=t1[25]-t0[56];
t0[58]=t1[26]-t0[57];
t0[59]=t1[27]-t0[58];
t0[60]=t1[28]-t0[59];
t0[61]=t1[29]-t0[60];
t0[62]=t1[30]-t0[61];
t0[63]=t1[31]-t0[62];
}
}
}
//04 乘完 A32 B32 的t0[0]--t0[63] 再乘以 P64
{
F=(double *)malloc(512);
F[0]=0.02209708691208*t0[0]; F[1]=0.03125000000000*t0[32];
F[2]=0.03125000000000*t0[1]; F[3]=0.03125000000000*t0[33];
F[4]=0.03125000000000*t0[2]; F[5]=0.03125000000000*t0[34];
F[6]=0.03125000000000*t0[3]; F[7]=0.03125000000000*t0[35];
F[8]=0.03125000000000*t0[4]; F[9]=0.03125000000000*t0[36];
F[10]=0.03125000000000*t0[5]; F[11]=0.03125000000000*t0[37];
F[12]=0.03125000000000*t0[6]; F[13]=0.03125000000000*t0[38];
F[14]=0.03125000000000*t0[7]; F[15]=0.03125000000000*t0[39];
F[16]=0.03125000000000*t0[8]; F[17]=0.03125000000000*t0[40];
F[18]=0.03125000000000*t0[9]; F[19]=0.03125000000000*t0[41];
F[20]=0.03125000000000*t0[10]; F[21]=0.03125000000000*t0[42];
F[22]=0.03125000000000*t0[11]; F[23]=0.03125000000000*t0[43];
F[24]=0.03125000000000*t0[12]; F[25]=0.03125000000000*t0[44];
F[26]=0.03125000000000*t0[13]; F[27]=0.03125000000000*t0[45];
F[28]=0.03125000000000*t0[14]; F[29]=0.03125000000000*t0[46];
F[30]=0.03125000000000*t0[15]; F[31]=0.03125000000000*t0[47];
F[32]=0.03125000000000*t0[16]; F[33]=0.03125000000000*t0[48];
F[34]=0.03125000000000*t0[17]; F[35]=0.03125000000000*t0[49];
F[36]=0.03125000000000*t0[18]; F[37]=0.03125000000000*t0[50];
F[38]=0.03125000000000*t0[19]; F[39]=0.03125000000000*t0[51];
F[40]=0.03125000000000*t0[20]; F[41]=0.03125000000000*t0[52];
F[42]=0.03125000000000*t0[21]; F[43]=0.03125000000000*t0[53];
F[44]=0.03125000000000*t0[22]; F[45]=0.03125000000000*t0[54];
F[46]=0.03125000000000*t0[23]; F[47]=0.03125000000000*t0[55];
F[48]=0.03125000000000*t0[24]; F[49]=0.03125000000000*t0[56];
F[50]=0.03125000000000*t0[25]; F[51]=0.03125000000000*t0[57];
F[52]=0.03125000000000*t0[26]; F[53]=0.03125000000000*t0[58];
F[54]=0.03125000000000*t0[27]; F[55]=0.03125000000000*t0[59];
F[56]=0.03125000000000*t0[28]; F[57]=0.03125000000000*t0[60];
F[58]=0.03125000000000*t0[29]; F[59]=0.03125000000000*t0[61];
F[60]=0.03125000000000*t0[30]; F[61]=0.03125000000000*t0[62];
F[62]=0.03125000000000*t0[31]; F[63]=0.03125000000000*t0[63];
}
free(t0); free(t1); free(t2); free(t3); free(t4); free(t5); free(t6); free(t7);
return(F);
}⌨️ 快捷键说明
复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?