/lib/cximage-6.0/raw/libdcr.c

Large files files are truncated, but you can click here to view the full file

/*
  libdcr version 0.1.8.89 11/Jan/2009

  libdcr : copyright (C) 2007-2009, Davide Pizzolato

  based on dcraw.c -- Dave Coffin's raw photo decoder
  Copyright 1997-2009 by Dave Coffin, dcoffin a cybercom o net

  Covered code is provided under this license on an "as is" basis, without warranty
  of any kind, either expressed or implied, including, without limitation, warranties
  that the covered code is free of defects, merchantable, fit for a particular purpose
  or non-infringing. The entire risk as to the quality and performance of the covered
  code is with you. Should any covered code prove defective in any respect, you (not
  the initial developer or any other contributor) assume the cost of any necessary
  servicing, repair or correction. This disclaimer of warranty constitutes an essential
  part of this license. No use of any covered code is authorized hereunder except under
  this disclaimer.

  No license is required to download and use libdcr.  However,
  to lawfully redistribute libdcr, you must either (a) offer, at
  no extra charge, full source code for all executable files
  containing RESTRICTED functions, (b) distribute this code under
  the GPL Version 2 or later, (c) remove all RESTRICTED functions,
  re-implement them, or copy them from an earlier, unrestricted
  revision of dcraw.c, or (d) purchase a license from the author
  of dcraw.c.

  --------------------------------------------------------------------------------

  dcraw.c home page: http://cybercom.net/~dcoffin/dcraw/
  libdcr  home page: http://www.xdp.it/libdcr/

 */

#define _GNU_SOURCE
#define _USE_MATH_DEFINES
#include <ctype.h>
#include <errno.h>
#include <fcntl.h>
#include <float.h>
#include <limits.h>
#include <math.h>
#include <setjmp.h>
#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <time.h>
#include <sys/types.h>

#include "libdcr.h"

// XYZ from RGB
const double xyz_rgb[3][3] = {
	{ 0.412453, 0.357580, 0.180423 },
	{ 0.212671, 0.715160, 0.072169 },
	{ 0.019334, 0.119193, 0.950227 } };

const float d65_white[3] = { 0.950456f, 1.0f, 1.088754f };

static int   dcr_sfile_read(dcr_stream_obj *obj, void *buf, int size, int cnt);
static int   dcr_sfile_write(dcr_stream_obj *obj, void *buf, int size, int cnt);
static long  dcr_sfile_seek(dcr_stream_obj *obj, long offset, int origin);
static int   dcr_sfile_close(dcr_stream_obj *obj);
static char* dcr_sfile_gets(dcr_stream_obj *obj, char *string, int n);
static int   dcr_sfile_eof(dcr_stream_obj *obj);
static long  dcr_sfile_tell(dcr_stream_obj *obj);
static int   dcr_sfile_getc(dcr_stream_obj *obj);
static int   dcr_sfile_scanf(dcr_stream_obj *obj,const char *format, void* output);

static dcr_stream_ops dcr_stream_fileops = {
	dcr_sfile_read,
	dcr_sfile_write,
	dcr_sfile_seek,
	dcr_sfile_close,
	dcr_sfile_gets,
	dcr_sfile_eof,
	dcr_sfile_tell,
	dcr_sfile_getc,
	dcr_sfile_scanf
};

	/*
	NO_JPEG disables decoding of compressed Kodak DC120 files.
	NO_LCMS disables the "-p" option.
	*/
#ifndef NO_JPEG
  #ifdef WIN32
    #include "../jpeg/jpeglib.h"
  #else
    #include <jpeglib.h>
  #endif
#endif

#ifndef NO_LCMS
  #ifdef WIN32
    #include "../lcms/lcms.h"
  #else
    #include <lcms.h>
  #endif
#endif

#ifdef LOCALEDIR
#include <libintl.h>
#define _(String) gettext(String)
#else
#define _(String) (String)
#endif

#ifdef __CYGWIN__
#include <io.h>
#endif
#ifdef WIN32
#include <sys/utime.h>
#include <winsock2.h>
#pragma comment(lib, "ws2_32.lib")
#define snprintf _snprintf
#define strcasecmp _stricmp
#define strncasecmp _strnicmp
	typedef __int64 INT64;
	typedef unsigned __int64 UINT64;
#else
#include <unistd.h>
#include <utime.h>
#include <netinet/in.h>
	typedef long long INT64;
	typedef unsigned long long UINT64;
	#define __int64		long long
#endif

#ifdef LJPEG_DECODE
#error Please compile dcraw.c by itself.
#error Do not link it with ljpeg_decode.
#endif


/*
   In order to inline this calculation, I make the risky
   assumption that all filter patterns can be described
   by a repeating pattern of eight rows and two columns

   Do not use the FC or BAYER macros with the Leaf CatchLight,
   because its pattern is 16x16, not 2x8.

   Return values are either 0/1/2/3 = G/M/C/Y or 0/1/2/3 = R/G1/B/G2

	PowerShot 600	PowerShot A50	PowerShot Pro70	Pro90 & G1
	0xe1e4e1e4:	0x1b4e4b1e:	0x1e4b4e1b:	0xb4b4b4b4:

	  0 1 2 3 4 5	  0 1 2 3 4 5	  0 1 2 3 4 5	  0 1 2 3 4 5
	0 G M G M G M	0 C Y C Y C Y	0 Y C Y C Y C	0 G M G M G M
	1 C Y C Y C Y	1 M G M G M G	1 M G M G M G	1 Y C Y C Y C
	2 M G M G M G	2 Y C Y C Y C	2 C Y C Y C Y
	3 C Y C Y C Y	3 G M G M G M	3 G M G M G M
			4 C Y C Y C Y	4 Y C Y C Y C
	PowerShot A5	5 G M G M G M	5 G M G M G M
	0x1e4e1e4e:	6 Y C Y C Y C	6 C Y C Y C Y
			7 M G M G M G	7 M G M G M G
	  0 1 2 3 4 5
	0 C Y C Y C Y
	1 G M G M G M
	2 C Y C Y C Y
	3 M G M G M G

   All RGB cameras use one of these Bayer grids:

	0x16161616:	0x61616161:	0x49494949:	0x94949494:

	  0 1 2 3 4 5	  0 1 2 3 4 5	  0 1 2 3 4 5	  0 1 2 3 4 5
	0 B G B G B G	0 G R G R G R	0 G B G B G B	0 R G R G R G
	1 G R G R G R	1 B G B G B G	1 R G R G R G	1 G B G B G B
	2 B G B G B G	2 G R G R G R	2 G B G B G B	2 R G R G R G
	3 G R G R G R	3 B G B G B G	3 R G R G R G	3 G B G B G B
 */

#define FC(row,col) \
	(p->filters >> ((((row) << 1 & 14) + ((col) & 1)) << 1) & 3)

#define BAYER(row,col) \
	p->image[((row) >> p->shrink)*p->iwidth + ((col) >> p->shrink)][FC(row,col)]

#define BAYER2(row,col) \
	p->image[((row) >> p->shrink)*p->iwidth + ((col) >> p->shrink)][dcr_fc(p,row,col)]

int DCR_CLASS dcr_fc (DCRAW* p, int row, int col)
{
	static const char filter[16][16] =
	{ { 2,1,1,3,2,3,2,0,3,2,3,0,1,2,1,0 },
	{ 0,3,0,2,0,1,3,1,0,1,1,2,0,3,3,2 },
	{ 2,3,3,2,3,1,1,3,3,1,2,1,2,0,0,3 },
	{ 0,1,0,1,0,2,0,2,2,0,3,0,1,3,2,1 },
	{ 3,1,1,2,0,1,0,2,1,3,1,3,0,1,3,0 },
	{ 2,0,0,3,3,2,3,1,2,0,2,0,3,2,2,1 },
	{ 2,3,3,1,2,1,2,1,2,1,1,2,3,0,0,1 },
	{ 1,0,0,2,3,0,0,3,0,3,0,3,2,1,2,3 },
	{ 2,3,3,1,1,2,1,0,3,2,3,0,2,3,1,3 },
	{ 1,0,2,0,3,0,3,2,0,1,1,2,0,1,0,2 },
	{ 0,1,1,3,3,2,2,1,1,3,3,0,2,1,3,2 },
	{ 2,3,2,0,0,1,3,0,2,0,1,2,3,0,1,0 },
	{ 1,3,1,2,3,2,3,2,0,2,0,1,1,0,3,0 },
	{ 0,2,0,3,1,0,0,1,1,3,3,2,3,2,2,1 },
	{ 2,1,3,2,3,1,2,1,0,3,0,2,0,2,0,2 },
	{ 0,3,1,0,0,2,0,3,2,1,3,1,1,3,1,3 } };

	if (p->filters != 1) return FC(row,col);
	return filter[(row+p->top_margin) & 15][(col+p->left_margin) & 15];
}

#ifndef __GLIBC__
char *dcr_memmem (char *haystack, size_t haystacklen,
	char *needle, size_t needlelen)
{
	char *c;
	for (c = haystack; c <= haystack + haystacklen - needlelen; c++)
		if (!memcmp (c, needle, needlelen))
			return c;
		return 0;
}
#define memmem dcr_memmem
#endif

