一天买了一些绿色的 LED，到手发现是黄绿色的，于是决定做一个可以更直观更方便地了解 LED 实际颜色的小工具，这是它的理论基础...

准备

在开始前，需要了解一点点简单的知识 ——

三刺激值函数

CIE 的特别委员会定义了 $\overline{x} (\lambda ), \overline{y} (\lambda ), \overline{z} (\lambda )$ 三个匹配函数来定义 CIE XYZ 色彩空间，引入了带有频谱功率分布 $P(\lambda)$ 的颜色的对应的 XYZ 三色刺激值：

\[ \begin{matrix} X = \int\limits_{\lambda} P(\lambda )\overline{x}(\lambda ) \mathrm{d}\lambda \\ Y = \int\limits_{\lambda} P(\lambda )\overline{y}(\lambda ) \mathrm{d}\lambda \\ Z = \int\limits_{\lambda} P(\lambda )\overline{z}(\lambda ) \mathrm{d}\lambda \end{matrix} \]

（积分在可见光谱范围内计算）

XYZ/RGB转换矩阵

对于一个给定的 RGB 色彩空间 $(x_{r}, y_{r}), (x_{g}, y_{g}), (x_{b}, y_{b})$ 和相同参考白色的 $(X_{W}, Y_{W}, Z_{W})$，可以利用 $3 \times 3$ 矩阵互相转换：

\[ \begin{bmatrix} X \\ Y \\ Z \end{bmatrix} = \begin{bmatrix} M \end{bmatrix} \begin{bmatrix} R \\ G \\ B \end{bmatrix} \\ [M] = \begin{bmatrix} S_r X_r & S_g X_g & S_b X_b \\ S_r Y_r & S_g Y_g & S_b Y_b \\ S_r Z_r & S_g Z_g & S_b Z_b \end{bmatrix} \\ X_r = x_r / y_r \\ Y_r = 1 \\ Z_r = (1 - x_r - y_r) / y_r \\ X_g = x_g / y_g \\ Y_g = 1 \\ Z_g = (1 - x_g - y_g) / y_g \\ X_b = x_b / y_b \\ Y_b = 1 \\ Z_b = (1 - x_b - y_b) / y_b \]

可以轻易求得将 XYZ 色彩空间转换为 RGB 色彩空间的矩阵，具体求法留给读者：

\[ \begin{bmatrix} S_r \\ S_g \\ S_b \end{bmatrix} = \begin{bmatrix} X_r & X_g & X_b \\ Y_r & Y_g & Y_b \\ Z_r & Z_g & Z_b\end{bmatrix} ^{-1}\begin{bmatrix}X_W \\Y_W \\Z_W\end{bmatrix} \]

从光谱转到 XYZ 色彩空间

显然我们可以利用三刺激值函数将光的波长映射到 XYZ 色彩空间，那么开始吧！

高斯光谱功率分布

在 CIE 定义的 XYZ 三色刺激值式子中，$P(\lambda)$ 就是大名鼎鼎的高斯分布的概率密度函数，通常它是长这个样的：

\[f(x) = \frac{1}{\sigma \sqrt{2\pi}} e^{-\frac{(x - \mu)^2}{2\sigma^2}}\]

我们可以利用给定的波长（中心波长）和半峰宽（FWHM）来计算光谱功率分布 $P(\lambda)$，稍稍变化一下：

\[ \begin{matrix}P(\lambda ) = \frac{1}{\sigma \sqrt{2\pi}} e^{-\frac{(\lambda - \lambda_0)^2}{2\sigma^2}} \\ \sigma = \frac{FWHM}{2\sqrt{2\ln{2} } } \end{matrix} \]

这样，对于给定的波长 $\lambda$，可以将 $\lambda$，$\lambda_0​$，和 $FWHM$的值代入上述高斯函数中，即可计算出 $P(\lambda)$。下面的代码使用 JS 计算了在给定中心波长（wavelength）和半峰宽（FWHM）下人眼可见光范围内的 $P(\lambda)$

let sigma = FWHM / 2.355;
let c = 1 / (sigma * math.sqrt(2 * math.pi));
let x = math.range(380, 750, 1).toArray();
let P = [];
x.forEach((value) => {
    P.push(c * math.exp(-math.pow(value - wavelength, 2) / (2 * math.pow(sigma, 2))));
});

现在我们已经解决了光谱功率分布这个难题，下面只需要求出对应 XYZ 的值即可

三色刺激值

$\overline{x} (\lambda ), \overline{y} (\lambda ), \overline{z} (\lambda )$ 三个匹配函数事实上并没有确切的函数来描述它们，它们通常来源于实验室测量。例如这是一个离散的列表 ColorMix/tristimulus.txt

因为它们是离散的，那么在实际应用中更多使用的是求和的形式：

