jpegdec.c

此是JPEG编码程序

// JPEG decoder module
// Copyright 1999 Cristi Cuturicu

#include "jpegdec.h"
// Used markers:
#define SOI 0xD8
#define EOI 0xD9
#define APP0 0xE0
#define SOF 0xC0
#define DQT 0xDB
#define DHT 0xC4
#define SOS 0xDA
#define DRI 0xDD
#define COM 0xFE

char error_string[90];
#define exit_func(err) { strcpy(error_string, err); return 0;}

static BYTE *buf; // the buffer we use for storing the entire JPG file

static BYTE bp; //current byte
static WORD wp; //current word

static DWORD byte_pos; // current byte position
#define BYTE_p(i) bp=buf[(i)++]
#define WORD_p(i) wp=(((WORD)(buf[(i)]))<<8) + buf[(i)+1]; (i)+=2

// WORD X_image_size,Y_image_size; // X,Y sizes of the image
static WORD X_round,Y_round; // The dimensions rounded to multiple of Hmax*8 (X_round)
			  // and Ymax*8 (Y_round)

static BYTE *im_buffer; // RGBA image buffer
static DWORD X_image_bytes; // size in bytes of 1 line of the image = X_round * 4
static DWORD y_inc_value ; // 32*X_round; // used by decode_MCU_1x2,2x1,2x2

BYTE YH,YV,CbH,CbV,CrH,CrV; // sampling factors (horizontal and vertical) for Y,Cb,Cr
static WORD Hmax,Vmax;


static BYTE zigzag[64]={ 0, 1, 5, 6,14,15,27,28,
				  2, 4, 7,13,16,26,29,42,
				  3, 8,12,17,25,30,41,43,
				  9,11,18,24,31,40,44,53,
				 10,19,23,32,39,45,52,54,
				 20,22,33,38,46,51,55,60,
				 21,34,37,47,50,56,59,61,
				 35,36,48,49,57,58,62,63 };
typedef struct {
   BYTE Length[17];  // k =1-16 ; L[k] indicates the number of Huffman codes of length k
   WORD minor_code[17];  // indicates the value of the smallest Huffman code of length k
   WORD major_code[17];  // similar, but the highest code
   BYTE V[65536];  // V[k][j] = Value associated to the j-th Huffman code of length k
	// High nibble = nr of previous 0 coefficients
	// Low nibble = size (in bits) of the coefficient which will be taken from the data stream
} Huffman_table;

static float *QT[4]; // quantization tables, no more than 4 quantization tables (QT[0..3])
static Huffman_table HTDC[4]; //DC huffman tables , no more than 4 (0..3)
static Huffman_table HTAC[4]; //AC huffman tables                  (0..3)

static BYTE YQ_nr,CbQ_nr,CrQ_nr; // quantization table number for Y, Cb, Cr
static BYTE YDC_nr,CbDC_nr,CrDC_nr; // DC Huffman table number for Y,Cb, Cr
static BYTE YAC_nr,CbAC_nr,CrAC_nr; // AC Huffman table number for Y,Cb, Cr

static BYTE Restart_markers; // if 1 => Restart markers on , 0 => no restart markers
static WORD MCU_restart; //Restart markers appears every MCU_restart MCU blocks
typedef void (*decode_MCU_func)(DWORD);


static SWORD DCY, DCCb, DCCr; // Coeficientii DC pentru Y,Cb,Cr
static SWORD DCT_coeff[64]; // Current DCT_coefficients
static BYTE Y[64],Cb[64],Cr[64]; // Y, Cb, Cr of the current 8x8 block for the 1x1 case
static BYTE Y_1[64],Y_2[64],Y_3[64],Y_4[64];
static BYTE tab_1[64],tab_2[64],tab_3[64],tab_4[64]; // tabelele de supraesantionare pt cele 4 blocuri

static SWORD Cr_tab[256],Cb_tab[256]; // Precalculated Cr, Cb tables
static SWORD Cr_Cb_green_tab[65536];

// Initial conditions:
// byte_pos = start position in the Huffman coded segment
// WORD_get(w1); WORD_get(w2);wordval=w1;

static BYTE d_k=0;  // Bit displacement in memory, relative to the offset of w1
			 // it's always <16
static WORD w1,w2; // w1 = First word in memory; w2 = Second word
static DWORD wordval ; // the actual used value in Huffman decoding.
static DWORD mask[17];
static SWORD neg_pow2[17]={0,-1,-3,-7,-15,-31,-63,-127,-255,-511,-1023,-2047,-4095,-8191,-16383,-32767};
static DWORD start_neg_pow2=(DWORD)neg_pow2;

