// LiftDwt.cpp: implementation of the CLiftDwt class.
// 作者：張增輝   國防科技大學理學院數學系
// QQ:15105473   email: zenghui1980@163.com
//
//  此程序為二進制9-7小波提升格式的變換，算法中只有整數加法和移位
//  程序中了使用的小波濾波器為2000年 IEEE Image Pro. 上的一篇介紹jpeg2000無損壓縮中的9-7小波濾波器，//  它其實是對cdf9-7提升格式中的幾個系數直接作二進制近似的到的
//
//  程序包括以下部分： 
//
//  Dwt2D(double buffer[], int height, int width, int lv) 二維提升小波變換，buffer進行變換的數據，he//  ight,width數據的尺寸，lv小波分解的層數
//  IDwt2D(double buffer[], int height, int width, int lv) 二維提升小波反變換
//  
//  Dwt1D(double *buffer, int buflen)  一維提升小波分解 buflen 數據的長度
//  IDwt1D(double *buffer, int buflen) 
//   
//
//  數據延拓方式有兩種： 周期延拓和對稱延拓
//////////////////////////////////////////////////////////////////////
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
#include "ImgPro.h"
#include "BinDwt97.h"
#include "math.h"

#ifdef _DEBUG
#undef THIS_FILE
static char THIS_FILE[]=__FILE__;
#define new DEBUG_NEW
#endif

//#define MAX(x,y) (x>y?x:y)
//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

CBinDwt97::CBinDwt97()
{
	shift = 4;
}

CBinDwt97::~CBinDwt97()
{

}

void CBinDwt97::Dwt2D(int **buffer, int height, int width, int lv)
{
	int i,j,k;
	int nheight,nwidth;
 
	for ( k=0; k<lv; k++){

       nheight=height>>k;
	   nwidth=width>>k;

       //行變換
	   int *pdata;
	   pdata = new int [nwidth];
		for (i=0; i<nheight; i++){
			for(j=0; j<nwidth; j++)
				*(pdata+j) = *(*(buffer+i)+j);
			Dwt1D(pdata,nwidth);
			for(j=0; j<nwidth; j++)
				*(*(buffer+i)+j) = *(pdata+j);
		}
		delete []pdata;
       
		//列變換
		int *p1data;
		p1data = new int [nheight];
		for(j=0; j<nwidth; j++){
			for (i=0; i<nheight; i++)
				*(p1data+i)=*(*(buffer+i)+j);
			Dwt1D(p1data,nheight);
			for(i=0; i<nwidth; i++)
				*(*(buffer+i)+j) = *(p1data+i);
		}
		delete []p1data;

	}

}

void CBinDwt97::IDwt2D(int **buffer, int height, int width, int lv)
{

	int i,j,k;
 	int nheight,nwidth;
  
	for (k=0; k<lv; k++){

 		nheight=height>>(lv-k-1);
 		nwidth=width>>(lv-k-1);

        //列反變換
		int *pdata;
		pdata = new int [nheight];
 		for(j=0; j<nwidth; j++){
 			for (i=0; i<nheight; i++)
 				*(pdata+i) = *(*(buffer+i)+j);
 			IDwt1D(pdata,nheight);
 			for(i=0; i<nwidth; i++)
 				*(*(buffer+i)+j) = *(pdata+i);
 		}
 		delete []pdata;
 
        //行反變換
		int *p1data;
		p1data = new int [nwidth];
 		for (i=0; i<nheight; i++){
 			for(j=0; j<nwidth; j++)
 				*(p1data+j)=*(*(buffer+i)+j);
 			IDwt1D(p1data,nwidth);
 			for(j=0; j<nwidth; j++)
 				*(*(buffer+i)+j) = *(p1data+j);
 		}
		delete []p1data;
 	}
	

}

