颜色迁移之三——Welsh经典算法

Welsh等人在Reinhard等人的彩色图像间色彩迁移算法的基础上，提出了灰度图像彩色化的思想，并提出了相应的算法。该算法主要利用查找匹配像素来实现灰度图像的色彩迁移，因为灰度图像只有亮度信息，所以该算法主要通过像素的亮度值匹配，来实现灰度图像的自动彩色化。

具体步骤如下：

(1)将参考图像和灰度图像分别由RGB空间转换到lαβ色彩空间。

(2)根据灰度图像的亮度及标准差，对参考图像进行亮度重映射。

由于参考图像和灰度图像的直方图取值不一定在同一范围，那么如果将一幅亮度值很低的参考图像与一幅亮度值很高的目标图像的直方图进行对比将会出现较大的误差。因此，需要对参考图像进行亮度重映射：

L = (nl’ / nl)* (l – ml) + ml’

其中，l是源图像l通道的数据，L是变换后得到新的源图像l通道的值，ml和ml’ 分别是源图像和着色图像的l通道的均值，nl和nl’表示它们的l通道标准方差。

(3)从参考图像中随机选取一批样本点，将像素点的亮度和邻域范围内亮度的标准差的线性组合值作为权值，计算公式如下：

W = l/2 + σ/2

其中，w为权值，l为像素点的亮度值，σ为该像素点周围某个邻域内亮度值的标准差。关于邻域的大小，Welsh指出一般取5X5，对个别图像取更大的邻域效果会更好。

(4)按行扫描灰度图像，对于每个像素点，由公式计算其权值，并在参考图像的样本点中找到一个最接近的样本点，即最佳匹配点，将匹配点的仅和口通道的值赋给灰度图像的对应像素点，从而实现色彩的传递。

(5)将参考图像和灰度图像从lαβ空间转换到RGB空间。

                                     

图 welsh效果图

Welsh等人的算法能有效地实现灰度图像的着色，但仍存在一些不足之处。首先，最佳匹配样本点的查找过程所花费时间较长，从而降低了算法的效率；其次，算法仅利用亮度和标准差来进行最佳匹配，对亮度与颜色对应的一致性要求很高。

代码在下面，执行起来花费时间较长，测试时请做好心理准备~~~

[cpp]  view plain  copy

 