static int shift_temp;
#define RIGHT_SHIFT(x,shft)  \
	((shift_temp = (x)) < 0 ? \
	 (shift_temp >> (shft)) | ((~(0L)) << (32-(shft))) : \
	 (shift_temp >> (shft)))
#define DESCALE(x,n)  RIGHT_SHIFT((x) + (1L << ((n)-1)), n)
#define RANGE_MASK 1023L
static BYTE *rlimit_table;
void prepare_range_limit_table()
/* Allocate and fill in the sample_range_limit table */
{
  int j;
  rlimit_table = (BYTE *)malloc(5 * 256L + 128) ;
  /* First segment of "simple" table: limit[x] = 0 for x < 0 */
  memset((void *)rlimit_table,0,256);
  rlimit_table += 256;	/* allow negative subscripts of simple table */
  /* Main part of "simple" table: limit[x] = x */
  for (j = 0; j < 256; j++) rlimit_table[j] = j;
  /* End of simple table, rest of first half of post-IDCT table */
  for (j = 256; j < 640; j++) rlimit_table[j] = 255;
  /* Second half of post-IDCT table */
  memset((void *)(rlimit_table + 640),0,384);
  for (j = 0; j < 128 ; j++) rlimit_table[j+1024] = j;
}

#ifdef _MSC_VER
WORD lookKbits(BYTE k)
{
 _asm {
	 mov dl, k
	 mov cl, 16
	 sub cl, dl
	 mov eax, [wordval]
	 shr eax, cl
 }
}

WORD WORD_hi_lo(BYTE byte_high,BYTE byte_low)
{
	_asm {
		  mov ah,byte_high
		  mov al,byte_low
		 }
}
SWORD get_svalue(BYTE k)
// k>0 always
// Takes k bits out of the BIT stream (wordval), and makes them a signed value
{
	_asm {
		   xor ecx, ecx
		   mov cl,k
		   mov eax,[wordval]
		   shl eax,cl
		   shr eax, 16
		   dec cl
		   bt eax,ecx
		   jc end_macro
	signed_value:inc cl
		   mov ebx,[start_neg_pow2]
		   add ax,word ptr [ebx+ecx*2]
		 end_macro:
	}
}
#endif

#ifdef __WATCOMC__

WORD lookKbits(BYTE k);
#pragma aux lookKbits=\
					  "mov eax,[wordval]"\
					  "mov cl, 16"\
					  "sub cl, dl"\
					  "shr eax, cl"\
					   parm [dl] \
					   value [ax] \
					   modify [eax cl];
WORD WORD_hi_lo(BYTE byte_high,BYTE BYTE_low);
#pragma aux WORD_hi_lo=\
			   parm [ah] [al]\
			   value [ax] \
			   modify [ax];

SWORD get_svalue(BYTE k);
// k>0 always
// Takes k bits out of the BIT stream (wordval), and makes them a signed value
#pragma aux get_svalue=\
			"xor ecx, ecx"\
			"mov cl, al"\
			"mov eax,[wordval]"\
			"shl eax, cl"\
			"shr eax, 16"\
			"dec cl"\
			"bt eax,ecx"\
			"jc end_macro"\
			"signed_value:inc cl"\
			"mov ebx,[start_neg_pow2]"\
			"add ax,word ptr [ebx+ecx*2]"\
			"end_macro:"\
			parm [al]\
			modify [eax ebx ecx]\
			value [ax];
#endif

void skipKbits(BYTE k)
{
 BYTE b_high,b_low;
 d_k+=k;
 if (d_k>=16) { d_k-=16;
		w1=w2;
		// Get the next word in w2
		BYTE_p(byte_pos);
		if (bp!=0xFF) b_high=bp;
		else {
			  if (buf[byte_pos]==0) byte_pos++; //skip 00
			  else byte_pos--; // stop byte_pos pe restart marker
			  b_high=0xFF;
		}
		BYTE_p(byte_pos);
		if (bp!=0xFF) b_low=bp;
		else {
			  if (buf[byte_pos]==0) byte_pos++; //skip 00
			  else byte_pos--; // stop byte_pos pe restart marker
			  b_low=0xFF;
		}
		w2=WORD_hi_lo(b_high,b_low);
 }

 wordval = ((DWORD)(w1)<<16) + w2;
 wordval <<= d_k;
 wordval >>= 16;
}

SWORD getKbits(BYTE k)
{
 SWORD signed_wordvalue;
 signed_wordvalue=get_svalue(k);
 skipKbits(k);
 return signed_wordvalue;
}

void calculate_mask()
{
  BYTE k;
  DWORD tmpdv;
  for (k=0;k<=16;k++) { tmpdv=0x10000;mask[k]=(tmpdv>>k)-1;} //precalculated bit mask
}

void init_QT()
{
 BYTE i;
 for (i=0;i<=3;i++) QT[i]=(float *)malloc(sizeof(float)*64);
}

