vcglib/vcg/space/colorspace.h

/****************************************************************************
* VCGLib                                                            o o     *
* Visual and Computer Graphics Library                            o     o   *
*                                                                _   O  _   *
* Copyright(C) 2006                                                \/)\/    *
* Visual Computing Lab                                            /\/|      *
* ISTI - Italian National Research Council                           |      *
*                                                                    \      *
* All rights reserved.                                                      *
*                                                                           *
* 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 (http://www.gnu.org/licenses/gpl.txt)          *
* for more details.                                                         *
*                                                                           *
****************************************************************************/

#ifndef __VCGLIB_COLORSPACE
#define __VCGLIB_COLORSPACE

#include <vcg/space/color4.h>
#include<algorithm>

/** 
 * EVERYTHING YOU ALWAYS WANTED TO KNOW ABOUT COLORS BUT WERE AFRAID TO ASK
 ***************************************************************************


 UNDERSTANDING ILLUMINANTS
 -------------------------


 ILLUMINANT A - INCANDESCENT

 The CIE standard illuminant A was first recommended in 1931 to represent incandescent light
 sources with an approximate color temperature of 2856 degrees Kelvin.

 ILLUMINANT C - SUNLIGHT

 CIE standard illuminant C was recommended in 1931 to represent average daylight with a
 color temperature of 6774 degrees Kelvin.

 ILLUMINANT D50 - DAYLIGHT AT 5000K - (DAYLIGHT - RED SHADE)

 CIE standard illuminant D50 is correlated color temperature D illuminant which is calculated
 from the standard illuminant D65. This standard illuminant represents Daylight at 5000
 degrees Kelvin and is the ANSI standard illuminant used in the graphic arts industry.

 ILLUMINANT D65 - DAYLIGHT AT 6500K (DAYLIGHT - NEUTRAL)

 The CIE standard illuminant D65 represents daylight at approximately 6500 degrees Kelvin.
 This standard and the method for calculating different correlating color temperatures was
 introduced in 1964.

 ILLUMINANT D75 - DAYLIGHT - BLUE SHADE

 ...

 ILLUMINANT F2 - COOL WHITE (COOL WHITE FLUORESCENT)

 The CIE standard illuminant F2 represents the typical cool white fluorescent source. The
 spikes which protrude above the continuous graph in the fluorescent illuminants represent the
 power measurements at the principle mercury emission lines.

 ILLUMINANT F7 - BROAD BAND DAYLIGHT (BROAD BAND WHITE FLUORESCENT)

 The CIE standard Illuminant F7 represents a broad band daylight fluorescent source.
 Examples of this type of source are the GE<EFBFBD> and Philips<EFBFBD> Daylight fluorescent sources.

 ILLUMINANT F11 - NARROW BAND WHITE (TL84 FLUORESCENT)

 The CIE standard illuminant F11 represents a narrow band white fluorescent source. 
 Sources similar to this illuminant are the Philips<EFBFBD> F40AX41 and TL841, as
 well as the GE<EFBFBD> SPX41.

 ILLUMINANT F12 - ULTRALUME 3000 FLUORESCENT

 ...


 
 References
 ----------

 * Danny Pascale, "RGB coordinates of the Macbeth ColorChecker".

 * X-Rite Incorporated, "Understanding Illuminants", Whitepaper,
   http://www.xrite.com/documents/apps/public/whitepapers/Ca00002a.pdf

 * http://www.brucelindbloom.com - a web site with a huge amount of information about color theory.


*/


namespace vcg 
{

/** \addtogroup space */
/*@{*/
    /**
        This (empty) class is used to convert colors between different color space.
				\note All computations are in double precision independently of color data type.
				\note All color components are assumed in the range [0.0, 1.0].
     */
template <typename T>
class ColorSpace
{

// definitions
public:

	enum ConeResponse	{	BRADFORD = 0,	VON_KRIES = 1, XYZ_SCALING = 2};

	enum Illuminant
	{
		VCG_ILLUMINANT_A = 0,	
		VCG_ILLUMINANT_B = 1,	
		VCG_ILLUMINANT_C = 2,
		VCG_ILLUMINANT_D50 = 3,	
		VCG_ILLUMINANT_D55 = 4,	
		VCG_ILLUMINANT_D65 = 5,	
		VCG_ILLUMINANT_D75 = 6,
		VCG_ILLUMINANT_E = 7,
		VCG_ILLUMINANT_INVALID = 8
	};

	enum RGBSpaces
	{
		ADOBE_RGB = 0,
		APPLE_RGB = 1,
		BEST_RGB = 2,
		BETA_RGB = 3,
		BRUCE_RGB = 4,
		CIE_RGB = 5,
		COLOR_MATCH = 6,
		DON_RGB4 = 7,
		ECI_RGB = 8,
		EKTA_SPACE = 9,
		NTSC_RGB = 10,
		PAL_RGB = 11,
		PROPHOTO = 12,
		SMPTE_C = 13,
		SRGB = 14,
		WIDE_GAMUT = 15
	};

// private methods
private:

	static double CIE_EPSILON()
	{
		return 0.008856;  // Intent of the CIE Standard --> 216/24389
	}

	static double CIE_KI()
	{
		return 903.3;     // Intent of the CIE Standard --> 24389/27
	}

	static double RGB2XYZ(int index)
	{
		// RGB WORKING SPACE DATA (RGB TO XYZ CONVERSION MATRIX)
		//
		//  144 elements (16 RGB spaces x matrices of 3 x 3)
		//
		static double RGB_TO_XYZ[]=
		{
			// Adobe RGB (1998)  (ref. D65)
			0.576700, 0.297361,  0.0270328,
			0.185556, 0.627355,  0.0706879,
			0.188212, 0.0752847, 0.991248,

			// Apple RGB  (ref. D65)
			0.449695, 0.244634,  0.0251829,
			0.316251, 0.672034,  0.141184,
			0.18452,  0.0833318, 0.922602,

			// BestRGB  (ref. D50)
			0.632670, 0.228457,  0.000000,
			0.204556, 0.737352,  0.00951424,
			0.126995, 0.0341908, 0.815696,

			// Beta RGB  (ref. D50)
			0.671254, 0.303273,  0.000000,
			0.174583, 0.663786,  0.040701,
			0.118383, 0.0329413, 0.784509,

			// BruceRGB  (ref. D65)
			0.467384, 0.240995,  0.0219086,
			0.294454, 0.683554,  0.0736135,
			0.188629, 0.0754517, 0.993447,

			// CIE  (ref. E)
			0.488718, 0.176204,  0.000000,
			0.310680, 0.812985,  0.0102048,
			0.200602, 0.0108109, 0.989795,

			// ColorMatch  (ref. D50)
			0.509344, 0.274884,  0.0242545,
			0.320907, 0.658132,  0.108782,
			0.133969, 0.0669845, 0.692174,

			// DonRGB4  (ref. D50)
			0.645771,   0.278350,   0.00371134,
			0.193351,   0.687970,   0.0179862,
			0.125098,   0.0336802,  0.803513,

			// ECI  (ref. D50)
			0.650204, 0.320250,  0.000000,
			0.178077, 0.602071,  0.0678390,
			0.135938, 0.0776791, 0.757371,

			// Ekta Space PS5  (ref. D50)
			0.593891,  0.260629,   0.000000,
			0.272980,  0.734946,   0.0419970,
			0.0973486, 0.00442493, 0.783213,

			// NTSC  (ref. C)
			0.606734, 0.298839, 0.000000,
			0.173564, 0.586811, 0.0661196,
			0.200112, 0.114350, 1.11491,

			// PAL / SECAM  (ref. D65)
			0.430587, 0.222021,  0.0201837,
			0.341545, 0.706645,  0.129551,
			0.178336, 0.0713342, 0.939234,

			// ProPhoto  (ref. D50)
			0.797675,  0.288040, 0.000000,
			0.135192,  0.711874, 0.000000,
			0.0313534, 0.000086, 0.825210,

			// SMPTE-C  (ref. D65)
			0.393555, 0.212395,  0.0187407,
			0.365253, 0.701049,  0.111932,
			0.191659, 0.0865558, 0.958297,

			// sRGB  (ref. D65)
			0.412424, 0.212656,  0.0193324,
			0.357579, 0.715158,  0.119193,  
			0.180464, 0.0721856, 0.950444,

			// WideGamut  (ref. D50)
			0.716105, 0.258187,  0.000000,
			0.100930, 0.724938,  0.0517813,
			0.147186, 0.0168748, 0.773429

		};

		return RGB_TO_XYZ[index];
	}

