options.done();	BinaryInputOpenerFromOptions input_opener(		options,input_options.filename,cin);	BinaryOutputOpenerFromOptions output_opener(		options,output_options.filename,cout);	cerr << input_options.errors();	cerr << output_options.errors();	cerr << options.errors();	cerr << input_opener.errors();	cerr << output_opener.errors();	if (!input_options.good()	 || !output_options.good()	 || !options.good()	 || !input_opener.good()	 || !output_opener.good()	 || !options	 || bad) {		cerr 	<< MMsgDC(Usage) << ": " << options.command()			<< input_options.usage()			<< output_options.usage()			<< " [-d|decompress]"			<< " -rows|height|h n"			<< " -columns|width|w n"			<< " -bits|depth n"			<< " [-near n]"			<< " [-T1|Ta n]"			<< " [-T2|Tb n]"			<< " [-T3|Tc n]"			<< " [-reset n]"			<< " [-nomarkers]"			<< " [-noruns]"			<< " [-v|verbose]"			<< " [" << MMsgDC(InputFile)				<< "[" << MMsgDC(OutputFile) << "]]"			<< " <" << MMsgDC(InputFile)			<< " >" << MMsgDC(OutputFile)			<< endl;		exit(1);	}	BinaryInputStream in(*(istream *)input_opener,input_options.byteorder);	BinaryOutputStream out(*(ostream *)output_opener,output_options.byteorder);	bool success=true;	Uint32 ROWS;	Uint32 COLUMNS;	Uint16 P;		// Sample precision	Uint32 MAXVAL;	bool haveLSE1=false;	if (decompressing && useJPEGmarkers) {		bool readrequiredmarkers=false;		Uint16 marker;		if (readJPEGMarker(in,marker) && readSOI(in,marker)) {			if (readJPEGMarker(in,marker) && readSOF55(in,marker,P,ROWS,COLUMNS)			 || (haveLSE1=readLSE1(in,marker,MAXVAL,T1,T2,T3,RESET)) && readJPEGMarker(in,marker) && readSOF55(in,marker,P,ROWS,COLUMNS)) {				if (readJPEGMarker(in,marker) && readSOS(in,marker,NEAR)				 || !haveLSE1 && (haveLSE1=readLSE1(in,marker,MAXVAL,T1,T2,T3,RESET)) && readJPEGMarker(in,marker) && readSOS(in,marker,NEAR)) {					readrequiredmarkers=true;				}				else {					cerr << "Corrupt JPEG stream ... expected SOS Marker" << endl;				}			}			else {				cerr << "Corrupt JPEG stream ... expected SOF55 Marker" << endl;			}		}		else {			cerr << "Corrupt JPEG stream ... expected SOI Marker" << endl;		}		if (!readrequiredmarkers) {			exit(1);		}	}	else {		P=bits;		ROWS=rows;		COLUMNS=cols;	}	if (!decompressing || !useJPEGmarkers || !haveLSE1) {	// Note that the LSE ID 1 marker is optional		MAXVAL=(1ul<<P)-1;		if (!RESET) RESET=64;		// May have been set on command line		// Initialization of default parameters as per A.1 reference to C.2.4.1.1.1		// Thresholds for context gradients ...		const Uint32 BASIC_T1 = 3;		const Uint32 BASIC_T2 = 7;		const Uint32 BASIC_T3 = 21;#define CLAMP_1(i)	((i > MAXVAL || i < NEAR+1) ? NEAR+1 : i)#define CLAMP_2(i)	((i > MAXVAL || i < T1) ? T1 : i)#define CLAMP_3(i)	((i > MAXVAL || i < T2) ? T2 : i)		// Only replace T1, T2, T3 if not set on command line ...		if (MAXVAL >= 128) {			Uint32 FACTOR=FloorDivision(Minimum(MAXVAL,4095)+128,256);			if (verbose) cerr << "MAXVAL >= 128" << endl;			if (verbose) cerr << "FACTOR = " << dec << FACTOR << endl;			if (!T1) T1=CLAMP_1(FACTOR*(BASIC_T1-2)+2+3*NEAR);			if (!T2) T2=CLAMP_2(FACTOR*(BASIC_T2-3)+3+5*NEAR);			if (!T3) T3=CLAMP_3(FACTOR*(BASIC_T3-4)+4+7*NEAR);		}		else {			Uint32 FACTOR=FloorDivision(256,MAXVAL+1);			if (verbose) cerr << "MAXVAL < 128" << endl;			if (verbose) cerr << "FACTOR = " << dec << FACTOR << endl;			if (!T1) T1=CLAMP_1(Maximum(2,BASIC_T1/FACTOR+3*NEAR));	// ? should these calculations be float since we are dividing ? :(			if (!T2) T2=CLAMP_2(Maximum(3,BASIC_T2/FACTOR+5*NEAR));			if (!T3) T3=CLAMP_3(Maximum(4,BASIC_T3/FACTOR+7*NEAR));		}	}	if (verbose) cerr << "NEAR = " << NEAR << endl;	if (verbose) cerr << "ROWS = " << ROWS << endl;	if (verbose) cerr << "COLUMNS = " << COLUMNS << endl;	if (verbose) cerr << "P = " << P << endl;	if (verbose) cerr << "MAXVAL = " << MAXVAL << endl;	if (verbose) cerr << "RESET = " << RESET << endl;	if (verbose) cerr << "T1 = " << T1 << endl;	if (verbose) cerr << "T2 = " << T2 << endl;	if (verbose) cerr << "T3 = " << T3 << endl;	Assert(ROWS);	Assert(COLUMNS);	Assert(P);	Assert(T1 == 0 || (NEAR+1 <= T1 && T1 <= MAXVAL));	Assert(T2 == 0 || (T1 <= T2 && T2 <= MAXVAL));	Assert(T3 == 0 || (T2 <= T3 && T3 <= MAXVAL));	// Initialization as per Annex A.2.1	Int32 RANGE = FloorDivision(MAXVAL+2*NEAR,2*NEAR+1)+1;	// int not unsigned to avoid need for cast when used	if (verbose) cerr << "RANGE = " << RANGE << endl;	Assert(MAXVAL == 0 || (1 <= MAXVAL && MAXVAL < 1ul<<P));	if (NEAR == 0) Assert(RANGE == MAXVAL+1);	Uint16 bpp = Maximum(2,Ceiling(Log(MAXVAL+1)));	// Number of bits needed to represent MAXVAL with a minumum of 2	Uint16 qbpp = Ceiling(Log(RANGE));		// Number of bits needed to represent a mapped error value	Uint16 LIMIT = 2*(bpp+Maximum(8,bpp));		// the value of glimit for a sample encoded in regular mode	if (verbose) cerr << "bpp = " << bpp << endl;	if (verbose) cerr << "qbpp = " << qbpp << endl;	if (verbose) cerr << "LIMIT = " << LIMIT << endl;	Assert(bpp >= 2);	Assert(LIMIT > qbpp);			// Else LIMIT-qbpp-1 will fail (see A.5.3)	if (!decompressing && useJPEGmarkers)  {		writeSOI(out);		writeSOF55(out,P,ROWS,COLUMNS);		writeLSE1(out,MAXVAL,T1,T2,T3,RESET);		writeSOS(out,NEAR);	}	// Fixed constants	const Uint16 nContexts = 365;	// plus two more run mode interruption contexts	const Int32 MIN_C = -128;	// Limits on values in bias correction array C	const Int32 MAX_C =  127;	// Initialization of variables ...	Int32	*N = new Int32[nContexts+2];	// counters for context type occurence [0..nContexts+2-1]						// [nContexts],[nContexts+1] for run mode interruption	Uint32	*A = new Uint32[nContexts+2];	// accumulated prediction error magnitude [0..nContexts-1]						// [nContexts],[nContexts+1] for run mode interruption	Int32	*B = new Int32[nContexts];	// auxilliary counters for bias cancellation [0..nContexts-1]	Int32	*C = new Int32[nContexts];	// counters indicating bias correction value [0..nContexts-1]						// (never -ve but often used as -N[Q] so int not unsigned saves cast)	Int32	*Nn = new Int32[2];		// negative prediction error for run interruption [365..366]	Assert(N);	Assert(A);	Assert(B);	Assert(C);	Assert(Nn);	{		Uint32 A_Init_Value=Maximum(2,FloorDivision(RANGE+(1lu<<5),(1lu<<6)));		if (verbose) cerr << "A_Init_Value = " << A_Init_Value << endl;		unsigned i;		for (i=0; i<nContexts; ++i) {			N[i]=1;			A[i]=A_Init_Value;			B[i]=C[i]=0;		}		N[nContexts]=1;		N[nContexts+1]=1;		A[nContexts]=A_Init_Value;		A[nContexts+1]=A_Init_Value;	}	Nn[365-365]=Nn[366-365]=0;	// The run variables seem to need to live beyond a single run or row !!!	unsigned RUNIndex = 0;	Uint16 *rowA = new Uint16[COLUMNS];	Uint16 *rowB = new Uint16[COLUMNS];	Assert(rowA);	Assert(rowB);	Uint32 row=0;	Uint16 *thisRow=rowA;	Uint16 *prevRow=rowB;	for (row=0; row<ROWS; ++row) {//cerr << "Row " << row << endl;		if (!decompressing)  {			Uint32 n=readRow(in,thisRow,COLUMNS,bpp);//cerr << "Row " << row << " read returns " << n << endl;			Assert (n==COLUMNS);		}		Uint32 col=0;		Uint16 prevRa0;		for (col=0; col<COLUMNS; ++col) {//cerr << "\tcol = " << col << endl;			//	c b d .			//	a x . .			//	. . . .			Int32 Rx;	// Reconstructed value - not Uint16 to allow overrange before clamping			// value at edges (first row and first col is zero) ...			Uint16 Ra;			Uint16 Rb;			Uint16 Rc;			Uint16 Rd;			if (row > 0) {				Rb=prevRow[col];				Rc=(col > 0) ? prevRow[col-1] : prevRa0;				Ra=(col > 0) ? thisRow[col-1] : (prevRa0=Rb);				Rd=(col+1 < COLUMNS) ? prevRow[col+1] : Rb;			}			else {				Rb=Rc=Rd=0;				Ra=(col > 0) ? thisRow[col-1] : (prevRa0=0);			}				//cerr << "\t\tRa = " << Ra << endl;//cerr << "\t\tRb = " << Rb << endl;//cerr << "\t\tRc = " << Rc << endl;//cerr << "\t\tRd = " << Rd << endl;			// NB. We want the Reconstructed values, which are the same			// in lossless mode, but if NEAR != 0 take care to write back			// reconstructed values into the row buffers in previous positions			// Compute local gradient ...			Int32 D1=(Int32)Rd-Rb;			Int32 D2=(Int32)Rb-Rc;			Int32 D3=(Int32)Rc-Ra;//cerr << "\t\tD1 = " << D1 << endl;//cerr << "\t\tD2 = " << D2 << endl;//cerr << "\t\tD3 = " << D3 << endl;			// Check for run mode ... (should check Abs() works ok for Int32)			if (Abs(D1) <= NEAR && Abs(D2) <= NEAR && Abs(D3) <= NEAR && useRunMode) {				// Run mode//cerr << "Row at run start " << row << endl;//cerr << "\tcol at run start " << col << endl;				if (decompressing) {//dumpReadBitPosition();					// Why is RUNIndex not reset to 0 here ?					Uint32 R;					while (readBit(in,R)) {//cerr << "\tcol " << col << endl;						if (R == 1) {							// Fill image with 2^J[RUNIndex] samples of Ra or till EOL							Int32 rm=J_rm[RUNIndex];//cerr << "\tRUNIndex " << RUNIndex << endl;//cerr << "\tFilling with " << rm << " samples of Ra " << Ra << endl;							while (rm-- && col < COLUMNS) {								thisRow[col]=Ra;//cerr << "pixel[" << row << "," << col << "] = " << thisRow[col] << endl;								++col;							}							// This will match when exact count coincides with end of row ...							if (rm == -1 && RUNIndex < 31) {								++RUNIndex;//cerr << "\tRUNIndex incremented to " << RUNIndex << endl;							}							if (col >= COLUMNS) {//cerr << "\tFilled to end of row" << endl;//cerr << "\tAfter having found end of row " << endl;//dumpReadBitPosition();								break;							}						}						else {							// Read J[RUNIndex] bits and fill image with that number of samples of Ra							Uint16 bits=J[RUNIndex];//cerr << "\tRUNIndex " << RUNIndex << endl;//cerr << "\tReading bits " << bits << endl;							Uint32 nfill=0;							Uint32 bit;							// msb bit is read first							while (bits-- && readBit(in,bit)) {								nfill=(nfill<<1) | bit;							}//cerr << "\tFill with " << nfill << " samples of Ra " << Ra << endl;							// Fill with nfill values of Ra							while (nfill--) {//if (!(col<(COLUMNS-1))) {//	cerr << "Fail at line 367 ... !(col<(COLUMNS-1))" << endl;//	cerr << "\tstill to fill " << nfill+1 << endl;//	cerr << "\trow is " << row << endl;//	cerr << "\tcol is " << col << endl;//}								Assert(col<(COLUMNS-1));								thisRow[col]=Ra;//cerr << "pixel[" << row << "," << col << "] = " << thisRow[col] << endl;								++col;							}							// Decode the run interruption sample ...//cerr << "\tcol at end of run " << col << endl;//cerr << "\tBefore decoding value that ends run " << endl;//dumpReadBitPosition();							// First update local context for interrupting sample, since weren't kept updated during run							if (row > 0) {								Rb=prevRow[col];								Ra=(col > 0) ? thisRow[col-1] : Rb;							}							else {								Rb=0;								Ra=(col > 0) ? thisRow[col-1] : 0;							}//cerr << "\t\tRa = " << Ra << endl;//cerr << "\t\tRb = " << Rb << endl;							codecRunEndSample(thisRow[col],Ra,Rb,RANGE,NEAR,MAXVAL,RESET,LIMIT,qbpp,J[RUNIndex],A,N,Nn,in,out,decompressing);//cerr << "pixel[" << row << "," << col << "] = " << thisRow[col] << endl;//cerr << "\tValue that ends run " << thisRow[col] << endl;//dumpReadBitPosition();							if (RUNIndex > 0) {								--RUNIndex;	// NB. Do this AFTER J[RUNIndex] used in the limited length Golomb coding//cerr << "\tRUNIndex decremented to " << RUNIndex << endl;							}							break;						}					}				}				else {//dumpWriteBitPosition();					// Scan to determine length of run (A.7.1.1) ...//cerr << "\tScan to determine length of run" << endl;//cerr << "\tRa is " << Ra << endl;					Int32 RUNval=Ra;					Int32 RUNcnt=0;					while (col < COLUMNS && (thisRow[col] == RUNval || (NEAR > 0 && Abs(Int32(thisRow[col])-RUNval) <= NEAR) ) ) {//cerr << "\tpixel[" << row << "," << col << "] = " << thisRow[col] << endl;						++RUNcnt;						if (NEAR > 0) thisRow[col]=RUNval;	// Replace with "reconstructed value"						++col;					}//cerr << "\tAt end of run, thisRow[" << row << "," << col << "] = " << thisRow[col] << " and RUNval = " << RUNval << endl;//cerr << "\tAbs(Int32(thisRow[col])-RUNval) = " << Abs(Int32(thisRow[col])-RUNval) << endl;//cerr << "\tNEAR = " << NEAR << endl;