拼豆与Potts模型

• 拼豆方案图生成器

拼豆是一种将彩色塑料豆按照网格排列，再通过熨烫融合成图案的手工材料。

问题定义

给定一张原始图片，目标是生成一个固定尺寸的拼豆网格图纸。这个过程中主要面临三个挑战。

1. 原图通常是高分辨率图像，而拼豆图纸只有有限格数。例如一张照片可能包含数百万像素，但最终图纸可能只有 $64\times64$ 个格点。
2. 原图颜色可以近似看作连续 RGB 空间，而真实拼豆品牌只提供有限色卡。例如 MARD、Perler、Hama、Artkal 等品牌通常只有几十到数百种颜色。
3. 如果直接把每个格点独立匹配到最近色卡颜色，结果往往会出现大量孤立杂色点。这些噪点在屏幕上看似细节丰富，但在实际拼豆制作中会增加换色难度，并使成品显得杂乱。

因此，拼豆图纸生成算法需要同时满足以下目标：

1. 每个格子的颜色应尽量接近原图；
2. 相邻格子的颜色应尽量形成稳定色块；
3. 原图中的真实边缘和主体结构不能被过度抹平；
4. 只出现一两颗的稀有颜色应尽量合并，降低制作复杂度。

数学抽象

设最终拼豆图纸是一个二维网格：

\[\Omega=\{(x,y)\mid 0\le x<W,\ 0\le y<H\}\]

其中 $W$ 是图纸宽度，$H$ 是图纸高度。每个格点 $p=(x,y)$ 需要分配一个颜色标签：

\[s_p\in\{1,2,\ldots,K\}\]

其中 $K$ 是当前色卡中的颜色总数。例如某个品牌色卡有 221 种颜色，则 $K=221$。

整个拼豆图纸可以表示为一个标签场：

\[S=\{s_p\}_{p\in\Omega}\]

每个标签 $s_p$ 对应一个真实拼豆颜色，例如 “Perler White” 或 “MARD A1”。

算法的目标是在所有可能的标签场 $S$ 中，寻找一个既接近原图、又足够平滑、同时便于实际制作的配置。

图像离散化

原始图片尺寸通常较大，并且宽高比例不固定。拼豆图纸需要固定格数，因此首先要把原图缩放到目标网格尺寸。

用户指定图纸宽度 $W$，算法根据原图宽高比计算高度：

\[H=\operatorname{round}\left(\frac{W}{\text{image width}/\text{image height}}\right)\]

然后将原图缩放到 $W\times H$ 的小画布上。缩小后，每个格点 $p$ 对应一个原图颜色：

\[I_p=(R_p,G_p,B_p)\]

这一步完成了从连续图像到拼豆离散网格的转换。

有限色卡

拼豆颜色不能任意选择，而只能来自用户指定品牌的色卡。设色卡为：

\[\mathcal{C}=\{c_1,c_2,\ldots,c_K\}\]

其中每个颜色可以表示为：

\[c_k=(R_k,G_k,B_k)\]

不同品牌可能有不同色卡，例如：

• MARD 全色；
• MARD 标准色；
• Artkal S 5mm；
• Perler 标准色；
• Hama 标准色。

生成图纸时，每个格点都必须从这个有限集合中选择颜色。这使得拼豆生成问题变成一个有限状态空间上的离散优化问题。

颜色差异度量：从 RGB 到 CIE Lab

最直接的颜色匹配方法是比较原图颜色和色卡颜色的 RGB 欧氏距离：

\[d_{\text{RGB}}(I_p,c_k)= \sqrt{(R_p-R_k)^2+(G_p-G_k)^2+(B_p-B_k)^2}\]

但 RGB 空间并不适合直接衡量人眼感知到的颜色差异。RGB 是屏幕显示编码空间，每个通道的数值变化并不等价于真实光强变化，人眼对不同颜色方向的敏感度也不同。

因此，算法先将 RGB 转换到更接近人眼感知的 CIE Lab 色彩空间。