void DCR_CLASS dcr_merror (DCRAW* p, void *ptr, char *where)
{
	if (ptr) return;
	if (p->sz_error){
		sprintf (p->sz_error,_("%s: Out of memory in %s\n"), p->ifname, where);
	} else {
		fprintf (stderr,_("%s: Out of memory in %s\n"), p->ifname, where);
	}
	longjmp (p->failure, 1);
}

void DCR_CLASS dcr_derror(DCRAW* p)
{
	if (!p->data_error) {
		fprintf (stderr, "%s: ", p->ifname);
		if (dcr_feof(p->obj_))
			fprintf (stderr,_("Unexpected end of file\n"));
		else
			fprintf (stderr,_("Corrupt data near 0x%llx\n"), (INT64) dcr_ftell(p->obj_));
	}
	p->data_error = 1;
}

ushort DCR_CLASS dcr_sget2 (DCRAW* p,uchar *s)
{
	if (p->order == 0x4949)		/* "II" means little-endian */
		return s[0] | s[1] << 8;
	else				/* "MM" means big-endian */
		return s[0] << 8 | s[1];
}

ushort DCR_CLASS dcr_get2(DCRAW* p)
{
	uchar str[2] = { 0xff,0xff };
	dcr_fread(p->obj_, str, 1, 2);
	return dcr_sget2(p, str);
}

unsigned DCR_CLASS dcr_sget4 (DCRAW* p,uchar *s)
{
	if (p->order == 0x4949)
		return s[0] | s[1] << 8 | s[2] << 16 | s[3] << 24;
	else
		return s[0] << 24 | s[1] << 16 | s[2] << 8 | s[3];
}
#define dcr_sget4(p, s) dcr_sget4(p, (uchar *)s)

unsigned DCR_CLASS dcr_get4(DCRAW* p)
{
	uchar str[4] = { 0xff,0xff,0xff,0xff };
	dcr_fread(p->obj_, str, 1, 4);
	return dcr_sget4(p, str);
}

unsigned DCR_CLASS dcr_getint (DCRAW* p,int type)
{
	return type == 3 ? dcr_get2(p) : dcr_get4(p);
}

float DCR_CLASS dcr_int_to_float (int i)
{
	union { int i; float f; } u;
	u.i = i;
	return u.f;
}

double DCR_CLASS dcr_getreal (DCRAW* p,int type)
{
	union { char c[8]; double d; } u;
	int i, rev;

	switch (type) {
	case 3: return (unsigned short) dcr_get2(p);
	case 4: return (unsigned int) dcr_get4(p);
	case 5:  u.d = (unsigned int) dcr_get4(p);
		return u.d / (unsigned int) dcr_get4(p);
	case 8: return (signed short) dcr_get2(p);
	case 9: return (signed int) dcr_get4(p);
	case 10: u.d = (signed int) dcr_get4(p);
		return u.d / (signed int) dcr_get4(p);
	case 11: return dcr_int_to_float (dcr_get4(p));
	case 12:
		rev = 7 * ((p->order == 0x4949) == (ntohs(0x1234) == 0x1234));
		for (i=0; i < 8; i++)
			u.c[i ^ rev] = dcr_fgetc(p->obj_);
		return u.d;
	default: return dcr_fgetc(p->obj_);
	}
}

void DCR_CLASS dcr_read_shorts (DCRAW* p,ushort *pixel, int count)
{
	if ((int)dcr_fread(p->obj_, pixel, 2, count) < count) dcr_derror(p);
	if ((p->order == 0x4949) == (ntohs(0x1234) == 0x1234))
		_swab ((char*)pixel, (char*)pixel, count*2);
}

void DCR_CLASS dcr_canon_black (DCRAW* p, double dark[2])
{
  int c, diff, row, col;

  if (p->raw_width < p->width+4) return;
  FORC(2) dark[c] /= (p->raw_width-p->width-2) * p->height >> 1;
  if ((diff = (int)(dark[0] - dark[1])))
    for (row=0; row < p->height; row++)
      for (col=1; col < p->width; col+=2)
	BAYER(row,col) += diff;
  dark[1] += diff;
  p->black = (unsigned int)((dark[0] + dark[1] + 1) / 2);
}

void DCR_CLASS dcr_canon_600_fixed_wb (DCRAW* p,int temp)
{
	static const short mul[4][5] = {
		{  667, 358,397,565,452 },
		{  731, 390,367,499,517 },
		{ 1119, 396,348,448,537 },
		{ 1399, 485,431,508,688 } };
	int lo, hi, i;
	float frac=0;

	for (lo=4; --lo; )
		if (*mul[lo] <= temp) break;
	for (hi=0; hi < 3; hi++)
		if (*mul[hi] >= temp) break;
	if (lo != hi)
		frac = (float) (temp - *mul[lo]) / (*mul[hi] - *mul[lo]);
	for (i=1; i < 5; i++)
		p->pre_mul[i-1] = 1 / (frac * mul[hi][i] + (1-frac) * mul[lo][i]);
}

/* Return values:  0 = p->white  1 = near p->white  2 = not p->white */
int DCR_CLASS dcr_canon_600_color (DCRAW* p,int ratio[2], int mar)
{
	int clipped=0, target, miss;

	if (p->flash_used) {
		if (ratio[1] < -104)
		{ ratio[1] = -104; clipped = 1; }
		if (ratio[1] >   12)
		{ ratio[1] =   12; clipped = 1; }
	} else {
		if (ratio[1] < -264 || ratio[1] > 461) return 2;
		if (ratio[1] < -50)
		{ ratio[1] = -50; clipped = 1; }
		if (ratio[1] > 307)
		{ ratio[1] = 307; clipped = 1; }
	}
	target = p->flash_used || ratio[1] < 197
		? -38 - (398 * ratio[1] >> 10)
		: -123 + (48 * ratio[1] >> 10);
	if (target - mar <= ratio[0] &&
		target + 20  >= ratio[0] && !clipped) return 0;
	miss = target - ratio[0];
	if (abs(miss) >= mar*4) return 2;
	if (miss < -20) miss = -20;
	if (miss > mar) miss = mar;
	ratio[0] = target - miss;
	return 1;
}

void DCR_CLASS dcr_canon_600_auto_wb(DCRAW* p)
{
	int mar, row, col, i, j, st, count[] = { 0,0 };
	int test[8], total[2][8], ratio[2][2], stat[2];

	memset (&total, 0, sizeof total);
	i = (int)(p->canon_ev + 0.5);
	if      (i < 10) mar = 150;
	else if (i > 12) mar = 20;
	else mar = 280 - 20 * i;
	if (p->flash_used) mar = 80;
	for (row=14; row < p->height-14; row+=4)
		for (col=10; col < p->width; col+=2) {
			for (i=0; i < 8; i++)
				test[(i & 4) + FC(row+(i >> 1),col+(i & 1))] =
				BAYER(row+(i >> 1),col+(i & 1));
			for (i=0; i < 8; i++)
				if (test[i] < 150 || test[i] > 1500) goto next;
			for (i=0; i < 4; i++)
				if (abs(test[i] - test[i+4]) > 50) goto next;
			for (i=0; i < 2; i++) {
				for (j=0; j < 4; j+=2)
					ratio[i][j >> 1] = ((test[i*4+j+1]-test[i*4+j]) << 10) / test[i*4+j];
				stat[i] = dcr_canon_600_color (p,ratio[i], mar);
			}
			if ((st = stat[0] | stat[1]) > 1) goto next;
			for (i=0; i < 2; i++)
				if (stat[i])
					for (j=0; j < 2; j++)
						test[i*4+j*2+1] = test[i*4+j*2] * (0x400 + ratio[i][j]) >> 10;
			for (i=0; i < 8; i++)
				total[st][i] += test[i];
			count[st]++;
next: ;
		}
	if (count[0] | count[1]) {
		st = count[0]*200 < count[1];
		for (i=0; i < 4; i++)
			p->pre_mul[i] = 1.0f / (total[st][i] + total[st][i+4]);
	}
}

void DCR_CLASS dcr_canon_600_coeff(DCRAW* p)
{
	static const short table[6][12] = {
		{ -190,702,-1878,2390,   1861,-1349,905,-393, -432,944,2617,-2105  },
		{ -1203,1715,-1136,1648, 1388,-876,267,245,  -1641,2153,3921,-3409 },
		{ -615,1127,-1563,2075,  1437,-925,509,3,     -756,1268,2519,-2007 },
		{ -190,702,-1886,2398,   2153,-1641,763,-251, -452,964,3040,-2528  },
		{ -190,702,-1878,2390,   1861,-1349,905,-393, -432,944,2617,-2105  },
		{ -807,1319,-1785,2297,  1388,-876,769,-257,  -230,742,2067,-1555  } };
	int t=0, i, c;
	float mc, yc;

	mc = p->pre_mul[1] / p->pre_mul[2];
	yc = p->pre_mul[3] / p->pre_mul[2];
	if (mc > 1 && mc <= 1.28 && yc < 0.8789) t=1;
	if (mc > 1.28 && mc <= 2) {
		if  (yc < 0.8789) t=3;
		else if (yc <= 2) t=4;
	}
	if (p->flash_used) t=5;
	for (p->raw_color = i=0; i < 3; i++)
		FORCC(p) p->rgb_cam[i][c] = table[t][i*4 + c] / 1024.0f;
}