1. void C***View::OnMenuWelsh() //Welsh灰度图像迁移算法  
2. {  
3.     int x,y,r,t;  
4.     for(x=0;x<m_bmp.lWidth;x++)  
5.        for(y=0;y<m_bmp.lHeight;y++)  
6.        {  
7.             //Luminance Remapping:使得颜色图像的l分量（等同于灰度分量）与形状图像保持一致  
8.             Carray_L[y][x]=(SsigmaL/CsigmaL)*(Carray_L[y][x]-CavgL)+SavgL;            
9.        }  
10.   
11.     for(t=0;t<m_bmp.lWidth;t++)  
12.         for(r=0;r<m_bmp.lHeight;r++)  
13.        {  
14.             //计算颜色图像5×5邻域统计量  
15.             int ns,a,b;  
16.             if(r>1&&r<m_bmp.lHeight-2&&t>1&&t<m_bmp.lWidth-2)  
17.             {     
18.                 ns=25;  
19.                 for(a=0;a<5;a++)  
20.                     for(b=0;b<5;b++)  
21.                     {  
22.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a+r-2][b+t-2];                             
23.                     }  
24.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
25.             }  
26.   
27.             else if(r==0&&t>1&&t<m_bmp.lWidth-2)  
28.             {  
29.                 ns=0;  
30.                 for(a=r;a<=r+2;a++)  
31.                     for(b=t-2;b<=t+2;b++)  
32.                     {  
33.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
34.                         ns++;  
35.                     }  
36.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;           
37.             }  
38.             else if(r==1&&t>1&&t<m_bmp.lWidth-2)  
39.             {  
40.                 ns=0;  
41.                 for(a=r-1;a<=r+2;a++)  
42.                     for(b=t-2;b<=t+2;b++)  
43.                     {  
44.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
45.                         ns++;  
46.                     }  
47.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;           
48.             }  
49.             else if(r==m_bmp.lHeight-1&&t>1&&t<m_bmp.lWidth-2)  
50.             {  
51.                 ns=0;  
52.                 for(a=r-2;a<=r;a++)  
53.                     for(b=t-2;b<=t+2;b++)  
54.                     {  
55.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
56.                         ns++;  
57.                     }  
58.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
59.             }  
60.             else if(r==m_bmp.lHeight-2&&t>1&&t<m_bmp.lWidth-2)  
61.             {  
62.                 ns=0;  
63.                 for(a=r-2;a<=r+1;a++)  
64.                     for(b=t-2;b<=t+2;b++)  
65.                     {  
66.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
67.                         ns++;  
68.                     }  
69.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
70.             }  
71. //  
72.             else if(t==0&&r>1&&r<m_bmp.lHeight-2)  
73.             {  
74.                 ns=0;  
75.                 for(a=r-2;a<=r+2;a++)  
76.                     for(b=t;b<=t+2;b++)  
77.                     {  
78.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
79.                         ns++;  
80.                     }  
81.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
82.             }  
83.             else if(t==1&&r>1&&r<m_bmp.lHeight-2)  
84.             {  
85.                 ns=0;  
86.                 for(a=r-2;a<=r+2;a++)  
87.                     for(b=t-1;b<=t+2;b++)  
88.                     {  
89.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
90.                         ns++;  
91.                     }  
92.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
93.             }  
94.             else if(t==m_bmp.lWidth-1&&r>1&&r<m_bmp.lHeight-2)  
95.             {  
96.                 ns=0;  
97.                 for(a=r-2;a<=r+2;a++)  
98.                     for(b=t-2;b<=t;b++)  
99.                     {  
100.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
101.                         ns++;  
102.                     }  
103.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
104.             }  
105.             else if(t==m_bmp.lWidth-2&&r>1&&r<m_bmp.lHeight-2)  
106.             {  
107.                 ns=0;  
108.                 for(a=r-2;a<=r+2;a++)  
109.                     for(b=t-2;b<=t+1;b++)  
110.                     {  
111.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
112.                         ns++;  
113.                     }  
114.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
115.             }  
116. /  
117.             else if(r==0&&t==0)  
118.             {  
119.                 ns=0;  
120.                 for(a=r;a<=r+2;a++)  
121.                     for(b=t;b<=t+2;b++)  
122.                     {  
123.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
124.                         ns++;  
125.                     }  
126.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
127.             }  
128.             else if(r==1&&t==0)  
129.             {  
130.                 ns=0;  
131.                 for(a=r-1;a<=r+2;a++)  
132.                     for(b=t;b<=t+2;b++)  
133.                     {  
134.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
135.                         ns++;  
136.                     }  
137.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
138.             }  
139.             else if(r==m_bmp.lHeight-2&&t==0)  
140.             {  
141.                 ns=0;  
142.                 for(a=r-2;a<=r+1;a++)  
143.                     for(b=t;b<=t+2;b++)  
144.                     {  
145.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
146.                         ns++;  
147.                     }  
148.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
149.             }  
150.             else if(r==m_bmp.lHeight-1&&t==0)  
151.             {  
152.                 ns=0;  
153.                 for(a=r-2;a<=r;a++)  
154.                     for(b=t;b<=t+2;b++)  
155.                     {  
156.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
157.                         ns++;  
158.                     }  
159.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
160.             }  
161. /  
162.             else if(r==0&&t==1)  
163.             {  
164.                 ns=0;  
165.                 for(a=r;a<=r+2;a++)  
166.                     for(b=t-1;b<=t+2;b++)  
167.                     {  
168.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
169.                         ns++;  
170.                     }  
171.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
172.             }  
173.             else if(r==1&&t==1)  
174.             {  
175.                 ns=0;  
176.                 for(a=r-1;a<=r+2;a++)  
177.                     for(b=t-1;b<=t+2;b++)  
178.                     {  
179.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
180.                         ns++;  
181.                     }  
182.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
183.             }  
184.             else if(r==m_bmp.lHeight-2&&t==1)  
185.             {  
186.                 ns=0;  
187.                 for(a=r-2;a<=r+1;a++)  
188.                     for(b=t-1;b<=t+2;b++)  
189.                     {  
190.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
191.                         ns++;  
192.                     }  
193.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
194.             }  
195.             else if(r==m_bmp.lHeight-1&&t==1)  
196.             {  
197.                 ns=0;  
198.                 for(a=r-2;a<=r;a++)  
199.                     for(b=t-1;b<=t+2;b++)  
200.                     {  
201.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
202.                         ns++;  
203.                     }  
204.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
205.             }  
206. /  
207.             else if(r==0&&t==m_bmp.lWidth-2)  
208.             {  
209.                 ns=0;  
210.                 for(a=r;a<=r+2;a++)  
211.                     for(b=t-2;b<=t+1;b++)  
212.                     {  
213.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
214.                         ns++;  
215.                     }  
216.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
217.             }  
218.             else if(r==1&&t==m_bmp.lWidth-2)  
219.             {  
220.                 ns=0;  
221.                 for(a=r-1;a<=r+2;a++)  
222.                     for(b=t-2;b<=t+1;b++)  
223.                     {  
224.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
225.                         ns++;  
226.                     }  
227.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
228.             }  
229.             else if(r==m_bmp.lHeight-2&&t==m_bmp.lWidth-2)  
230.             {  
231.                 ns=0;  
232.                 for(a=r-2;a<=r+1;a++)  
233.                     for(b=t-2;b<=t+1;b++)  
234.                     {  
235.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
236.                         ns++;  
237.                     }  
238.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
239.             }  
240.             else if(r==m_bmp.lHeight-1&&t==m_bmp.lWidth-2)  
241.             {  
242.                 ns=0;  
243.                 for(a=r-2;a<=r;a++)  
244.                     for(b=t-2;b<=t+1;b++)  
245.                     {  
246.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
247.                         ns++;  
248.                     }  
249.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
250.             }  
251.               
252.             else if(r==0&&t==m_bmp.lWidth-1)  
253.             {  
254.                 ns=0;  
255.                 for(a=r;a<=r+2;a++)  
256.                     for(b=t-2;b<=t;b++)  
257.                     {  
258.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
259.                         ns++;  
260.                     }  
261.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
262.             }  
263.             else if(r==1&&t==m_bmp.lWidth-1)  
264.             {  
265.                 ns=0;  
266.                 for(a=r-1;a<=r+2;a++)  
267.                     for(b=t-2;b<=t;b++)  
268.                     {  
269.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
270.                         ns++;  
271.                     }  
272.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
273.             }  
274.             else if(r==m_bmp.lHeight-2&&t==m_bmp.lWidth-1)  
275.             {  
276.                 ns=0;  
277.                 for(a=r-2;a<=r+1;a++)  
278.                     for(b=t-2;b<=t;b++)  
279.                     {  
280.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
281.                         ns++;  
282.                     }  
283.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
284.             }  
285.             else if(r==m_bmp.lHeight-1&&t==m_bmp.lWidth-1)  
286.             {  
287.                 ns=0;  
288.                 for(a=r-2;a<=r;a++)  
289.                     for(b=t-2;b<=t;b++)  
290.                     {  
291.                         Cavg_L[r][t]=Cavg_L[r][t]+Carray_L[a][b];     
292.                         ns++;  
293.                     }  
294.                 CSigma_L[r][t]=Cavg_L[r][t]/ns;  
295.             }  
296.     }  
297.   
298. ///全局搜索匹配//  
299.     double appr,apprtemp;  
300.     int apprnum[2];  
301.     int flag=0;  
302.     CDC *dc=GetDC();  
303.     t_bmp=s_bmp;  
304.     apprnum[0]=0;  
305.     apprnum[1]=0;  
306.     for(int a=0;a<s_bmp.lWidth;a++)  
307.        for(int b=0;b<s_bmp.lHeight;b++)  
308.        {  
309.            appr=(CSigma_L[0][0]+Carray_L[0][0])-(SSigma_L[b][a]+Sarray_L[b][a]);  
310.            appr=appr*appr;  
311.            for(t=0;t<m_bmp.lWidth;t++)  
312.               for(r=0;r<m_bmp.lHeight;r++)  
313.                {  
314.                    apprtemp=(CSigma_L[r][t]+Carray_L[r][t])-(SSigma_L[b][a]+Sarray_L[b][a]);  
315.                    apprtemp=apprtemp*apprtemp;  
316.                    if(appr>=apprtemp)     
317.                    {  
318.                        flag++;  
319.                        appr=apprtemp;  
320.                        apprnum[0]=r;  
321.                        apprnum[1]=t;  
322.                    }  
323.                }  
324.             Tarray_L[b][a]=Sarray_L[b][a];  
325.             Tarray_A[b][a]=Carray_A[apprnum[0]][apprnum[1]];  
326.             Tarray_B[b][a]=Carray_B[apprnum[0]][apprnum[1]];  
327.   
328.   
329.            l=Tarray_L[b][a]/sqrt(3)+Tarray_A[b][a]/sqrt(6)+Tarray_B[b][a]/sqrt(2);  
330.            m=Tarray_L[b][a]/sqrt(3)+Tarray_A[b][a]/sqrt(6)-Tarray_B[b][a]/sqrt(2);  
331.            s=Tarray_L[b][a]/sqrt(3)-2*Tarray_A[b][a]/sqrt(6);  
332.              
333.            //取反对数，10次方    
334.            l=pow(10,l);  
335.            m=pow(10,m);  
336.            s=pow(10,s);  
337.                  
338.            //转换回RGB颜色空间  
339.            t_bmp.pArray_R[b][a]=4.4687*l-3.5887*m+0.1196*s;  
340.            t_bmp.pArray_G[b][a]=-1.2197*l+2.3831*m-0.1626*s;  
341.            t_bmp.pArray_B[b][a]=0.0585*l-0.2611*m+1.2057*s;    
342.           
343.            //控制RGB值>255溢出  
344.            if(t_bmp.pArray_R[b][a]>255) {t_bmp.pArray_R[b][a]=255;}  
345.            if(t_bmp.pArray_G[b][a]>255) {t_bmp.pArray_G[b][a]=255;}  
346.            if(t_bmp.pArray_B[b][a]>255) {t_bmp.pArray_B[b][a]=255;}  
347.            //控制RGB值<0溢出  
348.            if(t_bmp.pArray_R[b][a]<0) {t_bmp.pArray_R[b][a]=0;}               
349.            if(t_bmp.pArray_G[b][a]<0) {t_bmp.pArray_G[b][a]=0;}               
350.            if(t_bmp.pArray_B[b][a]<0) {t_bmp.pArray_B[b][a]=0;}  
351.   
352.   
353.            dc->SetPixel(a+t_bmp.lWidth+m_bmp.lWidth+20,b,RGB(t_bmp.pArray_R[b][a],t_bmp.pArray_G[b][a],t_bmp.pArray_B[b][a]));  
354.        }  
355.           
356. }  

附下载链接

点击打开链接