一个 Lab 颜色写作：

\[(L^*,a^*,b^*)\]

其中：

• $L^*$ 表示亮度；
• $a^*$ 表示红绿方向，正值偏红，负值偏绿；
• $b^*$ 表示黄蓝方向，正值偏黄，负值偏蓝。

RGB 到 Lab 的转换通常分为四步。

RGB 归一化

给定 RGB 整数值：

\[(R,G,B),\quad R,G,B\in[0,255]\]

先归一化为：

\[r=\frac{R}{255},\quad g=\frac{G}{255},\quad b=\frac{B}{255}\]

此时的 $(r,g,b)$ 仍然是 sRGB 编码值，还不是线性光强。

Gamma 反校正

对 sRGB 做 gamma 反校正，得到线性 RGB：

\[r_{\text{lin}}= \begin{cases} r/12.92, & r\le 0.04045\\ \left(\frac{r+0.055}{1.055}\right)^{2.4}, & r>0.04045 \end{cases}\]

$g_{\text{lin}}$ 和 $b_{\text{lin}}$ 使用相同公式计算。

这一步的作用是将屏幕编码值还原为近似物理光强值。

转换为 CIE XYZ

线性 RGB 可转换为 CIE XYZ：

\[\begin{aligned} X &= 0.4124r_{\text{lin}} + 0.3576g_{\text{lin}} + 0.1805b_{\text{lin}}\\ Y &= 0.2126r_{\text{lin}} + 0.7152g_{\text{lin}} + 0.0722b_{\text{lin}}\\ Z &= 0.0193r_{\text{lin}} + 0.1192g_{\text{lin}} + 0.9505b_{\text{lin}} \end{aligned}\]

其中 $Y$ 近似表示亮度。绿色通道在 $Y$ 中权重最大，这反映了人眼对绿色附近亮度变化更敏感。

转换为 CIE Lab

使用 D65 参考白点：

\[X_n=0.95047,\quad Y_n=1.00000,\quad Z_n=1.08883\]

先进行归一化：

\[x=\frac{X}{X_n},\quad y=\frac{Y}{Y_n},\quad z=\frac{Z}{Z_n}\]

再定义函数：

\[f(t)= \begin{cases} \sqrt[3]{t}, & t>0.008856\\ 7.787t+\frac{16}{116}, & t\le 0.008856 \end{cases}\]

最后得到：

\[L^*=116f(y)-16\] \[a^*=500[f(x)-f(y)]\] \[b^*=200[f(y)-f(z)]\]

Delta E：感知颜色差异

Delta E 76

在 Lab 空间中，颜色差异通常用 Delta E 表示。最简单的形式是 $\Delta E_{76}$，即 Lab 空间中的欧氏距离：

\[\Delta E_{76}= \sqrt{(L_1-L_2)^2+(a_1-a_2)^2+(b_1-b_2)^2}\]

Delta E 94

不过 $\Delta E_{76}$ 对不同亮度、色度和色相区域的感知修正不足。因此，算法中使用更符合人眼感知的 $\Delta E_{94}$。

给定两个 Lab 颜色：

\[A=(L_1,a_1,b_1),\quad B=(L_2,a_2,b_2)\]

定义：

\[\Delta L=L_1-L_2\] \[C_1=\sqrt{a_1^2+b_1^2},\quad C_2=\sqrt{a_2^2+b_2^2}\] \[\Delta C=C_1-C_2\] \[\Delta H^2=(a_1-a_2)^2+(b_1-b_2)^2-\Delta C^2\]

则 $\Delta E_{94}$ 的一种常见形式为：

\[\Delta E_{94}(A,B)= \sqrt{ \Delta L^2+ \left(\frac{\Delta C}{1+0.045C_1}\right)^2+ \left(\frac{\Delta H}{1+0.015C_1}\right)^2 }\]

这里第一个颜色 $A$ 通常被视为标准色，第二个颜色 $B$ 被视为样品色。由于分母使用了 $C_1$，严格来说该公式不满足对称性：