void DCR_CLASS dcr_canon_600_load_raw(DCRAW* p)
{
	uchar  data[1120], *dp;
	ushort pixel[896], *pix;
	int irow, row, col, val;
	static const short mul[4][2] =
	{ { 1141,1145 }, { 1128,1109 }, { 1178,1149 }, { 1128,1109 } };

	for (irow=row=0; irow < p->height; irow++) {
		if ((long)dcr_fread(p->obj_, data, 1, p->raw_width*5/4) < (long)(p->raw_width*5/4)) dcr_derror(p);
		for (dp=data, pix=pixel; dp < data+1120; dp+=10, pix+=8) {
			pix[0] = (dp[0] << 2) + (dp[1] >> 6    );
			pix[1] = (dp[2] << 2) + (dp[1] >> 4 & 3);
			pix[2] = (dp[3] << 2) + (dp[1] >> 2 & 3);
			pix[3] = (dp[4] << 2) + (dp[1]      & 3);
			pix[4] = (dp[5] << 2) + (dp[9]      & 3);
			pix[5] = (dp[6] << 2) + (dp[9] >> 2 & 3);
			pix[6] = (dp[7] << 2) + (dp[9] >> 4 & 3);
			pix[7] = (dp[8] << 2) + (dp[9] >> 6    );
		}
		for (col=0; col < p->width; col++)
			BAYER(row,col) = pixel[col];
		for (col=p->width; col < p->raw_width; col++)
			p->black += pixel[col];
		if ((row+=2) > p->height) row = 1;
	}
	if (p->raw_width > p->width)
		p->black = p->black / ((p->raw_width - p->width) * p->height) - 4;
	for (row=0; row < p->height; row++)
		for (col=0; col < p->width; col++) {
			if ((val = BAYER(row,col) - p->black) < 0) val = 0;
			val = val * mul[row & 3][col & 1] >> 9;
			BAYER(row,col) = val;
		}
	dcr_canon_600_fixed_wb(p,1311);
	dcr_canon_600_auto_wb(p);
	dcr_canon_600_coeff(p);
	p->maximum = (0x3ff - p->black) * 1109 >> 9;
	p->black = 0;
}

void DCR_CLASS dcr_remove_zeroes(DCRAW* p)
{
	unsigned row, col, tot, n, r, c;

	for (row=0; row < p->height; row++)
		for (col=0; col < p->width; col++)
			if (BAYER(row,col) == 0) {
				tot = n = 0;
				for (r = row-2; r <= row+2; r++)
					for (c = col-2; c <= col+2; c++)
						if (r < p->height && c < p->width &&
							FC(r,c) == FC(row,col) && BAYER(r,c))
							tot += (n++,BAYER(r,c));
						if (n) BAYER(row,col) = tot/n;
			}
}

int DCR_CLASS dcr_canon_s2is(DCRAW* p)
{
	unsigned row;

	for (row=0; row < 100; row++) {
		dcr_fseek(p->obj_, row*3340 + 3284, SEEK_SET);
		if (dcr_fgetc(p->obj_) > 15) return 1;
	}
	return 0;
}

void DCR_CLASS dcr_canon_a5_load_raw(DCRAW* p)
{
	ushort data[2565], *dp, pixel;
	int vbits=0, buf=0, row, col, bc=0;

	p->order = 0x4949;
	for (row=-p->top_margin; row < p->raw_height-p->top_margin; row++) {
		dcr_read_shorts (p,dp=data, p->raw_width * 10 / 16);
		for (col=-p->left_margin; col < p->raw_width-p->left_margin; col++) {
			if ((vbits -= 10) < 0)
				buf = (vbits += 16, (buf << 16) + *dp++);
			pixel = buf >> vbits & 0x3ff;
			if ((unsigned) row < p->height && (unsigned) col < p->width)
				BAYER(row,col) = pixel;
			else if (col > 1-p->left_margin && col != p->width)
				p->black += (bc++,pixel);
		}
	}
	if (bc) p->black /= bc;
	p->maximum = 0x3ff;
	if (p->raw_width > 1600) dcr_remove_zeroes(p);
}

/*
dcr_getbits(p, -1) initializes the buffer
dcr_getbits(p, n) where 0 <= n <= 25 returns an n-bit integer
*/
unsigned DCR_CLASS dcr_getbits (DCRAW* p, int nbits)
{
	unsigned c;

	if (nbits == -1)
		return p->getbits_bitbuf = p->getbits_vbits = p->getbits_reset = 0;
	if (nbits == 0 || p->getbits_reset) return 0;
	while (p->getbits_vbits < nbits) {
		if ((c = dcr_fgetc(p->obj_)) == EOF) dcr_derror(p);
		if ((p->getbits_reset = p->zero_after_ff && c == 0xff && dcr_fgetc(p->obj_))) return 0;
		p->getbits_bitbuf = (p->getbits_bitbuf << 8) + (uchar) c;
		p->getbits_vbits += 8;
	}
	p->getbits_vbits -= nbits;
	return p->getbits_bitbuf << (32-nbits-p->getbits_vbits) >> (32-nbits);
}

void DCR_CLASS dcr_init_decoder(DCRAW* p)
{
	memset (p->first_decode, 0, sizeof p->first_decode);
	p->free_decode = p->first_decode;
}

/*
   Construct a decode tree according the specification in *source.
   The first 16 bytes specify how many codes should be 1-bit, 2-bit
   3-bit, etc.  Bytes after that are the leaf values.

   For example, if the source is

    { 0,1,4,2,3,1,2,0,0,0,0,0,0,0,0,0,
      0x04,0x03,0x05,0x06,0x02,0x07,0x01,0x08,0x09,0x00,0x0a,0x0b,0xff  },

   then the code is

	00		0x04
	010		0x03
	011		0x05
	100		0x06
	101		0x02
	1100		0x07
	1101		0x01
	11100		0x08
	11101		0x09
	11110		0x00
	111110		0x0a
	1111110		0x0b
	1111111		0xff
 */
uchar * DCR_CLASS dcr_make_decoder (DCRAW* p, const uchar *source, int level)
{
	struct dcr_decode *cur;
	int i, next;

	if (level==0) p->make_decoder_leaf=0;
	cur = p->free_decode++;
	if (p->free_decode > p->first_decode+2048) {
		fprintf (stderr,_("%s: decoder table overflow\n"), p->ifname);
		longjmp (p->failure, 2);
	}
	for (i=next=0; i <= p->make_decoder_leaf && next < 16; )
		i += source[next++];
	if (i > p->make_decoder_leaf) {
		if (level < next) {
			cur->branch[0] = p->free_decode;
			dcr_make_decoder (p, source, level+1);
			cur->branch[1] = p->free_decode;
			dcr_make_decoder (p, source, level+1);
		} else
			cur->leaf = source[16 + p->make_decoder_leaf++];
	}
	return (uchar *) source + 16 + p->make_decoder_leaf;
}