void load_quant_table(float *quant_table)
{
 float scalefactor[8]={1.0f, 1.387039845f, 1.306562965f, 1.175875602f,
					   1.0f, 0.785694958f, 0.541196100f, 0.275899379f};
 BYTE j,row,col;
// Load quantization coefficients from JPG file, scale them for DCT and reorder
// from zig-zag order
 for (j=0;j<=63;j++) quant_table[j]=buf[byte_pos+zigzag[j]];
 j=0;
 for (row=0;row<=7;row++)
   for (col=0;col<=7;col++) {
		quant_table[j]*=scalefactor[row]*scalefactor[col];
		j++;
   }
 byte_pos+=64;
}

void load_Huffman_table(Huffman_table *HT)
{
  BYTE k,j;
  DWORD code;

  for (j=1;j<=16;j++) {
	BYTE_p(byte_pos);
	HT->Length[j]=bp;
  }
  for (k=1;k<=16;k++)
	for (j=0;j<HT->Length[k];j++) {
		BYTE_p(byte_pos);
		HT->V[WORD_hi_lo(k,j)]=bp;
	}

  code=0;
  for (k=1;k<=16;k++) {
	 HT->minor_code[k] = (WORD)code;
	 for (j=1;j<=HT->Length[k];j++) code++;
	 HT->major_code[k]=(WORD)(code-1);
	 code*=2;
	 if (HT->Length[k]==0) {
		HT->minor_code[k]=0xFFFF;
		HT->major_code[k]=0;
	 }
  }
}

void process_Huffman_data_unit(BYTE DC_nr, BYTE AC_nr,SWORD *previous_DC)
{
// Process one data unit. A data unit = 64 DCT coefficients
// Data is decompressed by Huffman decoding, then the array is dezigzag-ed
// The result is a 64 DCT coefficients array: DCT_coeff
   BYTE nr,k,j,EOB_found;
   register WORD tmp_Hcode;
   BYTE size_val,count_0;
   WORD *min_code,*maj_code; // min_code[k]=minimum code of length k, maj_code[k]=similar but maximum
   WORD *max_val, *min_val;
   BYTE *huff_values;
   SWORD DCT_tcoeff[64];
   BYTE byte_temp;

// Start Huffman decoding
// First the DC coefficient decoding
   min_code=HTDC[DC_nr].minor_code;
   maj_code=HTDC[DC_nr].major_code;
   huff_values=HTDC[DC_nr].V;

   for (nr = 0; nr < 64 ; nr++) DCT_tcoeff[nr] = 0; //Initialize DCT_tcoeff

   nr=0;// DC coefficient

   min_val = &min_code[1]; max_val = &maj_code[1];
   for (k=1;k<=16;k++) {
	 tmp_Hcode=lookKbits(k);
//	   max_val = &maj_code[k]; min_val = &min_code[k];
	 if ( (tmp_Hcode<=*max_val)&&(tmp_Hcode>=*min_val) ) { //Found a valid Huffman code
		skipKbits(k);
		size_val=huff_values[WORD_hi_lo(k,(BYTE)(tmp_Hcode-*min_val))];
		if (size_val==0) DCT_tcoeff[0]=*previous_DC;
		else {
			DCT_tcoeff[0]=*previous_DC+getKbits(size_val);
			*previous_DC=DCT_tcoeff[0];
		}
		break;
	 }
	 min_val++; max_val++;
   }

// Second, AC coefficient decoding
   min_code=HTAC[AC_nr].minor_code;
   maj_code=HTAC[AC_nr].major_code;
   huff_values=HTAC[AC_nr].V;

   nr=1; // AC coefficient
   EOB_found=0;
   while ( (nr<=63)&&(!EOB_found) )
	{
	 max_val = &maj_code[1]; min_val =&min_code[1];
	 for (k=1;k<=16;k++)
	 {
	   tmp_Hcode=lookKbits(k);
//	   max_val = &maj_code[k]; &min_val = min_code[k];
	   if ( (tmp_Hcode<=*max_val)&&(tmp_Hcode>=*min_val) )
	   {
		skipKbits(k);
		byte_temp=huff_values[WORD_hi_lo(k,(BYTE)(tmp_Hcode-*min_val))];
		size_val=byte_temp&0xF;
		count_0=byte_temp>>4;
		if (size_val==0) {if (count_0==0) EOB_found=1;
						  else if (count_0==0xF) nr+=16; //skip 16 zeroes
						}
		else
		{
		 nr+=count_0; //skip count_0 zeroes
		 DCT_tcoeff[nr++]=getKbits(size_val);
		}
		break;
	   }
	   min_val++; max_val++;
	 }
	 if (k>16) nr++;  // This should not occur
	}
  for (j=0;j<=63;j++) DCT_coeff[j]=DCT_tcoeff[zigzag[j]]; // Et, voila ... DCT_coeff
}