\[\Delta E_{94}(A,B)\ne \Delta E_{94}(B,A)\]

在拼豆算法中，可以根据使用场景分别处理。

原图颜色到色卡颜色的匹配

当计算某个格点与色卡颜色的匹配误差时，方向是固定的：

\[\text{原图颜色}\rightarrow\text{色卡颜色}\]

因此可以直接使用：

\[E_{\text{data}}(p,k)=\Delta E_{94}(I_p,c_k)\]

这表示如果格点 $p$ 使用色卡颜色 $k$，会产生多少颜色误差。

色卡颜色之间的距离

在平滑项中，需要比较两个相邻格点所选色卡颜色之间的差异。因为邻接关系是无向的，所以距离最好是对称的。

可以使用平均色度构造对称版 $\Delta E_{94}$。设：

\[\bar C=\frac{C_1+C_2}{2}\]

则：

\[\Delta E_{94,\text{sym}}(A,B)= \sqrt{ \Delta L^2+ \left(\frac{\Delta C}{1+0.045\bar C}\right)^2+ \left(\frac{\Delta H}{1+0.015\bar C}\right)^2 }\]

这样可以保证：

\[\Delta E_{94,\text{sym}}(A,B)=\Delta E_{94,\text{sym}}(B,A)\]

数据项：颜色拟合代价

对每个格点 $p$ 和每个色卡颜色 $k$，算法计算数据项：

\[E_{\text{data}}(p,k)=\Delta E_{94}(I_p,c_k)\]

它表示：如果在格点 $p$ 放置颜色 $k$，该颜色与原图颜色的感知误差有多大。

如果只考虑数据项，最优策略就是逐格选择最近颜色：

\[s_p=\arg\min_k E_{\text{data}}(p,k)\]

这就是朴素颜色量化。

但是这种方法只考虑单个格点，不考虑周围邻居，因此容易产生大量碎点和孤立杂色。

Potts 模型：引入空间一致性

为了让拼豆图纸形成更稳定的色块，算法引入 Potts 模型。

经典 Potts 模型可以写作：

\[E(S)=\sum_p D_p(s_p)+\lambda\sum_{(p,q)\in\mathcal{N}}V(s_p,s_q)\]

其中：

• $D_p(s_p)$：数据项，表示当前颜色与原图颜色的误差；
• $\mathcal{N}$：邻接关系；
• $V(s_p,s_q)$：相邻格点之间的平滑惩罚；
• $\lambda$：平滑强度。

传统 Potts 模型通常使用简单的“相同/不同”惩罚：

\[V(s_p,s_q)= \begin{cases} 0, & s_p=s_q\\ 1, & s_p\ne s_q \end{cases}\]

但拼豆颜色之间的差异并不是二值的。浅黄和米白虽然不同，但视觉差异较小；黑色和亮黄差异很大，不应被同等看待。

颜色距离平滑 边缘调制 惩罚

算法会预先计算色卡中任意两种颜色之间的感知距离，并构造颜色距离矩阵：

\[C_{ij}=\frac{\Delta E_{94,\mathrm{sym}}(c_i,c_j)}{30}+0.2\]

其中，$c_i$ 和 $c_j$ 表示色卡中的两种颜色，$\Delta E_{94,\mathrm{sym}}$ 表示对称化后的 $\Delta E_{94}$ 颜色差异。

当两种颜色完全相同时，令：

\[C_{ii}=0\]

因此，这个矩阵具有如下含义：相同颜色相邻时不产生平滑惩罚；相近颜色相邻时惩罚较小；差异较大的颜色相邻时惩罚较大。

对于两个相邻格点 $p$ 和 $q$，如果它们当前选择的拼豆颜色标签分别为 $s_p$ 和 $s_q$，则二者之间的颜色差异惩罚为：

\[C_{s_p,s_q}\]