void DCR_CLASS dcr_crw_init_tables (DCRAW* p, unsigned table)
{
	static const uchar first_tree[3][29] = {
		{ 0,1,4,2,3,1,2,0,0,0,0,0,0,0,0,0,
			0x04,0x03,0x05,0x06,0x02,0x07,0x01,0x08,0x09,0x00,0x0a,0x0b,0xff  },
		{ 0,2,2,3,1,1,1,1,2,0,0,0,0,0,0,0,
		0x03,0x02,0x04,0x01,0x05,0x00,0x06,0x07,0x09,0x08,0x0a,0x0b,0xff  },
		{ 0,0,6,3,1,1,2,0,0,0,0,0,0,0,0,0,
		0x06,0x05,0x07,0x04,0x08,0x03,0x09,0x02,0x00,0x0a,0x01,0x0b,0xff  },
	};
	static const uchar second_tree[3][180] = {
		{ 0,2,2,2,1,4,2,1,2,5,1,1,0,0,0,139,
		0x03,0x04,0x02,0x05,0x01,0x06,0x07,0x08,
		0x12,0x13,0x11,0x14,0x09,0x15,0x22,0x00,0x21,0x16,0x0a,0xf0,
		0x23,0x17,0x24,0x31,0x32,0x18,0x19,0x33,0x25,0x41,0x34,0x42,
		0x35,0x51,0x36,0x37,0x38,0x29,0x79,0x26,0x1a,0x39,0x56,0x57,
		0x28,0x27,0x52,0x55,0x58,0x43,0x76,0x59,0x77,0x54,0x61,0xf9,
		0x71,0x78,0x75,0x96,0x97,0x49,0xb7,0x53,0xd7,0x74,0xb6,0x98,
		0x47,0x48,0x95,0x69,0x99,0x91,0xfa,0xb8,0x68,0xb5,0xb9,0xd6,
		0xf7,0xd8,0x67,0x46,0x45,0x94,0x89,0xf8,0x81,0xd5,0xf6,0xb4,
		0x88,0xb1,0x2a,0x44,0x72,0xd9,0x87,0x66,0xd4,0xf5,0x3a,0xa7,
		0x73,0xa9,0xa8,0x86,0x62,0xc7,0x65,0xc8,0xc9,0xa1,0xf4,0xd1,
		0xe9,0x5a,0x92,0x85,0xa6,0xe7,0x93,0xe8,0xc1,0xc6,0x7a,0x64,
		0xe1,0x4a,0x6a,0xe6,0xb3,0xf1,0xd3,0xa5,0x8a,0xb2,0x9a,0xba,
		0x84,0xa4,0x63,0xe5,0xc5,0xf3,0xd2,0xc4,0x82,0xaa,0xda,0xe4,
		0xf2,0xca,0x83,0xa3,0xa2,0xc3,0xea,0xc2,0xe2,0xe3,0xff,0xff  },
		{ 0,2,2,1,4,1,4,1,3,3,1,0,0,0,0,140,
		0x02,0x03,0x01,0x04,0x05,0x12,0x11,0x06,
		0x13,0x07,0x08,0x14,0x22,0x09,0x21,0x00,0x23,0x15,0x31,0x32,
		0x0a,0x16,0xf0,0x24,0x33,0x41,0x42,0x19,0x17,0x25,0x18,0x51,
		0x34,0x43,0x52,0x29,0x35,0x61,0x39,0x71,0x62,0x36,0x53,0x26,
		0x38,0x1a,0x37,0x81,0x27,0x91,0x79,0x55,0x45,0x28,0x72,0x59,
		0xa1,0xb1,0x44,0x69,0x54,0x58,0xd1,0xfa,0x57,0xe1,0xf1,0xb9,
		0x49,0x47,0x63,0x6a,0xf9,0x56,0x46,0xa8,0x2a,0x4a,0x78,0x99,
		0x3a,0x75,0x74,0x86,0x65,0xc1,0x76,0xb6,0x96,0xd6,0x89,0x85,
		0xc9,0xf5,0x95,0xb4,0xc7,0xf7,0x8a,0x97,0xb8,0x73,0xb7,0xd8,
		0xd9,0x87,0xa7,0x7a,0x48,0x82,0x84,0xea,0xf4,0xa6,0xc5,0x5a,
		0x94,0xa4,0xc6,0x92,0xc3,0x68,0xb5,0xc8,0xe4,0xe5,0xe6,0xe9,
		0xa2,0xa3,0xe3,0xc2,0x66,0x67,0x93,0xaa,0xd4,0xd5,0xe7,0xf8,
		0x88,0x9a,0xd7,0x77,0xc4,0x64,0xe2,0x98,0xa5,0xca,0xda,0xe8,
		0xf3,0xf6,0xa9,0xb2,0xb3,0xf2,0xd2,0x83,0xba,0xd3,0xff,0xff  },
		{ 0,0,6,2,1,3,3,2,5,1,2,2,8,10,0,117,
		0x04,0x05,0x03,0x06,0x02,0x07,0x01,0x08,
		0x09,0x12,0x13,0x14,0x11,0x15,0x0a,0x16,0x17,0xf0,0x00,0x22,
		0x21,0x18,0x23,0x19,0x24,0x32,0x31,0x25,0x33,0x38,0x37,0x34,
		0x35,0x36,0x39,0x79,0x57,0x58,0x59,0x28,0x56,0x78,0x27,0x41,
		0x29,0x77,0x26,0x42,0x76,0x99,0x1a,0x55,0x98,0x97,0xf9,0x48,
		0x54,0x96,0x89,0x47,0xb7,0x49,0xfa,0x75,0x68,0xb6,0x67,0x69,
		0xb9,0xb8,0xd8,0x52,0xd7,0x88,0xb5,0x74,0x51,0x46,0xd9,0xf8,
		0x3a,0xd6,0x87,0x45,0x7a,0x95,0xd5,0xf6,0x86,0xb4,0xa9,0x94,
		0x53,0x2a,0xa8,0x43,0xf5,0xf7,0xd4,0x66,0xa7,0x5a,0x44,0x8a,
		0xc9,0xe8,0xc8,0xe7,0x9a,0x6a,0x73,0x4a,0x61,0xc7,0xf4,0xc6,
		0x65,0xe9,0x72,0xe6,0x71,0x91,0x93,0xa6,0xda,0x92,0x85,0x62,
		0xf3,0xc5,0xb2,0xa4,0x84,0xba,0x64,0xa5,0xb3,0xd2,0x81,0xe5,
		0xd3,0xaa,0xc4,0xca,0xf2,0xb1,0xe4,0xd1,0x83,0x63,0xea,0xc3,
		0xe2,0x82,0xf1,0xa3,0xc2,0xa1,0xc1,0xe3,0xa2,0xe1,0xff,0xff  }
	};
	if (table > 2) table = 2;
	dcr_init_decoder(p);
	dcr_make_decoder (p, first_tree[table], 0);
	p->second_decode = p->free_decode;
	dcr_make_decoder (p, second_tree[table], 0);
}

/*
Return 0 if the image starts with compressed data,
1 if it starts with uncompressed low-order bits.

  In Canon compressed data, 0xff is always followed by 0x00.
*/
int DCR_CLASS dcr_canon_has_lowbits(DCRAW* p)
{
	uchar test[0x4000];
	int ret=1, i;

	dcr_fseek(p->obj_, 0, SEEK_SET);
	dcr_fread(p->obj_, test, 1, sizeof test);
	for (i=540; i < sizeof test - 1; i++)
		if (test[i] == 0xff) {
			if (test[i+1]) return 1;
			ret=0;
		}
	return ret;
}

void DCR_CLASS dcr_canon_compressed_load_raw(DCRAW* p)
{
	ushort *pixel, *prow;
	int nblocks, lowbits, i, row, r, col, save, val;
	unsigned irow, icol;
	struct dcr_decode *decode, *dindex;
	int block, diffbuf[64], leaf, len, diff, carry=0, pnum=0, base[2];
	double dark[2] = { 0,0 };
	uchar c;

	dcr_crw_init_tables (p, p->tiff_compress);
	pixel = (ushort *) calloc (p->raw_width*8, sizeof *pixel);
	dcr_merror (p, pixel, "canon_compressed_load_raw()");
	lowbits = dcr_canon_has_lowbits(p);
	if (!lowbits) p->maximum = 0x3ff;
	dcr_fseek(p->obj_, 540 + lowbits*p->raw_height*p->raw_width/4, SEEK_SET);
	p->zero_after_ff = 1;
	dcr_getbits(p, -1);
	for (row=0; row < p->raw_height; row+=8) {
		nblocks = MIN (8, p->raw_height-row) * p->raw_width >> 6;
		for (block=0; block < nblocks; block++) {
			memset (diffbuf, 0, sizeof diffbuf);
			decode = p->first_decode;
			for (i=0; i < 64; i++ ) {
				for (dindex=decode; dindex->branch[0]; )
					dindex = dindex->branch[dcr_getbits(p, 1)];
				leaf = dindex->leaf;
				decode = p->second_decode;
				if (leaf == 0 && i) break;
				if (leaf == 0xff) continue;
				i  += leaf >> 4;
				len = leaf & 15;
				if (len == 0) continue;
				diff = dcr_getbits(p, len);
				if ((diff & (1 << (len-1))) == 0)
					diff -= (1 << len) - 1;
				if (i < 64) diffbuf[i] = diff;
			}
			diffbuf[0] += carry;
			carry = diffbuf[0];
			for (i=0; i < 64; i++ ) {
				if (pnum++ % p->raw_width == 0)
					base[0] = base[1] = 512;
				if ((pixel[(block << 6) + i] = base[i & 1] += diffbuf[i]) >> 10)
					dcr_derror(p);
			}
		}
		if (lowbits) {
			save = dcr_ftell(p->obj_);
			dcr_fseek(p->obj_, 26 + row*p->raw_width/4, SEEK_SET);
			for (prow=pixel, i=0; i < p->raw_width*2; i++) {
				c = dcr_fgetc(p->obj_);
				for (r=0; r < 8; r+=2, prow++) {
					val = (*prow << 2) + ((c >> r) & 3);
					if (p->raw_width == 2672 && val < 512) val += 2;
					*prow = val;
				}
			}
			dcr_fseek(p->obj_, save, SEEK_SET);
		}
		for (r=0; r < 8; r++) {
			irow = row - p->top_margin + r;
			if (irow >= p->height) continue;
			for (col=0; col < p->raw_width; col++) {
				icol = col - p->left_margin;
				if (icol < p->width)
					BAYER(irow,icol) = pixel[r*p->raw_width+col];
				else if (col > 1)
					dark[icol & 1] += pixel[r*p->raw_width+col];
			}
		}
	}
	free (pixel);
	dcr_canon_black (p, dark);
}

/*
Not a full implementation of Lossless JPEG, just
enough to decode Canon, Kodak and Adobe DNG images.
*/