	static double XYZ2RGB(int index)
	{
		// RGB WORKING SPACE DATA (XYZ TO RGB CONVERSION MATRIX)
		//
		//  144 elements (16 RGB spaces x matrices of 3 x 3)
		//
		static double XYZ_TO_RGB[]=
		{
			// Adobe RGB (1998)  (ref. D65)
			2.04148,   -0.969258,   0.0134455, 
			-0.564977,   1.87599,  -0.118373,
			-0.344713,   0.0415557,  1.01527,

			// Apple RGB  (ref. D65)
			2.95176, -1.0851,   0.0854804,
			-1.28951, 1.99084,  -0.269456, 
			-0.47388, 0.0372023,   1.09113,

			// BestRGB  (ref. D50)
			1.75526,  -0.544134,   0.00634681,
			-0.483679, 1.50688,  -0.0175762, 
			-0.253000, 0.0215528, 1.22570,

			// Beta RGB  (ref. D50)
			1.68323,   -0.771023, 0.0400012,
			-0.428236, 1.70656,  -0.0885376,
			-0.236018, 0.0446899, 1.27236,

			// BruceRGB  (ref. D65)
			2.74566,  -0.969257,   0.0112707,
			-1.13589,  1.87599,  -0.113959,
			-0.435057, 0.0415557, 1.01311,

			// CIE  (ref. E)
			2.37067,  -0.513885,  0.00529818,
			-0.900040, 1.42530,  -0.0146949,
			-0.470634, 0.0885814, 1.00940,

			// ColorMatch  (ref. D50)
			2.64229,  -1.11198,   0.0821698,
			-1.22343,  2.05902,  -0.280725,  
			-0.393014, 0.0159614, 1.45599,

			// DonRGB4  (ref. D50)
			1.76039,  -0.712629,   0.00782072,
			-0.488120, 1.65274,  -0.0347411,
			-0.253613, 0.0416715, 1.24477,

			// ECI  (ref. D50)
			1.78276,   -0.959362,  0.0859318,
			-0.496985,  1.94780,  -0.174467, 
			-0.269010, -0.0275807, 1.32283,

			// Ekta Space PS5  (ref. D50)
			2.00438,  -0.711029,   0.0381263,
			-0.730484, 1.62021,  -0.0868780,
			-0.245005,  0.0792227, 1.27254,

			// NTSC  (ref. C)
			1.91049,   -0.984310,  0.0583744,
			-0.532592,  1.99845,  -0.118518,
			-0.288284, -0.0282980, 0.898611,

			// PAL / SECAM  (ref. D65)
			3.06313,   -0.969258,   0.0678674,
			-1.39328,   1.87599,  -0.228821,
			-0.475788,  0.0415557, 1.06919,

			// ProPhoto  (ref. D50)
			1.34594,   -0.544599,   0.000000,
			-0.255608,  1.50817,   0.000000,
			-0.0511118, 0.0205351, 1.21181,

			// SMPTE-C  (ref. D65)
			3.50570,   -1.06906,   0.0563117,
			-1.73964,   1.97781,  -0.196994,
			-0.544011,  0.0351720, 1.05005,

			// sRGB  (ref. D65)
			3.24071,  -0.969258,   0.0556352,
			-1.53726,  1.87599,   -0.203996,
			-0.498571, 0.0415557,  1.05707,

			// WideGamut  (ref. D50)
			1.46281,  -0.521793,   0.0349342,
			-0.184062, 1.44724,   -0.0968931,
			-0.274361, 0.0677228,  1.28841

		};
	
		return XYZ_TO_RGB[index];
	}