这个设计可以看作传统 Potts 模型的颜色距离扩展。传统 Potts 模型通常只判断两个标签是否相同：相同则惩罚为 $0$，不同则给出固定惩罚。而这里的惩罚不是固定值，而是由两种拼豆颜色的感知差异决定。

边缘调制

如果只使用颜色距离型平滑惩罚，算法会倾向于让相邻格点选择相同或相近颜色。这样可以减少噪声和孤立色点，但也可能破坏人物五官、轮廓线、文字边界、阴影交界等重要结构。

因此，程序会根据原图中相邻格点的颜色差异，对平滑惩罚进行边缘调制。

设原图缩小后，格点 $p$ 和 $q$ 的颜色分别为 $I_p$ 和 $I_q$，定义原图相邻颜色差异为：

\[d_I(p,q)=\Delta E_{94,\mathrm{sym}}(I_p,I_q)\]

边缘保留权重定义为：

\[w_{pq} = \frac{1}{ 1+\exp\left( \frac{d_I(p,q)-\tau}{\sigma} \right) }\]

其中，$\tau$ 是边缘判断的中点阈值，$\sigma$ 控制从“强平滑”到“弱平滑”的过渡宽度。

当原图中两个相邻格点颜色接近时，$d_I(p,q)$ 较小，$w_{pq}$ 较大，平滑作用较强。算法会更倾向于让这两个位置选择相同或相近的拼豆颜色。

当原图中两个相邻格点颜色差异明显时，$d_I(p,q)$ 较大，$w_{pq}$ 较小，平滑作用被削弱。算法允许这条邻接边两侧出现颜色断裂，从而保留轮廓、五官、文字边界和阴影分界等结构。

最终，两个相邻格点之间的平滑代价写作：

\[E_{\mathrm{smooth}}(p,q) = J\cdot w_{pq}\cdot C_{s_p,s_q}\]

展开后为：

\[E_{\mathrm{smooth}}(p,q) = J \cdot \frac{ \frac{\Delta E_{94,\mathrm{sym}}(c_{s_p},c_{s_q})}{30}+0.2 }{ 1+\exp\left( \frac{\Delta E_{94,\mathrm{sym}}(I_p,I_q)-\tau}{\sigma} \right) }\]

其中，$J$ 表示色块聚集强度。

并且当 $s_p=s_q$ 时，仍然令：$ C_{s_p,s_q}=0 $ 相同拼豆颜色相邻不会产生平滑惩罚。

一维示意函数

为了理解这个平滑机制的形状，可以用一个简化的一维函数表示：

\[f(x)= \frac{\frac{x}{30}+0.2}{ 1+e^{\frac{x-8}{2}} }\]

这里，$x$ 可以理解为颜色差异的示意变量。

函数单调性：先略微上升，然后下降到0。

• 当 $x$ 很小时：相邻颜色几乎一致，$f(x)$ 不大，颜色差异惩罚较小。
• 当 $x$ 处于中等范围时：颜色差异惩罚开始增大，但边缘抑制还没有完全生效，因此 $f(x)$ 会达到较高值。这个区域对应“相邻颜色不同，但还不像明确边缘”的情况，算法会倾向于压制这种不必要的碎色差异。
• 当 $x$ 很大时：说明相邻位置很可能跨过了真实边缘。此时指数项迅速增大，分母变大，$f(x)$ 下降。也就是说，算法不再强行平滑这条边，而是允许颜色断开，以保留图像结构。

邻域设计：四邻域与对角邻域

图像网格中的邻居可以分为四邻域和对角邻域。

四邻域包括：

\[(x-1,y),\ (x+1,y),\ (x,y-1),\ (x,y+1)\]

对角邻域包括：

\[(x-1,y-1),\ (x+1,y+1),\ (x+1,y-1),\ (x-1,y+1)\]

算法同时考虑这两类邻域。四邻域的平滑系数为 $J$，对角邻域的平滑系数为：

\[J_{\text{diag}}=0.707J\]

其中 $0.707\approx1/\sqrt{2}$，用于反映对角距离比水平或垂直距离更远。

