fdwt.c

这是我刚刚完成的关于JPEG2000的C语言实现的部分程序。小波变换是采用97变换

/******************************************************************************
功能描述:
   Implements 9/7 小波变换的实现
******************************************************************************/
/*****************************************************************************/
// 程序: 小波变换
// 版本:  V0.0
// 作者: 胡运平
// 最后修改时间 : 12, 5, 2005
/****************************************************************************/
#include <stdio.h>
#include<math.h>
#include"fdwt.h"
#define debug
#define fdwt	
void dwt(float buffer[LENGTH+8][WIDTH+8],float interleave[LENGTH+8][WIDTH+8],
		 float bufferout[LENGTH+8][WIDTH+8])
{
	int i,j,level;	
	float lifting_factors[5];
	lifting_factors[0] = (float) -1.586134342;
   	lifting_factors[1] = (float) -0.052980118;
 	lifting_factors[2] = (float)  0.882911075;
   	lifting_factors[3] = (float)  0.443506852;//
	lifting_factors[4] = (float)  1.230174105;//k
/*****************************************************************************
以下完成列变换
******************************************************************************/	
	for(level = 0;level<LEVEL;level++)
	{			
		//长度扩展	
		for(i = 0;i<LENGTH/(1<<level);i++)		
		{
			buffer[i+4][0] = buffer[i+4][8];
			buffer[i+4][1] = buffer[i+4][7];
			buffer[i+4][2] = buffer[i+4][6];
			buffer[i+4][3] = buffer[i+4][5];
			//把第二次变换时的数取出防止覆盖

			//516 =514 517= 513 518 = 512	
			buffer[i+4][WIDTH/(1<<level)+4] = buffer[i+4][WIDTH/(1<<level)+2];
			buffer[i+4][WIDTH/(1<<level)+5] = buffer[i+4][WIDTH/(1<<level)+1];
			buffer[i+4][WIDTH/(1<<level)+6] = buffer[i+4][WIDTH/(1<<level)];
			buffer[i+4][WIDTH/(1<<level)+7] = buffer[i+4][WIDTH/(1<<level)-1];
			}
		for(i = 0;i<LENGTH/(1<<level);i++)//行
		{	
			for(j =2;j<WIDTH/(1<<level)+8;j+=2)//第一步列3 5 7 9  11 13 --511 513 515 517 519 
			{//j = 2 4 6 8 --510 512 514 516 518		
				bufferout[i+4][j] =  buffer[i+4][j]+
				lifting_factors[0]*(buffer[i+4][j-1]+buffer[i+4][j+1]);//dl1奇数
			}
			for(j =3;j<WIDTH+7;j+=2)//第二步列 4 6 8 10 -- 510 512 514 516 518
			{//j = 3 5 7 9 --509 511 513 515 517
				bufferout[i+4][j] = buffer[i+4][j]+
				lifting_factors[1]*(bufferout[i+4][j-1]+bufferout[i+4][j+1]);//sl1偶数	
			}
			for(j =4;j<WIDTH/(1<<level)+6;j+=2)//第三步列5 7 9 11 -- 509 511 513 515 517 
			{//j = 4 6 8 10  -- 508 510 512 514 516
				bufferout[i+4][j] = bufferout[i+4][j]+
				lifting_factors[2]*(bufferout[i+4][j-1]+bufferout[i+4][j+1]);//dl1奇数
			}

			for(j =5;j<WIDTH/(1<<level)+5;j+=2)//第四步列 6 8 10 12 －－508 510 512 514 516
			{//j = 5 7 9 11 --507 509 511 513  515
				bufferout[i+4][j] = bufferout[i+4][j]+
				lifting_factors[3]*(bufferout[i+4][j-1]+bufferout[i+4][j+1]);//sl1偶数
			}
			for(j =4;j<WIDTH/(1<<level)+4;j+=2)//第五步列 5 7 9 11 －－ 509 511 513 515 517
			{//j = 4 6 8 10 －－ 508 510 512 514 516
				bufferout[i+4][j] = -lifting_factors[4]*bufferout[i+4][j];
			}
			for(j =5;j<WIDTH/(1<<level)+5;j+=2)//第六步列 6 8 10 12 -- 508 510 512 514 516
			{//j = 5 7 9 11 -- 507 509 511 513 515
				bufferout[i+4][j] = bufferout[i+4][j]/lifting_factors[4];
			}
	
		}//行数从第5行开始
		 //每一列的数放入5 6 7 8 －－511 512 513 514中
		 //数组中的j =  4 5 6 7 - - 510 511 512 513
/***********************************************************************
以上完成列变换，以下完成行变换
************************************************************************/
	//宽度扩展对称延拓
		for(j = 0;j<WIDTH/(1<<level);j++)
		{
			bufferout[0][j+4] = bufferout[8][j+4];
			bufferout[1][j+4] = bufferout[7][j+4];
			bufferout[2][j+4] = bufferout[6][j+4];
			bufferout[3][j+4] = bufferout[5][j+4];
			
			//132 = 130 133=129 134 = 128
			bufferout[LENGTH/(1<<level)+4][j+4] = bufferout[LENGTH/(1<<level)+2][j+4];
			bufferout[LENGTH/(1<<level)+5][j+4] = bufferout[LENGTH/(1<<level)+1][j+4];
			bufferout[LENGTH/(1<<level)+6][j+4] = bufferout[LENGTH/(1<<level)][j+4];
			bufferout[LENGTH/(1<<level)+7][j+4] = bufferout[LENGTH/(1<<level)-1][j+4];
		}

	for(j = 0;j<WIDTH/(1<<level); j++)
	{
		for(i = 2;i<LENGTH/(1<<level)+8;i+=2)//第一步行3 5 7 9  11 13 --127 129 131 133 135
		{//j = 2 4 6 8 --126 128 130 132 134
		
			buffer[i][j+4] =  bufferout[i][j+4]+
			lifting_factors[0]*(bufferout[i-1][j+4]+bufferout[i+1][j+4]);//dl1奇数
		//	bufferout[i+4][j] =  buffer[i+4][j]+
		//	lifting_factors[0]*(buffer[i+4][j-1]+buffer[i+4][j+1]);//dl1奇数
		}


		for(i = 3;i<LENGTH/(1<<level)+7;i+=2)//第二步行 4 6 8 10 -- 126 128 130 132 134
		{//j = 3 5 7 9 --125 127 129 131 133
			buffer[i][j+4] = bufferout[i][j+4]+
			lifting_factors[1]*(buffer[i-1][j+4]+buffer[i+1][j+4]);//sl1偶数	
			//bufferout[i+4][j] = buffer[i+4][j]+
			//lifting_factors[1]*(bufferout[i+4][j-1]+bufferout[i+4][j+1]);//sl1偶数	
		}

		for(i = 4;i<LENGTH/(1<<level)+6;i+=2)//第三步行5 7 9 11 -- 125 127 129 131 133 
		{//j = 4 6 8 10  -- 124 126 128 130 132 
			buffer[i][j+4] = buffer[i][j+4]+
			lifting_factors[2]*(buffer[i-1][j+4]+buffer[i+1][j+4]);//dl1奇数
			//bufferout[i+4][j] = bufferout[i+4][j]+
			//lifting_factors[2]*(bufferout[i+4][j-1]+bufferout[i+4][j+1]);//dl1奇数
		}

		for(i = 5;i<LENGTH/(1<<level)+5;i+=2)//第四步行 6 8 10 12 －－124 126 128 130 132 
		{//j = 5 7 9 11 --123 125 127 129 131 
			buffer[i][j+4] = buffer[i][j+4]+
			lifting_factors[3]*(buffer[i-1][j+4]+buffer[i+1][j+4]);//sl1偶数
			//bufferout[i+4][j] = bufferout[i+4][j]+
			//lifting_factors[3]*(bufferout[i+4][j-1]+bufferout[i+4][j+1]);//sl1偶数
		}
		for(i = 4;i<LENGTH/(1<<level)+4;i+=2)//第五步行 5 7 9 11 －－ 125 127 129 131 133
		{//j = 4 6 8 10 －－ 124 126 128 130 132
			buffer[i][j+4] = -lifting_factors[4]*buffer[i][j+4];
		    //bufferout[i+4][j] = -lifting_factors[4]*bufferout[i+4][j];

		}
		for(i = 5;i<LENGTH/(1<<level)+5;i+=2)//第六步行 6 8 10 12 -- 124 126 128 130 132
		{//j = 5 7 9 11 -- 123 125 127 129 131 
			buffer[i][j+4] = buffer[i][j+4]/lifting_factors[4];
			//bufferout[i+4][j] = bufferout[i+4][j]/lifting_factors[4];
		}
	
	}//每一行的数放入 5 6 7 8 －－中127 128 129 130 131 132
	  //数组中的j =  4 5 6 7 - - 126 127 128 129 130 131	
/*****************************************************************************
解交织		    LL = bufferout[i+4][j+4]    //偶数行偶数列
				HL = bufferout[i+5][j+5]   //奇数行偶数列
				LH = bufferout[i+4][j+5]   //偶数行奇数列
				HH = bufferout[i+5][j+5]   //奇数行奇数列
*****************************************************************************/
		for(i = 0;i<LENGTH/(1<<(level+1));i++)//行128 0--64
		{	
			for(j = 0;j<WIDTH/(1<<(level+1));j++)//列512  00--256
			{
				interleave[i+4][j+4] = buffer[i*2+5][j*2+5];    //偶数行偶数列
				interleave[i+4][j+WIDTH/(1<<(level+1))+4] = buffer[i*2+4][j*2+5];  //奇数行偶数列
				interleave[i+LENGTH/(1<<(level+1))+4][j+4] = buffer[i*2+5][j*2+4];  //偶数行奇数列
				interleave[i+LENGTH/(1<<(level+1))+4][j+WIDTH/(1<<(level+1))+4] =
					buffer[i*2+4][j*2+4];// 奇数行奇数列				
			}		
		}
		for(i = 0;i<LENGTH/(1<<(level+1));i++)//行128 0--64
		{	
			for(j = 0;j<WIDTH/(1<<(level+1));j++)//列512  00--256
			{
				buffer[i+4][j+4] = interleave[i+4][j+4];     //偶数行偶数列	
			}		
		}
#ifdef debug
//printf("\nthe interleave number is:");
//		for(i = 0;i<LENGTH;i++)
//			for(j = 0;j<WIDTH;j++)
//			{
//				if(!j%8) printf("\n");
//				printf("%3d ",(int)(interleave[i+4][j+4]+0.5));
//			}
//		for(i = 0;i<LENGTH;i++)
//			for(j = 0;j<WIDTH;j++)
//			{