	static double CA(int index)
	{
		// CHROMATIC ADAPTATION MATRIX 
		//
		//  196 (8 x 8 illuminants x 3 cone response) matrices of 3 x 3 elements
		//
		static double CAmatrix[]=
		{
			// A -> A
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,

			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,

			// A -> B (Bradford)
			0.890600, -0.097037, 0.054012,
			-0.082873, 1.075262, -0.091063,
			0.268563, 0.088084, 2.486772,
		
			// A -> B (Von Kries)
			0.957515, -0.018042, 0.000000,
			-0.164247, 1.018522, 0.000000,
			0.290746, 0.003635, 2.397466,

			// A -> B (XYZ scaling)
			0.902065, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 2.397474,

			// A -> C (Bradford)
			0.853036, -0.123884, 0.091195,
			-0.113032, 1.085373, -0.155374,
			0.440458, 0.142587, 3.477780,

			// A -> C (Von Kries)
			0.941839, -0.024699, 0.000000,
			-0.224854, 1.025358,  0.000000,
			0.480676, 0.004974, 3.322435,

			// A -> C (XYZ scaling)
			0.892678, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 3.322447,

			// A -> D50 (Bradford)
			0.878034, -0.111621, 0.050313,
			-0.091487, 1.092265, -0.083964,
			0.257070, 0.085312, 2.402220,

			// A -> D50 (Von Kries)
			0.953214, -0.019868, 0.000000,
			-0.180875, 1.020397, 0.000000,
			0.276267, 0.004004, 2.321454,

			// A -> D50 (XYZ scaling)
			0.877936, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 2.321462,

			// A -> D55 (Bradford)
			0.864256, -0.122489, 0.061011,
			-0.102279, 1.098455, -0.102340,
			0.307550, 0.101476, 2.689914,

			// A -> D55 (Von Kries)
			0.947636, -0.022237, 0.000000,
			-0.202440, 1.022829, 0.000000,
			0.331718, 0.004481, 2.590505,
		
			// A -> D55 (XYZ Scaling)
			0.870670, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 2.590514,

			// A -> D65 (Bradford)
			0.844763, -0.136524, 0.080011,
			-0.117903, 1.103952, -0.135190,
			0.395490, 0.129376, 3.196569,

			// A -> D65 (Von Kries)
			0.939518, -0.025685, 0.000000,
			-0.233826, 1.026369, 0.000000,
			0.428852, 0.005174, 3.063452,

			// A -> D65 (XYZ Scaling)
			0.865414, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 3.063462,

			// A -> D75 (Bradford)
			0.830991, -0.145710,  0.095536,
			-0.129133,  1.106062, -0.162123,
			0.466601,  0.151821,  3.608682,

			// A -> D75 (Von Kries)
			0.933660, -0.028172, 0.000000,
			-0.256474, 1.028923, 0.000000,
			0.507632, 0.005674, 3.447762,

			// A -> D75 (XYZ Scaling)
			0.864434,  0.000000,  0.000000,
			0.000000,  1.000000,  0.000000,
			0.000000,  0.000000,  3.447773,

			// A -> E (Bradford)
			0.881682, -0.100560,  0.071051,
			-0.090783,  1.070733, -0.120888,
			0.344476,  0.111711,  2.933736,

			// A -> E (Von Kries)
			0.953314, -0.019826, 0.000000,
			-0.180491, 1.020355, 0.000000,
			0.375531, 0.003993, 2.813168,

			// A -> E (XYZ Scaling)
			0.910515, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 2.813177,

			// B ->	A (Bradford)
			1.138867,  0.104490, -0.020910,
			0.077128,  0.934300,  0.032538,
			-0.125725, -0.044379,  0.403234,

			// B -> A (Von Kries)
			1.047555, 0.018556, 0.000000,
			0.168937, 0.984806, 0.000000,
			-0.127296, -0.003743, 0.417104,

			// B -> A (XYZ Scaling)
			1.108567, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.417106,

			// B -> B
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,

			// B -> C (Bradford)
			0.950475, -0.030658,  0.014905,
			-0.025481,  1.009149, -0.024973,
			0.075375,  0.024904,  1.397787,

			// B -> C (Von Kries)
			0.982455, -0.006847, 0.000000,
			-0.062331, 1.005606, 0.000000,
			0.081442, 0.001380, 1.385807,

			// B -> C (XYZ Scaling)
			0.989593, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.385811,

			// B -> D50 (Bradford)
			0.985029, -0.014775, -0.001704,
			-0.009391,  1.014671,  0.003596,
			-0.002672, -0.000039,  0.966056,

			// B -> D50 (Von Kries)
			0.995186, -0.001879, 0.000000,
			-0.017098, 1.001537, 0.000000,
			-0.005430, 0.000380, 0.968292,

			// B -> D50 (XYZ Scaling)
			0.973252, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.968295,

			// B -> D55 (Bradford)
			0.967155, -0.026843,  0.002545,
			-0.018894,  1.020141, -0.003387,
			0.019895,  0.007571,  1.081535,

			// B -> D55 (Von Kries)
			0.988943, -0.004315,  0.000000,
			-0.039278, 1.003532,  0.000000,
			0.018490,  0.000871,  1.080514,

			// B -> D55 (XYZ Scaling)
			0.965197, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.080518,

			// B -> D65 (Bradford)
			0.941484, -0.042836,  0.010157,
			-0.032134,  1.025102, -0.016128,
			0.058499,  0.020341,  1.284904,

			// B -> D65 (Von Kries)
			0.979857, -0.007861, 0.000000,
			-0.071558, 1.006436, 0.000000,
			0.060156, 0.001586, 1.277783,

			// B -> D65 (XYZ Scaling)
			0.959370, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.277788,

			// B -> D75 (Bradford)
			0.923139, -0.053546, 0.016407,
			-0.041374, 1.027096, -0.026684,
			0.089403, 0.030453, 1.450325,

			// B -> D75 (Von Kries)
			0.973300, -0.010419, 0.000000,
			-0.094851, 1.008531, 0.000000,
			0.093846, 0.002101, 1.438081,

			// B -> D75 (XYZ Scaling)
			0.958283, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.438086,

			// B -> E (Bradford)
			0.987430, -0.004980, 0.006942,
			-0.005608, 0.996265, -0.012009,
			0.032084, 0.010171, 1.179413,

			// B -> E (Von kries)
			0.995299, -0.001835, 0.000000,
			-0.016702, 1.001503, 0.000000,
			0.035959, 0.000370, 1.173388,

			// B -> E (XYZ Scaling)
			1.009367, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.173392,

			// C ->	A (Bradford)
			1.203987, 0.140744, -0.025283,
			0.102952, 0.928000, 0.038760,
			-0.156705, -0.055873, 0.289153,

			// C -> A (Von Kries)
			1.067896, 0.025724, 0.000000,
			0.234191, 0.980909, 0.000000,
			-0.154849, -0.005190, 0.300982,
		
			// C -> A (XYZ Scaling)
			1.120225, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.300983,

			// C -> B (Bradford)
			1.053816,   0.032278, -0.010661,
			0.025192,   0.991268,  0.017441,
			-0.057275, -0.019402,  0.715681,

			// C -> B (Von Kries)
			1.018301,   0.006933, 0.000000,
			0.063126,   0.994853, 0.000000,
			-0.059908, -0.001398, 0.721597,

			// C -> B (XYZ Scaling)	
			1.010516, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.721599,

			// C -> C
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,

			// C -> D50 (Bradford)
			1.037765, 0.017182, -0.011978,
			0.015460, 1.005438, 0.020371,
			-0.058148, -0.018868, 0.691416,

			// C -> D50 (Von Kries)
			1.013279, 0.005031, 0.000000,
			0.045808, 0.996264, 0.000000,
			-0.063514, -0.001014, 0.698718,

			// C -> D50 (XYZ Scaling)
			0.983487, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.698721,

			// C -> D55 (Bradford)
			1.018382, 0.004560, -0.008957,
			0.005982, 1.010689, 0.015570,
			-0.040789, -0.012837, 0.773954,

			// C -> D55 (Von Kries)
			1.006768, 0.002564, 0.000000,
			0.023348, 0.998095, 0.000000,
			-0.045848, -0.000516, 0.779698,

			// C -> D55 (XYZ Scaling)
			0.975347, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.779701,

			// C -> D65 (Bradford)
			0.990490, -0.012269, -0.003515,
			-0.007115, 1.015427, 0.006679,
			-0.011434, -0.002878, 0.919313,

			// C -> D65 (Von Kries)
			0.997291, -0.001026, 0.000000,
			-0.009340, 1.000760, 0.000000,
			-0.015193, 0.000208, 0.922047,

			// C -> D65 (XYZ Scaling)
			0.969459, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.922050,

			// C -> D75 (Bradford)
			0.970531, -0.023599,  0.000967,
			-0.016198, 1.017309, -0.000743,
			0.011914,  0.004934,  1.037548,

			// C -> D75 (Von Kries)
			0.990453, -0.003617, 0.000000,
			-0.032927, 1.002683, 0.000000,
			0.009543, 0.000730, 1.037718,

			// C -> D75 (XYZ Scaling)
			0.968360, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.037721,

			// C -> E (Bradford)
			1.040046, 0.026802, -0.005645,
			0.019876, 0.987617, 0.008842,
			-0.033485, -0.011765, 0.843919,

			// C -> E (Von Kries)
			1.013396, 0.005075, 0.000000,
			0.046208, 0.996233, 0.000000,
			-0.033655, -0.001024, 0.846716,

			// C -> E (XYZ Scaling)
			1.019981, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.846719,

			// D50 ->	A (Bradford)
			1.157264,   0.119830, -0.020050,
			0.087174,   0.922061,  0.030403,
			-0.126938, -0.045569,  0.417348,

			// D50 -> A (Von Kries)
			1.052976, 0.020503, 0.000000,
			0.186659, 0.983644, 0.000000,
			-0.125633, -0.004137, 0.430762,

			// D50 -> A (XYZ Scaling)
			1.139034, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.430763,

			// D50 -> B (Bradford)
			1.015344, 0.014785, 0.001735,
			0.009388, 0.985677, -0.003652,
			0.002809, 0.000081, 1.035142,

			// D50 -> B (Von Kries)
			1.004871, 0.001885, 0.000000,
			0.017164, 0.998496, 0.000000,
			0.005627, -0.000381, 1.032740,

			// D50 -> B (XYZ Scaling)
			1.027483, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.032743,

			// D50 -> C (Bradford)
			0.964813, -0.016165,  0.017191,
			-0.016469,  0.994317, -0.029580,
			0.080692,  0.025774,  1.446947,

			// D50 -> C (Von Kries)
			0.987122, -0.004985, 0.000000,
			-0.045378, 1.003977, 0.000000,
			0.089662, 0.001004, 1.431183,

			// D50 -> C (XYZ Scaling)
			1.016791, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.431187,

			// D50 -> D50
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,

			// D50 -> D55 (Bradford)
			0.981751, -0.012159, 0.004411,
			-0.009617, 1.005251, -0.007265,
			0.023309, 0.007844, 1.119548,

			// D50 -> D55 (Von Kries)
			0.993685, -0.002444,  0.000000,
			-0.022249, 1.001949,  0.000000,
			0.024676,  0.000493,  1.115894,

			// D50 -> D55 (XYZ Scaling)
			0.991724, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.115898,

			// D50 -> D65 (Bradford)
			0.955556, -0.028302, 0.012305,
			-0.023049, 1.009944, -0.020494,
			0.063197, 0.021018, 1.330084,

			// D50 -> D65 (Von Kries)
			0.984494, -0.006002, 0.000000,
		-0.054637, 1.004788, 0.000000,
			0.067667, 0.001210, 1.319622,

			// D50 -> D65 (XYZ Scaling)
			0.985737, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.319627,

			// D50 -> D75 (Bradford)
			0.936847, -0.039129, 0.018781,
			-0.032442, 1.011771, -0.031445,
			0.095135, 0.031456, 1.501335,

			// D50 -> D75 (Von Kries)
			0.977862, -0.008569, 0.000000,
			-0.078007, 1.006836, 0.000000,
			0.102432, 0.001727, 1.485168,

			// D50 -> D75 (XYZ Scaling)
			0.984620, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.485174,

			// D50 -> E (Bradford)
			1.002553, 0.009691, 0.008918,
			0.003624, 0.981912, -0.016079,
			0.035984, 0.010595, 1.220878,

			// D50 -> E (Von Kries)
			1.000115, 0.000044, 0.000000,
			0.000401, 0.999965, 0.000000,
			0.042744, -0.000010, 1.211809,

			// D50 -> E (XYZ Scaling)
			1.037108, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.211813,

			// D55 -> A (Bradford)
			1.180601, 0.133653, -0.021693,
			0.097012, 0.918163, 0.032732,
			-0.138643, -0.049919, 0.373005,

			// D55 -> A (Von Kries)
			1.060184, 0.023049, 0.000000,
			0.209842, 0.982241, 0.000000,
			-0.136121, -0.004650, 0.386023,

			// D55 -> A (XYZ Scaling)
			1.148540, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.386024,

			// D55 -> B (Bradford)
			1.034540, 0.027239, -0.002349,
			0.019098, 0.980736, 0.003026,
			-0.019163, -0.007366, 0.924635,

			// D55 -> B (Von Kries)
			1.011356, 0.004348, 0.000000,
			0.039593, 0.996649, 0.000000,
			-0.017339, -0.000878, 0.925479,

			// D55 -> B (XYZ Scaling)
			1.036058, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.925482,

			// D55 -> C (Bradford)
			0.982433, -0.004287,  0.011457,
			-0.006610, 0.989199, -0.019977,
			0.051668,  0.016181,  1.292341,

			// D55 -> C (Von Kries)
			0.993339, -0.002552, 0.000000,
			-0.023228, 1.001966, 0.000000,
			0.058394, 0.000514, 1.282539,
		
			// D55 -> C (XYZ Scaling)
			1.025276, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.282543,

			// D55 -> D50 (Bradford)
			1.018803, 0.012353, -0.003934,
			0.009594, 0.994842, 0.006418,
			-0.021278, -0.007227, 0.893255,

			// D55 -> D50 (Von Kries)
			1.006412, 0.002455, 0.000000,
			0.022357, 0.998107, 0.000000,
			-0.022266, -0.000495, 0.896136,

			// D55 -> D50 (XYZ Scaling)
			1.008345, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.896140,

			// D55 -> D55
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,

			// D55 -> D65 (Bradford)
			0.972990, -0.016440, 0.007051,
			-0.013357, 1.004598, -0.011734,
			0.036285, 0.012078, 1.187991,

			// D55 -> D65 (Von Kries)
			0.990671, -0.003573, 0.000000,
			-0.032527, 1.002753, 0.000000,
			0.038746, 0.000720, 1.182565,

			// D55 -> D65 (XYZ Scaling)
			0.993963, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.182570,

			// D55 -> D75 (Bradford)
			0.953688, -0.027490, 0.012840,
			-0.022677, 1.006379, -0.021468,
			0.065280, 0.021618, 1.340903,

			// D55 -> D75 (Von Kries)
			0.983939, -0.006152, 0.000000,
			-0.056002, 1.004740, 0.000000,
			0.070060, 0.001240, 1.330918,

			// D55 -> D75 (XYZ Scaling)
			0.992836, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.330923,

			// D55 -> E (Bradford)
			1.021308, 0.021962,  0.004085,
			0.013455, 0.977009, -0.008075,
			0.010784, 0.002161,  1.090481,

			// D55 -> E (Von Kries)
			1.006527, 0.002499, 0.000000,
			0.022756, 0.998075, 0.000000,
			0.016037, -0.000505, 1.085950,

			// D55 -> E (XYZ Scaling)
			1.045763, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.085953,

			// D65 -> A (Bradford)
			1.216371, 0.153235, -0.023966,
			0.110931, 0.915343, 0.035935,
			-0.154983, -0.056006, 0.314346,

			// D65 -> A (Von kries)
			1.071049, 0.026803, 0.000000,
			0.244014, 0.980413, 0.000000,
			-0.150348, -0.005408, 0.326427,

			// D65 -> A (XYZ Scaling)
			1.155516, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.326428,

			// D65 -> B (Bradford)
			1.064164, 0.044624, -0.007852,
			0.032589, 0.976635, 0.012001,
			-0.048965, -0.017493, 0.778436,

			// D65 -> B (Von Kries)
			1.021142, 0.007975, 0.000000,
			0.072613, 0.994171, 0.000000,
			-0.048164, -0.001609, 0.782600,

			// D65 -> B (XYZ Scaling)
			1.042351, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.782603,

			// D65 -> C (Bradford)
			1.009732, 0.012211, 0.003772,
			0.006993, 0.984871, -0.007129,
			0.012581, 0.003235, 1.087795,

			// D65 -> C (Von Kries)
			1.002728, 0.001028, 0.000000,
			0.009367, 0.999248, 0.000000,
			0.016519, -0.000208, 1.084536,

			// D65 -> C (XYZ Scaling)
			1.031503, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.084540,

			// D65 -> D50 (Bradford)
			1.047835, 0.029556, -0.009238,
			0.022897, 0.990481, 0.015050,
			-0.050147, -0.017056, 0.752034,

			// D65 -> D50 (Von Kries)
			1.016089, 0.006069, 0.000000,
			0.055260, 0.995564, 0.000000,
			-0.052154, -0.001224, 0.757788,

			// D65 -> D50 (XYZ Scaling)
			1.014470, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.757790,

			// D65 -> D55 (Bradford)
			1.028213, 0.016898, -0.005936,
			0.013304, 0.995522, 0.009754,
			-0.031540, -0.010637, 0.841840,

			// D65 -> D55 (Von Kries)
			1.009537, 0.003597, 0.000000,
			0.032756, 0.997370, 0.000000,
			-0.033098, -0.000726, 0.845614,

			// D65 -> D55 (XYZ Scaling)
			1.006074, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.845616,

			// D65 -> D65
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,

			// D65 -> D75 (Bradford)
			0.979823, -0.011388,  0.004880,
			-0.009251, 1.001719, -0.008121,
			0.025117,  0.008361,  1.128649,

			// D65 -> D75 (Von Kries)
			0.993120, -0.002596, 0.000000,
			-0.023629, 1.001897, 0.000000,
			0.026719, 0.000523, 1.125446,

			// D65 -> D75 (XYZ Scaling)
			0.998867, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.125450,

			// D65 -> E (Bradford)
			1.050285, 0.039078, -0.002409,
			0.027087, 0.972947, 0.002652,
			-0.023276, -0.009266, 0.917968,

			// D65 -> E (Von kries)
			1.016207, 0.006113, 0.000000,
			0.055662, 0.995532, 0.000000,
			-0.019769, -0.001234, 0.918297,

			// D65 -> E (XYZ Scaling)
			1.052114, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.918300,

			// D75 -> A (Bradford)	
			1.243649, 0.167322, -0.025407,
			0.120882, 0.914830, 0.037899,
			-0.165889, -0.060122, 0.278800,

			// D75 -> A (Von Kries)
			1.079173, 0.029548, 0.000000,
			0.269008, 0.979254, 0.000000,
			-0.159335, -0.005962, 0.290041,

			// D75 -> A (XYZ Scaling)
			1.156827, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.290042,

			// D75 -> B (Bradford)
			1.086904, 0.056997, -0.011247,
			0.042021, 0.975291, 0.017469,
			-0.067883, -0.023992, 0.689828,

			// D75 -> B (Von Kries)
			1.028470, 0.010625, 0.000000,
			0.096735, 0.992539, 0.000000,
			-0.067258, -0.002144, 0.695366,

			// D75 -> B (XYZ Scaling)
			1.043533, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.695369,

			// D75 -> C (Bradford)
			1.030774,   0.023916, -0.000943,
			0.016405,   0.983362,  0.000688,
			-0.011913, -0.004951,  0.963819,

			// D75 -> C (Von Kries)
			1.009762, 0.003643, 0.000000,
			0.033168, 0.997442, 0.000000,
			-0.009311, -0.000735,  0.963647,

			// D75 -> C (XYZ Scaling)
			1.032673, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.963650,

			// D75 -> D50 (Bradford)
			1.070127, 0.041775, -0.012512,
			0.032186, 0.988978, 0.020311,
			-0.068484, -0.023368, 0.666442,

			// D75 -> D50 (Von Kries)
			1.023337, 0.008709, 0.000000,
			0.079295, 0.993884, 0.000000,
			-0.070673, -0.001757, 0.673320,

			// D75 -> D50 (XYZ Scaling)
			1.015620, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.673322,

			// D75 -> D55 (Bradford)
			1.049904, 0.028885, -0.009591,
			0.022560, 0.993940, 0.015697,
			-0.051476, -0.017431, 0.745981,

			// D75 -> D55 (Von Kries)
			1.016680, 0.006225, 0.000000,
			0.056676, 0.995628, 0.000000,
			-0.053572, -0.001255, 0.751356,

			// D75 -> D55 (XYZ Scaling)
			1.007215, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.751359,

			// D75 -> D65 (Bradford)
			1.020813, 0.011641, -0.004330,
			0.009243, 0.998329, 0.007144,
			-0.022785, -0.007655, 0.886059,

			// D75 -> D65 (Von Kries)
			1.006992, 0.002609, 0.000000,
			0.023758, 0.998167, 0.000000,
			-0.023919, -0.000526,  0.888531,

			// D75 -> D65 (XYZ Scaling)
			1.001134, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.888534,

			// D75 -> D75
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,

			// D75 -> E (Bradford)
			1.072560, 0.051258, -0.006403,
			0.036583, 0.971617, 0.009183,
			-0.044763, -0.016548, 0.813408,

			// D75 -> E (Von Kries)
			1.023456, 0.008754, 0.000000,
			0.079698, 0.993854, 0.000000,
		-0.041900, -0.001766, 0.815937,

			// D75 -> E (XYZ Scaling)
			1.053308, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.815940,

			// E ->	A (Bradford)
			1.154755, 0.110892, -0.023397,
			0.082246, 0.937839, 0.036653,
			-0.138722, -0.048732, 0.342214,

			// E -> A (Von Kries)
			1.052848, 0.020457, 0.000000,
			0.186247, 0.983669, 0.000000,
			-0.140810, -0.004127, 0.355469,

			// E -> A (XYZ Scaling)
			1.098280, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.355470,

			// E -> B (Bradford)
			1.012951, 0.005123, -0.005910,
			0.005370, 1.003671, 0.010188,
			-0.027601, -0.008795, 0.847953,
		
			// E -> B (Von Kries)
			1.004757, 0.001841, 0.000000,
			0.016766, 0.998529, 0.000000,
			-0.030798, -0.000371, 0.852227,

			// E -> B (XYZ Scaling)
			0.990720, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.852230,

			// E -> C (Bradford)
			0.962209, -0.026032, 0.006709,
			-0.019703, 1.012944, -0.010744,
			0.037905, 0.013088, 1.185066,
		
			// E -> C (Von Kries)
			0.987012, -0.005028, 0.000000,
			-0.045772, 1.004013, 0.000000,
			0.039174, 0.001014, 1.181026,

			// E -> C (XYZ Scaling)
			0.980410, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.181030,

			// E -> D50 (Bradford)
			0.997754, -0.009768, -0.007417,
			-0.004163, 1.018316, 0.013442,
			-0.029371, -0.008549, 0.819186,

			// E -> D50 (Von Kries)
			0.999887, -0.000044, 0.000000,
			-0.000393, 1.000033, 0.000000,
			-0.035270, 0.000010, 0.825207,

			// E -> D50 (XYZ Scaling)
			0.964220, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.825210,

			// E -> D55 (Bradford)
			0.979467, -0.022009, -0.003832,
			-0.013568, 1.023820, 0.007632,
			-0.009659, -0.001811, 0.917050,

			// E -> D55 (Von Kries)
			0.993574, -0.002488, 0.000000,
			-0.022644, 1.001984, 0.000000,
			-0.014684, 0.000503, 0.920847,

			// E -> D55 (XYZ Scaling)
			0.956240, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 0.920850,

			// E -> D65 (Bradford)
			0.953167, -0.038259, 0.002612,
			-0.026600, 1.028843, -0.003042,
			0.023901, 0.009415, 1.089399,

			// E -> D65 (Von Kries)
			0.984385, -0.006045, 0.000000,
			-0.055029, 1.004824, 0.000000,
			0.021116, 0.001220, 1.088965,

			// E -> D65 (XYZ Scaling)
			0.950467, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.088969,

			// E -> D75 (Bradford)
			0.934355, -0.049157, 0.007910,
			-0.035658, 1.030889, -0.011919,
			0.050694, 0.018268, 1.229589,

			// E -> D75 (Von Kries)
			0.977754, -0.008612, 0.000000,
			-0.078398, 1.006874, 0.000000,
			0.050039, 0.001738, 1.225576,

			// E -> D75 (XYZ Scaling)
			0.949390, 0.000000, 0.000000,
			0.000000, 1.000000, 0.000000,
			0.000000, 0.000000, 1.225580,

			// E -> E
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0,
			
			1.0, 0.0, 0.0,
			0.0, 1.0, 0.0,
			0.0, 0.0, 1.0

		};

		return CAmatrix[index];
	}

// public methods
public:

	static double gamma(RGBSpaces rgb_space)
	{
		if (rgb_space == ADOBE_RGB)
		{
			return 2.2;
		}
		else if (rgb_space == APPLE_RGB)
		{
			return 1.8;
		}
		else if (rgb_space == BEST_RGB)
		{
			return 2.2;
		}
		else if (rgb_space == BETA_RGB)
		{
			return 2.2;
		}
		else if (rgb_space == BRUCE_RGB)
		{
			return 2.2;
		}
		else if (rgb_space == CIE_RGB)
		{
			return 2.2;
		}
		else if (rgb_space == DON_RGB4)
		{
			return 1.8;
		}
		else if (rgb_space == ECI_RGB)
		{
			return 2.2;
		}
		else if (rgb_space == EKTA_SPACE)
		{
			return 2.2;
		}
		else if (rgb_space == NTSC_RGB)
		{
			return 2.2;
		}
		else if (rgb_space == PAL_RGB)
		{
			return 2.2;
		}
		else if (rgb_space == PROPHOTO)
		{
			return 1.8;
		}
		else if (rgb_space == SMPTE_C)
		{
			return 2.2;
		}
		else if (rgb_space == SRGB)
		{
			return 2.2; // about 2.2
		}
		else if (rgb_space == WIDE_GAMUT)
		{
			return 2.2;
		}
		else
		{
			assert(false);
			return 0.0;
		}
	}