void IDCT_transform(SWORD *incoeff,BYTE *outcoeff,BYTE Q_nr)
// Fast float IDCT transform
{
 BYTE x;
 SBYTE y;
 SWORD *inptr;
 BYTE *outptr;
 float workspace[64];
 float *wsptr;//Workspace pointer
 float *quantptr; // Quantization table pointer
 float dcval;
 float tmp0,tmp1,tmp2,tmp3,tmp4,tmp5,tmp6,tmp7;
 float tmp10,tmp11,tmp12,tmp13;
 float z5,z10,z11,z12,z13;
 BYTE *range_limit=rlimit_table+128;
 // Pass 1: process COLUMNS from input and store into work array.
 wsptr=workspace;
 inptr=incoeff;
 quantptr=QT[Q_nr];
 for (y=0;y<=7;y++)
  {
   if( (inptr[8]|inptr[16]|inptr[24]|inptr[32]|inptr[40]|inptr[48]|inptr[56])==0)
	{
	 // AC terms all zero (in a column)
	 dcval=inptr[0]*quantptr[0];
	 wsptr[0]  = dcval;
	 wsptr[8]  = dcval;
	 wsptr[16] = dcval;
	 wsptr[24] = dcval;
	 wsptr[32] = dcval;
	 wsptr[40] = dcval;
	 wsptr[48] = dcval;
	 wsptr[56] = dcval;
	 inptr++;quantptr++;wsptr++;//advance pointers to next column
	 continue ;
	}
  //Even part
	tmp0 = inptr[0] *quantptr[0];
	tmp1 = inptr[16]*quantptr[16];
	tmp2 = inptr[32]*quantptr[32];
	tmp3 = inptr[48]*quantptr[48];

	tmp10 = tmp0 + tmp2;// phase 3
	tmp11 = tmp0 - tmp2;

	tmp13 = tmp1 + tmp3;// phases 5-3
	tmp12 = (tmp1 - tmp3) * 1.414213562f - tmp13; // 2*c4

	tmp0 = tmp10 + tmp13;// phase 2
	tmp3 = tmp10 - tmp13;
	tmp1 = tmp11 + tmp12;
	tmp2 = tmp11 - tmp12;

	// Odd part
	tmp4 = inptr[8] *quantptr[8];
	tmp5 = inptr[24]*quantptr[24];
	tmp6 = inptr[40]*quantptr[40];
	tmp7 = inptr[56]*quantptr[56];

	z13 = tmp6 + tmp5;// phase 6
	z10 = tmp6 - tmp5;
	z11 = tmp4 + tmp7;
	z12 = tmp4 - tmp7;

	tmp7 = z11 + z13;// phase 5
	tmp11= (z11 - z13) * 1.414213562f; // 2*c4

	z5 = (z10 + z12) * 1.847759065f; // 2*c2
	tmp10 = 1.082392200f * z12 - z5; // 2*(c2-c6)
	tmp12 = -2.613125930f * z10 + z5;// -2*(c2+c6)

	tmp6 = tmp12 - tmp7;// phase 2
	tmp5 = tmp11 - tmp6;
	tmp4 = tmp10 + tmp5;

	wsptr[0]  = tmp0 + tmp7;
	wsptr[56] = tmp0 - tmp7;
	wsptr[8]  = tmp1 + tmp6;
	wsptr[48] = tmp1 - tmp6;
	wsptr[16] = tmp2 + tmp5;
	wsptr[40] = tmp2 - tmp5;
	wsptr[32] = tmp3 + tmp4;
	wsptr[24] = tmp3 - tmp4;
	inptr++;
	quantptr++;
	wsptr++;//advance pointers to the next column
  }

//  Pass 2: process ROWS from work array, store into output array.
// Note that we must descale the results by a factor of 8 = 2^3
 outptr=outcoeff;
 wsptr=workspace;
 for (x=0;x<=7;x++)
  {
   // Even part
	tmp10 = wsptr[0] + wsptr[4];
	tmp11 = wsptr[0] - wsptr[4];

	tmp13 = wsptr[2] + wsptr[6];
	tmp12 =(wsptr[2] - wsptr[6]) * 1.414213562f - tmp13;

	tmp0 = tmp10 + tmp13;
	tmp3 = tmp10 - tmp13;
	tmp1 = tmp11 + tmp12;
	tmp2 = tmp11 - tmp12;

   // Odd part
	z13 = wsptr[5] + wsptr[3];
	z10 = wsptr[5] - wsptr[3];
	z11 = wsptr[1] + wsptr[7];
	z12 = wsptr[1] - wsptr[7];