int DCR_CLASS dcr_ljpeg_start (DCRAW* p, struct dcr_jhead *jh, int info_only)
{
	int c, tag, len;
	uchar data[0x10000], *dp;

	dcr_init_decoder(p);
	memset (jh, 0, sizeof *jh);
	FORC(6) jh->huff[c] = p->free_decode;
	jh->restart = INT_MAX;
	dcr_fread(p->obj_, data, 2, 1);
	if (data[1] != 0xd8) return 0;
	do {
		dcr_fread(p->obj_, data, 2, 2);
		tag =  data[0] << 8 | data[1];
		len = (data[2] << 8 | data[3]) - 2;
		if (tag <= 0xff00) return 0;
		dcr_fread(p->obj_, data, 1, len);
		switch (tag) {
		case 0xffc3:
			jh->sraw = ((data[7] >> 4) * (data[7] & 15) - 1) & 3;
		case 0xffc0:
			jh->bits = data[0];
			jh->high = data[1] << 8 | data[2];
			jh->wide = data[3] << 8 | data[4];
			jh->clrs = data[5] + jh->sraw;
			if (len == 9 && !p->dng_version) dcr_fgetc(p->obj_);
			break;
		case 0xffc4:
			if (info_only) break;
			for (dp = data; dp < data+len && *dp < 4; ) {
				jh->huff[*dp] = p->free_decode;
				dp = dcr_make_decoder (p, ++dp, 0);
			}
			break;
		case 0xffda:
			jh->psv = data[1+data[0]*2];
			jh->bits -= data[3+data[0]*2] & 15;
			break;
		case 0xffdd:
			jh->restart = data[0] << 8 | data[1];
		}
	} while (tag != 0xffda);
	if (info_only) return 1;
	if (jh->sraw) {
		FORC(4)        jh->huff[2+c] = jh->huff[1];
		FORC(jh->sraw) jh->huff[1+c] = jh->huff[0];
	}
	jh->row = (ushort *) calloc (jh->wide*jh->clrs, 4);
	dcr_merror (p, jh->row, "dcr_ljpeg_start()");
	return p->zero_after_ff = 1;
}

int DCR_CLASS dcr_ljpeg_diff (DCRAW* p, struct dcr_decode *dindex)
{
	int len, diff;

	while (dindex->branch[0])
		dindex = dindex->branch[dcr_getbits(p, 1)];
	len = dindex->leaf;
	if (len == 16 && (!p->dng_version || p->dng_version >= 0x1010000))
		return -32768;
	diff = dcr_getbits(p, len);
	if ((diff & (1 << (len-1))) == 0)
		diff -= (1 << len) - 1;
	return diff;
}

ushort * DCR_CLASS dcr_ljpeg_row (DCRAW* p, int jrow, struct dcr_jhead *jh)
{
	int col, c, diff, pred, spred=0;
	ushort mark=0, *row[3];

	if (jrow * jh->wide % jh->restart == 0) {
		FORC(6) jh->vpred[c] = 1 << (jh->bits-1);
		if (jrow)
			do mark = (mark << 8) + (c = dcr_fgetc(p->obj_));
			while (c != EOF && mark >> 4 != 0xffd);
			dcr_getbits(p, -1);
	}
	FORC3 row[c] = jh->row + jh->wide*jh->clrs*((jrow+c) & 1);
	for (col=0; col < jh->wide; col++)
		FORC(jh->clrs) {
			diff = dcr_ljpeg_diff (p, jh->huff[c]);
			if (jh->sraw && c <= jh->sraw && (col | c))
				pred = spred;
			else if (col) pred = row[0][-jh->clrs];
			else	    pred = (jh->vpred[c] += diff) - diff;
			if (jrow && col) switch (jh->psv) {
				case 1:	break;
				case 2: pred = row[1][0];					break;
				case 3: pred = row[1][-jh->clrs];				break;
				case 4: pred = pred +   row[1][0] - row[1][-jh->clrs];		break;
				case 5: pred = pred + ((row[1][0] - row[1][-jh->clrs]) >> 1);	break;
				case 6: pred = row[1][0] + ((pred - row[1][-jh->clrs]) >> 1);	break;
				case 7: pred = (pred + row[1][0]) >> 1;				break;
				default: pred = 0;
			}
			if ((**row = pred + diff) >> jh->bits) dcr_derror(p);
			if (c <= jh->sraw) spred = **row;
			row[0]++; row[1]++;
		}
		return row[2];
}

void DCR_CLASS dcr_lossless_jpeg_load_raw(DCRAW* p)
{
	int jwide, jrow, jcol, val, jidx, i, j, row=0, col=0;
	double dark[2] = { 0,0 };
	struct dcr_jhead jh;
	int min=INT_MAX;
	ushort *rp;

	if (!dcr_ljpeg_start (p,&jh, 0)) return;
	jwide = jh.wide * jh.clrs;

	for (jrow=0; jrow < jh.high; jrow++) {
		rp = dcr_ljpeg_row (p, jrow, &jh);
		for (jcol=0; jcol < jwide; jcol++) {
			val = *rp++;
			if (jh.bits <= 12)
				val = p->curve[val & 0xfff];
			if (p->cr2_slice[0]) {
				jidx = jrow*jwide + jcol;
				i = jidx / (p->cr2_slice[1]*jh.high);
				if ((j = i >= p->cr2_slice[0]))
					i  = p->cr2_slice[0];
				jidx -= i * (p->cr2_slice[1]*jh.high);
				row = jidx / p->cr2_slice[1+j];
				col = jidx % p->cr2_slice[1+j] + i*p->cr2_slice[1];
			}
			if (p->raw_width == 3984 && (col -= 2) < 0)
				col += (row--,p->raw_width);
			if ((unsigned) (row-p->top_margin) < p->height) {
				if ((unsigned) (col-p->left_margin) < p->width) {
					BAYER(row-p->top_margin,col-p->left_margin) = val;
					if (min > val) min = val;
				} else if (col > 1)
					dark[(col-p->left_margin) & 1] += val;
			}
			if (++col >= p->raw_width)
				col = (row++,0);
		}
	}
	free (jh.row);
	dcr_canon_black (p, dark);
	if (!strcasecmp(p->make,"KODAK"))
		p->black = min;
}

void DCR_CLASS dcr_canon_sraw_load_raw(DCRAW* p)
{
	struct dcr_jhead jh;
	short *rp=0, (*ip)[4];
	int jwide, slice, scol, ecol, row, col, jrow=0, jcol=0, pix[3], c;
	
	if (!dcr_ljpeg_start (p, &jh, 0)) return;
	jwide = (jh.wide >>= 1) * jh.clrs;
	
	for (ecol=slice=0; slice <= p->cr2_slice[0]; slice++) {
		scol = ecol;
		ecol += p->cr2_slice[1] * 2 / jh.clrs;
		if (!p->cr2_slice[0] || ecol > p->raw_width-1) ecol = p->raw_width & -2;
		for (row=0; row < p->height; row += (jh.clrs >> 1) - 1) {
			ip = (short (*)[4]) p->image + row*p->width;
			for (col=scol; col < ecol; col+=2, jcol+=jh.clrs) {
				if ((jcol %= jwide) == 0)
					rp = (short *) dcr_ljpeg_row (p, jrow++, &jh);
				if (col >= p->width) continue;
				FORC (jh.clrs-2)
					ip[col + (c >> 1)*p->width + (c & 1)][0] = rp[jcol+c];
				ip[col][1] = rp[jcol+jh.clrs-2] - 16384;
				ip[col][2] = rp[jcol+jh.clrs-1] - 16384;
			}
		}
	}
	ip = (short (*)[4]) p->image;
	rp = ip[0];
	for (row=0; row < p->height; row++, ip+=p->width) {
		if (row & (jh.sraw >> 1))
			for (col=0; col < p->width; col+=2)
				for (c=1; c < 3; c++)
					if (row == p->height-1)
						ip[col][c] =  ip[col-p->width][c];
					else ip[col][c] = (ip[col-p->width][c] + ip[col+p->width][c] + 1) >> 1;
					for (col=1; col < p->width; col+=2)
						for (c=1; c < 3; c++)
							if (col == p->width-1)
								ip[col][c] =  ip[col-1][c];
							else ip[col][c] = (ip[col-1][c] + ip[col+1][c] + 1) >> 1;
	}
	for ( ; rp < ip[0]; rp+=4) {
		if (p->unique_id < 0x80000200) {
			pix[0] = rp[0] + rp[2] - 512;
			pix[2] = rp[0] + rp[1] - 512;
			pix[1] = rp[0] + ((-778*rp[1] - (rp[2] << 11)) >> 12) - 512;
		} else {
			rp[1] += jh.sraw+1;
			rp[2] += jh.sraw+1;
			pix[0] = rp[0] + ((  200*rp[1] + 22929*rp[2]) >> 12);
			pix[1] = rp[0] + ((-5640*rp[1] - 11751*rp[2]) >> 12);
			pix[2] = rp[0] + ((29040*rp[1] -   101*rp[2]) >> 12);
		}
		FORC3 rp[c] = CLIP(pix[c] * p->sraw_mul[c] >> 10);
	}
	free (jh.row);
	p->maximum = 0x3fff;
}

void DCR_CLASS dcr_adobe_copy_pixel (DCRAW* p, int row, int col, ushort **rp)
{
	unsigned r, c;

	r = row -= p->top_margin;
	c = col -= p->left_margin;
	if (p->is_raw == 2 && p->opt.shot_select) (*rp)++;
	if (p->filters) {
		if (p->fuji_width) {
			r = row + p->fuji_width - 1 - (col >> 1);
			c = row + ((col+1) >> 1);
		}
		if (r < p->height && c < p->width)
			BAYER(r,c) = **rp < 0x1000 ? p->curve[**rp] : **rp;
		*rp += p->is_raw;
	} else {
		if (r < p->height && c < p->width)
			FORC(p->tiff_samples)
				p->image[row*p->width+col][c] = (*rp)[c] < 0x1000 ? p->curve[(*rp)[c]]:(*rp)[c];
			*rp += p->tiff_samples;
	}
	if (p->is_raw == 2 && p->opt.shot_select) (*rp)--;
}

