coder.c

傅立叶变换和小波变换是图像压缩的重要工具。该代大戏是利用小波变换进行图像压缩。

#include "coder.h"
#include <crblib/arithc.h>
#include <crblib/codeutil.h>
#include <crblib/rungae.h>
#include "mathutil.h"
#include "dpcm.h"
#include "wave.h"

// <> could share the rung contexts from the previous bands

//#define DO_STATISTICS

#define USE_FAR_NEIGHBORS

//#define USE_COUNT_CNTX // non-count context is much faster & (barely) better cmpsion

//#define CODE_ROWON_RAW

#define ROW_ISON_LIKELY		(7) // give it its own rung; helps alot
#define ROW_ISON_UNLIKELY	(7)

/***************
*
	In the very-high compression scenarios we're interested in,
	there may be huge runs of all-zero values, and we could get a
	lot of speed from run-coding or something

<> to try : code a run-length of zeros every time you code a zero
	even on the non-zerorow rows, 60% of vals are zero		
did it : the average run len is 30 ! and we still get no speed help
	over the whole-zero-row coder

row zeros improves time like this: (CDFxHaar)
	0.039886 secs = 182.6 cycles / pel = 1643066 pels /sec
	0.030420 secs = 139.3 cycles / pel = 2154337 pels /sec

--------------
Far neighbors are a very even compression/speed tradeoff as expected:

No far neighbors:

lena.256.raw         : 65536 -> 4128 = 0.50 bpp
DecodePlane : 168.3 cycles / pel = 1782633 pels /sec
mse = 33.178757, rmse = 5.760100, psnr = 32.956398

With far neighbors:

lena.256.raw         : 65536 -> 4150 = 0.51 bpp
DecodePlane : 181.3 cycles / pel = 1654727 pels /sec
mse = 32.442520, rmse = 5.695833, psnr = 33.053852

*
*****************/

void encodeBands(SubBand *LH,SubBand *HL,SubBand *HH,arithInfo *ari);
void decodeBands(SubBand *LH,SubBand *HL,SubBand *HH,arithInfo *ari);

#ifdef DO_STATISTICS
static int Count_Pels=0,Count_Zeros=0,Count_ZeroRows=0;
#endif

/*}{*** subband loop ***********/

void codeWavelet(Wavelet *w,bool forward)
{
int l,b;
arithInfo *ari;
#ifdef _DEBUG
int lastPos;
#endif

	ari = w->ari;

#ifdef DO_STATISTICS
	Count_Pels=Count_Zeros=Count_ZeroRows=0;
#endif

	#ifdef _DEBUG
	if ( forward )
	{
		lastPos = arithTellEncPos(ari);
	}
	#endif

	for(b=0;b<3;b++)
	{
	SubBandTree * sbt;
	int levels;
		sbt = w->trees + b;

		levels = w->levels;
		if ( b < 2 )
			levels--;

		dpcmCode(&(sbt->LL),levels,w->dpcmLoss,ari,forward);

		if ( forward )	
		{
			for(l=0;l<levels;l++)
			{
				encodeBands(&(sbt->bands[0][l]), &(sbt->bands[1][l]), &(sbt->bands[2][l]), ari);
			}
		}
		else
		{
			for(l=0;l<levels;l++)
			{
				decodeBands(&(sbt->bands[0][l]), &(sbt->bands[1][l]), &(sbt->bands[2][l]), ari);
			}
		}

	#if 0 //{
		if ( forward )
		{
		int pos;
			pos = arithTellEncPos(ari);
			printf("plane %d coded to %d bytes\n",b,pos - lastPos);
			lastPos = pos;
		} //}
	#endif
	}

#ifdef DO_STATISTICS
	printf("zero = %d of %d = %2.1f %%\n",Count_Zeros,Count_Pels,Count_Zeros*100.0f/Count_Pels);
	printf("rows = %d of %d = %2.1f %%\n",Count_ZeroRows,Count_Pels,Count_ZeroRows*100.0f/Count_Pels);
#endif

}

/*}{***************** the functions ******************************/

int INLINE misoncontext(		uint *row,uint *prow,uint *pprow,subBandType sbType);
int INLINE munarycontext(	uint *row,uint *prow,uint *pprow);
int INLINE msigncontext(	uint *row,uint *prow,uint *pprow,subBandType sbType);

void LapInfo_Reset(LapInfo *li)
{
int x;
	for(x=0;x<NUM_CONTEXTS;x++)
		rungModelInit(&(li->isOn[x]));
	for(x=0;x<NUM_UNARY_CONTEXTS;x++)
		rungModelInit(&(li->unary[x]));
	for(x=0;x<NUM_SIGN_CONTEXTS;x++)
		rungModelInit(&(li->signs[x]));
	li->lastRowOn = true;
}