	static Illuminant refIlluminant(RGBSpaces rgb_space)
	{
		// RGB WORKING SPACE DATA
		// 
		// Illuminant Reference for each RGB Space.
		//
		
		if (rgb_space == ADOBE_RGB)
		{
			return VCG_ILLUMINANT_D65;
		}
		else if (rgb_space == APPLE_RGB)
		{
			return VCG_ILLUMINANT_D65;
		}
		else if (rgb_space == BEST_RGB)
		{
			return VCG_ILLUMINANT_D50;
		}
		else if (rgb_space == BETA_RGB)
		{
			return VCG_ILLUMINANT_D50;
		}
		else if (rgb_space == BRUCE_RGB)
		{
			return VCG_ILLUMINANT_D65;
		}
		else if (rgb_space == CIE_RGB)
		{
			return VCG_ILLUMINANT_E;
		}
		else if (rgb_space == DON_RGB4)
		{
			return VCG_ILLUMINANT_D50;
		}
		else if (rgb_space == ECI_RGB)
		{
			return VCG_ILLUMINANT_D50;
		}
		else if (rgb_space == EKTA_SPACE)
		{
			return VCG_ILLUMINANT_D50;
		}
		else if (rgb_space == NTSC_RGB)
		{
			return VCG_ILLUMINANT_C;
		}
		else if (rgb_space == PAL_RGB)
		{
			return VCG_ILLUMINANT_D65;
		}
		else if (rgb_space == PROPHOTO)
		{
			return VCG_ILLUMINANT_D50;
		}
		else if (rgb_space == SMPTE_C)
		{
			return VCG_ILLUMINANT_D65;
		}
		else if (rgb_space == SRGB)
		{
			return VCG_ILLUMINANT_D65;
		}
		else if (rgb_space == WIDE_GAMUT)
		{
			return VCG_ILLUMINANT_D50;
		}
		else
		{
			assert(false);
			return VCG_ILLUMINANT_INVALID;
		}
	}

