www.pudn.com > xiaobo.zip.zip > coder.c


 
/** 
 
<> todo : if (!coder->encodeBand) 
	provide a generic stub which calls encodeBandBP several times 
	(same for BandZT) 
	this is a pain in the ass, and perhaps obsolete; the calls to encodeBand 
		should probably all go through encodeImage, which is the master. 
 
***/ 
 
#define CHECK_EOF 
 
#include  
#include  
#include  
#include  
#include  
#include  
 
#include "coder.h" 
#include "image.h" 
#include "subbands.h" 
 
int tune_param = 0; 
 
extern coder coderNone, coderOzero, coderOone, coderSigMap,coderSM2,coderSMo0 , 
	coderBitPlane,coderBP,coderBPsorted,coderBPbin,coderBPB2,coderBPBF, coderVQ, coderFrac, coderO1_1,coderO2, coderO1_SB, 
	coderSH_O1,coderSH_O1SB,coderZF,coderBPZT,coderBPZT2,coderBPBFZT,coderNOP; 
 
const num_coders = 24;	// (sizeof(coder_list)/sizeofpointer)-1 
const coder * coder_list[] = { 
		&coderBP, 
		&coderNOP, 
		&coderBPsorted, 
		&coderBitPlane, 
		&coderBPZT, 
		&coderBPZT2, 
		&coderBPBF, 
		&coderBPBFZT, 
		&coderBPbin, 
		&coderBPB2, 
		&coderZF, 
		&coderSigMap, 
		&coderSM2, 
		&coderSMo0, 
		&coderO1_1, 
		&coderO1_SB, 
		&coderO2, 
		&coderOone, 
		&coderOzero, 
		&coderVQ, 
		&coderFrac, 
		&coderNone, 
		&coderSH_O1, 
		&coderSH_O1SB, 
		NULL 
}; 
 
#define SAFE_PAD (1<<16) 
 
extern wavelet * newWavelet(const image *im,int levels) 
{ 
wavelet *w; 
	if ( (w = new(wavelet)) == NULL ) return NULL; 
	w->width = im->width; 
	w->height = im->height; 
	w->planes = im->planes; 
	w->levels = levels; 
	w->complen = im->tot_bytes; 
	if ( (w->comp = malloc(w->complen)) == NULL ) { 
		freeWavelet(w); return NULL; 
	} 
 
	w->stopline = -1; 
	w->im = (image *)im; 
 
return w; 
} 
extern void freeWavelet(wavelet *w) 
{ 
	if ( w ) { 
		freeSubbands(w->subband_root); 
		smartfree(w->qi); 
		smartfree(w->comp); 
		free(w); 
	} 
} 
 
 
coder * coder_create_read(wavelet *w) 
{ 
coder *ret; 
 
	if ( (ret = new(coder)) == NULL ) return NULL; 
	if ( w->coder_template ) memcpy(ret,w->coder_template,sizeof(coder)); 
 
	ret->w = w; 
 
	w->stoplen = min(w->complen,w->stoplen); 
	 
	if ( (ret->arith = arithInit()) == NULL ){ 
		coder_destroy(ret); return NULL; 
	} 
	arithDecodeInit(ret->arith,ret->w->comp); 
 
return ret; 
} 
 
coder * coder_create_write(const coder *template,wavelet *w,int stoplen) 
{ 
coder *ret; 
 
	if ( (ret = new(coder)) == NULL ) return NULL; 
	if ( template ) memcpy(ret,template,sizeof(coder)); 
 
	w->coder_template = template; 
  
	ret->w = w; 
 
	w->stoplen = stoplen; 
 
	if ( (ret->arith = arithInit()) == NULL ){ 
		coder_destroy(ret); return NULL; 
	} 
	arithEncodeInit(ret->arith,ret->w->comp); 
 
return ret; 
} 
 
void coder_flush_read(coder *c) 
{ 
#ifdef CHECK_EOF 
int got; 
got = arithGet(c->arith,79); 
if ( got == 43 ) arithDecode(c->arith,43,44,79); 
else errputs("warning : didn't get EOF"); 
#endif // EOF 
 
arithDecodeDone(c->arith); 
} 
void coder_flush_write(coder *c) 
{ 
#ifdef CHECK_EOF 
arithEncode(c->arith,43,44,79); 
#endif //EOF 
 
c->w->complen = arithEncodeDone(c->arith); 
} 
 
void coder_destroy(coder *c) 
{ 
	if ( c ) { 
		if ( c->arith ) arithFree(c->arith); 
		free(c); 
	} 
} 
 
