Change colourspace handling round a bit:1.0-colour-cleanup

- move the essence of GammaLUT into TransferFunction and handle
  different bit depths more neatly
- add ColourConversion to describe input gamma correction, colour
  transformation and then output gamma correction in one class.
- add default ColourConversions for sRGB->XYZ, Rec709->XYZ and
  XYZ->RGB.

16 files changed, 526 insertions, 147 deletions

diff --git a/src/colour_conversion.cc b/src/colour_conversion.cc
new file mode 100644
index 00000000..28edf9bc
--- /dev/null
+++ b/src/colour_conversion.cc
@@ -0,0 +1,75 @@
+/*
+    Copyright (C) 2014 Carl Hetherington <cth@carlh.net>
+
+    This program is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    This program is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with this program; if not, write to the Free Software
+    Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#include "colour_conversion.h"
+#include "gamma_transfer_function.h"
+#include "modified_gamma_transfer_function.h"
+#include "colour_matrix.h"
+
+using boost::shared_ptr;
+using namespace dcp;
+
+ColourConversion ColourConversion::srgb_to_xyz (
+	shared_ptr<const TransferFunction> (new ModifiedGammaTransferFunction (2.4, 0.04045, 0.055, 12.92)),
+	dcp::colour_matrix::srgb_to_xyz,
+	shared_ptr<const TransferFunction> (new GammaTransferFunction (2.6))
+	);
+
+/* XXX: what sort of RGB is this...? */
+ColourConversion ColourConversion::xyz_to_rgb (
+	shared_ptr<const TransferFunction> (new GammaTransferFunction (2.6)),
+	dcp::colour_matrix::xyz_to_rgb,
+	shared_ptr<const TransferFunction> (new GammaTransferFunction (2.2))
+	);
+
+ColourConversion ColourConversion::rec709_to_xyz (
+	shared_ptr<const TransferFunction> (new ModifiedGammaTransferFunction (2.4, 0.081, 0.099, 4.5)),
+	dcp::colour_matrix::srgb_to_xyz,
+	shared_ptr<const TransferFunction> (new GammaTransferFunction (2.6))
+	);
+
+ColourConversion::ColourConversion (
+	shared_ptr<const TransferFunction> in,
+	double const matrix[3][3],
+	shared_ptr<const TransferFunction> out
+	)
+	: _in (in)
+	, _matrix (3, 3)
+	, _out (out)
+{
+	for (int i = 0; i < 3; ++i) {
+		for (int j = 0; j < 3; ++j) {
+			_matrix (i, j) = matrix[i][j];
+		}
+	}
+}
+
+bool
+ColourConversion::about_equal (ColourConversion const & other, float epsilon) const
+{
+	for (int i = 0; i < 3; ++i) {
+		for (int j = 0; j < 3; ++j) {
+			if (fabs (_matrix(i, j) - other._matrix(i, j)) > epsilon) {
+				return false;
+			}
+		}
+	}
+
+	return _in->about_equal (other._in, epsilon) && _out->about_equal (other._out, epsilon);
+}
diff --git a/src/colour_conversion.h b/src/colour_conversion.h
new file mode 100644
index 00000000..9bb7c6ee
--- /dev/null
+++ b/src/colour_conversion.h
@@ -0,0 +1,76 @@
+/*
+    Copyright (C) 2014 Carl Hetherington <cth@carlh.net>
+
+    This program is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    This program is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with this program; if not, write to the Free Software
+    Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#include <boost/shared_ptr.hpp>
+#include <boost/numeric/ublas/matrix.hpp>
+
+namespace dcp {
+
+class TransferFunction;
+
+class ColourConversion
+{
+public:
+	ColourConversion ()
+		: _matrix (3, 3)
+	{}
+	
+	ColourConversion (
+		boost::shared_ptr<const TransferFunction> in,
+		double const matrix[3][3],
+		boost::shared_ptr<const TransferFunction> out
+		);
+
+	boost::shared_ptr<const TransferFunction> in () const {
+		return _in;
+	}
+
+	boost::numeric::ublas::matrix<double> matrix () const {
+		return _matrix;
+	}
+
+	boost::shared_ptr<const TransferFunction> out () const {
+		return _out;
+	}
+
+	void set_in (boost::shared_ptr<const TransferFunction> f) {
+		_in = f;
+	}
+
+	void set_matrix (boost::numeric::ublas::matrix<double> m) {
+		_matrix = m;
+	}
+
+	void set_out (boost::shared_ptr<const TransferFunction> f) {
+		_out = f;
+	}
+
+	bool about_equal (ColourConversion const & other, float epsilon) const;
+
+	static ColourConversion srgb_to_xyz;
+	static ColourConversion xyz_to_rgb;
+	static ColourConversion rec709_to_xyz;
+
+protected:
+	boost::shared_ptr<const TransferFunction> _in;
+	boost::numeric::ublas::matrix<double> _matrix;
+	boost::shared_ptr<const TransferFunction> _out;
+};
+
+}
diff --git a/src/gamma_transfer_function.cc b/src/gamma_transfer_function.cc
new file mode 100644
index 00000000..b4aa20f1
--- /dev/null
+++ b/src/gamma_transfer_function.cc
@@ -0,0 +1,55 @@
+/*
+    Copyright (C) 2012-2014 Carl Hetherington <cth@carlh.net>
+
+    This program is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    This program is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with this program; if not, write to the Free Software
+    Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#include "gamma_transfer_function.h"
+#include <cmath>
+
+using std::pow;
+using boost::shared_ptr;
+using boost::dynamic_pointer_cast;
+using namespace dcp;
+
+GammaTransferFunction::GammaTransferFunction (float gamma)
+	: _gamma (gamma)
+{
+
+}
+
+float *
+GammaTransferFunction::make_lut (int bit_depth) const
+{
+	int const bit_length = pow (2, bit_depth);
+	float* lut = new float[int(std::pow(2.0f, bit_depth))];
+	for (int i = 0; i < bit_length; ++i) {
+		lut[i] = pow(float(i) / (bit_length - 1), _gamma);
+	}
+
+	return lut;
+}
+
+bool
+GammaTransferFunction::about_equal (shared_ptr<const TransferFunction> other, float epsilon) const
+{
+	shared_ptr<const GammaTransferFunction> o = dynamic_pointer_cast<const GammaTransferFunction> (other);
+	if (!o) {
+		return false;
+	}
+
+	return fabs (_gamma - o->_gamma) < epsilon;
+}
diff --git a/src/gamma_lut.h b/src/gamma_transfer_function.h
index 76a63ccf..85035f66 100644
--- a/src/gamma_lut.h
+++ b/src/gamma_transfer_function.h
@@ -17,47 +17,26 @@
 
 */
 
