#include <stdio.h>
#include <memory.h>
#include <math.h>
#include <malloc.h>
int row=8,col=8;
char **sTable,**mainTable,minorTable[64],MortonOrder[64];
//小波实验数据
const int example[64] =
{
63,-34, 49, 10, 7, 13,-12, 7 ,
-31, 23, 14,-13, 3, 4, 6, -1 ,
15, 14, 3,-12, 5, -7, 3, 9 ,
-9, -7,-14, 8, 4, -2, 3, 2 ,
-5, 9, -1, 47, 4, 6, -2, 2,
3, 0, -3, 2, 3, -2, 0, 4 ,
2, -3, 6, -4, 3, 6, 3, 6 ,
5, 11, 5, 6, 0, 3, -4, 4
};
//小波系数的处理----扫描顺序
void MortonOrderScan()
{
int i,j,k,minx,miny,maxx,maxy,x,y,order;
MortonOrder[0]=0;order=0;
for(k=0;k<16;k++)
{
y=MortonOrder[k]/8;x=MortonOrder[k]%8;
minx=x<<1;miny=y<<1;
maxx=(x+1)<<1;maxy=(y+1)<<1;
for(i=miny;i<maxy;i++)
{
for(j=minx;j<maxx;j++)
{
MortonOrder[order++]=i*col+j;
}
}
}
}
//Returns max(abs(matrix)), the absolute maximum value of a matrix.
int matrix_2d_abs_max(int *m)
{
int i,tMax=0;
for(i=0;i<row*col;i++)
if(m[i]>tMax)
tMax=m[i];
return tMax;
}
//Returns 1 if descendance tree is a zerotree.
int zeroTree(int *m,int x,int y,int threshold)
{
int i,j,tMax=0;
int minx,miny,maxx,maxy;
minx=x<<1;miny=y<<1;
maxx=(x+1)<<1;maxy=(y+1)<<1;
if(minx>=col || miny>=row)
return 0;
while(maxx<=col && maxy<=row)
{
for(i=miny;i<maxy;i++)
{
for(j=minx;j<maxx;j++)
{
if(m[i*col+j]>=threshold && minorTable[i*col+j]!=1)
return 0;
}
}
minx=minx<<1;miny=miny<<1;maxx=maxx<<1;maxy=maxy<<1;
}
return 1;
}
// Returns a dominant-pass code from the alphabet [P,N,T,Z].
void code(int *m,int x,int y,int threshold,int times)
{
int i,j,minx,miny,maxx,maxy;
if(minorTable[y*col+x]==1)
{
mainTable[times][y*col+x]=' ';
return ;
}
//如果该小波系数的幅值大于当前阈值,则是大系数(但是不知道正负)
if(abs(m[y*col+x])>=threshold)
{
//如果它大于0,则返回为P(即"正大系数"),否则返回为NEG("负大系数")
if(m[y*col+x]>=0)
mainTable[times][y*col+x]='P';
else
mainTable[times][y*col+x]='N';
}
//如果它的幅值小于当前阈值,则是小系数(但不知道是否"零树根")
else
{
//如果零树判断为真,则返回T("零树根"),否则返回Z("孤独零")
if (zeroTree(m,x,y,threshold))
{
mainTable[times][y*col+x]='T';
minx=x<<1;miny=y<<1;
maxx=(x+1)<<1;maxy=(y+1)<<1;
while(maxx<=col && maxy<=row)
{
for(i=miny;i<maxy;i++)
{
for(j=minx;j<maxx;j++)
{
if(minorTable[i*col+j]!=1)
mainTable[times][i*col+j]='*';
}
}
minx=minx<<1;miny=miny<<1;maxx=maxx<<1;maxy=maxy<<1;
}
}
else
mainTable[times][y*col+x]='Z';
}
}
//编码
void EZW(int *m,int row,int col,int threshold,int times,FILE *fw)
{
printf("第%d次主扫描结果:\n",times+1);
fprintf(fw,"第%d次主扫描结果:\n",times+1);
int x,y,temp;
for(y=0;y<row;y++)
{
for(x=0;x<col;x++)
{
temp=y*col+x;
if(mainTable[times][temp]!='*')
code(m,x,y,threshold,times);
if(minorTable[temp]!=1)
{
printf("%c\t",mainTable[times][temp]);fprintf(fw,"%c\t",mainTable[times][temp]);
}
else
{
printf("%c\t",32);fprintf(fw,"%c\t",32);
}
}
printf("\n");fprintf(fw,"\n");
}
printf("为解码器提供的信息:\n");fprintf(fw,"为解码器提供的信息:\n");
printf("T%d=%d,",times,threshold);fprintf(fw,"T%d=%d,",times,threshold);
printf("D%d:",times+1);fprintf(fw,"D%d:",times+1);
for(y=0;y<row;y++)
{
for(x=0;x<col;x++)
{
temp=MortonOrder[y*col+x];
if(minorTable[temp]!=1 && mainTable[times][temp]!='*')
{
printf("%c",mainTable[times][temp]);fprintf(fw,"%c",mainTable[times][temp]);
if(mainTable[times][temp]=='P' || mainTable[times][temp]=='N')
minorTable[temp]=1;
}
}
}
printf(";S%d:",times+1);fprintf(fw,";S%d:",times+1);
for(y=0;y<row;y++)
{
for(x=0;x<col;x++)
{
temp=MortonOrder[y*col+x];
if(minorTable[temp]==1)
{
if(abs(m[temp])<abs(m[temp]/threshold*threshold)+threshold/2)
{
sTable[times][temp]=0;
printf("0");fprintf(fw,"0");
}
else
{
sTable[times][temp]=1;
printf("1");fprintf(fw,"1");
}
}
}
}
printf("\n\n");fprintf(fw,"\n\n");
}
//解码
void IEZW(int *iEZWData,int row,int col,int threshold,int times,FILE *fw)
{
printf("第%d次解码结果:\n",times+1);fprintf(fw,"第%d次解码结果:\n",times+1);
int i,j,k,x,y,temp,minx,miny,maxx,maxy;
for(i=0;i<64;i++)
{
if(iEZWData[i]!=0 && iEZWData[i]!=64)
{
if(sTable[times][i]==1)
{
if(iEZWData[i]<0)
{
iEZWData[i]=(int)(abs(iEZWData[i])/threshold)*threshold+(int)(threshold*0.75);
iEZWData[i]*=-1;
}
else
iEZWData[i]=(int)(iEZWData[i]/threshold)*threshold+(int)(threshold*0.75);
}
else if(sTable[times][i]==0)
{
if(iEZWData[i]<0)
{
iEZWData[i]=(int)(abs(iEZWData[i])/threshold)*threshold+(int)(threshold*0.25);
iEZWData[i]*=-1;
}
else
iEZWData[i]=(int)(iEZWData[i]/threshold)*threshold+(int)(threshold*0.25);
}
}
else
{
if(mainTable[times][i]=='P')
{
if(sTable[times][i]==0)
iEZWData[i]=(int)(threshold*1.25);
else
iEZWData[i]=(int)(threshold*1.75);
}
else if(mainTable[times][i]=='N')
{
if(sTable[times][i]==0)
iEZWData[i]=-1*(int)(threshold*1.25);
else
iEZWData[i]=-1*(int)(threshold*1.75);
}
else if(mainTable[times][i]=='Z')
{
iEZWData[i]=0;
}
else if(mainTable[times][i]=='T')
{
iEZWData[i]=0;
x=i%8;y=i/8;
minx=x<<1;miny=y<<1;
maxx=(x+1)<<1;maxy=(y+1)<<1;
while(maxx<=col && maxy<=row)
{
for(k=miny;k<maxy;k++)
{
for(j=minx;j<maxx;j++)
{
// if(minorTable[k*col+j]!=1)
if(iEZWData[k*col+j]==0)
iEZWData[k*col+j]=64;
}
}
minx=minx<<1;miny=miny<<1;maxx=maxx<<1;maxy=maxy<<1;
}
}
}
}
for(y=0;y<row;y++)
{
for(x=0;x<col;x++)
{
temp=y*col+x;
if(iEZWData[temp]==64)
{
printf("*\t");fprintf(fw,"*\t");
}
else
{
printf("%d\t",iEZWData[temp]);fprintf(fw,"%d\t",iEZWData[temp]);
}
}
printf("\n");fprintf(fw,"\n");
}
printf("\n");fprintf(fw,"\n");
}
void main()
{
int times,Num=6,*threshold;
//创建主表与辅助表
threshold=(int *)malloc(sizeof(int)*Num);
sTable = (char **)malloc(sizeof(char *)*Num);
mainTable = (char **)malloc(sizeof(char *)*Num);
for( int i=0; i<Num; i++ )
{
sTable[i] = (char *)malloc(sizeof(char)*64);
mainTable[i] = (char *)malloc(sizeof(char)*64);
memset(sTable[i],0,sizeof(char)*64);
memset(mainTable[i],0,sizeof(char)*64);
}
//准备编码后的数据流输出的文件
FILE *fw;
fw=fopen("result.txt","w");
int tMax=0;
int *imageData;
imageData=(int *)malloc(sizeof(int)*col*row);
memset(imageData,0,sizeof(int)*col*row);
memcpy(imageData,example,sizeof(int)*col*row);
//求矩阵的最大值
tMax=matrix_2d_abs_max(imageData);
//设置初始阈值
threshold[0] = 1 << (int)(floor(log10(tMax)/log10(2)));
for(times=1;times<Num;times++)
threshold[times]=threshold[times-1]>>1;
//小波系数扫描顺序
MortonOrderScan();
//编码
for(times=0;times<Num;times++)
{
EZW(imageData,row,col,threshold[times],times,fw);
}
//解码
int iEZWData[64];
for(i=0;i<64;i++)
iEZWData[i]=64;
for(times=0;times<Num;times++)
{
IEZW(iEZWData,row,col,threshold[times],times,fw);
}
//释放内存空间
fclose(fw);
for(i=0; i<Num; i++ )
{
free(mainTable[i]);free(sTable[i]);
}
free(sTable);free(mainTable);free(imageData);
}