void CBinDwt97::Dwt1D(int *buffer, int buflen)
{
	int i;
	int itemp;
	int data_tmp;

	int *d, *s, *p;
	p = new int [buflen+(shift<<1)];
	d = new int [(buflen>>1)+shift];
	s = new int [(buflen>>1)+shift];

/*    使用周期延拓方式	*/
/*	for (i=0; i<shift; i++){
		*(p+i) = *(buffer+buflen-shift+i);
		*(p+buflen+shift+i) = *(buffer+i);
	}
	for (i=0; i<buflen; i++) *(p+i+shift) = *(buffer+i);
/* 數據周期延拓結束  計算p1 = d1          p2 = s1 */

	/*  對稱延拓 */
	for (i=0; i<shift; i++){
		*(p+i) = *(buffer+shift-i);
		*(p+buflen+shift+i) = *(buffer+buflen-i-2);
	}
	for (i=0; i<buflen; i++) *(p+i+shift) = *(buffer+i);
    /*  結束  */



	for (i=0; i<(buflen>>1)+shift; i++){
		itemp = i-(shift>>1);
		data_tmp = - (*(p+shift+(itemp<<1)+2) + *(p+shift+(itemp<<1)));
		*(d+i) = *(p+shift+(itemp<<1)+1) + data_tmp + (data_tmp>>1) + (data_tmp>>4) + (data_tmp >>6) + (data_tmp >>7);
	}

	for (i=0; i<(buflen>>1)+shift-1; i++){
		itemp = i+1-(shift>>1);
		data_tmp = -(*(d+i+1) + *(d+i));
        *(s+i+1) = *(p+shift+(itemp<<1)) + (data_tmp >>5) +(data_tmp >>6) + (data_tmp >>8) + (data_tmp >>9) + (data_tmp >>12);
	}
   /*  d1 s1 計算結束 下計算 d2  */
	for (i=0; i<(buflen>>1)+shift-1; i++){
		data_tmp = *(s+i) + *(s+i+1);
        *(p+shift+(buflen>>1)+i) = *(d+i) + (data_tmp>>1) + (data_tmp >>2) + (data_tmp >>3) + (data_tmp >>7);
	}

	/*  s2 */
    for (i=0; i<(buflen>>1)+shift-1; i++){
		data_tmp = *(p+(buflen>>1)+shift+i+1) + *(p+(buflen>>1)+shift+i);
		*(p+i+1) = *(s+i+1) + (data_tmp >>2) + (data_tmp >>3) + (data_tmp >>4) + (data_tmp >>8) + (data_tmp >>9) + (data_tmp >>12);
	}

	///  存在冗余計算
	for (i=0; i<(buflen>>1); i++)
		*(buffer+i+(buflen>>1)) = *(p+shift+(buflen>>1) + (shift>>1) +i);

	for (i=0; i<(buflen>>1); i++)
		*(buffer+i) = *(p + (shift>>1) + i);

    delete []p;
	delete []d;
	delete []s;
    
}

void CBinDwt97::IDwt1D(int *buffer, int buflen)
{
	int i;
    int *p1, *p2, *s, *d;
	int data_tmp;

	p1 = new int [(buflen>>1)+shift];
	p2 = new int [(buflen>>1)+shift];
	s  = new int [(buflen>>1)+shift];
	d  = new int [(buflen>>1)+shift];

/*    延拓方式為周期延拓     /*
	for (i=0; i<(shift>>1); i++){
		p1[i] = buffer[i+(buflen>>1)-(shift>>1)];
		p1[i+(shift>>1)+(buflen>>1)] = buffer[i];
		p2[i] = buffer[i+buflen-(shift>>1)];
		p2[i+(shift>>1)+(buflen>>1)] = buffer[i+(buflen>>1)];
	}
	for (i=0; i<(buflen>>1); i++){
		p1[i+(shift>>1)] = buffer[i];
		p2[i+(shift>>1)] = buffer[i+(buflen>>1)];
	}

/*  延拓結束 */	

	/*  對稱延拓	*/
   	for (i=0; i<(shift>>1); i++){
		p1[i] = buffer[(shift>>1)-i];
		p1[i+(shift>>1)+(buflen>>1)] = buffer[(buflen>>1)-i-1];
		p2[i] = buffer[(buflen>>1)+(shift>>1)-i-1];
		p2[i+(shift>>1)+(buflen>>1)] = buffer[buflen-i-2];
	}
	for (i=0; i<(buflen>>1); i++){
		p1[i+(shift>>1)] = buffer[i];
		p2[i+(shift>>1)] = buffer[i+(buflen>>1)];
	}
/*   對稱延拓結束  */

 
	
	/*  s1  */
	for (i=0; i<(buflen>>1)+shift-1; i++){
	   data_tmp = *(p2+i+1) + *(p2+i/*+1-1*/);
       *(s+i+1) = *(p1+i+1) - ((data_tmp >>2) + (data_tmp >>3) + (data_tmp >>4) + (data_tmp >>8) + (data_tmp >>9) + (data_tmp >>12));
	}
   /*  d1 */
	for (i=0; i<(buflen>>1)+shift-1; i++){
		data_tmp = *(s+i) + *(s+i+1);
        *(d+i) = *(p2+i) - ((data_tmp>>1) + (data_tmp >>2) + (data_tmp >>3) + (data_tmp >>7));
	}

    /* p1 = s0 */
	for (i=0; i<(buflen>>1)+shift-1; i++){
		data_tmp = -(*(d+i+1)+*(d+i/*+1-1*/));
		*(p1+i+1) = *(s+i+1) - ((data_tmp >>5) +(data_tmp >>6) + (data_tmp >>8) + (data_tmp >>9) + (data_tmp >>12));
	}

   /*  p2 = d0 */
	for (i=0; i<(buflen>>1)+shift-1; i++){
		data_tmp = - (*(p1+i) + *(p1+i+1));
		*(p2+i) = *(d+i) - (data_tmp + (data_tmp>>1) + (data_tmp>>4) + (data_tmp >>6) + (data_tmp >>7));
	}

   for (i=0; i<(buflen>>1); i++){
	   *(buffer+2*i) = *(p1+i+(shift>>1)) ;
	   *(buffer+2*i+1) = *(p2+i+(shift>>1));
   }


	delete []p1;
	delete []p2;
	delete []s;
	delete []d;
}