	static Color4<T> RGBtoHSV(const Color4<T> & color)
	{
		double h,s,v;
		RGBtoHSV(static_cast<double>(color[0]), static_cast<double>(color[1]),
      static_cast<double>(color[2]), h,s,v);
		Color4<T> c(h,s,v,color[3]);
		return c;
	}

	static void RGBtoHSV(double R, double G, double B, double &H, double &S, double &V)
	{
    double v_min = std::min(std::min(R, G), B);  // Min value of RGB
    double v_max = std::max(std::max(R, G), B);  // Max value of RGB

		double delta = v_max - v_min;   //Delta RGB value

		V = v_max;

		if (delta < 0.00000000001)   // This is a gray, no chroma...
		{
			H = 0.0;                   // HSV results = 0 <20> 1
			S = 0.0;
		}
		else                         // Chromatic data...
		{
			S = delta / v_max;

			double deltaR = (((v_max - R) / 6) + (delta/2.0) ) / delta;
			double deltaG = (((v_max - G) / 6) + (delta/2.0) ) / delta;
			double deltaB = (((v_max - B) / 6) + (delta/2.0) ) / delta;

			if ( R == v_max ) 
				H = deltaB - deltaG;
			else if (G == v_max) 
				H = (1.0 / 3.0) + deltaR - deltaB;
			else if (B == v_max) 
				H = (2.0 / 3.0) + deltaG - deltaR;
      else H = 0; // this final branch should never happen...

      if ( H < 0 ) H += 1.0;
      if ( H > 1 ) H -= 1.0;
		}
	}

	static Color4<T> HSVtoRGB(const Color4<T> & color)
	{
		double r,g,b;
		HSVtoRGB(static_cast<double>(color[0]), static_cast<double>(color[1]),
      static_cast<double>(color[2]), r,g,b);
		Color4<T> c(r,g,b,color[3]);
		return c;
	}

	static void HSVtoRGB(double H, double S, double V, double &R, double &G, double &B)
	{
		if (S == 0)
		{
			R = V;
			G = V;
			B = V;
		}
		else
		{
			double var_h = H * 6.0;

			if (var_h == 6.0)
				var_h = 0.0;     // H must be < 1

			int var_i = static_cast<int>(var_h);
			double var_1 = V * (1.0 - S);
			double var_2 = V * (1.0 - S * (var_h - var_i ));
			double var_3 = V * (1.0 - S * (1.0 - (var_h - var_i)));

			double var_r, var_g, var_b;

			if (var_i == 0) 
			{ 
				var_r = V; 
				var_g = var_3; 
				var_b = var_1;
			}
			else if (var_i == 1) 
			{
				var_r = var_2;
				var_g = V; 
				var_b = var_1;
			}
			else if (var_i == 2) 
			{ 
				var_r = var_1;
				var_g = V;
				var_b = var_3;
			}
			else if (var_i == 3)
			{ 
				var_r = var_1;
				var_g = var_2;
				var_b = V;
			}
			else if (var_i == 4)
			{
				var_r = var_3;
				var_g = var_1;
				var_b = V;
			}
			else
			{
				var_r = V;
				var_g = var_1;
				var_b = var_2;
			}

			R = var_r;
			G = var_g;
			B = var_b;
		}

	}