-#ifndef LIBDCP_GAMMA_LUT_H
-#define LIBDCP_GAMMA_LUT_H
-
-#include "lut_cache.h"
+#include "transfer_function.h"
 
 namespace dcp {
 
-class GammaLUT
+class GammaTransferFunction : public TransferFunction
 {
 public:
-	GammaLUT (int bit_depth, float gamma, bool linearised);
-
-	~GammaLUT () {
-		delete[] _lut;
-	}
-	
-	float const * lut () const {
-		return _lut;
-	}
-
-	int bit_depth () const {
-		return _bit_depth;
-	}
+	GammaTransferFunction (float gamma);
 
 	float gamma () const {
 		return _gamma;
 	}
 
-	bool linearised () const {
-		return _linearised;
-	}
+	bool about_equal (boost::shared_ptr<const TransferFunction> other, float epsilon) const;
+
+protected:
+	float * make_lut (int bit_depth) const;
 
-	static LUTCache<GammaLUT> cache;
-	
 private:
-	float* _lut;
-	int _bit_depth;
 	float _gamma;
-	bool _linearised;
 };
 
 }
-
-#endif
diff --git a/src/modified_gamma_transfer_function.cc b/src/modified_gamma_transfer_function.cc
new file mode 100644
index 00000000..7878773b
--- /dev/null
+++ b/src/modified_gamma_transfer_function.cc
@@ -0,0 +1,68 @@
+/*
+    Copyright (C) 2012-2014 Carl Hetherington <cth@carlh.net>
+
+    This program is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    This program is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with this program; if not, write to the Free Software
+    Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#include "modified_gamma_transfer_function.h"
+#include <cmath>
+
+using std::pow;
+using boost::shared_ptr;
+using boost::dynamic_pointer_cast;
+using namespace dcp;
+
+ModifiedGammaTransferFunction::ModifiedGammaTransferFunction (float power, float threshold, float A, float B)
+	: _power (power)
+	, _threshold (threshold)
+	, _A (A)
+	, _B (B)
+{
+
+}
+
+float *
+ModifiedGammaTransferFunction::make_lut (int bit_depth) const
+{
+	int const bit_length = pow (2, bit_depth);
+	float* lut = new float[int(std::pow(2.0f, bit_depth))];
+	for (int i = 0; i < bit_length; ++i) {
+		float const p = static_cast<float> (i) / (bit_length - 1);
+		if (p > _threshold) {
+			lut[i] = pow ((p + _A) / (1 + _A), _power);
+		} else {
+			lut[i] = p / _B;
+		}
+	}
+
+	return lut;
+}
+
+bool
+ModifiedGammaTransferFunction::about_equal (shared_ptr<const TransferFunction> other, float epsilon) const
+{
+	shared_ptr<const ModifiedGammaTransferFunction> o = dynamic_pointer_cast<const ModifiedGammaTransferFunction> (other);
+	if (!o) {
+		return false;
+	}
+
+	return (
+		fabs (_power - o->_power) < epsilon &&
+	        fabs (_threshold - o->_threshold) < epsilon &&
+		fabs (_A - o->_A) < epsilon &&
+		fabs (_B - o->_B) < epsilon
+		);
+}
diff --git a/src/modified_gamma_transfer_function.h b/src/modified_gamma_transfer_function.h
new file mode 100644
index 00000000..0c6a6b84
--- /dev/null
+++ b/src/modified_gamma_transfer_function.h
@@ -0,0 +1,62 @@
+/*
+    Copyright (C) 2012-2014 Carl Hetherington <cth@carlh.net>
+
+    This program is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    This program is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with this program; if not, write to the Free Software
+    Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#include "transfer_function.h"
+
+namespace dcp {
+
+/** A transfer function which for an input x gives an output y where
+ *
+ *  y = x / B                      for x <= threshold
+ *  y = ((x + A) / (1 + A))^power  for x >  threshold
+ */
+class ModifiedGammaTransferFunction : public TransferFunction
+{
+public:
+	ModifiedGammaTransferFunction (float power, float threshold, float A, float B);
+
+	float power () const {
+		return _power;
+	}
+
+	float threshold () const {
+		return _threshold;
+	}
+
+	float A () const {
+		return _A;
+	}
+
+	float B () const {
+		return _B;
+	}
+
+	bool about_equal (boost::shared_ptr<const TransferFunction>, float epsilon) const;
+	
+protected:
+	float * make_lut (int bit_depth) const;
+	
+private:
+	float _power;
+	float _threshold;
+	float _A;
+	float _B;
+};
+
+}
diff --git a/src/mono_picture_frame.cc b/src/mono_picture_frame.cc
index f3063c0c..31499c64 100644
--- a/src/mono_picture_frame.cc
+++ b/src/mono_picture_frame.cc
@@ -24,9 +24,9 @@
 #include "mono_picture_frame.h"
 #include "exceptions.h"
 #include "argb_frame.h"
