advancejpeg-lscode.cpp

改进的JPEG-LS算法

#include "stdafx.h"
#include <math.h>
#include "DoubleCode.h"
#include "AdvanceJPEG-LSCode.h"
//图像行编码器
CAdvanceJepgLSOneLineCode::CAdvanceJepgLSOneLineCode()
{
	lpImage=NULL;
	lpCodeStream=NULL;
	CodeBitcp=0;
	OneSampleBits=0;
	qNEAR=0;//误差界
	lpPrecedingLine=lpCurrentLine=lpThirdLine=lpBuf1=lpBuf2=lpBuf3=NULL;//存放编码上下文象素数据
	lpReconstructPrecedingLine=lpReconstructCurrentLine=lpBuf4=lpBuf5=NULL;//存放重建象素数据
	//用于大误差编码
	lpDoubleErrList=NULL;
	DoubleErrN=0;
	DErrN=0;
	DDErrN=0;
	lpDErrNumberFastCoding=NULL;
	lpDErrRunFastCoding=NULL;
	lpDDErrNumberFastCoding=NULL;
	lpDDErrRunFastCoding=NULL;

	AdjustWinLength=BakAdjustWinLength=6;//调节窗口长度。
	AdjustNumber=1;
	AllTestNumber=2;
	lpTwoAdjustNumberFastCoding=NULL;
	lpTwoAdjustRunFastCoding=NULL;
	lpAdjustNumberFastCoding=NULL;
	lpAdjustRunFastCoding=NULL;
	//用于多波段编码
	AllowBandPrediction=FALSE;//波段预测允许开关
}
CAdvanceJepgLSOneLineCode::~CAdvanceJepgLSOneLineCode()
{
}
CString CAdvanceJepgLSOneLineCode::GetString()
{
	CString str,s;
	str.Empty();
	s.Format(_T("\r\n%-4d"),LineN);str+=s;
	s.Format(_T(" %-8d"),CodeBitcp);str+=s;
	s.Format(_T(" %08x"),AndCode);str+=s;
	s.Format(_T(" %-4d"),RUNindex);str+=s;
	s.Format(_T(" %-3d"),DoubleMaxError);str+=s;
	s.Format(_T(" %-3d"),RunContext);str+=s;
	s.Format(_T(" %1d"),RunTestThresholdStep);str+=s;
	s.Format(_T(" %-2d"),MaxErrAdjust1);str+=s;
	s.Format(_T(" %-2d"),MaxErrAdjust2);str+=s;
	s.Format(_T(" %-3d"),DoubleErrN);str+=s;
	s.Format(_T(" %-3d"),TestErrorMeLimen);str+=s;
	s.Format(_T(" %-3d"),AdjustErrorMe);str+=s;
	s.Format(_T(" %-3d"),BakAdjustWinLength);str+=s;
	s.Format(_T(" %-3d"),AdjustWinLength);str+=s;
	s.Format(_T(" %-3d"),AdjustNumber);str+=s;
	s.Format(_T(" %-3d"),AllTestNumber);str+=s;
	s.Format(_T(" %-8.5f"),((double)A_Mul)/UNIFORM_FLOAT_TO_INT_MUL);str+=s;
	s.Format(_T(" %-5d"),B_Add);str+=s;
	s.Format(_T(" %-8.5f"),((double)p)/UNIFORM_FLOAT_TO_INT_MUL);str+=s;
	s.Format(_T(" %-8.4f"),((double)PredictErr_D)/UNIFORM_FLOAT_TO_INT_MUL);str+=s;
	return str;
}
void CAdvanceJepgLSOneLineCode::PrepareCode(int LineL,CAdvanceJepgLSOneLineCode *lpPredictionBand)//为行编码器准备缓冲区
{
	X=LineL;//行长度
	lpImage=new int[X];
	lpCurrentLine=lpBuf1=new int[X];
	ZeroMemory(lpCurrentLine,X*sizeof(int));
	lpPrecedingLine=lpBuf2=new int[X];
	ZeroMemory(lpPrecedingLine,X*sizeof(int));
	lpThirdLine=lpBuf3=new int[X];
	ZeroMemory(lpThirdLine,X*sizeof(int));
	if((lpPredictionBandCode=lpPredictionBand)!=NULL)
	{//允许波段预测
		AllowBandPrediction=TRUE;
		lpReconstructCurrentLine=lpBuf4=new int[X];
		lpReconstructPrecedingLine=lpBuf5=new int[X];
	}
	else
	{//不作波段预测时
		AllowBandPrediction=FALSE;
		lpReconstructPrecedingLine=lpPrecedingLine;
		lpReconstructCurrentLine=lpCurrentLine;
	}
	lpDoubleErrList=new AdvanceJPEGLSMaxErrList[X];

	lpDErrNumberFastCoding=new CDoubleErrorNumberFastCoding;
	lpDErrRunFastCoding=new CDoubleErrorRunFastCoding;
	lpDDErrNumberFastCoding=new CDoubleErrorNumberFastCoding;
	lpDDErrRunFastCoding=new CDoubleErrorRunFastCoding;

	lpTwoAdjustRunFastCoding=new CDoubleErrorRunFastCoding;
	lpTwoAdjustNumberFastCoding=new CDoubleErrorNumberFastCoding;
	lpAdjustRunFastCoding=new CDoubleErrorRunFastCoding;
	lpAdjustNumberFastCoding=new CDoubleErrorNumberFastCoding;
}
void CAdvanceJepgLSOneLineCode::EndCode()//主要是清除申请的缓冲区,要在编码数据处理完后方可结束
{
	if(lpImage!=NULL)delete lpImage;lpImage=NULL;
	if(lpBuf1!=NULL)delete lpBuf1;lpBuf1=NULL;
	if(lpBuf2!=NULL)delete lpBuf2;lpBuf2=NULL;
	if(lpBuf3!=NULL)delete lpBuf3;lpBuf3=NULL;
	if(lpPredictionBandCode!=NULL)
	{
		if(lpBuf4!=NULL)delete lpBuf4;
		lpBuf4=NULL;
		if(lpBuf5!=NULL)delete lpBuf5;
		lpBuf5=NULL;
	}
	if(lpDoubleErrList!=NULL)delete lpDoubleErrList;lpDoubleErrList=NULL;

	if(lpDErrNumberFastCoding!=NULL)delete lpDErrNumberFastCoding;lpDErrNumberFastCoding=NULL;
	if(lpDErrRunFastCoding!=NULL)delete lpDErrRunFastCoding;lpDErrRunFastCoding=NULL;
	if(lpDDErrNumberFastCoding!=NULL)delete lpDDErrNumberFastCoding;lpDDErrNumberFastCoding=NULL;
	if(lpDDErrRunFastCoding!=NULL)delete lpDDErrRunFastCoding;lpDDErrRunFastCoding=NULL;
	if(lpTwoAdjustNumberFastCoding!=NULL)delete lpTwoAdjustNumberFastCoding;lpTwoAdjustNumberFastCoding=NULL;
	if(lpTwoAdjustRunFastCoding!=NULL)delete lpTwoAdjustRunFastCoding;lpTwoAdjustRunFastCoding=NULL;
	if(lpAdjustNumberFastCoding!=NULL)delete lpAdjustNumberFastCoding;lpAdjustNumberFastCoding=NULL;
	if(lpAdjustRunFastCoding!=NULL)delete lpAdjustRunFastCoding;lpAdjustRunFastCoding=NULL;

	lpImage=NULL;
	lpCodeStream=NULL;
	CodeBitcp=0;
	OneSampleBits=0;
}
void CAdvanceJepgLSOneLineCode::InitBandPrediction()
{
	int i;
	OpenOffBandPrediction=TRUE;
	Last_A_Mul=A_Mul=UNIFORM_FLOAT_TO_INT_MUL;
	Last_B_Add=B_Add=0;
	p=UNIFORM_FLOAT_TO_INT_MUL;
	D_Cun=D_Pre=E_Cun=E_Pre=R_CunPre=0;
	StatisticLinN=0;
	for(i=0;i<X;i++)
	{
		lpCurrentLine[i]=lpPrecedingLine[i]=lpThirdLine[i]=Pixel_ADD;
	}
	AndCode<<=2;AndCode|=0x3;
	MAXVAL=AndCode;
	RANGE=(MAXVAL+2*qNEAR)/(2*qNEAR+1)+1;
	for(qbpp=15;(RANGE&(1<<qbpp))==0;qbpp--);if(RANGE!=(1<<qbpp))qbpp++;
	for(bpp=15;((MAXVAL+1)&(1<<bpp))==0;bpp--);if((MAXVAL+1)!=(1<<bpp))bpp++;if(bpp<2)bpp=2;
	LIMIT=2*(bpp+((bpp>8)?bpp:8));//为规则方式中的GOLOMB编码用的最大码长限度
}
void CAdvanceJepgLSOneLineCode::OpenBandPrediction()
{
	OpenOffBandPrediction=TRUE;//打开波段预测
	int i,pp,ppx;
	//准备预测误差的上下文
	for(i=0;i<X;i++)
	{//将上下文变成预测误差上下文.
		pp=lpPredictionBandCode->lpReconstructCurrentLine[i];//预测波段重建值
		if((ppx=(A_Mul*pp/UNIFORM_FLOAT_TO_INT_MUL+B_Add))>Pixel_ADD)ppx=pp;//预测值
		lpCurrentLine[i]=lpCurrentLine[i]-ppx+PIX_ADD;//当前上下文由重建值变成预测误差值
		pp=lpPredictionBandCode->lpReconstructPrecedingLine[i];//前行处理(同上)
		if((ppx=(Last_A_Mul*pp/UNIFORM_FLOAT_TO_INT_MUL+Last_B_Add))>PIX_ADD)ppx=pp;
		lpPrecedingLine[i]=lpPrecedingLine[i]-ppx+PIX_ADD;
	}
	AndCode<<=2;AndCode|=0x3;
	MAXVAL=AndCode;
	RANGE=(MAXVAL+2*qNEAR)/(2*qNEAR+1)+1;
	for(qbpp=15;(RANGE&(1<<qbpp))==0;qbpp--);if(RANGE!=(1<<qbpp))qbpp++;
	for(bpp=15;((MAXVAL+1)&(1<<bpp))==0;bpp--);if((MAXVAL+1)!=(1<<bpp))bpp++;if(bpp<2)bpp=2;
	LIMIT=2*(bpp+((bpp>8)?bpp:8));//为规则方式中的GOLOMB编码用的最大码长限度
}
void CAdvanceJepgLSOneLineCode::CloseBandPrediction()
{
	int i,pp,ppx;
	OpenOffBandPrediction=FALSE;//关闭波段预测
	//准备预测误差的上下文
	AndCode>>=2;//预测误差的范围要扩大两比特
	MAXVAL=AndCode;
	RANGE=(MAXVAL+2*qNEAR)/(2*qNEAR+1)+1;
	for(qbpp=15;(RANGE&(1<<qbpp))==0;qbpp--);if(RANGE!=(1<<qbpp))qbpp++;
	for(bpp=15;((MAXVAL+1)&(1<<bpp))==0;bpp--);if((MAXVAL+1)!=(1<<bpp))bpp++;if(bpp<2)bpp=2;
	LIMIT=2*(bpp+((bpp>8)?bpp:8));//为规则方式中的GOLOMB编码用的最大码长限度
	for(i=0;i<X;i++)
	{//将当前预测误差上下文变成重建图像上下文.
		pp=lpPredictionBandCode->lpReconstructCurrentLine[i];
		if((ppx=(A_Mul*pp/UNIFORM_FLOAT_TO_INT_MUL+B_Add))>PIX_ADD)ppx=pp;//预测值
		lpCurrentLine[i]=lpCurrentLine[i]+ppx-PIX_ADD;//还原预测误差为重建成值
		if(lpCurrentLine[i]<0)lpCurrentLine[i]=0;
		else if(lpCurrentLine[i]>MAXVAL)lpCurrentLine[i]=MAXVAL;
		pp=lpPredictionBandCode->lpReconstructPrecedingLine[i];
		if((ppx=(Last_A_Mul*pp/UNIFORM_FLOAT_TO_INT_MUL+Last_B_Add))>PIX_ADD)ppx=pp;
		lpPrecedingLine[i]=lpPrecedingLine[i]+ppx-PIX_ADD;
		if(lpPrecedingLine[i]<0)lpPrecedingLine[i]=0;
		else if(lpPrecedingLine[i]>MAXVAL)lpPrecedingLine[i]=MAXVAL;
	}
}
void CAdvanceJepgLSOneLineCode::UpdateBandPrediction()//更新波段预测参数，更新变量A_Mul,B_Add,p。
{
	//时序要求: 参与统计的重建数据必须是当前象素之前的.
	int i,pp,pc,num,dp;
	__int64 Rxy,Sx,Mx,My,Sy;
	int aMul,bAdd;
	Rxy=Sx=Mx=My=Sy=0;
	for(i=0,num=0;i<X;i++)
	{
		pc=lpReconstructCurrentLine[i];//当前波段当前行象素重建值
		pp=lpPredictionBandCode->lpReconstructCurrentLine[i];//预测波段当前行象素重建值
		if((dp=pc-pp)<0)dp=-dp;//取得误差绝对值
		if(dp<128||(((pp+pc)*4)/dp)>6/*和与差之比大于3/2时*/)
		{//当两波段因特殊象素造成过大的误异时不作为统计量
			num++;//统计量加1
			Rxy+=pc*pp;//相关统计量
			My+=pc;//当前波段和
			Sy+=pc*pc;//当前波段平方和
			Sx+=pp*pp;//预测波段平方和
			Mx+=pp;//预测波段和
		}
	}
	if(num<=0)return;
	//预测误差方差的平方
	if(LineN<=1)PredictErr_D=(int)((Sy+Sx-Rxy-Rxy-(My-Mx)*(My-Mx)/num)/num);
	else PredictErr_D=(int)((PredictErr_D+(Sy+Sx-Rxy-Rxy-(My-Mx)*(My-Mx)/num)/num+1)/2);
	//累计统计量
	E_Cun+=My;
	E_Pre+=Mx;
	D_Cun+=Sy;
	D_Pre+=Sx;
	R_CunPre+=Rxy;
	StatisticLinN+=num;//相关系数统计的行计数
	if(StatisticLinN>=(UPDATE_RELATIVE_LINE_NUMBER*X-X/2)/*必须大于等于2*/)
	{//用即时相关系数控制是否预测
		//计算相关系数
		D_Cun=(D_Cun-(E_Cun*E_Cun/StatisticLinN))/StatisticLinN;//当前波段方差平方
		D_Pre=(D_Pre-(E_Pre*E_Pre/StatisticLinN))/StatisticLinN;//预测波段方差平方
		R_CunPre=(R_CunPre-E_Cun*E_Pre/StatisticLinN)/StatisticLinN;
		if((D_Cun*D_Pre)==0)p=UNIFORM_FLOAT_TO_INT_MUL;
		else p=(int)(R_CunPre*R_CunPre*UNIFORM_FLOAT_TO_INT_MUL/(D_Cun*D_Pre));//整数化相关系数平方
			
		//清0累计统量
		D_Cun=D_Pre=E_Cun=E_Pre=R_CunPre=0;
		StatisticLinN=0;
			
		//判断是否预测
		if(p>ALLOW_PREDICTION_MIN_RELATIVE&&OpenOffBandPrediction==FALSE)
		{//相关系数大于最小允许预测的相关系数值，同时当前处在关闭波段预测的状态时
			OpenBandPrediction();//打开波段预测允许
		}
		else if(p<NOT_PREDICTION_MAX_RELATIVE&&OpenOffBandPrediction==TRUE)
		{//相关系数小于最大不允许预测的相关系数值，同时当前处在打开波段预测的状态时
			CloseBandPrediction();//关闭波段预测允许
		}
	}
	if(OpenOffBandPrediction==TRUE)
	{//允许波段预测时，用预测误差的方差控制平坦区象素调节参数。
		if(PredictErr_D<(qNEAR*qNEAR)/2)//相当于预测误差的方差处在误差界2/3时
		{
			AdjustWinLength=BakAdjustWinLength+BakAdjustWinLength/2;//平坦区调节窗口加长
			if(LineN>2)
			{//将游程检测及游程编码入口门限设成最小值。
				RunTestThresholdStep=0;
				RunContext=SetRunContext[RunTestThresholdStep];
				MaxErrAdjust1=SetMaxErrAdjust1[RunTestThresholdStep];
				MaxErrAdjust2=SetMaxErrAdjust2[RunTestThresholdStep];
				DoubleMaxError=SetDoubleMaxError[RunTestThresholdStep];
				if(TestErrorMeLimen==0)
				{
					RunContext=DoubleMaxError=qNEAR;
				}
			}
		}
		else 
		{
			AdjustWinLength=BakAdjustWinLength;
			if(LineN>2)
			{//将游程检测及游程编码入口门限设成中间值。
				RunTestThresholdStep=1;
				RunContext=SetRunContext[RunTestThresholdStep];
				MaxErrAdjust1=SetMaxErrAdjust1[RunTestThresholdStep];
				MaxErrAdjust2=SetMaxErrAdjust2[RunTestThresholdStep];
				DoubleMaxError=SetDoubleMaxError[RunTestThresholdStep];
				if(TestErrorMeLimen==0)
				{
					RunContext=DoubleMaxError=qNEAR;
				}
			}
		}
	}
	else AdjustWinLength=BakAdjustWinLength;
	//更新预测系数
	Last_A_Mul=A_Mul;
	Last_B_Add=B_Add;
	if((Sx-Mx*Mx/num)==0)aMul=MAX_A_MUL;
	else aMul=(int)((Rxy-My*Mx/num)*UNIFORM_FLOAT_TO_INT_MUL/(Sx-Mx*Mx/num));
	bAdd=(int)((My-aMul*Mx/UNIFORM_FLOAT_TO_INT_MUL)/num);
	A_Mul=(A_Mul+aMul+1)/2;//新老系数取平均
	if(A_Mul>MAX_A_MUL)A_Mul=MAX_A_MUL;
	if(A_Mul<MIN_A_MUL)A_Mul=MIN_A_MUL;
	B_Add=(B_Add+bAdd+1)/2;//新老系数取平均
}
void CAdvanceJepgLSOneLineCode::InitialisationsCode(int OnePixelBit,int AllowMaxError)
{//初始化编码器,在编码器缓冲区准备如后
	int FACTOR;
	LineN=0;//行计数
	OneSampleBits=OnePixelBit;
	Pixel_ADD=1<<(OnePixelBit+1);
	AndCode=(1<<OnePixelBit)-1;
	MAXVAL=AndCode;
	qNEAR=AllowMaxError;
	if(MAXVAL>=128)
	{
		FACTOR=(((MAXVAL>4095)?4095:MAXVAL) + 128)/256;

		T1 = FACTOR * (BASIC_T1 - 2) + 2 + 3*qNEAR;
		if(T1>MAXVAL||T1<(qNEAR+1))T1=qNEAR+1;//CLAMP_1()

		T2 = FACTOR * (BASIC_T2 - 3) + 3 + 5*qNEAR;
		if(T2>MAXVAL||T2<T1)T2=T1;

		T3 = FACTOR * (BASIC_T3 - 4) + 4 + 7*qNEAR;
		if(T3>MAXVAL||T3<T2)T3=T2;
	}
	else
	{
		FACTOR = 256/(MAXVAL + 1);

		T1=BASIC_T1/FACTOR+3*qNEAR;
		if(T1<2)T1=2;
		if(T1>MAXVAL||T1<(qNEAR+1))T1=qNEAR+1;//CLAMP_1()

		T2=BASIC_T2/FACTOR+5*qNEAR;
		if(T2<3)T1=3;
		if(T2>MAXVAL||T2<T1)T2=T1;

		T3=BASIC_T3/FACTOR+7*qNEAR;
		if(T3<4)T1=4;
		if(T3>MAXVAL||T3<T2)T3=T2;
	}
/*1.	Compute the parameter RANGE，Compute the parameters 
	qbpp = [log RANGE], bpp = max(2, [log(MAXVAL+1)]), and LIMIT = 2 * (bpp + max(8,bpp))

qbpp 是待编码的误差动态范围所需的比特数表示。
bpp  原图像的比特数。
*/

	RANGE=(MAXVAL+2*qNEAR)/(2*qNEAR+1)+1;
	for(qbpp=15;(RANGE&(1<<qbpp))==0;qbpp--);if(RANGE!=(1<<qbpp))qbpp++;
	for(bpp=15;((MAXVAL+1)&(1<<bpp))==0;bpp--);if((MAXVAL+1)!=(1<<bpp))bpp++;if(bpp<2)bpp=2;
	LIMIT=2*(bpp+((bpp>8)?bpp:8));//为规则方式中的GOLOMB编码用的最大码长限度

	lpCurrentLine=lpBuf1;
	ZeroMemory(lpCurrentLine,X*sizeof(int));
	lpPrecedingLine=lpBuf2;
	ZeroMemory(lpPrecedingLine,X*sizeof(int));
	lpThirdLine=lpBuf3;
	ZeroMemory(lpThirdLine,X*sizeof(int));
	//附加编码器的复位
	lpDErrNumberFastCoding->Reset();
	lpDErrRunFastCoding->Reset();
	lpDDErrNumberFastCoding->Reset();
	lpDDErrRunFastCoding->Reset();
	lpTwoAdjustNumberFastCoding->Reset();
	lpTwoAdjustRunFastCoding->Reset();
	lpAdjustNumberFastCoding->Reset();
	lpAdjustRunFastCoding->Reset();

	//初始化游程检测及游程编码入口门限
	SetDoubleMaxError[0]=MIN_MDoubleMaxErr(qNEAR);
	SetRunContext[0]=MIN_MINCONTEXT(qNEAR);
	SetMaxErrAdjust1[0]=MIN_MDoubleMaxErrAdjust(qNEAR);
	SetMaxErrAdjust2[0]=MIN_MDoubleMaxErrAdjust2(qNEAR);

	SetDoubleMaxError[1]=MID_MDoubleMaxErr(qNEAR);
	SetRunContext[1]=MID_MINCONTEXT(qNEAR);
	SetMaxErrAdjust1[1]=MID_MDoubleMaxErrAdjust(qNEAR);
	SetMaxErrAdjust2[1]=MID_MDoubleMaxErrAdjust2(qNEAR);

	SetDoubleMaxError[2]=MAX_MDoubleMaxErr(qNEAR);
	SetRunContext[2]=MAX_MINCONTEXT(qNEAR);
	SetMaxErrAdjust1[2]=MAX_MDoubleMaxErrAdjust(qNEAR);
	SetMaxErrAdjust2[2]=MAX_MDoubleMaxErrAdjust2(qNEAR);

	RunTestThresholdStep=1;
	RunContext=SetRunContext[RunTestThresholdStep];
	MaxErrAdjust1=SetMaxErrAdjust1[RunTestThresholdStep];
	MaxErrAdjust2=SetMaxErrAdjust2[RunTestThresholdStep];
	DoubleMaxError=SetDoubleMaxError[RunTestThresholdStep];

	//初始化平坦区象素调节编码控制参数
	TestErrorMeLimen=MTestErrMeLimen(qNEAR);
	if(TestErrorMeLimen==0)
	{
		RunContext=qNEAR;
		DoubleMaxError=qNEAR;
	}
	else
	{
		RunContext=SetRunContext[RunTestThresholdStep];
		DoubleMaxError=SetDoubleMaxError[RunTestThresholdStep];
	}
	AdjustNumber=1000;
	AllTestNumber=20;
	AdjustWinLength=BakAdjustWinLength=DEFAULT_ADJUST_WIN_LENGTH;
	//初始化波段预测参数
	if(AllowBandPrediction==TRUE)
	{
		InitBandPrediction();
	}
	else OpenOffBandPrediction=FALSE;
	//初始口化上下文
	for(Q=0;Q<365;Q++)
	{
		B[Q]=C[Q]=0;
	}
	for(Q=0;Q<367;Q++)
	{
		A[Q]=(RANGE+32)/64;
		if(A[Q]<2)A[Q]=2;
		N[Q]=1;
	}
	//初始化游程编码码表
	RUNindex=0;
	J[0]=J[1]=J[2]=J[3] = 0;
	J[4]=J[5]=J[6]=J[7] = 1;
	J[8]=J[9]=J[10]=J[11]=2;
	J[12]=J[13]=J[14]=J[15]=3;
	J[16]=J[17]=4;
	J[18]=J[19]=5;