也谈CIE-RGB到CIE-XYZ的转换(附MATLAB程序及数据)

也谈CIERGB到CIEXYZ的转换SMT傅亚炜CIE1931XYZ标准色度观察者色度空间（也有文献中叫颜色模型、色彩空间、色度学系统的，都是一个意思），以下简称CIEXYZ色度空间，是国际上色度计算、颜色测量和颜色表征的统一标准，是几乎所有测色仪器的设计与制造依据。CIEXYZ色度空间的三原色X、Y、Z是假想的，实际上并不存在，因此也无法进行实验测试。XYZ系统所使用的数据，都是根据CIE1931RGB色度空间的实验数据（颜色匹配函数或称光谱三刺激值）、、rgb转换而来的：2.76891.75171.13024.59070.06010.05595.5949XRGBYRGBZGB(1)0.490.310.20.6671.13231.20070.1770.81240.01060.6671.13231.20070.00990.99010.6671.13231.2007rgbxrgbrgbyrgbgbzrgb(2)由于计算精度问题，不同文献中的数值会有微小差异（参见《附录一》）。在一般的文献和教材中都是直接引用上述转换公式，很少有给出推导过程的。少数给出过正统推导过程的[1][2][3]，都是从色度图坐标直接先推导出三刺激值，其过程繁琐，且不直观，不利于理解（参见《附录二》）。所以我根据我在MATLAB中绘制色度图的过程，给出了概念清晰的简明推导过程，先从rgb色度坐标导出xyz色度坐标，这样能很直观地看到转换过程对色度图所带来的影响。CIE在推出XYZ色度空间时，提了以下条件：(1)所有三刺激值都是正的，即所有真实色都要落在新的麦克斯韦三角形之内。(2)假想色空间尽量小，即真实色在新麦克斯韦三角形内占的面积尽可能地大。(3)X轴和Z轴上亮度为零，即Y轴代表所有的亮度（光通量）。(4)X、Y、Z相等时同样要合成等能白光E，即白光E的xyz坐标和rgb坐标相同，都是(0.333,0333,0.333)。根据上述要求(1)、(2)、(3)，在rg色度图上画出了如下三角形：其中XZ这条线为无亮度线。CIERGB色度空间的亮度方程为：4.59070.0601CRGB(3)其中R、G、B是三刺激值。而色度坐标r、g、b与三刺激值的关系为：RrRGBGgRGBBbRGB(4)这个变换，好多文献称之为“归一化”，实际上是不准确的，我觉得应该称之为“降维”（想起了《三体》中的二向箔），就是把三维图拍扁，全部挤进一个平面里。不多解释，直接上图：降维--上面右图中的三角形就是rgb空间的麦克斯韦颜色三角，可以看到有很多颜色落在麦克斯韦三角之外；三角形三个顶点就是三原色700nm、546.1nm和435.8nm；而其背后是它的投影，也就是rg色度图。把(4)式代入(3)式，无亮度的情况即ΦC等于0，此时分母中的(R+G+B)可以约掉，原式变为：4.59070.06010rgb(5)考虑到1rgb，将1brg代入上式，即得上式在r-g平面内的切线：0.93994.58060.06010rg(6)XY线设为过700nm的切线，对应方程为：0.9910rg(7)因为从546.1nm到700nm虽然近似为一条直线，但多少有点弧度，所以这条直线比135度线略微有点上翘，与g轴的交点是（0，1.01）。YZ线，所有文献中都告诉我们说这是过503nm的一条切线，原因是条件(2)，就是令假想色空间最小，但没人给出过严格的证明，估计多少有点拍脑袋决定的成分：1.450.5510rg(8)两两联立上述方程(6)、(7)、(8)，利用小学水平的解二元一次方程组可以解得（r，g，b）的三组坐标值：（1.2750，-0.2778，0.0028）、（-1.7391，2.7668，-0.0277）、（-0.7431，0.1409，1.6022）。数值都经过了四舍五入，所以在不同文献中也会略有不同。这三组坐标值的点，我们的目标是要将它们转换为（x，y，z）为（1，0，0）、（0，1，0）、（0，0，1）的点，所以：1.27501.73910.74310.27782.76680.14090.00280.02771.6022rxgybz(9)将（1，0，0）、（0，1，0）、（0，0，1）这三组值代入方程右边，就可以很清楚地得到原先的三组（r，g，b）坐标。对上述矩阵求逆即可得：1110.90880.57490.37090.09120.41880.005500.00620.6236xrygzb(10)之所以加了下标1，是因为这并非我们最终所要的（x，y，z）。如下图所示，图中红点代表等能白光E，其rgb坐标是（0.333，0.333，0.333），但是在上述变换后变成了（0.6182，1.1718，0.2099）。从下图中也可以直观地看出，由于rg色度图中在r为负值的左边区域还有着很大一块，所以在将左图的三角形均匀地变形为右图的三角形之后，白点必然是偏向右下的。所以，通过三个顶点坐标求逆得到的图形，三个顶点坐标是正确的；由于其线性性质，三条边的坐标也都是正确的；但是内部是变形的，不能满足条件(4)。为了将白点E变换到（0.333，0.333，0.333），可以再做一个比例变换，即让(10)式中的矩阵的每一行都乘以一个常数，使其满足：0.3330.90880.57490.38090.3330.3330.09120.41880.00550.3330.33300.00620.62360.333xxxyyyzzCCCCCCCC(11)方程(11)两边都乘以3，即可轻松求得：10.53920.90880.57490.370911.93970.09120.41880.005511.58770.00620.6236xyzCCC(12)将(12)式中的系数乘到(10)式，得到一组新的（x，y，z）：2220.49000.31000.20000.17700.81240.010600.00990.9901xrygzb(13)这个变换，才是真正意义上的“归一化”，即系数矩阵的每一行加起来都等于1。但是（2x，2y，2z）还不是我们所要的（x，y，z）。因为如上图所示，虽然白点的坐标如愿变成了（0.333，0.333，0.333），但是三个顶点的坐标又不对了，不再是（1，0，0）、（0，1，0）、（0，0，1）。所以还得对这个看起来已经是扁扁的图形再做一次“降维”打击，即令：2222222222220.490.310.20.6671.13231.20070.1770.81240.01060.6671.13231.20070.00990.99010.6671.13231.2007xrgbxxyzrgbyrgbyxyzrgbzgbzxyzrgb(14)这也就是一开始所给出的(2)式。把光谱色的rgb数据代进去，画出图形来一看，正是我们熟悉的CIE1931标准观察者色度图，也叫“马蹄图”或“舌形图”：如果将直接求逆所得的图形与归一化再降维之后的图形相比较，可以发现它们外面的麦克斯韦三角形是一致的，但是从绿线的对应关系可以直观地看出，降维之后，三角形内部发生了疏密度的变化。我们所要的xy色度图，其绿色端被进行了压缩，而红色端被进行了放大，从而使得等能白光E点往左上移动，到了（0.333，0.333，0.333）的位置。所以XYZ系统与RGB系统相比，不仅是取消了负值，而且颜色分布也更加均匀了（虽然以后会发现这还不够均匀。而均匀，一直是CIE努力的方向，所以才会有后续的CIE1960、1964和1976色度空间）。如下图所示，标上颜色后比较起来更加直观，可明显看到红色区的扩大：以上就是从CIErgb色度坐标到CIExyz色度坐标的转换过程。直观地来说，就是把左边的三角形拉成一个等边直角三角形，然后往绿色端挤一挤。这个过程，并非CIE的原始推理过程，但是数学意义和物理意义很清晰，用到的数学知识很浅显，也很适合于用图形来演示，有助于帮助人们对色度空间的转换有一个直观的感性认识。再来考虑三刺激值的转换。已知CIERGB色度空间的三刺激值R、G、B与色度坐标r、g、b的关系如(4)式所示：RrRGBGgRGBBbRGB(4)在CIEXYZ色度空间中我们同样定义假想三原色的三刺激值X、Y、Z为：XxXYZYyXYZZzXYZ(15)将(4)式和(15)式代入(14)式，约去分母(R+G+B)，可得：0.490.310.20.6671.13231.20070.1770.81240.01060.6671.13231.20070.00990.99010.6671.13231.2007XRGBxXYZRGBXRGByXYZRGBXGBzXYZRGB(16)联立(14)式和(16)式，因为对于所有的（x，y，z）这两式都成立，不做严格的数学证明，可以直观地看出至少存在如下解：0.49000.31000.2000K0.17700.81240.010600.00990.9901XRYGZB(17)其中K是一个常数。由于上式中的系数矩阵是归一化的，所以当R=G=B=1时，我们有X=Y=Z=K。再利用一次条件(3)，即仅有Y代表亮度，则根据亮度方程(3)，有：4.59070.0601YRGB(18)将R=G=B=1代入上式，可得：K=14.59070.06015.6508Y(19)代回(17)式，即得如(1)式所示的三刺激值转换公式。这个系数K是亮度方程中三个系数之和。它的物理意义是，当R、G、B三原色的光通量都为1光瓦（或1流明）时，混合而得的白光的光通量即为5.6508光瓦（流明）。其数学意义是，亮度方程(18)中，当匹配等能光谱色时，在555nm（黄绿光）处Y达到最大值为1，起到归一化的作用。也就是说在实验的精度范围之内，亮度方程Y与明视觉相对光敏函数（光谱光视效率）V应该是一致的。参考文献[1]周恕义,王少华.1931CIE-RGB系统向XYZ系统的转换[J].纺织基础科学学报,1990,(3):72-75.[2]周太明.色度学讲座[M].上海:复旦大学,2004.[3]薛兵,张晓军.高等学校“十二五”规划教材：应用物理[M].西安:西安电子科技大学出版社,2012.附录一：计算精度问题在一般的文献中，三刺激值转换公式均如下所示：2.76891.75171.13024.59070.06010.05655.5943XRGBYRGBZGB(20)而我利用MATLAB求得的公式是：2.76891.75171.13024.59070.06010.05595.5949XRGBYRGBZGB(1)最大误差出现在Z表达式的G分量系数上，其误差超过1%。究其原因，在于1931年的时候，电子计算机或计算器还没有出现，还是用手工和计算尺进行计算的。在联立解方程(6)、(7)、(8)时，如果将结果四舍五入保留到小数点后4位，即：1.27501.73910.74310.27782.76680.14090.00280.02771.6022rxgybz(9)并且计算过程中的中间数据全部保