-#include "gamma_lut.h"
 #include "util.h"
 #include "rgb_xyz.h"
+#include "colour_conversion.h"
 #include "KM_fileio.h"
 #include "AS_DCP.h"
 #include <openjpeg.h>
@@ -123,7 +123,7 @@ MonoPictureFrame::argb_frame (int reduce, float srgb_gamma) const
 {
 	return xyz_to_rgba (
 		decompress_j2k (const_cast<uint8_t*> (_buffer->RoData()), _buffer->Size(), reduce),
-		GammaLUT::cache.get (12, DCI_GAMMA, false), GammaLUT::cache.get (12, 1 / srgb_gamma, false)
+		ColourConversion::xyz_to_rgb
 		);
 }
 
@@ -132,7 +132,7 @@ MonoPictureFrame::rgb_frame (uint8_t* buffer) const
 {
 	xyz_to_rgb (
 		decompress_j2k (const_cast<uint8_t*> (_buffer->RoData()), _buffer->Size(), 0),
-		GammaLUT::cache.get (12, DCI_GAMMA, false), GammaLUT::cache.get (12, 1 / 2.4, false),
+		ColourConversion::xyz_to_rgb,
 		buffer
 		);
 }
diff --git a/src/rgb_xyz.cc b/src/rgb_xyz.cc
index f53a51ea..de8d484a 100644
--- a/src/rgb_xyz.cc
+++ b/src/rgb_xyz.cc
@@ -20,9 +20,10 @@
 #include "rgb_xyz.h"
 #include "argb_frame.h"
 #include "xyz_frame.h"
-#include "gamma_lut.h"
 #include "image.h"
 #include "colour_matrix.h"
+#include "colour_conversion.h"
+#include "transfer_function.h"
 #include <cmath>
 
 using std::min;