void encodeRow(uint *inrow,SubBand *sb,arithInfo *ari,LapInfo *li)
{
uint *row,*prow,*pprow,*rowbase;
int w,cnt;
subBandType sbType;
RowBuffer * buf;
bool rowOn;

	sbType = sb->type;
	w = sb->width;
	buf = sb->wavelet->CodeBuffer[ sbType - 1 ];

	row   = (uint *) RowBuffer_Row(buf, 0);
	prow  = (uint *) RowBuffer_Row(buf,-1);
	pprow = (uint *) RowBuffer_Row(buf,-2);
	rowbase = row;

	memcpy(row, inrow,4*w);

#ifdef DO_STATISTICS
	Count_Pels += w;
#endif
	
	{
	uint *ptr;
		// is row all zero?
		rowOn = false;
		ptr = row;
		while(w--)
		{
			if ( *ptr++ )
			{
				rowOn = true;
				break;
			}
		}

	#ifdef CODE_ROWON_RAW
		arithEncBitRaw(ari,rowOn);
	#else
		if ( li->lastRowOn )	rungModelEncBit(ari,rowOn,li->isOn + ROW_ISON_LIKELY);
		else					rungModelEncBit(ari,rowOn,li->isOn + ROW_ISON_UNLIKELY);
	#endif
	}

	li->lastRowOn = rowOn;

	if ( ! rowOn )
	{
	#ifdef DO_STATISTICS
		Count_ZeroRows += sb->width;
	#endif
		goto earlyDone;
	}

	w = sb->width;
	
	while(w--)
	{
		cnt = misoncontext(row,prow,pprow,sbType);

		assert( *row != 1);

		if ( ! *row )
		{
#ifdef DO_STATISTICS
	Count_Zeros ++;
#endif
			rungModelEncBit(ari,0,li->isOn + cnt);
		}
		else
		{
		int val,sign;

			rungModelEncBit(ari,1,li->isOn + cnt);

			val = *row;
			sign = val&1;
			val = (val>>1);
			assert(val);

			cnt = munarycontext(row,prow,pprow);

			rungEncodeUnary(ari,val-1,li->unary + cnt);

			cnt = msigncontext(row,prow,pprow,sbType);

			rungModelEncBit(ari,sign,li->signs + cnt);
		}

		row++; prow++; pprow++;
	}
	
earlyDone:

	row = rowbase + sb->width;

	// this row will be the prev for the next row ; mirror into the pad zones
	rowbase[-2] = rowbase[-1] = rowbase[0];
	row[1] = row[0] = row[-1];
}

uint * decodeRow(SubBand *sb,arithInfo *ari,LapInfo *li)
{
uint *row,*prow,*pprow,*rowbase;
int w;
subBandType sbType;
RowBuffer * buf;
rung_t *isOn,*unary,*signs;
bool rowOn;

	sbType = sb->type;
	w = sb->width;
	buf = sb->wavelet->CodeBuffer[ sbType - 1 ];

	row   = (uint *) RowBuffer_Row(buf, 0);
	prow  = (uint *) RowBuffer_Row(buf,-1);
	pprow = (uint *) RowBuffer_Row(buf,-2);
	rowbase = row;

	isOn  = li->isOn;
	unary = li->unary;
	signs = li->signs;	

	#ifdef CODE_ROWON_RAW
		rowOn = arithDecBitRaw(ari);
	#else
		if ( li->lastRowOn )	rowOn = rungModelDecBit(ari,li->isOn + ROW_ISON_LIKELY);
		else					rowOn = rungModelDecBit(ari,li->isOn + ROW_ISON_UNLIKELY);
		li->lastRowOn = rowOn;
	#endif

	if ( rowOn )
	{

		/******
		<> todos :

			make 3 different decodeRow funcs for each sbtype
				(not until all else is stable)

			inline rungModelDecBit ?

		*******/

		while(w--)
		{
			if ( rungModelDecBit(ari,isOn + misoncontext(row,prow,pprow,sbType)) )
			{
			int val;
			
				val = rungDecodeUnary(ari,unary + munarycontext(row,prow,pprow)) + 1;

				val <<= 1;

				val += rungModelDecBit(ari,signs + msigncontext(row,prow,pprow,sbType) );

				*row = val;
				row++; prow++; pprow++;
			}
			else
			{
				*row = 0;
				row++; prow++; pprow++;
			}
		}
	}
	else
	{
		memset(row,0,w<<2);
		row += w;
	}
	
	// this row will be the prev for the next row ; mirror into the pad zones
	rowbase[-2] = rowbase[-1] = rowbase[0];
	row[1] = row[0] = row[-1];

return rowbase;
}

