//cerr << "\t\tdetermineGolombParameter: a = " << a << endl;	Assert(n);			// Make sure we don't get out of control	for (k=0;(n<<k) < a; ++k)	// Number of occurrences vs accumulated error magnitude		Assert(k<31);		// ... internal limit ... don't exceed width of Uint32//cerr << "\t\tdetermineGolombParameter: k = " << k << endl;	return k;}static BinaryInputStream &decodeMappedErrvalWithGolomb(Uint16 k,Uint16 glimit,Uint16 qbpp,Uint32 &value,BinaryInputStream &in){//cerr << "\t\tdecodeMappedErrvalWithGolomb: k = " << k << endl;//cerr << "\t\tdecodeMappedErrvalWithGolomb: glimit = " << glimit << endl;//cerr << "\t\tdecodeMappedErrvalWithGolomb: qbpp = " << qbpp << endl;	// Read unary representation of remaining most significant bits	Uint32 bit;	Uint32 unarycode=0;	while (readBit(in,bit) && !bit) ++unarycode;	// stops after bit is 1 (having read and discared trailing 1 bit)//cerr << endl;	Uint32 offset;	Uint16 bitstoread;	Assert(glimit > qbpp+1);	Uint16 limit=glimit-qbpp-1;//cerr << "\t\tdecodeMappedErrvalWithGolomb: unarycode = " << unarycode << endl;//cerr << "\t\tdecodeMappedErrvalWithGolomb: limit = " << limit << endl;	if (unarycode < limit) {		// use it to form most significant bits//cerr << "\t\tdecodeMappedErrvalWithGolomb: not limited, read " << unarycode << " zero bits (as value) followed by 1 then will read remaining " << k << " bits" << endl;		value=unarycode;		// will later get shifted into ms bits		bitstoread=k;		offset=0;	}	else {//cerr << "\t\tdecodeMappedErrvalWithGolomb: limited, read " << unarycode << " zero bits followed by 1 then will read remaining " << qbpp << " bits of value-1" << endl;		value=0;			// no contribution from unary code ... whole value is next		bitstoread=qbpp;		offset=1;	}	// Read least significant k bits	while (bitstoread-- && readBit(in,bit)) value=(value<<1) | bit;	// msb bit is read first//cerr << endl;	value+=offset;				// correct for limited case //cerr << "\t\tdecodeMappedErrvalWithGolomb: value = " << value << endl;	return in;}static BinaryOutputStream &encodeMappedErrvalWithGolomb(Uint16 k,Uint16 glimit,Uint16 qbpp,Uint32 value,BinaryOutputStream &out){//cerr << "\t\tencodeMappedErrvalWithGolomb: k = " << k << " value = " << value << endl;//cerr << "\t\tencodeMappedErrvalWithGolomb: glimit = " << glimit << endl;//cerr << "\t\tencodeMappedErrvalWithGolomb: qbpp = " << qbpp << endl;	// A.5.3 Mapped-error encoding	Uint32 unarycode=value>>k;					// Most significant bits go into unary code	Assert(glimit > qbpp+1);	Uint16 limit=glimit-qbpp-1;//cerr << "\t\tencodeMappedErrvalWithGolomb: unarycode = " << unarycode << endl;//cerr << "\t\tencodeMappedErrvalWithGolomb: limit = " << limit << endl;	if (unarycode < limit) {//cerr << "\t\tencodeMappedErrvalWithGolomb: not limited, writing " << unarycode << " zero bits followed by 1 then remaining " << k << " bits" << endl;		while (unarycode--) writeBit(out,0);			// Append unary representation of remaining most significant bits		writeBit(out,1);					// Flag the end of the unary code		Uint16 bits=k;						// Append least significant k bits		while (bits--) { writeBit(out,(value>>bits)&1); } 	// msb bit is written first & use the decremented bits as shift	}	else {//cerr << "\t\tencodeMappedErrvalWithGolomb: limited, writing " << limit << " zero bits followed by 1 then remaining " << qbpp << " bits of value-1" << endl;		while (limit--) writeBit(out,0);			// Append limit 0 bits		writeBit(out,1);					// Flag the end of the unary code		value-=1;		while (qbpp--) { writeBit(out,(value>>qbpp)&1); } 	// write whole value (always of length qbpp)	}//cerr << endl;	return out;}static voidquantizeErrval(Uint16 NEAR,Int32 &Errval){//cerr << "\t\tquantizeErrval: before Errval = " << Errval << endl;	if (NEAR) {		if (Errval > 0)			Errval=(Errval+NEAR)/(2*NEAR+1);		else			Errval=(Errval-NEAR)/(2*NEAR+1);		// in A.4.4 it is actually -(NEAR-Errval)/(2*NEAR+1)	}	// else leave Errval as it is for lossless mode//cerr << "\t\tquantizeErrval: after Errval = " << Errval << endl;}static voiddeQuantizeErrval(Uint16 NEAR,Int32 &Errval){//cerr << "\t\tdeQuantizeErrval: before Errval = " << Errval << endl;	if (NEAR) Errval=Errval*(2*NEAR+1);//cerr << "\t\tdeQuantizeErrval: after Errval = " << Errval << endl;}static inline voidclampPredictedValue(Int32 &X,Int32 MAXVAL){//cerr << "\t\tclampPredictedValue: before value = " << X << endl;	if      (X > MAXVAL)	X=MAXVAL;	else if (X < 0)		X=0;//cerr << "\t\tclampPredictedValue: after value = " << X << endl;}static voidcodecRunEndSample(Uint16 &Ix,Int32 Ra,Int32 Rb,Int32 RANGE,Uint16 NEAR,Uint32 MAXVAL,Uint16 RESET,		Uint16 LIMIT,Uint16 qbpp,Uint16 rk,		Uint32 *A,Int32 *N,Int32 *Nn,		BinaryInputStream &in,BinaryOutputStream &out,bool decompressing){//cerr << "\t\tcodecRunEndSample: " << (decompressing ? "decoding" : "encoding") << endl;//if (!decompressing) cerr << "\t\tcodecRunEndSample: value = " << Ix << endl;	bool RItype = (Ra == Rb || Abs(Ra-Rb) <= NEAR);	Int16 SIGN = (!RItype && Ra > Rb) ? -1 : 1;	Int32 Px = RItype ? Ra : Rb;//cerr << "\t\tcodecRunEndSample: Ra = " << Ra << endl;//cerr << "\t\tcodecRunEndSample: Rb = " << Rb << endl;//cerr << "\t\tcodecRunEndSample: RItype = " << (RItype ? "1":"0") << endl;//cerr << "\t\tcodecRunEndSample: SIGN = " << SIGN << endl;//cerr << "\t\tcodecRunEndSample: Px = " << Px << endl;	Uint32 TEMP = RItype ? A[366]+(N[366]>>1) : A[365];	Uint16 Q = 365 + (RItype ? 1 : 0);//cerr << "\t\tcodecRunEndSample: TEMP = " << TEMP << endl;//cerr << "\t\tcodecRunEndSample: Q = " << Q << endl;	Uint16 k = determineGolombParameter(N[Q],TEMP);//cerr << "\t\tcodecRunEndSample: k = " << k << endl;	Int32  Errval;	Int32  updateErrval;	Uint32 EMErrval;	if (decompressing) {		decodeMappedErrvalWithGolomb(k,LIMIT-rk-1,qbpp,EMErrval,in);	// needs work :(//cerr << "\t\tcodecRunEndSample: EMErrval = " << EMErrval << endl;		Uint32 tEMErrval = EMErrval + (RItype ? 1 : 0);		// use local copy to leave original for parameter update later//cerr << "\t\tcodecRunEndSample: tEMErrval = " << tEMErrval << endl;		if (tEMErrval == 0) {			Errval = 0;		}		else if (k == 0) {			if (2*Nn[Q-365] < N[Q]) {				if (tEMErrval%2 == 0) {					Errval = -Int32(tEMErrval>>1);		// "map = 0"	2 becomes -1, 4 becomes -2, 6 becomes -3				}				else {					Errval = (tEMErrval+1)>>1;		// "map = 1"	1 becomes 1, 3 becomes 2, 5 becomes 3				}			}			else {	// 2*Nn[Q-365] >= N[Q]				if (tEMErrval%2 == 0) {					Errval = tEMErrval>>1;			// "map = 0"	2 becomes 1, 4 becomes 2, 6 becomes 3				}				else {					Errval = -Int32((tEMErrval+1)>>1);	// "map = 1"	1 becomes -1, 3 becomes -2, 5 becomes -3				}			}		}		else {			if (tEMErrval%2 == 0) {				Errval = tEMErrval>>1;				// "map = 0"	2 becomes  1, 4 becomes  2, 6 becomes 3			}			else {				Errval = -Int32((tEMErrval+1)>>1);		// "map = 1"	1 becomes -1, 3 becomes -2, 5 becomes -3			}		}//cerr << "\t\tcodecRunEndSample: Errval after sign unmapping = " << Errval << endl;		updateErrval=Errval;		if (NEAR > 0) deQuantizeErrval(NEAR,Errval);//cerr << "\t\tcodecRunEndSample: Errval SIGN uncorrected = " << Errval << endl;		if (SIGN < 0) Errval=-Errval;		// if "context type" was negative//cerr << "\t\tcodecRunEndSample: Errval result = " << Errval << endl;		Int32 Rx = Px+Errval;		// modulo(RANGE*(2*NEAR+1)) as per F.1 Item 14		// (NB. Is this really the reverse of the encoding procedure ???)		if (Rx < -NEAR)			Rx+=RANGE*(2*NEAR+1);		else if (Rx > MAXVAL+NEAR)			Rx-=RANGE*(2*NEAR+1);		clampPredictedValue(Rx,MAXVAL);		// Apply inverse point transform and mapping table when implemented		Ix=(Uint16)Rx;	}	else {//cerr << "\t\tIx " << Ix << endl;		Errval = Int32(Ix) - Px;//cerr << "\t\tErrval " << Errval << endl;		if (SIGN < 0) Errval=-Errval;		// if "context type" was negative//cerr << "\t\tErrval sign corrected " << Errval << endl;		// Figure out sign to later correct Errval (Figure A.19) ...		if (NEAR > 0) {	// For near-lossless, quantize Errval and derive reconstructed value (A.4.4)			quantizeErrval(NEAR,Errval);//cerr << "\t\tErrval quantized " << Errval << endl;			// Replace with the reconstructed value the decoder will have			// (obviously if in lossless mode there will be no difference)			Int32 Rx=Px+SIGN*Errval*(2*NEAR+1);			clampPredictedValue(Rx,MAXVAL);								Ix=(Uint16)Rx;//cerr << "\t\tReplaced Rx " << Rx << endl;		}		// Modulo reduction of the prediction error (A.4.5)		if (Errval < 0)			Errval=Errval+RANGE;		if (Errval >= (RANGE+1)/2)	Errval=Errval-RANGE;//cerr << "\t\tErrval modulo " << RANGE << " = " << Errval << endl;		updateErrval=Errval;		// Golomb stuff is outside decompress/compress decision since same		// Map error to non-negative ...		// Int16 map = ((k == 0 && Errval > 0 && 2*Nn[Q-365] < N[Q]) || (Errval < 0 && (2*Nn[Q-365] >= N[Q] || k != 0)) ? 1 : 0;		// EMErrval = 2*Abs(Errval) - RItype - map;		if (k == 0) {			if (Errval > 0) {				if (2*Nn[Q-365] < N[Q]) {					EMErrval = 2*Errval - 1;		// "map = 1"	1 becomes 1, 2 becomes 3, 3 becomes 5				}				else {	// 2*Nn[Q-365] >= N[Q]					EMErrval = 2*Errval;			// "map = 0"	1 becomes 2, 2 becomes 4, 3 becomes 6				}			}			else if (Errval < 0) {				if (2*Nn[Q-365] < N[Q]) {					EMErrval = -2*Errval;			// "map = 0"	-1 becomes 2, -2 becomes 4, -3 becomes 6				}				else {	// 2*Nn[Q-365] >= N[Q]					EMErrval = -2*Errval - 1;		// "map = 1"	-1 becomes 1, -2 becomes 3, -3 becomes 5				}			}			else { // Errval == 0				EMErrval = 0;					// "map = 0"	0 stays 0			}		}		else {	// k != 0			if (Errval > 0) {				EMErrval = 2*Errval;				// "map = 0"	1 becomes 2, 2 becomes 4, 3 becomes 6			}			else if (Errval < 0) {				EMErrval = -2*Errval - 1;			// "map = 1"	-1 becomes 1, -2 becomes 3, -3 becomes 5			}			else { // Errval == 0				EMErrval = 0;					// "map = 0"	0 stays 0			}		}//cerr << "\t\tcodecRunEndSample: EMErrval before subtraction of RItype = " << EMErrval << endl;		EMErrval-=(RItype ? 1 : 0);//cerr << "\t\tcodecRunEndSample: EMErrval after subtraction of RItype = " << EMErrval << endl;		encodeMappedErrvalWithGolomb(k,LIMIT-rk-1,qbpp,EMErrval,out);	}	// Update parameters ...//cerr << "\t\tcodecRunEndSample: Update parameters ... updateErrval used = " << updateErrval << endl;//cerr << "\t\tcodecRunEndSample: Update parameters ... EMErrval used = " << EMErrval << endl;//cerr << "\t\tcodecRunEndSample: A[" << Q << "]  before = " << A[Q]  << endl;//cerr << "\t\tcodecRunEndSample: N[" << Q << "]  before = " << N[Q]  << endl;//cerr << "\t\tcodecRunEndSample: Nn[" << (Q-365) << "] before = " << Nn[Q-365] << endl;	if (updateErrval < 0) ++Nn[Q-365];	A[Q]+=(EMErrval+1-(RItype ? 1 : 0))>>1;	if (N[Q] == RESET) {		A[Q]=A[Q]>>1;		N[Q]=N[Q]>>1;		Nn[Q-365]=Nn[Q-365]>>1;	}	++N[Q];//cerr << "\t\tcodecRunEndSample: A[" << Q << "]  updated = " << A[Q]  << endl;//cerr << "\t\tcodecRunEndSample: N[" << Q << "]  updated = " << N[Q]  << endl;//cerr << "\t\tcodecRunEndSample: Nn[" << (Q-365) << "] updated = " << Nn[Q-365] << endl;//if (decompressing) cerr << "\t\tcodecRunEndSample: value = " << Ix << endl;}intmain(int argc,char **argv){	bool bad=false;	GetNamedOptions				options(argc,argv);	BinaryInputOptionsWithByteOrder		input_options(options);	BinaryOutputOptionsWithByteOrder	output_options(options);	bool verbose=options.get("v") || options.get("verbose");	bool decompressing=options.get("d") || options.get("decompress");	bool useJPEGmarkers=!options.get("nomarkers");	bool useRunMode=!options.get("noruns");	unsigned rows=0;	if ((!decompressing || !useJPEGmarkers) && !options.get("rows",rows) && !options.get("height",rows) && !options.get("h",rows)) {		cerr << EMsgDC(NeedOption) << " - rows" << endl;		bad=true;	}	unsigned cols=0;	if ((!decompressing || !useJPEGmarkers) && !options.get("columns",cols) && !options.get("width",cols) && !options.get("w",cols)) {		cerr << EMsgDC(NeedOption) << " - columns" << endl;		bad=true;	}	unsigned bits=0;	if ((!decompressing || !useJPEGmarkers) && !options.get("bits",bits) && !options.get("depth",bits)) {		cerr << EMsgDC(NeedOption) << " - bits" << endl;		bad=true;	}	Assert(bits <= 16);	Uint16 NEAR=0;		// Lossless if zero	Uint16 T1=0;	Uint16 T2=0;	Uint16 T3=0;	Uint16 RESET=0;	if ((!decompressing || !useJPEGmarkers)) {		unsigned near;	if (options.get("near",near))				NEAR=near;		unsigned t1;	if (options.get("T1",t1) || options.get("Ta",t1))	T1=t1;		unsigned t2;	if (options.get("T2",t2) || options.get("Tb",t2))	T2=t2;		unsigned t3;	if (options.get("T3",t3) || options.get("Tc",t3))	T3=t3;		unsigned reset;	if (options.get("reset",reset))				RESET=reset;	}	input_options.done();	output_options.done();