@@ -34,18 +35,15 @@ using namespace dcp;
 
 /** Convert an openjpeg XYZ image to RGBA.
  *  @param xyz_frame Frame in XYZ.
- *  @param lut_in Input Gamma LUT to use.
- *  @param lut_out Output Gamma LUT to use.
  *  @return RGB image.
  */
 shared_ptr<ARGBFrame>
 dcp::xyz_to_rgba (
 	boost::shared_ptr<const XYZFrame> xyz_frame,
-	boost::shared_ptr<const GammaLUT> lut_in,
-	boost::shared_ptr<const GammaLUT> lut_out
+	ColourConversion const & conversion
 	)
 {
-	int const max_colour = pow (2, lut_out->bit_depth()) - 1;
+	int const max_colour = pow (2, 12) - 1;
 
 	struct {
 		double x, y, z;
@@ -60,9 +58,12 @@ dcp::xyz_to_rgba (
 	int* xyz_z = xyz_frame->data (2);
 
 	shared_ptr<ARGBFrame> argb_frame (new ARGBFrame (xyz_frame->size ()));
-
 	uint8_t* argb = argb_frame->data ();
 	
+	float const * lut_in = conversion.in()->lut (16);
+	float const * lut_out = conversion.out()->lut (12);
+	boost::numeric::ublas::matrix<double> matrix = conversion.matrix ();
+
 	for (int y = 0; y < xyz_frame->size().height; ++y) {
 		uint8_t* argb_line = argb;
 		for (int x = 0; x < xyz_frame->size().width; ++x) {
@@ -70,9 +71,9 @@ dcp::xyz_to_rgba (
 			assert (*xyz_x >= 0 && *xyz_y >= 0 && *xyz_z >= 0 && *xyz_x < 4096 && *xyz_y < 4096 && *xyz_z < 4096);
 			
 			/* In gamma LUT */
-			s.x = lut_in->lut()[*xyz_x++];
-			s.y = lut_in->lut()[*xyz_y++];
-			s.z = lut_in->lut()[*xyz_z++];
+			s.x = lut_in[*xyz_x++];
+			s.y = lut_in[*xyz_y++];
+			s.z = lut_in[*xyz_z++];
 
 			/* DCI companding */
 			s.x /= DCI_COEFFICIENT;
@@ -80,9 +81,9 @@ dcp::xyz_to_rgba (
 			s.z /= DCI_COEFFICIENT;
 
 			/* XYZ to RGB */
-			d.r = ((s.x * colour_matrix::xyz_to_rgb[0][0]) + (s.y * colour_matrix::xyz_to_rgb[0][1]) + (s.z * colour_matrix::xyz_to_rgb[0][2]));
-			d.g = ((s.x * colour_matrix::xyz_to_rgb[1][0]) + (s.y * colour_matrix::xyz_to_rgb[1][1]) + (s.z * colour_matrix::xyz_to_rgb[1][2]));
-			d.b = ((s.x * colour_matrix::xyz_to_rgb[2][0]) + (s.y * colour_matrix::xyz_to_rgb[2][1]) + (s.z * colour_matrix::xyz_to_rgb[2][2]));
+			d.r = ((s.x * matrix(0, 0)) + (s.y * matrix(0, 1)) + (s.z * matrix(0, 2)));
+			d.g = ((s.x * matrix(1, 0)) + (s.y * matrix(1, 1)) + (s.z * matrix(1, 2)));
+			d.b = ((s.x * matrix(2, 0)) + (s.y * matrix(2, 1)) + (s.z * matrix(2, 2)));
 			
 			d.r = min (d.r, 1.0);
 			d.r = max (d.r, 0.0);
@@ -94,9 +95,9 @@ dcp::xyz_to_rgba (
 			d.b = max (d.b, 0.0);
 			
 			/* Out gamma LUT */
-			*argb_line++ = lut_out->lut()[(int) (d.b * max_colour)] * 0xff;
-			*argb_line++ = lut_out->lut()[(int) (d.g * max_colour)] * 0xff;
-			*argb_line++ = lut_out->lut()[(int) (d.r * max_colour)] * 0xff;
+			*argb_line++ = lut_out[(int) (d.b * max_colour)] * 0xff;
+			*argb_line++ = lut_out[(int) (d.g * max_colour)] * 0xff;
+			*argb_line++ = lut_out[(int) (d.r * max_colour)] * 0xff;
 			*argb_line++ = 0xff;
 		}
 		
@@ -117,12 +118,11 @@ dcp::xyz_to_rgba (
 void
 dcp::xyz_to_rgb (
 	boost::shared_ptr<const XYZFrame> xyz_frame,
-	boost::shared_ptr<const GammaLUT> lut_in,
-	boost::shared_ptr<const GammaLUT> lut_out,
+	ColourConversion const & conversion,
 	uint8_t* buffer
 	)
 {
-	int const max_colour = pow (2, lut_out->bit_depth()) - 1;
+	int const max_colour = pow (2, 12) - 1;
 
 	struct {
 		double x, y, z;
@@ -136,6 +136,10 @@ dcp::xyz_to_rgb (
 	int* xyz_y = xyz_frame->data (1);
 	int* xyz_z = xyz_frame->data (2);
 
+	float const * lut_in = conversion.in()->lut (16);
+	float const * lut_out = conversion.out()->lut (12);
+	boost::numeric::ublas::matrix<double> matrix = conversion.matrix ();
+	
 	for (int y = 0; y < xyz_frame->size().height; ++y) {
 		uint8_t* buffer_line = buffer;
 		for (int x = 0; x < xyz_frame->size().width; ++x) {
@@ -143,9 +147,9 @@ dcp::xyz_to_rgb (
 			assert (*xyz_x >= 0 && *xyz_y >= 0 && *xyz_z >= 0 && *xyz_x < 4096 && *xyz_y < 4096 && *xyz_z < 4096);
 			
 			/* In gamma LUT */
-			s.x = lut_in->lut()[*xyz_x++];
-			s.y = lut_in->lut()[*xyz_y++];
-			s.z = lut_in->lut()[*xyz_z++];
+			s.x = lut_in[*xyz_x++];
+			s.y = lut_in[*xyz_y++];
+			s.z = lut_in[*xyz_z++];
 
 			/* DCI companding */
 			s.x /= DCI_COEFFICIENT;
@@ -153,9 +157,9 @@ dcp::xyz_to_rgb (
 			s.z /= DCI_COEFFICIENT;
 
 			/* XYZ to RGB */
-			d.r = ((s.x * colour_matrix::xyz_to_rgb[0][0]) + (s.y * colour_matrix::xyz_to_rgb[0][1]) + (s.z * colour_matrix::xyz_to_rgb[0][2]));
-			d.g = ((s.x * colour_matrix::xyz_to_rgb[1][0]) + (s.y * colour_matrix::xyz_to_rgb[1][1]) + (s.z * colour_matrix::xyz_to_rgb[1][2]));
-			d.b = ((s.x * colour_matrix::xyz_to_rgb[2][0]) + (s.y * colour_matrix::xyz_to_rgb[2][1]) + (s.z * colour_matrix::xyz_to_rgb[2][2]));
+			d.r = ((s.x * matrix(0, 0)) + (s.y * matrix(0, 1)) + (s.z * matrix(0, 2)));
+			d.g = ((s.x * matrix(1, 0)) + (s.y * matrix(1, 1)) + (s.z * matrix(1, 2)));
+			d.b = ((s.x * matrix(2, 0)) + (s.y * matrix(2, 1)) + (s.z * matrix(2, 2)));
 			
 			d.r = min (d.r, 1.0);
 			d.r = max (d.r, 0.0);
@@ -167,9 +171,9 @@ dcp::xyz_to_rgb (
 			d.b = max (d.b, 0.0);
 			
 			/* Out gamma LUT */
-			*buffer_line++ = lut_out->lut()[(int) (d.r * max_colour)] * 0xff;
-			*buffer_line++ = lut_out->lut()[(int) (d.g * max_colour)] * 0xff;
-			*buffer_line++ = lut_out->lut()[(int) (d.b * max_colour)] * 0xff;
+			*buffer_line++ = lut_out[(int) (d.r * max_colour)] * 0xff;
+			*buffer_line++ = lut_out[(int) (d.g * max_colour)] * 0xff;
+			*buffer_line++ = lut_out[(int) (d.b * max_colour)] * 0xff;
 		}
 		
 		buffer += xyz_frame->size().width * 3;
@@ -182,14 +186,9 @@ dcp::xyz_to_rgb (
 shared_ptr<dcp::XYZFrame>
 dcp::rgb_to_xyz (
 	boost::shared_ptr<const Image> rgb,
-	boost::shared_ptr<const GammaLUT> lut_in,
-	boost::shared_ptr<const GammaLUT> lut_out,
-	double const colour_matrix[3][3]
+	ColourConversion const & conversion
 	)
 {
-	assert (lut_in->bit_depth() == 12);
-	assert (lut_out->bit_depth() == 16);
-	
 	shared_ptr<XYZFrame> xyz (new XYZFrame (rgb->size ()));
 
 	struct {
@@ -200,28 +199,24 @@ dcp::rgb_to_xyz (
 		double x, y, z;
 	} d;
 
+	float const * lut_in = conversion.in()->lut (12);
+	float const * lut_out = conversion.out()->lut (16);
+	boost::numeric::ublas::matrix<double> matrix = conversion.matrix ();
+
 	int jn = 0;
 	for (int y = 0; y < rgb->size().height; ++y) {
 		uint16_t* p = reinterpret_cast<uint16_t *> (rgb->data()[0] + y * rgb->stride()[0]);
 		for (int x = 0; x < rgb->size().width; ++x) {
 
 			/* In gamma LUT (converting 16-bit to 12-bit) */
-			s.r = lut_in->lut()[*p++ >> 4];
-			s.g = lut_in->lut()[*p++ >> 4];
-			s.b = lut_in->lut()[*p++ >> 4];
+			s.r = lut_in[*p++ >> 4];
+			s.g = lut_in[*p++ >> 4];
+			s.b = lut_in[*p++ >> 4];
 			
 			/* RGB to XYZ Matrix */
-			d.x = ((s.r * colour_matrix[0][0]) +
-			       (s.g * colour_matrix[0][1]) +
-			       (s.b * colour_matrix[0][2]));
-			
-			d.y = ((s.r * colour_matrix[1][0]) +
-			       (s.g * colour_matrix[1][1]) +
-			       (s.b * colour_matrix[1][2]));
-			
-			d.z = ((s.r * colour_matrix[2][0]) +
-			       (s.g * colour_matrix[2][1]) +
-			       (s.b * colour_matrix[2][2]));
+			d.x = ((s.r * matrix(0, 0)) + (s.g * matrix(0, 1)) + (s.b * matrix(0, 2)));
+			d.y = ((s.r * matrix(1, 0)) + (s.g * matrix(1, 1)) + (s.b * matrix(1, 2)));
+			d.z = ((s.r * matrix(2, 0)) + (s.g * matrix(2, 1)) + (s.b * matrix(2, 2)));
 			
 			/* DCI companding */
 			d.x = d.x * DCI_COEFFICIENT * 65535;
@@ -233,9 +228,9 @@ dcp::rgb_to_xyz (
 			assert (d.z >= 0 && d.z < 65536);
 			
 			/* Out gamma LUT */
-			xyz->data(0)[jn] = lut_out->lut()[(int) d.x] * 4096;
-			xyz->data(1)[jn] = lut_out->lut()[(int) d.y] * 4096;
-			xyz->data(2)[jn] = lut_out->lut()[(int) d.z] * 4096;
+			xyz->data(0)[jn] = lut_out[(int) d.x] * 4096;
+			xyz->data(1)[jn] = lut_out[(int) d.y] * 4096;
+			xyz->data(2)[jn] = lut_out[(int) d.z] * 4096;
 
 			++jn;
 		}
diff --git a/src/rgb_xyz.h b/src/rgb_xyz.h
index 871473c0..a6463b0d 100644
--- a/src/rgb_xyz.h
+++ b/src/rgb_xyz.h
@@ -24,21 +24,12 @@ namespace dcp {
 
 class ARGBFrame;	
 class XYZFrame;
-class GammaLUT;
 class Image;
+class ColourConversion;
 	
-extern boost::shared_ptr<ARGBFrame> xyz_to_rgba (
-	boost::shared_ptr<const XYZFrame>, boost::shared_ptr<const GammaLUT>, boost::shared_ptr<const GammaLUT>
-	);
-
-extern void xyz_to_rgb (
-	boost::shared_ptr<const XYZFrame>, boost::shared_ptr<const GammaLUT>, boost::shared_ptr<const GammaLUT>, uint8_t* buffer
-	);
-	
-extern boost::shared_ptr<XYZFrame> rgb_to_xyz (
-	boost::shared_ptr<const Image>, boost::shared_ptr<const GammaLUT>, boost::shared_ptr<const GammaLUT>, double const colour_matrix[3][3]
-	);
-
+extern boost::shared_ptr<ARGBFrame> xyz_to_rgba (boost::shared_ptr<const XYZFrame>, ColourConversion const & conversion);
+extern void xyz_to_rgb (boost::shared_ptr<const XYZFrame>, ColourConversion const & conversion, uint8_t* buffer);
+extern boost::shared_ptr<XYZFrame> rgb_to_xyz (boost::shared_ptr<const Image>, ColourConversion const & conversion);
 extern boost::shared_ptr<XYZFrame> xyz_to_xyz (boost::shared_ptr<const Image>);
 	
 }
diff --git a/src/stereo_picture_frame.cc b/src/stereo_picture_frame.cc
index bcc17e3a..7d989777 100644
--- a/src/stereo_picture_frame.cc
+++ b/src/stereo_picture_frame.cc
@@ -21,8 +21,8 @@
 #include "exceptions.h"
 #include "argb_frame.h"
 #include "util.h"
-#include "gamma_lut.h"
 #include "rgb_xyz.h"
+#include "colour_conversion.h"
 #include "AS_DCP.h"
 #include "KM_fileio.h"
 #include <openjpeg.h>
@@ -89,7 +89,7 @@ StereoPictureFrame::argb_frame (Eye eye, int reduce, float srgb_gamma) const
 		break;
 	}
 	
-	return xyz_to_rgba (xyz_frame, GammaLUT::cache.get (12, DCI_GAMMA, false), GammaLUT::cache.get (12, 1 / srgb_gamma, false));
+	return xyz_to_rgba (xyz_frame, ColourConversion::xyz_to_rgb);
 }
 
 void
@@ -105,7 +105,7 @@ StereoPictureFrame::rgb_frame (Eye eye, uint8_t* buffer) const
 		break;
 	}
 	
-	return xyz_to_rgb (xyz_frame, GammaLUT::cache.get (12, DCI_GAMMA, false), GammaLUT::cache.get (12, 1 / 2.4, false), buffer);
+	return xyz_to_rgb (xyz_frame, ColourConversion::xyz_to_rgb, buffer);
 }
 
 uint8_t const *
diff --git a/src/gamma_lut.cc b/src/transfer_function.cc
index abe35b7c..2c7e3e52 100644
--- a/src/gamma_lut.cc
+++ b/src/transfer_function.cc
@@ -17,34 +17,30 @@
 
 */
 
-#include "gamma_lut.h"
-#include "lut_cache.h"
+#include "transfer_function.h"
 #include <cmath>
 
+using std::pow;
+using std::map;
 using namespace dcp;
 
-LUTCache<GammaLUT> GammaLUT::cache;
+TransferFunction::~TransferFunction ()
+{
+	for (map<int, float*>::const_iterator i = _luts.begin(); i != _luts.end(); ++i) {
+		delete[] i->second;
+	}
 
-GammaLUT::GammaLUT (int bit_depth, float gamma, bool linearised)
-	: _bit_depth (bit_depth)
-	, _gamma (gamma)
-	, _linearised (linearised)
+	_luts.clear ();
+}
+
+float const *
+TransferFunction::lut (int bit_depth) const
 {
-	_lut = new float[int(std::pow(2.0f, _bit_depth))];
-	int const bit_length = pow (2, _bit_depth);
-
-	if (_linearised) {
-		for (int i = 0; i < bit_length; ++i) {
-			float const p = static_cast<float> (i) / (bit_length - 1);
-			if (p > 0.04045) {
-				_lut[i] = pow ((p + 0.055) / 1.055, gamma);
-			} else {
-				_lut[i] = p / 12.92;
-			}
-		}
-	} else {
-		for (int i = 0; i < bit_length; ++i) {
-			_lut[i] = pow(float(i) / (bit_length - 1), gamma);
-		}
+	map<int, float*>::const_iterator i = _luts.find (bit_depth);
+	if (i != _luts.end ()) {
+		return i->second;
 	}
+
+	_luts[bit_depth] = make_lut (bit_depth);
+	return _luts[bit_depth];
 }
diff --git a/src/lut_cache.h b/src/transfer_function.h
index 9bffa7b0..31a794e2 100644
--- a/src/lut_cache.h
+++ b/src/transfer_function.h
@@ -1,5 +1,5 @@
 /*
-    Copyright (C) 2013-2014 Carl Hetherington <cth@carlh.net>
+    Copyright (C) 2012-2014 Carl Hetherington <cth@carlh.net>
 
     This program is free software; you can redistribute it and/or modify
     it under the terms of the GNU General Public License as published by
@@ -17,34 +17,29 @@
 
 */
 
-#ifndef LIBDCP_LUT_CACHE_H
-#define LIBDCP_LUT_CACHE_H
+#ifndef LIBDCP_TRANSFER_FUNCTION_H
+#define LIBDCP_TRANSFER_FUNCTION_H
 
-#include <list>
-#include <boost/shared_ptr.hpp>
 #include <boost/noncopyable.hpp>
+#include <boost/shared_ptr.hpp>
+#include <map>
 
 namespace dcp {
 
-template<class T>
-class LUTCache : public boost::noncopyable
+class TransferFunction : public boost::noncopyable
 {
 public:
-	boost::shared_ptr<T> get (int bit_depth, float gamma, bool linearised)
-	{
-		for (typename std::list<boost::shared_ptr<T> >::iterator i = _cache.begin(); i != _cache.end(); ++i) {
-			if ((*i)->bit_depth() == bit_depth && (*i)->gamma() == gamma) {
-				return *i;
-			}
-		}
-
-		boost::shared_ptr<T> lut (new T (bit_depth, gamma, linearised));
-		_cache.push_back (lut);
-		return lut;
-	}
-
-private:
-	std::list<boost::shared_ptr<T> > _cache;
+	virtual ~TransferFunction ();
+
+	float const * lut (int bit_depth) const;
+
+	virtual bool about_equal (boost::shared_ptr<const TransferFunction> other, float epsilon) const = 0;
+
+protected:
+	virtual float * make_lut (int bit_depth) const = 0;
+
+private:	
+	mutable std::map<int, float*> _luts;
 };
 
 }
diff --git a/src/util.cc b/src/util.cc
index 1334b21d..58eacbfb 100644
--- a/src/util.cc
+++ b/src/util.cc
@@ -26,7 +26,6 @@
 #include "types.h"
 #include "argb_frame.h"
 #include "certificates.h"
-#include "gamma_lut.h"
 #include "xyz_frame.h"
 #include "compose.hpp"
 #include "KM_util.h"
diff --git a/src/wscript b/src/wscript
index 583c3608..4b08da49 100644
--- a/src/wscript
+++ b/src/wscript
@@ -6,6 +6,7 @@ def build(bld):
              asset.cc
              certificate_chain.cc
              certificates.cc
+             colour_conversion.cc
              colour_matrix.cc
              content.cc
              cpl.cc
@@ -17,13 +18,14 @@ def build(bld):
              exceptions.cc
              file.cc
              font.cc
-             gamma_lut.cc
+             gamma_transfer_function.cc
              image.cc
              interop_load_font.cc
              interop_subtitle_content.cc
              key.cc
              local_time.cc
              metadata.cc
+             modified_gamma_transfer_function.cc
              mono_picture_mxf.cc
              mono_picture_mxf_writer.cc
              mono_picture_frame.cc
@@ -54,6 +56,7 @@ def build(bld):
              subtitle_content.cc
              subtitle_string.cc
              text.cc
+             transfer_function.cc
              types.cc
              util.cc
              version.cc
@@ -64,6 +67,7 @@ def build(bld):
               asset.h
               certificate_chain.h
               certificates.h
+              colour_conversion.h
               colour_matrix.h
               cpl.h
               content.h
@@ -73,16 +77,16 @@ def build(bld):
               decrypted_kdm_key.h
               encrypted_kdm.h
               exceptions.h
-              gamma_lut.h
+              gamma_transfer_function.h
               image.h
               interop_load_font.h
               interop_subtitle_content.h
               key.h
               local_time.h
-              lut_cache.h
               metadata.h
               mono_picture_mxf.h
               mono_picture_frame.h
+              modified_gamma_transfer_function.h
               mxf.h
               mxf_writer.h
               object.h
@@ -111,6 +115,7 @@ def build(bld):
               subtitle.h
               subtitle_content.h
               subtitle_string.h
+              transfer_function.h
               types.h
               util.h
               version.h
diff --git a/test/colour_conversion_test.cc b/test/colour_conversion_test.cc
new file mode 100644
index 00000000..060148a6
--- /dev/null
+++ b/test/colour_conversion_test.cc
@@ -0,0 +1,82 @@
+/*
+    Copyright (C) 2014 Carl Hetherington <cth@carlh.net>
+
+    This program is free software; you can redistribute it and/or modify
+    it under the terms of the GNU General Public License as published by
+    the Free Software Foundation; either version 2 of the License, or
+    (at your option) any later version.
+
+    This program is distributed in the hope that it will be useful,
+    but WITHOUT ANY WARRANTY; without even the implied warranty of
+    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+    GNU General Public License for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with this program; if not, write to the Free Software
+    Foundation, Inc., 675 Mass Ave, Cambridge, MA 02139, USA.
+
+*/
+
+#include "gamma_transfer_function.h"
+#include "colour_conversion.h"
+#include "modified_gamma_transfer_function.h"
+#include <boost/test/unit_test.hpp>
+#include <cmath>
+
+using std::pow;
+using boost::shared_ptr;
+using namespace dcp;
+
+static void
+check_gamma (shared_ptr<const TransferFunction> tf, int bit_depth, float gamma)
+{
+	float const * lut = tf->lut (bit_depth);
+	int const count = pow (2, bit_depth);
+
+	for (int i = 0; i < count; ++i) {
+		BOOST_CHECK_CLOSE (lut[i], pow (float(i) / (count - 1), gamma), 0.001);
+	}
+}
+
+static void
+check_modified_gamma (shared_ptr<const TransferFunction> tf, int bit_depth, float power, float threshold, float A, float B)
+{
+	float const * lut = tf->lut (bit_depth);
+	int const count = pow (2, bit_depth);
+
+	for (int i = 0; i < count; ++i) {
+		float const x = float(i) / (count - 1);
+		if (x > threshold) {
+			BOOST_CHECK_CLOSE (lut[i], pow ((x + A) / (1 + A), power), 0.001);
+		} else {
+			BOOST_CHECK_CLOSE (lut[i], (x / B), 0.001);
+		}
+	}
+}
+
+BOOST_AUTO_TEST_CASE (colour_conversion_test1)
+{
+	ColourConversion cc = ColourConversion::srgb_to_xyz;
+
+	check_modified_gamma (cc.in(), 8, 2.4, 0.04045, 0.055, 12.92);
+	check_modified_gamma (cc.in(), 12, 2.4, 0.04045, 0.055, 12.92);
+	check_modified_gamma (cc.in(), 16, 2.4, 0.04045, 0.055, 12.92);
+
+	check_gamma (cc.out(), 8, 2.6);
+	check_gamma (cc.out(), 12, 2.6);
+	check_gamma (cc.out(), 16, 2.6);
+}
+
+BOOST_AUTO_TEST_CASE (colour_conversion_test2)
+{
+	ColourConversion cc = ColourConversion::rec709_to_xyz;
+
+	check_modified_gamma (cc.in(), 8, 2.4, 0.081, 0.099, 4.5);
+	check_modified_gamma (cc.in(), 12, 2.4, 0.081, 0.099, 4.5);
+	check_modified_gamma (cc.in(), 16, 2.4, 0.081, 0.099, 4.5);
+
+	check_gamma (cc.out(), 8, 2.6);
+	check_gamma (cc.out(), 12, 2.6);
+	check_gamma (cc.out(), 16, 2.6);
+}
+
diff --git a/test/wscript b/test/wscript
index 0ac8b800..dbbb51d5 100644
--- a/test/wscript
+++ b/test/wscript
@@ -28,6 +28,7 @@ def build(bld):
                  argb_frame_test.cc
                  certificates_test.cc
                  colour_test.cc
+                 colour_conversion_test.cc
                  cpl_sar_test.cc
                  dcp_test.cc
                  dcp_time_test.cc