OKLAB颜色空间详解 2 ——颜色空间转换、oklab空间中的色域确定

本文主要介绍以下两个方面：

介绍并用python实现了oklab空间与XYZ空间的转换
介绍并用python实现了在oklab空间如何确定一个色域的范围

1. oklab与XYZ互相转换

（1）XYZ向oklab空间转换总共三步，两步矩阵乘法+一步非线性变换，计算简单是oklab的一大优点；
第一步是将XYZ转换到一个近似的锥体细胞反应 lms:
在这里插入图片描述
第二步是一个非线性变换：

第三步是转换至Lab的坐标下：

上式中的M1、M2的数值如下：

（2）从oklab到XYZ的变换是上述过程的逆过程，如下所示：

python实现：

MATRIX_1_XYZ_TO_LMS = np.array([
        [0.8189330101, 0.3618667424, -0.1288597137],
        [0.0329845436, 0.9293118715, 0.0361456387],
        [0.0482003018, 0.2643662691, 0.6338517070],
    ])
MATRIX_2_LMS_TO_LAB = np.array([
    [0.2104542553, 0.7936177850, -0.0040720468],
    [1.9779984951, -2.4285922050, 0.4505937099],
    [0.0259040371, 0.7827717662, -0.8086757660],
])
MATRIX_1_LMS_TO_XYZ = np.linalg.inv(MATRIX_1_XYZ_TO_LMS)
MATRIX_2_LAB_TO_LMS = np.linalg.inv(MATRIX_2_LMS_TO_LAB)

def XYZ_to_oklab(XYZ): 
    LMS = np.einsum('...ij,...j->...i', MATRIX_1_XYZ_TO_LMS, XYZ)
    LMS_prime = np.power(LMS, 1/3)
    lab = np.einsum('...ij,...j->...i', MATRIX_2_LMS_TO_LAB, LMS_prime)
    return lab

def oklab_to_XYZ(oklab):
    LMS = np.einsum('...ij,...j->...i', MATRIX_2_LAB_TO_LMS, oklab)
    LMS_prime = np.power(LMS, 3)
    XYZ = np.einsum('...ij,...j->...i', MATRIX_1_LMS_TO_XYZ, LMS_prime)
    return XYZ

2. oklab颜色空间的色域范围确定

oklab空间拥有众多优秀的品质，我们可以在此空间进行很多图像处理操作，当进行色彩相关操作时，往往我们需要确定色域(gamut)的范围。下面讲一下在oklab空间确定gamut的方法。

首先简单介绍下色域的概念。虽然人眼能够感知的颜色很丰富，但由于技术的限制，拍摄或显示终端所能产生的颜色是有限的，色域标准则定义了颜色的子集合。其经常在CIE 1931 xy Color Coordinates中进行可视化，经过从三维到二维平面的映射后，呈现为三角形。
在这里插入图片描述

比较常见的三个色域为BT709、DCI-P3、BT2020。

在oklab颜色空间，由于其分别建模了亮度、色度和色相，因此可以固定色相，在单一色相上确定色域范围。
方法是固定一个色相后，我们分别对亮度和色度在[0,1]上采样，采样后进行oklch -> oklab -> XYZ -> rgb的转换，
r/g/b在[0,1]范围外，说明对应的oklch也在色域外。

下面为确认bt2020色域的python代码

N = 1000 # N越大色域图画得越细腻，如果想看一个点是不是在色域内，N过大的话，一个点显示不出来。可以根据需求适当调整
def oklab2rgb_wo_clip(oklab):
    XYZ = oklab_to_XYZ(oklab)
    XYZ_to_rgb_matrix_hdr = colour.models.BT2020_COLOURSPACE.XYZ_to_RGB_matrix
    rgb_out_lin = colour.utilities.dot_vector(XYZ_to_rgb_matrix_hdr, XYZ)
    return rgb_out_lin

def draw_gamut_in_oklab(input_h):
    c, l = np.meshgrid(np.linspace(0,1,N), np.linspace(0, 1, N))
    h = np.ones_like(l)*input_h
    a = c*np.cos(h)
    b = c*np.sin(h)
    oklab_30 = np.stack((l, a, b), axis=2)
    rgb_lin = oklab2rgb_wo_clip(oklab_30)
    r = rgb_lin[:,:,0]
    g = rgb_lin[:,:,1]
    b = rgb_lin[:,:,2]
    gamut = np.minimum(r, np.minimum(g, b))
    gamut_2 = np.maximum(r, np.maximum(g, b))
    gamut[gamut_2>1]=-1
    gamut_flip = np.flip(gamut, 0)
    return gamut_flip

def find_cusp(gamut_flip):
    gamut_map = np.sign(gamut_flip)
    gamut_point = np.argwhere(gamut_map==1)
    cusp_c = np.max(gamut_point[:,1])
    cusp_l = gamut_point[:,0][np.argmax(gamut_point[:,1])]
    return np.float(cusp_c)/N, np.float(N-cusp_l)/N


# 比如一个点的lch的值如下：
l_ex = 0.15417300711845494
c_ex = 0.09907472692884031
h_ex = 0.40446415294979943

c_index = np.int(c_ex*N)
l_index = np.int(l_ex*N)
# print(c_index, l_index)

gamut_flip = draw_gamut_in_oklab(h_ex)
cusp_c, cusp_l = find_cusp(gamut_flip)
print("hue: {}, cusp_c:{}, cusp_l:{}".format(str(h_ex), str(cusp_c), str(cusp_l)))


guamut_vis = np.sign(gamut_flip)
guamut_vis[guamut_vis==-1]=64
guamut_vis[guamut_vis==0]=128
guamut_vis[guamut_vis==1]=128
guamut_vis[N-l_index, c_index] = 255 # 因为进行了上下flip，所以在高度上要进行翻转
guamut_vis_c3 = np.stack((guamut_vis, guamut_vis, guamut_vis), 2)
plt.figure(figsize=(20, 20))
plt.imshow(guamut_vis_c3.astype(np.uint8))
plt.show()

画出的图为：
在这里插入图片描述
浅灰色为该色相上的bt2020色域范围；oklab上的色域范围均为三角形形状，作者的blog中有一些例子如下：