经典论文回顾：Palette-based Photo Recoloring

Palette-based Photo Recoloring 是图像重着色领域的经典论文，发表于SIGGRAPH 2015，目前引用量180+。因为跟我目前的工作高度相关，所以今天简要回顾一下这篇论文。

先上一张论文的teaser，下图第1列是输入图像和提取的调色板，第2、3、4列展示了修改后的调色板以及原始图像的重着色效果，看上去效果很不错。基于调色板的图像重着色主要解决两个问题：1）调色板提取，2）基于调色板的图像重着色。针对调色板的提取，目前主要存在两类方法，分别是基于聚类的方法和基于几何凸包的方法。这篇论文是基于聚类的方法提取调色板。

1、调色板提取

作者使用改进的K-Means算法提取图像的调色板。直接在所有图像像素上聚类复杂度较高，为了加速聚类过程，作者对图像空间进行采样，并在这些采样点上进行聚类。具体而言，作者将整个RGB颜色空间归一化到$[0,1{]}^3$，并将RGB颜色空间平均划分为$16\times 16\times 16$个网格。然后，统计每一个网格内部的像素数目 $n_i$，并计算网格内像素的平均颜色 $c_i$ 代表该网格的颜色。需要注意的是，这里的网格颜色需要从RGB空间转化到Lab空间，一方面作者认为Lab空间比RGB空间感知更加均匀，另一方面后续对调色板的操作会涉及到Lab中的L通道。

采样点搞定后，接下来就是选取初始的聚类中心。由于K-Means算法对初值较为敏感，为了选取差异较大的初始颜色。作者设计了一种类似于最远点采样的方式选取初始聚类中心。首先给每个采样颜色$C_i$赋权值$n_i$，接下来迭代的选取权重最大的颜色加入到初始聚类中心集合，因此，第一次就选取包含最多像素的网格的平均颜色作为第1个聚类中心 $C_1$；接下来，对其他候选颜色权重进行更新：
$$n_j=(1-exp(-d_{ij}^2/{\sigma }_{\alpha }^2))n_j\text{\hspace{1em}(1)}$$
其中，$d_{ij}$表示Lab空间中$C_i$与$C_j$的距离，$\sigma$是一个经验值，这里取值为80。容易看出来，权重跟$d_{ij}$正相关，颜色差异越大，分数越高，后续越有可能被选入初始聚类中心集合，反之则不太可能选入初始聚类中心。重复该过程，直到所有K个颜色都加入到初始值集合。

最后一步就是使用K-Means对所有采样颜色进行聚类。聚类算法仅针对有效的采样点（内部没有像素的网格被丢弃）实施，聚类速度跟图像像素数目无关，只跟采样点数目有关，因此可以快速收敛。最终收敛的聚类中心组成图像调色板。需要注意的是，图像中的黑色一般难以调整，因此一般默认先将黑色放入调色板内，待聚类收敛后，再将黑色剔除。

2、图像重着色

基于调色板的图像编辑操作简单，用户只需要对调色板的颜色进行简单的修改就可以驱动图像重新着色。经过第一步，我们已经提取到输入图像对应的调色板。接下来探讨如何将调色板的变化映射到图像的颜色上。