/************* routines for DPCM coding the top plane  
 
as of v1.5 this coder is quite bad-ass , even beats things like CALIC on your teeny 
images (because context coders take too long to rev up) : 
 
we do something cheezy like all the standard image coders.  Since the LL band 
is so small, we don't have time to adapt a context coder.  Also, the different LL 
bands of images have widely varying statistics, so we can't just use the plain 
order -1 coder. 
 
instead we try to make a local estimate of the mean prediction error (MPE).  We then 
code using a static coder : 
 
	0  		+ (1 < MPE) 
	01 		+ (MPE < 3*MPE ) 
	001 	+ (3*MPE < 7*MPE ) 
	0001	+ (7*MPE < 15 * MPE) 
	... 
 
it would make more sense to use a Guassian probability distribution 
sent to the arithcoder, with MPE as the standard dev, but the cumprobs 
for a gaussian are the error function.  you would think we could 
handle that by tabulation, but there are intricacies with coding 
(sym/mpe) as a real number vs. an integer. (eg. hard to define the 
"next" and "prev" symbols to get the neighboring cumprobs) 
 
------ 
 
our sign coder is even cheezier.  We build the N+W sign neighbors context. 
if N == W we have confidence, and use a binary adaptive coder on (sign == neighbor) 
else we just send sign raw 
 
***************/ 
 
#define EST_ERROR_SHIFT			3 
#define EST_ERROR_PAD			0 
#define EST_ERROR(grad,prev)	(((grad + grad + prev)>>EST_ERROR_SHIFT) + 2 + EST_ERROR_PAD) 
#define ERROR_CODER_MULT		2		// tunes to 1 !!! almost an 0.2 b gain (on the LL) 
 
#define INC 5 
 
#define bitModel(bit,P0,PT)		do { PT += INC; if (!(bit)) P0 += INC;  if ( PT > 1000 ) { PT >>= 1; P0 >>= 1; P0++; PT += 2; } } while(0) 
#define bitEnc(bit,ari,P0,PT)	do { arithEncBit(ari,P0,PT,bit);	bitModel(bit,P0,PT); } while(0) 
#define bitDec(bit,ari,P0,PT)	do { bit = arithDecBit(ari,P0,PT);	bitModel(bit,P0,PT); } while(0) 
 
static int signs_p0,signs_pt; 
 
void coder_encodesign(arithInfo *ari,bool sign,bool W,bool N) 
{ 
	if ( N == W ) { 
		if ( ! N ) sign = !sign; 
		bitEnc(sign,ari,signs_p0,signs_pt); 
	} else { 
		arithEncBitRaw(ari,sign); 
	} 
} 
bool coder_decodesign(arithInfo *ari,bool W,bool N) 
{ 
	if ( N == W ) { 
		bool sign; 
		bitDec(sign,ari,signs_p0,signs_pt); 
		return sign ? N : !N; 
	} else { 
		return arithDecBitRaw(ari); 
	} 
} 
 
void coder_encodeDPCM(coder *c,int *plane,int width,int height,int rowpad) 
{ 
int x,y,pred,grad,val,sign,prev_err,est,fullw; 
arithInfo * ari = c->arith; 
int *ptr,*pline; 
 
	fullw = width + rowpad; 
	signs_p0 = 10; signs_pt = 20; 
	ptr = plane; prev_err = 99; 
	for(y=0;y>1; grad = abs(ptr[-1] - ptr[-2]); } 
			} else if ( x == 0 ) { pred = (pline[0] + pline[1])>>1; grad = abs(pline[0] - pline[1]);  
			} else if ( x == width-1 ) { 
				pred = (ptr[-1] + pline[0])>>1; 
				grad = abs(ptr[-1] - pline[0]); 
			} else { 
				pred = (ptr[-1]*3 + pline[0]*3 + ptr[-1] + pline[1])>>3; 
				grad = max( abs(ptr[-1] - ptr[-1]) , abs( pline[0] - pline[1]) ); 
			} 
 
			val = (*ptr) - pred; 
 
			if ( val < 0 ) { sign = 1; val = -val; } 
			else sign = 0; 
	 
			est = EST_ERROR(grad,prev_err); 
			cu_putMulting_ari(val,ari,est,ERROR_CODER_MULT); 
 
			if ( val > 0 ) { 
				if ( x == 0 || y == 0 ) { 
					arithEncBitRaw(ari,sign); 
				} else { 
					coder_encodesign(ari,sign,isneg(ptr[-1]),isneg(pline[0])); 
				} 
			} 
 
			ptr++; 
			prev_err = val; 
		} 
		ptr += rowpad; 
	} 
} 
 