void encodeBands(SubBand *LH,SubBand *HL,SubBand *HH,arithInfo *ari)
{
int y,s,h,i;
LapInfo liLH,liHL,liHH;
RowBuffer ** bufs;

	s = LH->stride;
	h = LH->height;
	bufs = LH->wavelet->CodeBuffer;
	
	// start up the codebuffer:
	for(i=0;i<3;i++)
		RowBuffer_Clear(bufs[i],LH->width);

	LapInfo_Reset(&liLH);
	LapInfo_Reset(&liHL);
	LapInfo_Reset(&liHH);

	i = 0;
	for(y=0;y<h;y++)
	{
		encodeRow( (uint *)(LH->band) + i,LH,ari,&liLH);
		encodeRow( (uint *)(HL->band) + i,HL,ari,&liHL);
		encodeRow( (uint *)(HH->band) + i,HH,ari,&liHH);

		bufs[0]->rowCenter ++;
		bufs[1]->rowCenter ++;
		bufs[2]->rowCenter ++;
		i += s;
	}
}

void decodeBands(SubBand *LH,SubBand *HL,SubBand *HH,arithInfo *ari)
{
int y,w,s,h,i;
LapInfo liLH,liHL,liHH;
RowBuffer ** bufs;

	s = LH->stride;
	w = LH->width;
	h = LH->height;
	bufs = LH->wavelet->CodeBuffer;
	
	// start up the codebuffer:
	// start up the codebuffer:
	for(i=0;i<3;i++)
		RowBuffer_Clear(bufs[i],LH->width);

	LapInfo_Reset(&liLH);
	LapInfo_Reset(&liHL);
	LapInfo_Reset(&liHH);

	i = 0;
	for(y=0;y<h;y++)
	{
	uint * dPtr[3];

		dPtr[0] = decodeRow(LH,ari,&liLH);
		dPtr[1] = decodeRow(HL,ari,&liHL);
		dPtr[2] = decodeRow(HH,ari,&liHH);

		memcpy((uint *)(LH->band) + i,dPtr[0],4*w);
		memcpy((uint *)(HL->band) + i,dPtr[1],4*w);
		memcpy((uint *)(HH->band) + i,dPtr[2],4*w);

		bufs[0]->rowCenter ++;
		bufs[1]->rowCenter ++;
		bufs[2]->rowCenter ++;
		i += s;
	}
}

/*}{***************** the context ******************************/

#define	N	prow[0]
#define	W	row[-1]
#define	NW	prow[1]
#define	NE	prow[-1]
#define	WW	row[-2]
#define	NN	pprow[0]

int INLINE msigncontext(uint *row,uint *prow,uint *pprow,subBandType sbType)
{
	switch(sbType)
	{
		case SBT_HH:
			return (NW &1) + (NE &1);
		case SBT_HL: // horizontal
			return (W&1) + (WW&1);
		case SBT_LH: // vertical
			return (N&1) + (NN&1);

	#ifdef _DEBUG
		default:
			assert(0);
			return 0;
	#endif

	}
}

int INLINE munarycontext(uint *row,uint *prow,uint *pprow)
{
uint count;

	count = N + W + ((NW + NE)>>2);

	// basically a log2 :

	if ( count > 10 )
	{
		return 4;
	}
	else
	{
	static const int utable[] = { 0, 1, 2, 2, 3, 3, 3, 3, 3, 3, 3 };
		return utable[count];
	}
}

int INLINE misoncontext(uint *row,uint *prow,uint *pprow,subBandType sbType)
{

	// <> most of our time is in here
	
	switch(sbType)
	{
		case SBT_HL:
		{
			switch( ( W ? 1 : 0 ) | ( WW ? 2 : 0 ) )
			{
				case 0:
					if ( NW || NE ) return 2;
					if ( N  || NN ) return 1;

				#ifdef USE_FAR_NEIGHBORS
					if(    pprow[-1] || pprow[-2] 
						|| pprow[+1] || pprow[+2] ) return 1;
				#endif

					return 0;

				case 1:
				case 2:
					if ( NW || NE ) return 5;
					if ( N  || NN ) return 4;
					return 3;

				case 3:
					return 6; // horizontal active

			}
		}

		case SBT_LH:
		{
			switch( ( N ? 1 : 0 ) | ( NN ? 2 : 0 ) )
			{
				case 0:
					if ( NW || NE ) return 2;
					if ( W  || WW ) return 1;
					return 0;

				case 1:
				case 2:
					if ( NW || NE ) return 5;
					if ( W  || WW ) return 4;
					return 3;

				case 3:
					return 6; // vertical active
			}
		}

		case SBT_HH:
		{
			switch( ( NE  ? 1 : 0 ) | ( NW ? 2 : 0 ) )
			{
				case 0:
					if ( N  || W ) return 3;			
					if ( NN || WW ) return 2;

				#ifdef USE_FAR_NEIGHBORS
					if(    pprow[-1] || pprow[-2] 
						|| pprow[+1] || pprow[+2] ) return 1;
				#endif

					return 0;

				case 1:
				case 2:
					if ( N || W ) return 5;
					return 4;

				case 3:
					return 6;	// diagonals active
			}
		}
	}
	
	#ifdef _DEBUG
	assert(0);
	return 0;
	#endif
}

/*}{********* EOF *******************************************/