因此，局部平滑项可以概括为：

\[\sum_{q\in\mathcal{N}_4(p)} J\cdot w_{pq}\cdot C_{s_p,s_q} + \sum_{q\in\mathcal{N}_{\text{diag}}(p)} J_{\text{diag}}\cdot w_{pq}\cdot C_{s_p,s_q}\]

去孤立噪点

拼豆图纸中常见的问题是“孤豆”：某种颜色只在一个位置出现，周围全是其他颜色。这类孤立点通常会增加制作难度，也可能只是颜色量化造成的噪声。

为减少孤豆，算法加入孤立惩罚。对于格点 $p$ 的候选颜色 $k$，统计周围邻居中与它颜色相同的数量：

\[n_{\text{same}}(p,k)=|\{q\in\mathcal{N}(p):s_q=k\}|\]

如果相同邻居太少，就增加惩罚：

\[E_{\text{curve}}(p,k)= \begin{cases} 1.5K_{\text{curve}}, & n_{\text{same}}\le1\\ 0.5K_{\text{curve}}, & n_{\text{same}}=2\\ 0, & n_{\text{same}}>2 \end{cases}\]

其中 $K_{\text{curve}}$ 对应界面中的“去孤立噪点”参数。

总能量函数

综合颜色拟合、空间平滑和孤立惩罚，某个格点 $p$ 选择颜色 $k$ 的局部能量可以写作：

\[H_p(k)= E_{\text{data}}(p,k) + \sum_{q\in\mathcal{N}(p)} J_{pq}w_{pq}C_{k,s_q} + E_{\text{curve}}(p,k)\]

其中：

• $E_{\text{data}}(p,k)$：颜色拟合项；
• $J_{pq}$：方向相关的平滑强度，四邻域为 $J$，对角邻域为 $0.707J$；
• $w_{pq}$：边缘保留权重；
• $C_{k,s_q}$：候选颜色与邻居颜色之间的色卡距离；
• $E_{\text{curve}}(p,k)$：孤立颜色惩罚。

整个图纸的总能量可以抽象为：

\[E(S)= \sum_pE_{\text{data}}(p,s_p) + \sum_{(p,q)\in\mathcal{N}} J_{pq}w_{pq}C_{s_p,s_q} + \sum_pE_{\text{curve}}(p,s_p)\]

目标是寻找低能量配置：

\[S^*=\arg\min_SE(S)\]

这就是拼豆图纸生成问题的 Potts 模型表述。

为什么不能穷举求解

如果图纸有 $N=W\times H$ 个格点，每个格点有 $K$ 种颜色，那么所有可能图纸数量为：

\[K^N\]

例如 $64\times64$ 的图纸有 4096 个格点。如果色卡有 221 种颜色，组合数量为：

\[221^{4096}\]

这是天文数字，不可能穷举搜索。因此，程序采用局部迭代优化方法，而不是直接求全局最优解。

Top-K 候选颜色：降低计算量

如果每次更新都遍历全部 $K$ 个色卡颜色，计算量会很大。为提高效率，算法对每个格点只保留一小部分候选颜色。

首先，按照数据项误差排序，取最接近原图颜色的前 $K_{\text{top}}$ 个颜色：

\[\mathcal{T}_p=\operatorname{TopK}_{k}E_{\text{data}}(p,k)\]

程序中通常取 $K_{\text{top}}=7$。

在实际更新格点 $p$ 时，候选集合包括：

1. 当前颜色；
2. 原图颜色最接近的 Top-K 颜色；
3. 周围邻居已经使用的颜色。

即：

\[\mathcal{A}_p=\{s_p\}\cup\mathcal{T}_p\cup\{s_q:q\in\mathcal{N}(p)\}\]

这样做有两个好处。

一方面，Top-K 颜色保证局部颜色不会偏离原图太远。另一方面，邻居颜色允许已有色块向周围扩张，从而形成更稳定的区域。

四相位更新

