#include <lcms.h>
#include <gsl/gsl_spline.h>

#include "pnm.h"

static int
pam_transform(FILE *in_fp, FILE *out_fp,
	const struct pnm *in_pnm, cmsHTRANSFORM *hTransform);

static LPGAMMATABLE
buildGammaCurve(const char *ctrl_points, int resolution);

int
main(int argc, char **argv) {
	int c, i, nProf = 0;
	char *in_prof = NULL, *out_prof = NULL;
	char *curve = NULL;
	double gamma = 0.0, bright = 0.0, contrast = 1.0, hue = 0.0;
	double saturation = 0.0;
	int tempSrc = 5000, tempDst = 5000;
	cmsHTRANSFORM hTransform;
	cmsHPROFILE profiles[3];
	DWORD type;
	struct pnm in_pnm;

	while ((c = getopt(argc, argv, "i:o:g:b:c:h:s:S:D:C:")) != EOF) {
		switch (c) {  
			case 'i':
				in_prof = optarg;
				break;
			case 'o':
				out_prof = optarg;
				break;
			case 'g':
				gamma = atof(optarg);
				break;
			case 'b':
				bright = atof(optarg);
				break;
			case 'c':
				contrast = atof(optarg);
				break;
			case 'h':
				hue = atof(optarg);
				break;
			case 's':
				saturation = atof(optarg);
				break;
			case 'S':
				tempSrc = atoi(optarg);
				break;
			case 'D':
				tempDst = atoi(optarg);
			case 'C':
				curve = optarg;
				break;

			default:
				break;
		}
	}
#if 0
	if (in_prof) {
		profiles[nProf++] = cmsOpenProfileFromFile(in_prof, "r");
	} else {
		profiles[nProf++] = cmsCreate_sRGBProfile();
	}
#endif

	if (gamma) {
		LPGAMMATABLE gTable[3];

		gTable[0] = cmsBuildGamma(256, gamma);
		gTable[1] = gTable[2] = cmsBuildGamma(256, 1.0);

		profiles[nProf++] = cmsCreateLinearizationDeviceLink(
			icSigLabData,
			gTable);
	}

	if (curve) {
		LPGAMMATABLE gTable[3];

		gTable[0] = gTable[1] = gTable[2] = buildGammaCurve(curve, 256);

		profiles[nProf++] = cmsCreateLinearizationDeviceLink(
			icSigLabData,
			gTable);
	}

#if 0
	profiles[nProf++] = cmsCreateBCHSWabstractProfile(
		5,
		bright,
		contrast,
		hue,
		saturation,
		tempSrc,
		tempDst);

	if (out_prof) {
		profiles[nProf++] = cmsOpenProfileFromFile(out_prof, "r");
	} else {
		profiles[nProf++] = cmsCreate_sRGBProfile();
	}
#endif

	if (readPnmHeader(stdin, &in_pnm) != 0) {
		exit(1);
	}

	if (in_pnm.maxval == 255) {
		type = TYPE_RGB_8;
	} else if (in_pnm.maxval == 65535) {
		type = TYPE_RGB_16_SE;
	} else {
		fprintf(stderr, "unsupported PNM maxval %d\n", in_pnm.maxval);
		exit(1);
	}

	hTransform = cmsCreateMultiprofileTransform(
		profiles,
		nProf,
		type,
		type,
		INTENT_PERCEPTUAL,
		cmsFLAGS_BLACKPOINTCOMPENSATION);
	
	writePnmHeader(stdout, &in_pnm);

	pam_transform(stdin, stdout, &in_pnm, &hTransform);

	cmsDeleteTransform(hTransform);
	
	for (i = 0; i < nProf; i++) {
	//	cmsCloseProfile(profiles[i]);
	}

	return 0;
}

static int
pam_transform(FILE *in_fp, FILE *out_fp,
	const struct pnm *in_pnm, cmsHTRANSFORM *hTransform)
{
	int row, col;
	int nbytes = in_pnm->maxval == 65535?2:1;
	char *buf = malloc(in_pnm->width * nbytes * 3);

	for (row = 0; row < in_pnm->height; row++) {
		fread(buf, in_pnm->width, nbytes * 3, in_fp);

		cmsDoTransform(*hTransform, buf, buf, in_pnm->width);

		fwrite(buf, in_pnm->width, nbytes * 3, out_fp);
	}

	free(buf);
	return 0;
}

#define MAX_CTRL 256

static LPGAMMATABLE
buildGammaCurve(const char *ctrl_points, int resolution) {
	LPGAMMATABLE table;
	char *pstr, *buf = strdup(ctrl_points);
	int i, n = 0;
	double X[MAX_CTRL], Y[MAX_CTRL];
	gsl_interp_accel *acc;
	gsl_spline *spline;

	while (pstr = strsep(&buf, ";")) {
		if (n == MAX_CTRL) {
			fprintf(stderr, "Maximum number of control points (%d) reached.\n",
				MAX_CTRL);
			break;
		}

		if (sscanf(pstr, "%lf,%lf", &X[n], &Y[n]) != 2 ||
			X[n] < 0.0 || X[n] > 1.0 || Y[n] < 0.0 || Y[n] > 1.0) {
			fprintf(stderr, "Could not parse control point %s.\n", pstr);
			free(buf);
			return NULL;
		}

		n++;
	}

	free(buf);

	table = (LPGAMMATABLE) calloc(resolution, sizeof(GAMMATABLE));
	table->nEntries = resolution;

	acc = gsl_interp_accel_alloc();
	spline = gsl_spline_alloc(gsl_interp_cspline, n);
	gsl_spline_init (spline, X, Y, n);

	for (i = 0; i < resolution; i++) {
		double _x = (double) i / (resolution - 1);
		double _y = gsl_spline_eval(spline, _x, acc);
        double dval = _y * 65535.0 + .5;
        if (dval > 65535.) dval = 65535.0;
        if (dval < 0) dval = 0;

        table->GammaTable[i] = (WORD) floor(dval);
	}

	gsl_spline_free(spline);
	gsl_interp_accel_free(acc);
	
	return table;
}