\[ \begin{matrix} X = \sum_{i}\overline{x}_{i} P_{i} \Delta \lambda \\ Y = \sum_{i}\overline{y}_{i} P_{i} \Delta \lambda \\ Z = \sum_{i}\overline{z}_{i} P_{i} \Delta \lambda \end{matrix} \]

在刚刚计算光谱功率分布时我们也是按照离散的来计算的，那么只需要将 $\overline{x} (\lambda ), \overline{y} (\lambda ), \overline{z} (\lambda )$ 三个匹配函数在 $\lambda$ 处的取值与其功率概率相乘求和叠加即可！

const tristimulus = [...];  // 这里导入了 https://github.com/HuangdaxianZaixian/ColorMix/blob/master/tristimulus.txt
let T = [];
tristimulus.forEach((value) => {
    T.push(value.slice(1));
});
let XYZ = [];
T.forEach((value, index) => {
    XYZ.push(math.multiply(value, P[index]));
});
let XYZSum = math.sum(XYZ, 0);

现在已经获得了在给定中心波长（wavelength）和半峰宽（FWHM）下所有波长叠加的 XYZ 色彩空间的具体值，最后只需要通过前面的逆矩阵 $[M]^{-1}$ 即可将其映射到 RGB 色彩空间

从 XYZ 色彩空间到 RGB 色彩空间

事实上通常并不需要计算矩阵 $[M]$ 和 $[M]^{-1}$，已经有人完成了这些计算过程（RGB/XYZ Matrices）。这里以 sRGB 色彩空间为例：

\[ [M]^{-1} = \begin{bmatrix} 2.0413690 & -0.5649464 & -0.3446944 \\ -0.9692660 & 1.8760108 & 0.0415560 \\ 0.0134474 & -0.1183897 & 1.0154096\end{bmatrix} \]

看起来直接计算这个矩阵即可：

\[ \begin{bmatrix}r \\ g \\ b \end{bmatrix} = [M]^{-1}\begin{bmatrix}X \\ Y \\ Z\end{bmatrix} \]

但在进行计算前，还需要将 $[X,Y,Z]$ 三个分量进行归一化到 $[0, 1]$

let XYZMax = math.max(...XYZSum);
XYZSum.forEach((value, index) => {
    XYZSum[index] = value / XYZMax;
});

在完成归一化后进行矩阵乘法得到的参量属于线性 RGB 通道（记为 $(r,g,b)$ 和 $v$），需要变换到非线性 RGB（记为 $(R,G,B)$ 和 $V$）

\[ \begin{matrix} v \in {r, g, b} \\ V \in {R, G, B} \end{matrix} \]

对三个通道进行相同的变换即可，这里以 sRGB 为例：

\[ V = \begin{cases} 12.92v & v \le 0.0031308 \\ 1.055v^{1/2.4 } - 0.055 & else \end{cases} \]

RGB = RGB.map((value) => {
    if (value <= 0.0031308) {
        return 12.92 * value;
    } else {
        return 1.055 * math.pow(value, 1 / 2.4) - 0.055;
    }
});

经过变换后（即使不经过变换）输出的 RGB 值是在 $[0,1]$ 范围内的，将其乘以 255 即可得到 $[0, 255]$ 范围内的 RGB 值

但有可能输出的值超出了 $[0,1]$ 范围，小于 0 的部分近似为 0，大于 1 的部分近似于 1 即可

单一波长情况

若只给出了中心波长，则可以将其看作冲击函数 $\delta(\lambda-\lambda_0)$：

\[ \begin{matrix} X = \int\limits_{\lambda} \delta(\lambda-\lambda_0)\overline{x}(\lambda ) \mathrm{d}\lambda \\ Y = \int\limits_{\lambda} \delta(\lambda-\lambda_0)\overline{y}(\lambda ) \mathrm{d}\lambda \\ Z = \int\limits_{\lambda} \delta(\lambda-\lambda_0)\overline{z}(\lambda ) \mathrm{d}\lambda \end{matrix} \]

根据冲激函数的性质 $\int_{-\infty }^{+\infty } \delta(t)\mathrm{d}t = 1$，所以式子变成了：

\[ \begin{matrix} X = \overline{x}(\lambda_0 ) \\ Y = \overline{y}(\lambda_0 ) \\ Z = \overline{z}(\lambda_0 ) \end{matrix} \]

按照上面的过程使用 $[M]^{-1}$ 计算即可

尝试一下！

嘉立创EDA & EasyEDA 专业版扩展

参考

[1] Wikipedia. Cie 1931色彩空间[EB/OL]. 2023[2024-11-22]. https://zh.wikipedia.org/wiki/CIE_1931%E8%89%B2%E5%BD%A9%E7%A9%BA%E9%97%B4.
[2] HuangdaxianZaixian. ColorMix: RGB Leds color mixture[EB/OL]. [2024-11-22]. https://github.com/HuangdaxianZaixian/ColorMix/.
[3] Lindbloom, B.J. Brucelindbloom[EB/OL]. [2024-11-22]. http://brucelindbloom.com/index.html.