	static Color4<T> XYZtoRGB(const Color4<T> & color, Illuminant src, 
		RGBSpaces dest, ConeResponse response = BRADFORD)
	{
		double r,g,b;
		XYZtoRGB(static_cast<double>(color[0]), static_cast<double>(color[1]),
      static_cast<double>(color[2]), src, r,g,b, dest, response);
		Color4<T> c(r,g,b,color[3]);
		return c;
	}

	static void XYZtoRGB(double X, double Y, double Z, Illuminant src, 
		double &R, double &G, double &B, RGBSpaces space, ConeResponse response = BRADFORD)
	{
		double Xp, Yp, Zp;
		chromaticAdaptation(X, Y, Z, src, Xp, Yp, Zp, refIlluminant(space), response);

		int index = static_cast<int>(space) * 3 * 3;
		double r = Xp * XYZ2RGB(index)   + Yp * XYZ2RGB(index+3) + Zp * XYZ2RGB(index+6);
		double g = Xp * XYZ2RGB(index+1) + Yp * XYZ2RGB(index+4) + Zp * XYZ2RGB(index+7);
		double b = Xp * XYZ2RGB(index+2) + Yp * XYZ2RGB(index+5) + Zp * XYZ2RGB(index+8);

		// Account for gamma correction
		if (space == SRGB)
		{
			if (r > 0.0031308)
				R = 1.055 * pow(r, 1.0/2.4) - 0.055;
			else
				R = 12.92 * r;

			if (g > 0.0031308)
				G = 1.055 * pow(g, 1.0/2.4) - 0.055;
			else
				G = 12.92 * g;

			if (b > 0.0031308)
				B = 1.055 * pow(b, 1.0/2.4) - 0.055;
			else
				B = 12.92 * b;
		}
		else
		{
			double lambda = gamma(space);
			R = pow(r, 1.0/lambda);
			G = pow(g, 1.0/lambda);
			B = pow(b, 1.0/lambda);
		}
	}

	static Color4<T> RGBtoXYZ(const Color4<T> & color, RGBSpaces space, 
		Illuminant dest, ConeResponse response = BRADFORD)
	{
		double x,y,z;
		RGBtoXYZ(static_cast<double>(color[0]), static_cast<double>(color[1]), 
			static_cast<double>(color[2]), space, x,y,z, dest, response);
		Color4<T> c(x,y,z,color[3]);
		return c;
	}

	static void RGBtoXYZ(double R, double G, double B, RGBSpaces space, 
		double &X, double &Y, double &Z, Illuminant dest, ConeResponse response = BRADFORD)
	{
		// Account for Gamma Correction
		double r,g,b;
		if (space == SRGB)
		{
			if (R <= 0.04045)
				r = R / 12.92;
			else
				r = pow((R + 0.055) / 1.055, 2.4);

			if (G <= 0.04045)
				g = G / 12.92;
			else
				g = pow((G + 0.055) / 1.055, 2.4);

			if (B <= 0.04045)
				b = B / 12.92;
			else
				b = pow((B + 0.055) / 1.055, 2.4);
		}
		else
		{
			double lambda = gamma(space);
			r = pow(R, lambda);
			g = pow(G, lambda);
			b = pow(B, lambda);
		}

		int index = static_cast<int>(space) * 3 * 3;
		double Xt = r * RGB2XYZ(index) + g * RGB2XYZ(index+3) + b * RGB2XYZ(index+6);
		double Yt = r * RGB2XYZ(index+1) + g * RGB2XYZ(index+4) + b * RGB2XYZ(index+7);
		double Zt = r * RGB2XYZ(index+2) + g * RGB2XYZ(index+5) + b * RGB2XYZ(index+8);

		// Convert from reference illuminant to the desired one
		chromaticAdaptation(Xt, Yt, Zt, refIlluminant(space), X, Y, Z, dest, response);
	}

	static Color4<T> XYZtoCIELab(const Color4<T> & color, Illuminant ref)
	{
		double L,a,b;
		XYZtoCIELab(static_cast<double>(color[0]), static_cast<double>(color[1]),
			static_cast<double>(color[2]), L,a,b, ref);
		Color4<T> c(L,a,b,color[3]);
		return c;
	}

	static void XYZtoCIELab(double X, double Y, double Z, double &L, double &a, double &b,
		Illuminant whiteRef)
	{	
		// Reference white computation 
		// Note: The reference white depends only on the illuminant
		//       (it is independent on the RGB Working Space)
		double Xr, Yr, Zr; 
		RGBtoXYZ(1.0, 1.0, 1.0, SRGB, Xr, Yr, Zr, whiteRef);

		double xr = X / Xr;
		double yr = Y / Yr;
		double zr = Z / Zr;

		double cieEpsilon = CIE_EPSILON();
		double cieKi = CIE_KI();

		double fx;
		if (xr > cieEpsilon)
			fx = pow(xr, 1.0/3.0);
		else
			fx = (cieKi * xr + 16.0) / 116.0;

		double fy;
		if (yr > cieEpsilon)
			fy = pow(yr, 1.0/3.0);
		else
			fy = (cieKi * yr + 16.0) / 116.0;

		double fz;
		if (zr > cieEpsilon)
			fz = pow(zr, 1.0/3.0);
		else fz = (cieKi * zr + 16.0) / 116.0;

		L = 116.0 * fy - 16.0;
		a = 500.0 * (fx - fy);
		b = 200.0 * (fy - fz);
	}

	static Color4<T> CIELabtoXYZ(const Color4<T> & color, Illuminant ref)
	{
		double x,y,z;
		CIELabtoXYZ(static_cast<double>(color[0]), static_cast<double>(color[1]),
			static_cast<double>(color[2]), x,y,z, ref);
		Color4<T> c(x,y,z,color[3]);
		return c;
	}

	static void CIELabtoXYZ(double L, double a, double b, double &X, double &Y, double &Z,
		Illuminant whiteRef)
	{
		double cieEpsilon = CIE_EPSILON();
		double cieKi = CIE_KI();

		double value;

		double yr;
		if (L > cieKi * cieEpsilon)
		{
			value = (L + 16.0) / 116.0;
			yr = value * value * value;
		}
		else
			yr = L / cieKi;

		double fy;
		if (yr > cieEpsilon)
			fy = (L + 16.0) / 116.0;
		else
			fy = (cieKi * yr + 16.0) / 116.0;

		double fz = fy - (b / 200.0);
		double fz_cubed = fz * fz * fz;

		double fx = (a / 500.0) + fy;
		double fx_cubed = fx * fx * fx;

		double xr;
		if (fx_cubed > cieEpsilon)
			xr = fx_cubed;
		else
			xr = (116.0 * fx - 16.0) / cieKi;

		double zr;
		if (fz_cubed > cieEpsilon)
			zr = fz_cubed;
		else 
			zr = (116.0 * fz - 16.0) / cieKi;

		// Reference white computation 
		// Note: The reference white depends only on the illuminant
		//       (it is independent on the RGB Working Space)
		double Xr, Yr, Zr; 
		RGBtoXYZ(1.0, 1.0, 1.0, SRGB, Xr, Yr, Zr, whiteRef);

		X = xr * Xr;
		Y = yr * Yr;
		Z = zr * Zr;
	}