void DCR_CLASS dcr_adobe_dng_load_raw_lj(DCRAW* p)
{
	unsigned save, trow=0, tcol=0, jwide, jrow, jcol, row, col;
	struct dcr_jhead jh;
	ushort *rp;

	while (trow < p->raw_height) {
		save = dcr_ftell(p->obj_);
		if (p->tile_length < INT_MAX)
			dcr_fseek(p->obj_, dcr_get4(p), SEEK_SET);
		if (!dcr_ljpeg_start (p,&jh, 0)) break;
		jwide = jh.wide;
		if (p->filters) jwide *= jh.clrs;
		jwide /= p->is_raw;
		for (row=col=jrow=0; (int)jrow < jh.high; jrow++) {
			rp = dcr_ljpeg_row (p, jrow, &jh);
			for (jcol=0; jcol < jwide; jcol++) {
				dcr_adobe_copy_pixel (p,trow+row, tcol+col, &rp);
				if (++col >= p->tile_width || col >= p->raw_width)
					row += 1 + (col = 0);
			}
		}
		dcr_fseek(p->obj_, save+4, SEEK_SET);
		if ((tcol += p->tile_width) >= p->raw_width)
			trow += p->tile_length + (tcol = 0);
		free (jh.row);
	}
}

void DCR_CLASS dcr_adobe_dng_load_raw_nc(DCRAW* p)
{
	ushort *pixel, *rp;
	unsigned int row, col;

	pixel = (ushort *) calloc (p->raw_width * p->tiff_samples, sizeof *pixel);
	dcr_merror (p, pixel, "adobe_dng_load_raw_nc()");
	for (row=0; row < p->raw_height; row++) {
		if (p->tiff_bps == 16)
			dcr_read_shorts (p, pixel, p->raw_width * p->tiff_samples);
		else {
			dcr_getbits(p, -1);
			for (col=0; col < p->raw_width * p->tiff_samples; col++)
				pixel[col] = dcr_getbits(p, p->tiff_bps);
		}
		for (rp=pixel, col=0; col < p->raw_width; col++)
			dcr_adobe_copy_pixel (p,row, col, &rp);
	}
	free (pixel);
}

void DCR_CLASS dcr_pentax_k10_load_raw(DCRAW* p)
{
	static const uchar pentax_tree[] =
	{ 0,2,3,1,1,1,1,1,1,2,0,0,0,0,0,0,
	3,4,2,5,1,6,0,7,8,9,10,11,12 };
	int row, col, diff;
	ushort vpred[2][2] = {{0,0},{0,0}}, hpred[2];

	dcr_init_decoder(p);
	dcr_make_decoder (p, pentax_tree, 0);
	dcr_getbits(p, -1);
	for (row=0; row < p->height; row++)
		for (col=0; col < p->raw_width; col++) {
			diff = dcr_ljpeg_diff (p, p->first_decode);
			if (col < 2) hpred[col] = vpred[row & 1][col] += diff;
			else	   hpred[col & 1] += diff;
			if (col < p->width)
				BAYER(row,col) = hpred[col & 1];
			if (hpred[col & 1] >> 12) dcr_derror(p);
		}
}

void DCR_CLASS dcr_nikon_compressed_load_raw(DCRAW* p)
{
	static const uchar nikon_tree[][32] = {
		{ 0,1,5,1,1,1,1,1,1,2,0,0,0,0,0,0,	/* 12-bit lossy */
			5,4,3,6,2,7,1,0,8,9,11,10,12 },
		{ 0,1,5,1,1,1,1,1,1,2,0,0,0,0,0,0,	/* 12-bit lossy after split */
		0x39,0x5a,0x38,0x27,0x16,5,4,3,2,1,0,11,12,12 },
		{ 0,1,4,2,3,1,2,0,0,0,0,0,0,0,0,0,  /* 12-bit lossless */
		5,4,6,3,7,2,8,1,9,0,10,11,12 },
		{ 0,1,4,3,1,1,1,1,1,2,0,0,0,0,0,0,	/* 14-bit lossy */
		5,6,4,7,8,3,9,2,1,0,10,11,12,13,14 },
		{ 0,1,5,1,1,1,1,1,1,1,2,0,0,0,0,0,	/* 14-bit lossy after split */
		8,0x5c,0x4b,0x3a,0x29,7,6,5,4,3,2,1,0,13,14 },
		{ 0,1,4,2,2,3,1,2,0,0,0,0,0,0,0,0,	/* 14-bit lossless */
		7,6,8,5,9,4,10,3,11,12,2,0,1,13,14 } };
	struct dcr_decode *dindex;
	ushort ver0, ver1, vpred[2][2], hpred[2], csize;
	int i, min, max, step=0, huff=0, split=0, row, col, len, shl, diff;

	dcr_fseek(p->obj_, p->meta_offset, SEEK_SET);
	ver0 = dcr_fgetc(p->obj_);
	ver1 = dcr_fgetc(p->obj_);
	if (ver0 == 0x49 || ver1 == 0x58)
		dcr_fseek(p->obj_, 2110, SEEK_CUR);
	if (ver0 == 0x46) huff = 2;
	if (p->tiff_bps == 14) huff += 3;
	dcr_read_shorts (p, vpred[0], 4);
	max = 1 << p->tiff_bps & 0x7fff;
	if ((csize = dcr_get2(p)) > 1)
		step = max / (csize-1);
	if (ver0 == 0x44 && ver1 == 0x20 && step > 0) {
		for (i=0; i < csize; i++)
			p->curve[i*step] = dcr_get2(p);
		for (i=0; i < max; i++)
			p->curve[i] = ( p->curve[i-i%step]*(step-i%step) +
			p->curve[i-i%step+step]*(i%step) ) / step;
		dcr_fseek(p->obj_, p->meta_offset+562, SEEK_SET);
		split = dcr_get2(p);
	} else if (ver0 != 0x46 && csize <= 0x4001)
		dcr_read_shorts (p, p->curve, max=csize);
	while (p->curve[max-2] == p->curve[max-1]) max--;
	dcr_init_decoder(p);
	dcr_make_decoder (p, nikon_tree[huff], 0);
	dcr_fseek(p->obj_, p->data_offset, SEEK_SET);
	dcr_getbits(p, -1);
	for (min=row=0; row < p->height; row++) {
		if (split && row == split) {
			dcr_init_decoder(p);
			dcr_make_decoder (p, nikon_tree[huff+1], 0);
			max += (min = 16) << 1;
		}
		for (col=0; col < p->raw_width; col++) {
			for (dindex=p->first_decode; dindex->branch[0]; )
				dindex = dindex->branch[dcr_getbits(p, 1)];
			len = dindex->leaf & 15;
			shl = dindex->leaf >> 4;
			diff = ((dcr_getbits(p, len-shl) << 1) + 1) << shl >> 1;
			if ((diff & (1 << (len-1))) == 0)
				diff -= (1 << len) - !shl;
			if (col < 2) hpred[col] = vpred[row & 1][col] += diff;
			else	   hpred[col & 1] += diff;
			if ((ushort)(hpred[col & 1] + min) >= max) dcr_derror(p);
			if ((unsigned) (col-p->left_margin) < p->width)
				BAYER(row,col-p->left_margin) = p->curve[LIM((short)hpred[col & 1],0,0x3fff)];
		}
	}
}

/*
Figure out if a NEF file is compressed.  These fancy heuristics
are only needed for the D100, thanks to a bug in some cameras
that tags all images as "compressed".
*/
int DCR_CLASS dcr_nikon_is_compressed(DCRAW* p)
{
	uchar test[256];
	int i;

	dcr_fseek(p->obj_, p->data_offset, SEEK_SET);
	dcr_fread(p->obj_, test, 1, 256);
	for (i=15; i < 256; i+=16)
		if (test[i]) return 1;
		return 0;
}

/*
Returns 1 for a Coolpix 995, 0 for anything else.
*/
int DCR_CLASS dcr_nikon_e995(DCRAW* p)
{
	int i, histo[256];
	const uchar often[] = { 0x00, 0x55, 0xaa, 0xff };

	memset (histo, 0, sizeof histo);
	dcr_fseek(p->obj_, -2000, SEEK_END);
	for (i=0; i < 2000; i++)
		histo[dcr_fgetc(p->obj_)]++;
	for (i=0; i < 4; i++)
		if (histo[often[i]] < 200)
			return 0;
		return 1;
}

/*
Returns 1 for a Coolpix 2100, 0 for anything else.
*/
int DCR_CLASS dcr_nikon_e2100(DCRAW* p)
{
	uchar t[12];
	int i;

	dcr_fseek(p->obj_, 0, SEEK_SET);
	for (i=0; i < 1024; i++) {
		dcr_fread(p->obj_, t, 1, 12);
		if (((t[2] & t[4] & t[7] & t[9]) >> 4
			& t[1] & t[6] & t[8] & t[11] & 3) != 3)
			return 0;
	}
	return 1;
}