void coder_decodeDPCM(coder *c,int *plane,int width,int height,int rowpad) 
{ 
int x,y,pred,grad,val,prev_err,est,fullw; 
arithInfo * ari = c->arith; 
int *ptr,*pline; 
 
	fullw = width + rowpad; 
	signs_p0 = 10; signs_pt = 20; 
	ptr = plane; prev_err = 99; 
	for(y=0;y>1; grad = abs(ptr[-1] - ptr[-2]); } 
			} else if ( x == 0 ) { pred = (pline[0] + pline[1])>>1; grad = abs(pline[0] - pline[1]);  
			} else if ( x == width-1 ) { 
				pred = (ptr[-1] + pline[0])>>1; 
				grad = abs(ptr[-1] - pline[0]); 
			} else { 
				pred = (ptr[-1]*3 + pline[0]*3 + ptr[-1] + pline[1])>>3; 
				grad = max( abs(ptr[-1] - ptr[-1]) , abs( pline[0] - pline[1]) ); 
			} 
 
			est = EST_ERROR(grad,prev_err); 
			val = cu_getMulting_ari(ari,est,ERROR_CODER_MULT); 
 
			prev_err = val; 
 
			if ( val > 0 ) { 
				if ( x == 0 || y == 0 ) { 
					if ( arithDecBitRaw(ari) ) val = -val; 
				} else { 
					if ( coder_decodesign(ari,isneg(ptr[-1]),isneg(pline[0])) ) 
						val = -val; 
				} 
			} 
 
			*ptr++ = val + pred; 
		} 
		ptr += rowpad; 
	} 
} 
 
/***** 
 
the LL DPCM no longer uses the _m1 routines , but some of 
the coders use them to do order (-1) coding. 
 
******/ 
 
#define M1_THRESH_1 3 
#define M1_THRESH_2 8 
#define M1_THRESH_3 30 
#define M1_THRESH_4 128 
 
void encode_m1(arithInfo * ari,int sym) 
{ 
	if ( sym < M1_THRESH_1 ) { 
		arithBit(ari,0); 
		arithEncode(ari,sym,sym+1,M1_THRESH_1); 
	} else { 
		arithBit(ari,1); sym -= M1_THRESH_1; 
		if ( sym < M1_THRESH_2 ) { 
			arithBit(ari,0); 
			arithEncode(ari,sym,sym+1,M1_THRESH_2); 
		} else { 
			arithBit(ari,1); sym -= M1_THRESH_2; 
			if ( sym < M1_THRESH_3 ) { 
				arithBit(ari,0); 
				arithEncode(ari,sym,sym+1,M1_THRESH_3); 
			} else { 
				int escape = M1_THRESH_4; 
 
				arithBit(ari,1); sym -= M1_THRESH_3; 
 
				while( sym >= escape ) { 
					arithEncode(ari,escape,escape+1,escape+1); 
					sym -= escape;	escape += escape; 
					if ( escape > ari->safeProbMax ) escape = ari->safeProbMax; 
				} 
				arithEncode(ari,sym,sym+1,escape+1); 
			} 
		} 
	} 
} 
 
int decode_m1(arithInfo *ari) 
{ 
int sym,cur; 
sym=0; 
	if ( arithGetBit(ari) == 0 ) { 
		cur = arithGet(ari,M1_THRESH_1); 
		arithDecode(ari,cur,cur+1,M1_THRESH_1); 
		sym+=cur; 
	} else { 
		sym += M1_THRESH_1; 
		if ( arithGetBit(ari) == 0 ) { 
			cur = arithGet(ari,M1_THRESH_2); 
			arithDecode(ari,cur,cur+1,M1_THRESH_2); 
			sym+=cur; 
		} else { 
			sym += M1_THRESH_2; 
			if ( arithGetBit(ari) == 0 ) { 
				cur = arithGet(ari,M1_THRESH_3); 
				arithDecode(ari,cur,cur+1,M1_THRESH_3); 
				sym+=cur; 
			} else { 
				int escape = M1_THRESH_4>>1; 
 
				sym += M1_THRESH_3; 
				do { 
					escape += escape; 
					if ( escape > ari->safeProbMax ) escape = ari->safeProbMax; 
					cur = arithGet(ari,escape+1); 
					arithDecode(ari,cur,cur+1,escape+1); 
					sym += cur; 
				} while( cur == escape ); 
			} 
		} 
	} 
return sym; 
}