	// RGB --> HSL
	static Color4<T> RGBtoHSL(const Color4<T> & color)
	{
		double h,s,l;
		RGBtoHSL(static_cast<double>(color[0]), static_cast<double>(color[1]),
			static_cast<double>(color[2]), h,s,l);
		Color4<T> c(h,s,l,color[3]);
		return c;
	}

	// RGB --> HSL
	static void RGBtoHSL(double R, double G, double B, double &H, double &S, double &L)
	{
		double v = std::min(R,G);
		double v_min = std::min(v, B);    // Min value of RGB

		v = std::max(R,G);
		double v_max = std::max(v, B);    // Max value of RGB

		double delta = v_max - v_min;         // Delta RGB value

		L = (v_max + v_min) / 2.0;

		if (delta == 0.0)                     // This is a gray, no chroma...
		{
			H = 0.0;
			S = 0.0;
		}
		else                                  // Chromatic data...
		{
			if ( L < 0.5 ) 
				S = delta / (v_max + v_min);
			else
				S = delta / (2.0 - v_max - v_min);

			double deltaR = (((v_max - R) / 6.0) + (delta/2.0)) / delta;
			double deltaG = (((v_max - G) / 6.0) + (delta/2.0)) / delta;
			double deltaB = (((v_max - B) / 6.0) + (delta/2.0)) / delta;

			if (R == v_max) 
				H = deltaB - deltaG;
			else if (G == v_max) 
				H = (1.0 / 3.0) + deltaR - deltaB;
			else if (B == v_max) 
				H = (2.0 / 3.0) + deltaG - deltaR;

			if ( H < 0.0 )
				H += 1.0;

			if ( H > 1.0 )
				H -= 1.0;
		}
	}

	static Color4<T> HSLtoRGB(const Color4<T> & color)
	{
		double r,g,b;
		HSLtoRGB(static_cast<double>(color[0]), static_cast<double>(color[1]),
			static_cast<double>(color[2]), r,g,b);
		Color4<T> c(r,g,b,color[3]);
		return c;
	}
	
	// HSL --> RGB
	static void HSLtoRGB(double H, double S, double L, double &R, double &G, double &B)
	{
		if (S == 0.0)
		{
			R = L;
			G = L;
			B = L;
		}
		else
		{
			double var_1, var_2;

			if (L < 0.5) 
				var_2 = L * ( 1 + S );
			else
				var_2 = ( L + S ) - ( S * L );

			var_1 = 2 * L - var_2;

			R = Hue2RGB(var_1, var_2, H + (1.0/3.0));
			G = Hue2RGB(var_1, var_2, H );
			B = Hue2RGB(var_1, var_2, H - (1.0/3.0));
		}

	}

	static double Hue2RGB(double v1, double v2, double vH)
	{
		if ( vH < 0 ) 
			vH += 1.0;

		if ( vH > 1 ) 
			vH -= 1.0;

		if ( (6.0 * vH) < 1.0) 
			return (v1 + ( v2 - v1 ) * 6.0 * vH);
		
		if ( (2.0 * vH) < 1.0) 
			return v2;

		if ( (3.0 * vH) < 2.0)
			return ( v1 + ( v2 - v1 ) * ((2.0/3.0) - vH) * 6.0);
		
		return v1;
	}

	static Color4<T> xyYtoXYZ(const Color4<T> & color)
	{
		double X,Y,Z;
		xyYtoXYZ(static_cast<double>(color[0]), static_cast<double>(color[1]),
			static_cast<double>(color[2]), X,Y,Z);

		Color4<T> c(X,Y,Z,color[3]);
		return c;
	}

	// CIE xyY --> CIE XYZ
	static void xyYtoXYZ(double x, double y, double _Y, double &X, double &Y, double &Z)
	{
		X = x * (Y / y);
		Y = _Y;
		Z = (1.0 - x - y) * (Y / y);
	}

	static Color4<T> XYZtoxyY(const Color4<T> & color)
	{
		double x,y,Y;
		XYZtoxyY(static_cast<double>(color[0]), static_cast<double>(color[1]),
			static_cast<double>(color[2]), x, y, Y);

		Color4<T> c(x,y,Y,color[3]);
		return c;
	}

	// CIE XYZ --> CIE xyY
	static void XYZtoxyY(double X, double _Y, double Z, double &x, double &y, double &Y)
	{
		x = X / (X + Y + Z);
		y = Y / (X + Y + Z);
		Y = _Y;
	}

	// XYZ (Illuminant src) --> XYZ (Illuminant dest) - [ILLUMINANT CHANGE]
	static Color4<T> chromaticAdaptation(const Color4<T> & color, Illuminant src, 
		Illuminant dst, ConeResponse response = BRADFORD)
	{
		double X,Y,Z;
    chromaticAdaptation(static_cast<double>(color[0]), static_cast<double>(color[1]),
			static_cast<double>(color[2]), src, X, Y, Z, dst, response);

		Color4<T> c(X,Y,Z,color[3]);
		return c;
	}

	// XYZ (Illuminant src) --> XYZ (Illuminant dest) - [ILLUMINANT CHANGE]
	static void chromaticAdaptation(double Xi, double Yi, double Zi, Illuminant src,
		double &X, double &Y, double &Z, Illuminant dst, ConeResponse response = BRADFORD)
	{
		int index = static_cast<int>(src) * 8 * 3 * 3 * 3 + static_cast<int>(dst) * 3 * 3 * 3 +
			static_cast<int>(response) * 3 * 3;

		X = CA(index)   * Xi + CA(index+3) * Yi + CA(index+6) * Zi;
		Y = CA(index+1) * Xi + CA(index+4) * Yi + CA(index+7) * Zi;
		Z = CA(index+2) * Xi + CA(index+5) * Yi + CA(index+8) * Zi;
	}

	// RGB (working space src) --> RGB (working space dest) - [RGB WORKING SPACE CHANGE]
	static Color4<T> RGBtoRGB(const Color4<T> & color, RGBSpaces rgbsrc, RGBSpaces rgbdest, 
		ConeResponse response = BRADFORD)
	{
		double R,G,B;
		RGBtoRGB(static_cast<double>(color[0]), static_cast<double>(color[1]),
			static_cast<double>(color[2]), rgbsrc, R,G,B, rgbdest, response);

		Color4<T> c(R,G,B,color[3]);
		return c;
	}

	// RGB (working space src) --> RGB (working space dest) - [RGB WORKING SPACE CHANGE]
	static void RGBtoRGB(double Ri, double Gi, double Bi, RGBSpaces rgbsrc,
		double &R, double &G, double &B, RGBSpaces rgbdest, ConeResponse response = BRADFORD)
	{
		double X,Y,Z;
		RGBtoXYZ(Ri, Gi, Bi, rgbsrc, X,Y,Z, refIlluminant(rgbsrc));
		XYZtoRGB(X,Y,Z, refIlluminant(rgbsrc), R,G,B, rgbdest);
	}
	
};


} // end of NameSpace VCG

#endif /*  __VCGLIB_COLORSPACE */