在局部能量计算中，程序不仅考虑上下左右四邻域，也考虑四个对角邻域。因此，一个格点 $p$ 的邻域可以写成：

\[\mathcal{N}_8(p) = \mathcal{N}_4(p) \cup \mathcal{N}_{\mathrm{diag}}(p)\]

其中四邻域为：

\[\mathcal{N}_4(p)= \{(x-1,y),(x+1,y),(x,y-1),(x,y+1)\}\]

对角邻域为：

\[\mathcal{N}_{\mathrm{diag}}(p)= \{(x-1,y-1),(x+1,y+1),(x+1,y-1),(x-1,y+1)\}\]

四相位划分策略

将网格划分为四个相位：

\[\mathrm{phase}(x,y) = 2(x\bmod 2)+(y\bmod 2)\]

也就是根据 $x$ 和 $y$ 的奇偶性组合，将格点分为四类：

\[(0,0),\quad (0,1),\quad (1,0),\quad (1,1)\]

这四类格点分别对应四个 phase：

\[0,\quad 1,\quad 2,\quad 3\]

每轮迭代依次更新这四类格点。由于同一相位中的任意两个格点在横向、纵向和对角方向上都至少间隔两个格子，它们不会互为八邻域邻居。因此，四相位更新更适合当前这种包含对角平滑项的 Potts/MRF 模型。

代码实现

程序中的相位划分逻辑如下：

\[\mathrm{phase}(x,y) = 2(x\bmod 2)+(y\bmod 2)\]

对应 JavaScript 代码为：

const cellPhase = 2 * (x & 1) + (y & 1);

完整更新框架为：

for (let phase = 0; phase < 4; phase++) {
  for (let y = 0; y < h; y++) {
    for (let x = 0; x < w; x++) {
      const cellPhase = 2 * (x & 1) + (y & 1);
      if (cellPhase !== phase) continue;

      // 更新当前格点
    }
  }
}

模拟退火

如果每一步都选择局部能量最低的颜色，算法很容易陷入局部最优。例如某个区域一旦形成错误色块，后续可能难以修正。

因此，算法引入温度参数 $T$。初始温度较高，例如：$ T_0=4.0 $

每轮迭代后逐渐降低：$ T\leftarrow0.987T $

在温度较高时，程序不会总是选择最低能量颜色，而是根据 Boltzmann 分布随机选择：

\[P(k)= \frac{ \exp\left(-\frac{H_p(k)-H_{\min}}{T}\right) }{ \sum_{j\in\mathcal{A}_p} \exp\left(-\frac{H_p(j)-H_{\min}}{T}\right) }\]

其中：

\[H_{\min}=\min_{j\in\mathcal{A}_p}H_p(j)\]

当温度较高时，能量稍高的候选颜色仍有一定概率被选中，算法可以探索更多配置。

当温度降低后，概率分布逐渐集中到最低能量颜色，结果趋于稳定。

当温度低于阈值，例如 $T<0.2$ 时，程序可以改为直接选择局部能量最低的颜色。

模拟退火的作用可以概括为：

• 前期保留随机性，避免过早陷入局部最优；
• 后期降低随机性，使图纸逐渐收敛稳定。

程序实际使用轮盘赌采样实现：先计算所有候选颜色的权重总和 $Z$，再从 $[0,Z)$ 中随机取一个数，根据累计权重落点选择颜色。

稀有颜色合并

Potts 优化结束后，图纸中仍可能存在只使用一两颗的颜色。虽然这些颜色可能在数学上有一定拟合意义，但从手工制作角度看，它们会增加购买和摆放复杂度。

因此，程序会检测使用量极少的颜色。若颜色 $k$ 的使用量 $ n_k\le2 $ 则该颜色被视为稀有颜色。

如果用户选择自动合并，程序会将这些稀有颜色替换为非稀有颜色集合中最适合该位置的颜色。

这一步不是 Potts 模型本身的一部分，而是面向实际制作的后处理策略。它的目的不是追求理论上的最优能量，而是减少实际制作成本。