#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <math.h>

unsigned char buf[6];

long getlong(FILE *fp)
{
	fread(buf, 4, 1, fp);
	return ((buf[0] << 24) | (buf[1] << 16) | (buf[2] << 8) | buf[3]);
}

#define EOTF(x) ((x < 0.08145)? (x / 4.5) : (pow((x + 0.0993) / 1.0993, 1.0 / 0.45)))
static double lut_eotf[256];  /* LUT for Electro-Optical Transfer Function */
static double lut_oetf[256];  /* LUT for Opto-Electronic Transfer Function */

void init_LUTs()
{
	int idx;
	double e;

	for (idx = 0; idx < 256; idx++) {
		e = (idx - 16) / 219.0;
		e = EOTF(e);
		lut_eotf[idx] = e;
		e = (idx - 16.5) / 219.0;
		e = EOTF(e);
		lut_oetf[idx] = e;
	}
	return;
}

/* Inverse Quantization and Electro-Optical Transfer Function */
double iqt_eotf(int idx)
{
	return (lut_eotf[idx]);
}

/*  Opto-Electronic Transfer Function and Quantization */
int oetf_qt(double e)
{
	int step, idx = 128;

	/* Inverse table look-up by binary search */
	for (step = 64; step > 0; step >>= 1) {
		if (e < lut_oetf[idx]) {
			idx -= step;
		}
		else {
			idx += step;
		}
	}
	if (e < lut_oetf[idx]) idx--;
	return (idx);
}

/* Color conversion */
void BT2020toBT709(int *rgb)
{
	const double mat[3][3] = {
		{ 1.6605, -0.5876, -0.0729 },
		{ -0.1244,  1.1328, -0.0084 },
		{ -0.0181, -0.1006,  1.1188 }
	};
	int i, j;
	double rgb2020[3], rgb709[3];

	for (i = 0; i < 3; i++) {
		rgb2020[i] = iqt_eotf(rgb[i]);
	}
	for (i = 0; i < 3; i++) {
		rgb709[i] = 0.0;
		for (j = 0; j < 3; j++) {
			rgb709[i] += mat[i][j] * rgb2020[j];
		}
	}
	for (i = 0; i < 3; i++) {
		rgb[i] = oetf_qt(rgb709[i]);
	}
	return;
}

int main(int argc, char **argv)
{
	int x, y, i, width, height, rgb[3];
	int colorconv = 0;
	long offset;
	char *infile, *outfile;
	FILE *fpi, *fpo;

	if (argc < 3) {
		printf("Usage: dpx2ppm (-c) infile outfile\n");
		return (0);
	}
	if (strcmp(argv[1], "-c") == 0) {
		colorconv = 1;
		argv++;
	}
	infile = *++argv;
	outfile = *++argv;
	if ((fpi = fopen(infile, "rb")) == NULL) {
		fprintf(stderr, "Can't open %s !\n", infile);
		return (1);
	}
	if ((fpo = fopen(outfile, "wb")) == NULL) {
		fprintf(stderr, "Can't open %s !\n", outfile);
		return (1);
	}
	if (getlong(fpi) != 0x53445058) {
		fprintf(stderr, "Not a DPX file !\n");
		return (1);
	}
	offset = getlong(fpi);
	fseek(fpi, 772L, SEEK_SET);
	width = getlong(fpi);
	height = getlong(fpi);
	fseek(fpi, offset, SEEK_SET);
	fprintf(fpo, "P6\n%d %d\n255\n", width, height);
	init_LUTs();
	for (y = 0; y < height; y++) {
		for (x = 0; x < width; x++) {
			fread(buf, 6, 1, fpi);
			for (i = 0; i < 3; i++) {
				rgb[i] = buf[i << 1] & 0xff;
			}
			if (colorconv) {
				BT2020toBT709(rgb);
			}
			for (i = 0; i < 3; i++) {
				fputc(rgb[i], fpo);
			}
		}
	}
	fclose(fpi);
	fclose(fpo);
	return 0;
}