void DCR_CLASS dcr_nikon_3700(DCRAW* p)
{
	int bits, i;
	uchar dp[24];
	const struct {
		int bits;
		char make[12], model[15];
	} table[] = {
		{ 0x00, "PENTAX",  "Optio 33WR" },
		{ 0x03, "NIKON",   "E3200" },
		{ 0x32, "NIKON",   "E3700" },
		{ 0x33, "OLYMPUS", "C740UZ" } };

	dcr_fseek(p->obj_, 3072, SEEK_SET);
	dcr_fread(p->obj_, dp, 1, 24);
	bits = (dp[8] & 3) << 4 | (dp[20] & 3);
	for (i=0; i < sizeof table / sizeof *table; i++)
		if (bits == table[i].bits) {
			strcpy (p->make,  table[i].make );
			strcpy (p->model, table[i].model);
		}
}

/*
Separates a Minolta DiMAGE Z2 from a Nikon E4300.
*/
int DCR_CLASS dcr_minolta_z2(DCRAW* p)
{
	int i, nz;
	char tail[424];

	dcr_fseek(p->obj_, -(long)(sizeof tail), SEEK_END);
	dcr_fread(p->obj_, tail, 1, sizeof tail);
	for (nz=i=0; i < sizeof tail; i++)
		if (tail[i]) nz++;
	return nz > 20;
}

/* Here p->raw_width is in bytes, not pixels. */
void DCR_CLASS dcr_nikon_e900_load_raw(DCRAW* p)
{
	int offset=0, irow, row, col;

	for (irow=0; irow < p->height; irow++) {
		row = irow * 2 % p->height;
		if (row == 1)
			offset = - (-offset & -4096);
		dcr_fseek(p->obj_, offset, SEEK_SET);
		offset += p->raw_width;
		dcr_getbits(p, -1);
		for (col=0; col < p->width; col++)
			BAYER(row,col) = dcr_getbits(p, 10);
	}
}

/*
The Fuji Super CCD is just a Bayer grid rotated 45 degrees.
*/
void DCR_CLASS dcr_fuji_load_raw(DCRAW* p)
{
	ushort *pixel;
	int wide, row, col, r, c;

	dcr_fseek(p->obj_, (p->top_margin*p->raw_width + p->left_margin) * 2, SEEK_CUR);
	wide = p->fuji_width << !p->fuji_layout;
	pixel = (ushort *) calloc (wide, sizeof *pixel);
	dcr_merror (p, pixel, "fuji_load_raw()");
	for (row=0; row < p->raw_height; row++) {
		dcr_read_shorts (p, pixel, wide);
		dcr_fseek(p->obj_, 2*(p->raw_width - wide), SEEK_CUR);
		for (col=0; col < wide; col++) {
			if (p->fuji_layout) {
				r = p->fuji_width - 1 - col + (row >> 1);
				c = col + ((row+1) >> 1);
			} else {
				r = p->fuji_width - 1 + row - (col >> 1);
				c = row + ((col+1) >> 1);
			}
			BAYER(r,c) = pixel[col];
		}
	}
	free (pixel);
}

void DCR_CLASS dcr_jpeg_thumb (DCRAW* p, FILE *tfp);

void DCR_CLASS dcr_ppm_thumb (DCRAW* p, FILE *tfp)
{
	char *thumb;
	p->thumb_length = p->thumb_width*p->thumb_height*3;
	thumb = (char *) malloc (p->thumb_length);
	dcr_merror (p, thumb, "ppm_thumb()");
	fprintf (tfp, "P6\n%d %d\n255\n", p->thumb_width, p->thumb_height);
	dcr_fread(p->obj_, thumb, 1, p->thumb_length);
	fwrite (thumb, 1, p->thumb_length, tfp);
	free (thumb);
}

void DCR_CLASS dcr_layer_thumb (DCRAW* p, FILE *tfp)
{
	int i, c;
	char *thumb, map[][4] = { "012","102" };

	p->colors = p->thumb_misc >> 5 & 7;
	p->thumb_length = p->thumb_width*p->thumb_height;
	thumb = (char *) calloc (p->colors, p->thumb_length);
	dcr_merror (p, thumb, "layer_thumb()");
	fprintf (tfp, "P%d\n%d %d\n255\n",
		5 + (p->colors >> 1), p->thumb_width, p->thumb_height);
	dcr_fread(p->obj_, thumb, p->thumb_length, p->colors);
	for (i=0; i < (int)p->thumb_length; i++)
		FORCC(p) putc (thumb[i+p->thumb_length*(map[p->thumb_misc >> 8][c]-'0')], tfp);
	free (thumb);
}

void DCR_CLASS dcr_rollei_thumb (DCRAW* p, FILE *tfp)
{
	unsigned i;
	ushort *thumb;

	p->thumb_length = p->thumb_width * p->thumb_height;
	thumb = (ushort *) calloc (p->thumb_length, 2);
	dcr_merror (p, thumb, "rollei_thumb()");
	fprintf (tfp, "P6\n%d %d\n255\n", p->thumb_width, p->thumb_height);
	dcr_read_shorts (p, thumb, p->thumb_length);
	for (i=0; i < p->thumb_length; i++) {
		putc (thumb[i] << 3, tfp);
		putc (thumb[i] >> 5  << 2, tfp);
		putc (thumb[i] >> 11 << 3, tfp);
	}
	free (thumb);
}

void DCR_CLASS dcr_rollei_load_raw(DCRAW* p)
{
	uchar pixel[10];
	unsigned iten=0, isix, i, buffer=0, row, col, todo[16];

	isix = p->raw_width * p->raw_height * 5 / 8;
	while (dcr_fread(p->obj_, pixel, 1, 10) == 10) {
		for (i=0; i < 10; i+=2) {
			todo[i]   = iten++;
			todo[i+1] = pixel[i] << 8 | pixel[i+1];
			buffer    = pixel[i] >> 2 | buffer << 6;
		}
		for (   ; i < 16; i+=2) {
			todo[i]   = isix++;
			todo[i+1] = buffer >> (14-i)*5;
		}
		for (i=0; i < 16; i+=2) {
			row = todo[i] / p->raw_width - p->top_margin;
			col = todo[i] % p->raw_width - p->left_margin;
			if (row < p->height && col < p->width)
				BAYER(row,col) = (todo[i+1] & 0x3ff);
		}
	}
	p->maximum = 0x3ff;
}

int DCR_CLASS dcr_bayer (DCRAW* p, unsigned row, unsigned col)
{
	return (row < p->height && col < p->width) ? BAYER(row,col) : 0;
}

void DCR_CLASS dcr_phase_one_flat_field (DCRAW* p, int is_float, int nc)
{
	ushort head[8];
	unsigned wide, y, x, rend, cend, row, col;
	int c;
	float *mrow, num, mult[4];

	dcr_read_shorts (p, head, 8);
	wide = head[2] / head[4];
	mrow = (float *) calloc (nc*wide, sizeof *mrow);
	dcr_merror (p, mrow, "phase_one_flat_field()");
	for (y=0; y < (unsigned int)(head[3] / head[5]); y++) {
		for (x=0; x < wide; x++)
			for (c=0; c < nc; c+=2) {
				num = is_float ? (float)dcr_getreal(p, 11) : dcr_get2(p)/32768.0f;
				if (y==0) mrow[c*wide+x] = num;
				else mrow[(c+1)*wide+x] = (num - mrow[c*wide+x]) / head[5];
			}
		if (y==0) continue;
		rend = head[1]-p->top_margin + y*head[5];
		for (row = rend-head[5]; row < p->height && row < rend; row++) {
			for (x=1; x < wide; x++) {
				for (c=0; c < nc; c+=2) {
					mult[c] = mrow[c*wide+x-1];
					mult[c+1] = (mrow[c*wide+x] - mult[c]) / head[4];
				}
				cend = head[0]-p->left_margin + x*head[4];
				for (col = cend-head[4]; col < p->width && col < cend; col++) {
					c = nc > 2 ? FC(row,col) : 0;
					if (!(c & 1)) {
						c = (int)(BAYER(row,col) * mult[c]);
						BAYER(row,col) = LIM(c,0,65535);
					}
					for (c=0; c < nc; c+=2)
						mult[c] += mult[c+1];
				}
			}
			for (x=0; x < wide; x++)
				for (c=0; c < nc; c+=2)
					mrow[c*wide+x] += mrow[(c+1)*wide+x];
		}
	}
	free (mrow);
}