输入：原始调色板： { C i } \{C_i\} {Ci​}， 修改后的调色板： { C i ′ } \{C'_i\} {Ci′​}
输出：重着色的图像，对任意像素$I_i$ : $I_i^{\prime }=f(I_i)$。 （$I_i^{\prime }$ 即是重着色后的颜色，$f$ 是要求的映射函数）

基于调色板的图像编辑本质上是一个插值问题，文中列出的一些重要的插值性质包括：

1. Interpolation： 图像中的颜色如果跟初始调色板的颜色相同，则能精确的映射到修改后的调色板颜色：$f(C_i)=C_i^{\prime }$。
2. In Gamut: 输出的颜色位于色域内 （比如RGB的话，像素的三个通道均位于[0,255]区间）。
3. Pixel Continuity： 相似的颜色映射后也要相似：$li{m}_{q\to p}f(q)=f(p)$。
4. One-to-One：满足映射$f$是一个单射，不同的输入对应不同的输出：$f(p)=f(q)=>p=q$
5. Monotonicity in L：相对亮度保持：$L(p)<L(q)=>L(f(p))<=L(f(q))$

为了得到满足以上所有性质的映射函数，作者将L通道跟ab通道分开处理。

2.1 亮度通道L变换

之前提到，调色板的所有颜色位于Lab空间，这里首先对调色板的颜色根据亮度排序使得：$L(C_{i<j})<L(C_j)$.
为了满足插值性质5，作者提出了两种方案：

1）当用户修改调色板的某些颜色后，保持整个调色板的相对亮度不变，从而使得图像重着色后的像素保持相对亮度。

1. 修改后的调色板中某个颜色$C_i^{\prime }$亮度增加。 $C_i^{\prime }$之前的所有颜色亮度依然保持递增，对于后续的颜色的亮度L修改为：
   $$L(C_{j>i}^{\prime })=max(L(C_j^{\prime }),L(C_{j-1}^{\prime }))$$
2. 修改后的调色板中某个颜色$C_i^{\prime }$亮度降低。 $C_i^{\prime }$之后的所有颜色亮度依然保持递增，对于之前的颜色的亮度L修改为：
   $$L(C_{j<i}^{\prime })=min(L(C_j^{\prime }),L(C_{j+1}^{\prime }))$$

2）对于任意像素颜色映射时，其亮度可以通过最接近的两个调色板的亮度加权平均取得，比如在原始图像中像素$x$的亮度可以通过两个临近的调色板颜色插值得到：
$$L(x)=\lambda L(C_j)+(1-\lambda )L(C_k)$$
那么当调色板的亮度变化后依然可以保持这种线性权重不变即：
$$L(x^{\prime })=\lambda L(C_j^{\prime })+(1-\lambda )L(C_k^{\prime })$$

2.2 颜色通道ab插值

2.2.1 单个调色板颜色修改

首先讲简单的情况，即调色板只包含一种颜色，且变化为$C\to C^{\prime }$, 那么理想的情况是所有像素都经过一个类似的变化：
$$f_1(x)=x+C^{\prime }-C$$
但可能导致某些位于Lab空间边界上的像素可能越界，超过色域范围。如下图(b)中下方的情况：$x_0=x+C^{\prime }-C$. 为此，作者提出如下的解决方案。

首先找到射线 $C{C}^{\prime }$ 与ab空间边界的交点 $C_b$；
对于任意像素点$x$, 我们令$x_0=x+C^{\prime }-C$， 根据得到的值是否越过色域分为两种情况：

1. 如果 $x_0$ 位于色域内，这里称为far case。 那么通过二分搜索确定射线 $x{x}_0$ 与ab空间边界的交点 $x_b$；
2. 如果 $x_0$ 位于色域外，这里称为near case。 那么通过二分搜索确定射线 $C^{\prime }{x}_0$ 与ab空间边界的交点 $x_b$；

针对任意的像素点$x$, 重着色后的颜色$x^{\prime }$满足：
(1) $x^{\prime }$ 位于线段$x{x}_b$上；（2）如下等式：
$$\frac{||x^{\prime }-x||}{||C^{\prime }-C||}=min(1,\frac{||x_b-x||}{||C_b-C||})\text{\hspace{1em}(2)}$$
因此，进一步分两种情况：
(1) 如果$||x_b-x||>||C_b-C||$, 则：
$$\frac{||x^{\prime }-x||}{||C^{\prime }-C||}=min(1,\frac{||x_b-x||}{||C_b-C||})=1$$
$$x^{\prime }=x+C^{\prime }-C\text{\hspace{1em}(3)}$$
(2) 如果$||x_b-x||<||C_b-C||$, 则：
$$\frac{||x^{\prime }-x||}{||C^{\prime }-C||}=min(1,\frac{||x_b-x||}{||C_b-C||})=\frac{||x_b-x||}{||C_b-C||}$$
$$||x^{\prime }-x||=||x_b-x||\frac{||C^{\prime }-C||}{||C_b-C||}$$
$$x^{\prime }=x+(x_b-x)\frac{||C^{\prime }-C||}{||C_b-C||}\text{\hspace{1em}(4)}$$
其实我没有完全理解这样做的深意，我不知道对于near case，$x_b$ 取为 $x{x}_0$ 与ab空间边界的交点 会如何。

需要注意的是：第2部分计算颜色转换时，作者的实现跟论文描述有较大差异，最主要的区别在于：寻找边界的时候，文中是在LAB空间边界寻找$x_b,C_b$，而具体实现确是将LAB转换到RGB空间，然后看对应的RGB是否越界。

2.2.2 多个调色板颜色修改

如上图 （c）所示，当有多个调色板的颜色改变时，某个点最终的颜色应当受到到多个调色板颜色影响。
$$\Delta x=\sum _i^k{w}_i(x)f_i(x),and\sum _i^k{w}_i(x)=1\text{\hspace{1em}(5)}$$
改变后的颜色为：(这个很重要！！)
$$x^{\prime }=x+\Delta x$$

其中$k$表示调色板的颜色数目。$w_i(x)$表示混合权重，这里采用RBF插值的方法求这些权重:
$$w_i(x)=\sum _j^k{\lambda }_{ij}\varphi (x,C_j)\text{\hspace{1em}(6)}$$
其中,
$$\varphi (x,C_j)=exp(-(x-C_j{)}^2/2{\sigma }_r^2)\text{\hspace{1em}(7)}$$

可以看到每个$w_i$关联到$k$个$\lambda$: ${\lambda }_{i1},{\lambda }_{i2}，....{\lambda }_{ik}$. 因此，所有的权重关联到$k^2$个$\lambda$。因此要求$w$，实际需要求解一个$k\times k$的$\lambda$矩阵。

给定原始调色板: $C_1,C_2,...C_k$，以及修改后的调色板： $C_1^{\prime },C_2^{\prime },...C_k^{\prime }$.
那么对于任意的颜色 $x=C_i$,则一定能精确的转变为 $x^{\prime }=C_i^{\prime }$ (插值性质1定义)。这样的话， $x=C_i$将只受到调色板$C_i\to C_i^{\prime }$的影响，而不会受到其他调色板的颜色变化影响，因此$w_i=1,w_{j\ne i}=0$. 为此，我们可以将公式（6）展开成矩阵乘积的形式，给出如下的线性方程组求解$\lambda$矩阵：

$$\left[\begin{array}{cccc}{\lambda }_{11}& {\lambda }_{12}& ...& {\lambda }_{1k}\\ {\lambda }_{21}& {\lambda }_{22}& ...& {\lambda }_{2k}\\ .& .& .& .\\ {\lambda }_{k1}& {\lambda }_{k2}& ...& {\lambda }_{kk}\end{array}\right]\left[\begin{array}{cccc}{d}_{11}& d_{21}& ...& d_{k1}\\ d_{12}& d_{22}& ...& d_{k2}\\ .& .& .& .\\ d_{1k}& d_{2k}& ...& d_{kk}\end{array}\right]=\left[\begin{array}{cccc}1& 0& ...& 0\\ 0& 1& ...& 0\\ .& .& 1& .\\ 0& 0& ...& 1\end{array}\right]\text{\hspace{1em}(8)}$$

其中，$d_{ij}=\varphi (C_i,C_j)$, 上式右侧为单位矩阵，实际要求的$\lambda$矩阵等于上述第2个矩阵的逆。

3、一些加速方法

文章并不是直接计算所有像素的插值，而是借助之前的网格，先计算所有网格点的颜色，然后再通过三线性插值的方法对网格内部的像素点插值。

实现效果和对比

注：OpenCV中RGB2LUV的方法：

Vec3d RGB2LUV(Vec3d c) {
	Vec3f c_f(c[0], c[1], c[2]);
	Mat3f rgb(c_f), luv;
	cvtColor(rgb, luv, cv::COLOR_RGB2Luv);
	Vec3f LUV(luv.at<Vec3f>(0, 0)[0], luv.at<Vec3f>(0, 0)[1], luv.at<Vec3f>(0, 0)[2]);
	return Vec3d(LUV[0], LUV[1], LUV[2]);
}