void DCR_CLASS dcr_phase_one_correct(DCRAW* p)
{
	unsigned entries, tag, data, save, col, row, type;
	int len, i, j, k, cip, val[4], dev[4], sum, max;
	int head[9], diff, mindiff=INT_MAX, off_412=0;
	static const signed char dir[12][2] =
	{ {-1,-1}, {-1,1}, {1,-1}, {1,1}, {-2,0}, {0,-2}, {0,2}, {2,0},
	{-2,-2}, {-2,2}, {2,-2}, {2,2} };
	float poly[8], num, cfrac, frac, mult[2], *yval[2];
	ushort curve[0x10000], *xval[2];

	if (p->opt.half_size || !p->meta_length) return;
	if (p->opt.verbose) fprintf (stderr,_("Phase One correction...\n"));
	dcr_fseek(p->obj_, p->meta_offset, SEEK_SET);
	p->order = dcr_get2(p);
	dcr_fseek(p->obj_, 6, SEEK_CUR);
	dcr_fseek(p->obj_, p->meta_offset+dcr_get4(p), SEEK_SET);
	entries = dcr_get4(p);  dcr_get4(p);
	while (entries--) {
		tag  = dcr_get4(p);
		len  = dcr_get4(p);
		data = dcr_get4(p);
		save = dcr_ftell(p->obj_);
		dcr_fseek(p->obj_, p->meta_offset+data, SEEK_SET);
		if (tag == 0x419) {				/* Polynomial curve */
			for (dcr_get4(p), i=0; i < 8; i++)
				poly[i] = (float)dcr_getreal(p, 11);
			poly[3] += (p->ph1.tag_210 - poly[7]) * poly[6] + 1;
			for (i=0; i < 0x10000; i++) {
				num = (poly[5]*i + poly[3])*i + poly[1];
				curve[i] = (unsigned short)LIM(num,0,65535);
			} goto apply;				/* apply to right half */
		} else if (tag == 0x41a) {			/* Polynomial curve */
			for (i=0; i < 4; i++)
				poly[i] = (float)dcr_getreal(p, 11);
			for (i=0; i < 0x10000; i++) {
				for (num=0, j=4; j--; )
					num = num * i + poly[j];
				curve[i] = (unsigned short)LIM(num+i,0,65535);
			} apply:					/* apply to whole image */
			for (row=0; row < p->height; row++)
				for (col = (tag & 1)*p->ph1.split_col; col < p->width; col++)
					BAYER(row,col) = curve[BAYER(row,col)];
		} else if (tag == 0x400) {			/* Sensor defects */
			while ((len -= 8) >= 0) {
				col  = dcr_get2(p) - p->left_margin;
				row  = dcr_get2(p) - p->top_margin;
				type = dcr_get2(p); dcr_get2(p);
				if (col >= p->width) continue;
				if (type == 131)			/* Bad column */
					for (row=0; row < p->height; row++)
						if (FC(row,col) == 1) {
							for (sum=i=0; i < 4; i++)
								sum += val[i] = dcr_bayer (p, row+dir[i][0], col+dir[i][1]);
							for (max=i=0; i < 4; i++) {
								dev[i] = abs((val[i] << 2) - sum);
								if (dev[max] < dev[i]) max = i;
							}
							BAYER(row,col) = (unsigned short)((sum - val[max])/3.0 + 0.5);
						} else {
							for (sum=0, i=8; i < 12; i++)
								sum += dcr_bayer (p, row+dir[i][0], col+dir[i][1]);
							BAYER(row,col) = (unsigned short)(0.5 + sum * 0.0732233 +
								(dcr_bayer(p, row,col-2) + dcr_bayer(p, row,col+2)) * 0.3535534);
						}
						else if (type == 129) {			/* Bad pixel */
							if (row >= p->height) continue;
							j = (FC(row,col) != 1) * 4;
							for (sum=0, i=j; i < j+8; i++)
								sum += dcr_bayer (p, row+dir[i][0], col+dir[i][1]);
							BAYER(row,col) = (sum + 4) >> 3;
						}
			}
		} else if (tag == 0x401) {			/* All-color flat fields */
			dcr_phase_one_flat_field (p, 1, 2);
		} else if (tag == 0x416 || tag == 0x410) {
			dcr_phase_one_flat_field (p, 0, 2);
		} else if (tag == 0x40b) {			/* Red+blue flat field */
			dcr_phase_one_flat_field (p, 0, 4);
		} else if (tag == 0x412) {
			dcr_fseek(p->obj_, 36, SEEK_CUR);
			diff = abs (dcr_get2(p) - p->ph1.tag_21a);
			if (mindiff > diff) {
				mindiff = diff;
				off_412 = dcr_ftell(p->obj_) - 38;
			}
		}
		dcr_fseek(p->obj_, save, SEEK_SET);
	}
	if (off_412) {
		dcr_fseek(p->obj_, off_412, SEEK_SET);
		for (i=0; i < 9; i++) head[i] = dcr_get4(p) & 0x7fff;
		yval[0] = (float *) calloc (head[1]*head[3] + head[2]*head[4], 6);
		dcr_merror (p, yval[0], "phase_one_correct()");
		yval[1] = (float  *) (yval[0] + head[1]*head[3]);
		xval[0] = (ushort *) (yval[1] + head[2]*head[4]);
		xval[1] = (ushort *) (xval[0] + head[1]*head[3]);
		dcr_get2(p);
		for (i=0; i < 2; i++)
			for (j=0; j < head[i+1]*head[i+3]; j++)
				yval[i][j] = (float)dcr_getreal(p, 11);
		for (i=0; i < 2; i++)
			for (j=0; j < head[i+1]*head[i+3]; j++)
				xval[i][j] = dcr_get2(p);
		for (row=0; row < p->height; row++)
			for (col=0; col < p->width; col++) {
				cfrac = (float) col * head[3] / p->raw_width;
				cfrac -= cip = (int)cfrac;
				num = (float)(BAYER(row,col) * 0.5);
				for (i=cip; i < cip+2; i++) {
					for (k=j=0; j < head[1]; j++)
						if (num < xval[0][k = head[1]*i+j]) break;
						frac = (j == 0 || j == head[1]) ? 0 :
						(xval[0][k] - num) / (xval[0][k] - xval[0][k-1]);
						mult[i-cip] = yval[0][k-1] * frac + yval[0][k] * (1-frac);
				}
				i = (int)(((mult[0] * (1-cfrac) + mult[1] * cfrac)
					* (row + p->top_margin) + num) * 2);
				BAYER(row,col) = LIM(i,0,65535);
			}
		free (yval[0]);
	}
}

void DCR_CLASS dcr_phase_one_load_raw(DCRAW* p)
{
	int row, col, a, b;
	ushort *pixel, akey, bkey, mask;

	dcr_fseek(p->obj_, p->ph1.key_off, SEEK_SET);
	akey = dcr_get2(p);
	bkey = dcr_get2(p);
	mask = p->ph1.format == 1 ? 0x5555:0x1354;
	dcr_fseek(p->obj_, p->data_offset + p->top_margin*p->raw_width*2, SEEK_SET);
	pixel = (ushort *) calloc (p->raw_width, sizeof *pixel);
	dcr_merror (p, pixel, "phase_one_load_raw()");
	for (row=0; row < p->height; row++) {
		dcr_read_shorts (p, pixel, p->raw_width);
		for (col=0; col < p->raw_width; col+=2) {
			a = pixel[col+0] ^ akey;
			b = pixel[col+1] ^ bkey;
			pixel[col+0] = (a & mask) | (b & ~mask);
			pixel[col+1] = (b & mask) | (a & ~mask);
		}
		for (col=0; col < p->width; col++)
			BAYER(row,col) = pixel[col+p->left_margin];
	}
	free (pixel);
	dcr_phase_one_correct(p);
}

unsigned DCR_CLASS dcr_ph1_bits (DCRAW* p,int nbits)
{
	if (nbits == -1)
		return (unsigned int)(p->ph1_bits_bitbuf = p->ph1_bits_vbits = 0);
	if (nbits == 0) return 0;
	if ((p->ph1_bits_vbits -= nbits) < 0) {
		p->ph1_bits_bitbuf = p->ph1_bits_bitbuf << 32 | dcr_get4(p);
		p->ph1_bits_vbits += 32;
	}
	return (unsigned int)(p->ph1_bits_bitbuf << (64-nbits-p->ph1_bits_vbits) >> (64-nbits));
}

void DCR_CLASS dcr_phase_one_load_raw_c(DCRAW* p)
{
	static const int length[] = { 8,7,6,9,11,10,5,12,14,13 };
	int *offset, len[2], pred[2], row, col, i, j;
	ushort *pixel;
	short (*black)[2];

	pixel = (ushort *) calloc (p->raw_width + p->raw_height*4, 2);
	dcr_merror (p, pixel, "phase_one_load_raw_c()");
	offset = (int *) (pixel + p->raw_width);
	dcr_fseek(p->obj_, p->strip_offset, SEEK_SET);
	for (row=0; row < p->raw_height; row++)
		offset[row] = dcr_get4(p);
	black = (short (*)[2]) offset + p->raw_height;
	dcr_fseek(p->obj_, p->ph1.black_off, SEEK_SET);
	if (p->ph1.black_off)
		dcr_read_shorts (p, (ushort *) black[0], p->raw_height*2);
	for (i=0; i < 256; i++)
		p->curve[i] = (unsigned short)(i*i / 3.969 + 0.5);
	for (row=0; row < p->raw_height; row++) {
		dcr_fseek(p->obj_, p->data_offset + offset[row], SEEK_SET);
		dcr_ph1_bits(p,-1);
		pred[0] = pred[1] = 0;
		for (col=0; col < p->raw_width; col++) {
			if (col >= (p->raw_width & -8))
				len[0] = len[1] = 14;
			else if ((col & 7) == 0)
				for (i=0; i < 2; i++) {
					for (j=0; j < 5 && !dcr_ph1_bits(p,1); j++);
					if (j--) len[i] = length[j*2 + dcr_ph1_bits(p,1)];
				}
				if ((i = len[col & 1]) == 14)
					pixel[col] = pred[col & 1] = dcr_ph1_bits(p,16);
				else
					pixel[col] = pred[col & 1] += dcr_ph1_bits(p,i) + 1 - (1 << (i - 1));
				if (pred[col & 1] >> 16) dcr_derror(p);
				if (p->ph1.format == 5 && pixel[col] < 256)
					pixel[…