二值图像的距离变换研究

桃扇骨 2022-09-27 15:22 202阅读 0赞

[原文地址][Link 1]

\[研究内容\]

二值图像距离变换

\[正文\]

二值图像距离变换的概念由Rosenfeld和Pfaltz于1966年在论文\[1\]中提出，目前广泛应用于计算机图形学，目标识别及GIS空间分析等领域，其主要思想是通过表识空间点(目标点与背景点)距离的过程，最终将二值图像转换为灰度图像。

距离变换按照距离的类型可以分为欧式距离变换(Eudlidean Distance Transfrom)和非欧式距离变换两种，其中，非欧式距离变换又包括棋盘距离变换(Chessboard Distance Transform)，城市街区距离变换(Cityblock Distance Transform)，倒角距离变换(Chamfer Distance Transform)等；

距离变换的主要过程：

假设一幅二值图像I，包含一个连通区域S，其中有目标集O和背景集B，距离图为D，则距离变换的定义如公式1-(1)：

![SouthEast][]                                        1-(1)

具体步骤如下：

1，将图像中的目标像素点分类，分为内部点，外部点和孤立点。

以中心像素的四邻域为例，如果中心像素为目标像素(值为1)且四邻域都为目标像素(值为1)，则该点为内部点。如果该中心像素为目标像素，四邻域为背景像素(值为0)，则该中心点为孤立点，如下图Fig.1所示。除了内部点和孤立点之外的目标区域点为边界点。

![SouthEast 1][]

2，计算图像中所有的内部点和非内部点，点集分别为S1，S2。

3，对于S1中的每一个内部点(x,y)，使用距离公式disf()计算骑在S2中的最小距离，这些最小距离构成集合S3。

4，计算S3中的最大最小值Max,Min。

5，对于每一个内部点，转换后的灰度值G计算如下公式1-(2)所示：

![SouthEast 2][]    1-(2)

其中，S3(x,y)表示S1中的点(x,y)在S2中的最短距离

6，对于孤立点保持不变。

在以上距离变换的过程中，距离函数disf()的选取如果是欧式距离，则该距离变换称为欧式距离变换，依次类推。对于距离的求取，目前主要的距离公式如下：

欧式距离：

![SouthEast 3][]     1-(3)

棋盘距离：

![SouthEast 4][]               1-(4)

城市街区距离：

![SouthEast 5][]  1-(5)

对于欧式距离变换，由于其结果准确，而计算相比非欧式距离变换较为复杂，因此，出现了较多的快速欧式距离变换[算法][Link 2]，这里笔者介绍一种基于3\*3模板的快速欧式距离变换算法(文献2)，具体过程如下：

1，按照从上到下，从左到右的顺序，使用模板如图Fig.2，依次循环遍历图像I，此过程称为前向循环。

![SouthEast 6][]

对于p对应的像素(x,y)，我们计算五个距离：d0,d1,d2,d3,d4：

d0=p(x,y)

d1=p(x-1,y)+disf((x-1,y),(x,y))

d2=p(x-1,y-1)+disf((x-1,y-1),(x,y))

d3=p(x,y-1)+disf((x,y-1),(x,y))

d4=p(x+1,y-1)+disf((x-1,y-1),(x,y))

则p(x,y)变换后的像素值为：

p(x,y)=Min(d0,d1,d2,d3,d4);

使用上述算法得到图像I'。

2，按照从下到上，从右到左的顺序，使用Fig.2所示模板依次循环遍历图像I’，此过程称为后向循环。

对于p对应的像素(x,y)，我们计算五个距离：d0,d5,d6,d7,d8：

d0=p(x,y)

d5=p(x+1,y)+disf((x+1,y),(x,y))

d6=p(x+1,y+1)+disf((x+1,y+1),(x,y))

d7=p(x,y+1)+disf((x,y+1),(x,y))

d8=p(x-1,y+1)+disf((x-1,y+1),(x,y))

则p(x,y)后向变换后的像素值为：

p(x,y)=Min(d0,d5,d6,d7,d8);

使用上述算法得到的图像即为距离变换得到的灰度图像。

以上过程即文献2所示快速欧式距离变换算法。如果我们需要非欧氏距离变换的快速算法，只需要修改文献2算法中的欧式距离公式disf()为非欧式距离公式，如棋盘距离，城市街区距离等，过程依次类推。

对于欧式距离变换算法，相关学者研究了速度更快的倒角距离变换算法，来近似欧式距离变换的效果。具体过程如下：

1，使用前向模板如图Fig.3中左边3\*3模板，对图像从上到下，从左到右进行扫描，模板中心0点对应的像素值如果为0则跳过，如果为1则计算模板中每个元素与其对应的像素值的和，分别为Sum1,Sum2,Sum3,Sum4，Sum5，而中心像素值为这五个和值中的最小值。

2，使用后向模板如图Fig.3中右边的3\*3模板，对图像从下到上，从右到左进行扫描，方法同上。

3，一般我们使用的倒角距离变换模板为3\*3和5\*5，分别如下图所示：

![SouthEast 7][]

\[实验结果\]

实验采用512\*512大小的图像进行[测试][Link 3]，测试PC为64位，Intel(R) Core(TM) i5-3230 CPU, 2.6GHz, 4G RAM，运行环境为VS2008，C\#。

实验结果如下：

![SouthEast 8][]

(a)原图

![Center][]

(b)Euclidean Distance Transfrom

![SouthEast 9][]

![SouthEast 10][]

(d) Chessboard Distance Transform

![SouthEast 11][]

(e) Chamfer Distance Transform

对于以上欧式距离变换与非欧式距离变换，我们做了时间分析，结果如下：

![SouthEast 12][]

对于Table 1的数据，是通过计算50张512\*512大小的图像得到的平均结果，代码未曾优化，距离变换结果均做了均衡化处理，对于不同配置，不同程序语言可能存在一定差异，总体而言，基于3\*3模板的倒角距离变换速度最快，大概是欧氏距离快速算法的一半。

\[参考文献\]

\[1\] Rosenfeld A,PfaltzJ.L, Sequential operations in digital pic ture processing. Journal of ACM,1966, 13(4):471-494.

\[2\] Frank Y.Shih,Yi-Ta Wu, Fast Euclidean distance transformation in two scans using a 3\*3 neighborhood. Journal of Computer Vision and Image Understanding 2004,195–205.

附录

本人使用C\#编写的代码如下：

本人的完整demo，包含参考文献，测试图像，地址：[http://download.csdn.net/detail/trent1985/6841125][http_download.csdn.net_detail_trent1985_6841125]

**\[csharp\]**  [view plain][]   [copy][view plain]

[![在CODE上查看代码片][CODE]][CODE_CODE] [![派生到我的代码片][ico_fork.svg]][ico_fork.svg 1]

1.  /// <summary>  
2.  /// Distance transform for binary image.  
3.  /// </summary>  
4.  /// <param name="src">The source image.</param>  
5.  /// <param name="distanceMode">One parameter to choose the mode of distance transform from 1 to 3, 1 means Euclidean Distance, 2 means CityBlock Distance, 3 means ChessBoard Distance.</param>  
6.  /// <returns>The result image.</returns>  
7.  **public** Bitmap DistanceTransformer(Bitmap src, **int** distanceMode)  
8.  \{  
9.      **int** w = src.Width;  
10.     **int** h = src.Height;  
11.     **double** p0 = 0, p1 = 0, p2 = 0, p3 = 0, p4 = 0, p5 = 0, p6 = 0, p7 = 0, p8 = 0, min = 0;  
12.     **int** mMax = 0, mMin = 255;  
13.     System.Drawing.Imaging.BitmapData srcData = src.LockBits(**new** Rectangle(0, 0, w, h), System.Drawing.Imaging.ImageLockMode.ReadOnly, System.Drawing.Imaging.PixelFormat.Format24bppRgb);  
14.     **unsafe**  
15.     \{  
16.         **byte**\* p = (**byte**\*)srcData.Scan0.ToPointer();  
17.         **int** stride = srcData.Stride;  
18.         **byte**\* pTemp;  
19.         **for** (**int** y = 0; y < h; y++)  
20.         \{  
21.             **for** (**int** x = 0; x < w; x++)  
22.             \{  
23.                 **if** (x > 0 && x < w - 1 && y > 0 && y < h - 1)  
24.                 \{  
25.                     p0 = (p + x \* 3 + y \* stride)\[0\];  
26.                     **if** (p0 != 0)  
27.                     \{  
28.                         pTemp = p + (x - 1) \* 3 + (y - 1) \* stride;  
29.                         p2 = pTemp\[0\] + GetDistance(x, y, x - 1, y - 1, distanceMode);  
30.                         pTemp = p + x \* 3 + (y - 1) \* stride;  
31.                         p3 = pTemp\[0\] + GetDistance(x, y, x, y - 1, distanceMode);  
32.                         pTemp = p + (x + 1) \* 3 + (y - 1) \* stride;  
33.                         p4 = pTemp\[0\] + GetDistance(x, y, x + 1, y - 1, distanceMode);  
34.                         pTemp = p + (x - 1) \* 3 + y \* stride;  
35.                         p1 = pTemp\[0\] + GetDistance(x, y, x - 1, y, distanceMode);  
36.                         min = GetMin(p0, p1, p2, p3, p4);  
37.                         pTemp = p + x \* 3 + y \* stride;  
38.                         pTemp\[0\] = (**byte**)Math.Min(min, 255);  
39.                         pTemp\[1\] = (**byte**)Math.Min(min, 255);  
40.                         pTemp\[2\] = (**byte**)Math.Min(min, 255);  
41.                     \}  
42.                 \}  
43.                 **else**  
44.                 \{  
45.                     pTemp = p + x \* 3 + y \* stride;  
46.                     pTemp\[0\] = 0;  
47.                     pTemp\[1\] = 0;  
48.                     pTemp\[2\] = 0;  
49.                 \}  
50.   
51.             \}  
52.   
53.         \}  
54.         **for** (**int** y = h - 1; y > 0; y--)  
55.         \{  
56.             **for** (**int** x = w - 1; x > 0; x--)  
57.             \{  
58.                 **if** (x > 0 && x < w - 1 && y > 0 && y < h - 1)  
59.                 \{  
60.                     p0 = (p + x \* 3 + y \* stride)\[0\];  
61.                     **if** (p0 != 0)  
62.                     \{  
63.                         pTemp = p + (x + 1) \* 3 + y \* stride;  
64.                         p5 = pTemp\[0\] + GetDistance(x, y, x + 1, y, distanceMode);  
65.                         pTemp = p + (x + 1) \* 3 + (y + 1) \* stride;  
66.                         p6 = pTemp\[0\] + GetDistance(x, y, x + 1, y + 1, distanceMode);  
67.                         pTemp = p + x \* 3 + (y + 1) \* stride;  
68.                         p7 = pTemp\[0\] + GetDistance(x, y, x, y + 1, distanceMode);  
69.                         pTemp = p + (x - 1) \* 3 + (y + 1) \* stride;  
70.                         p8 = pTemp\[0\] + GetDistance(x, y, x - 1, y + 1, distanceMode);  
71.                         min = GetMin(p0, p5, p6, p7, p8);  
72.                         pTemp = p + x \* 3 + y \* stride;  
73.                         pTemp\[0\] = (**byte**)Math.Min(min, 255);  
74.                         pTemp\[1\] = (**byte**)Math.Min(min, 255);  
75.                         pTemp\[2\] = (**byte**)Math.Min(min, 255);  
76.                         mMax = (**int**)Math.Max(min, mMax);  
77.                     \}  
78.                 \}  
79.                 **else**  
80.                 \{  
81.                     pTemp = p + x \* 3 + y \* stride;  
82.                     pTemp\[0\] = 0;  
83.                     pTemp\[1\] = 0;  
84.                     pTemp\[2\] = 0;  
85.                 \}  
86.   
87.             \}  
88.   
89.         \}  
90.         mMin = 0;  
91.         **for** (**int** y = 0; y < h; y++)  
92.         \{  
93.             **for** (**int** x = 0; x < w; x++)  
94.             \{  
95.                 pTemp = p + x \* 3 + y \* stride;  
96.                 **if** (pTemp\[0\] != 0)  
97.                 \{  
98.                     **int** temp = pTemp\[0\];  
99.                     pTemp\[0\] = (**byte**)((temp - mMin) \* 255 / (mMax - mMin));  
100.                     pTemp\[1\] = (**byte**)((temp - mMin) \* 255 / (mMax - mMin));  
101.                     pTemp\[2\] = (**byte**)((temp - mMin) \* 255 / (mMax - mMin));  
102.                 \}  
103.             \}  
104.   
105.         \}  
106.         src.UnlockBits(srcData);  
107.         **return** src;  
108.     \}  
109. \}  
110. /// <summary>  
111. /// Chamfer distance transform(using two 3\*3 windows:forward window434 300 000;backward window 000 003 434).  
112. /// </summary>  
113. /// <param name="src">The source image.</param>  
114. /// <returns>The result image.</returns>  
115. **public** Bitmap ChamferDistancetransfrom(Bitmap src)  
116. \{  
117.     **int** w = src.Width;  
118.     **int** h = src.Height;  
119.     **double** p0 = 0, p1 = 0, p2 = 0, p3 = 0, p4 = 0, p5 = 0, p6 = 0, p7 = 0, p8 = 0, min = 0;  
120.     **int** mMax = 0, mMin = 255;  
121.     System.Drawing.Imaging.BitmapData srcData = src.LockBits(**new** Rectangle(0, 0, w, h), System.Drawing.Imaging.ImageLockMode.ReadOnly, System.Drawing.Imaging.PixelFormat.Format24bppRgb);  
122.     **unsafe**  
123.     \{  
124.         **byte**\* p = (**byte**\*)srcData.Scan0.ToPointer();  
125.         **int** stride = srcData.Stride;  
126.         **byte**\* pTemp;  
127.         **for** (**int** y = 0; y < h; y++)  
128.         \{  
129.             **for** (**int** x = 0; x < w; x++)  
130.             \{  
131.                 **if** (x > 0 && x < w - 1 && y > 0 && y < h - 1)  
132.                 \{  
133.                     p0 = (p + x \* 3 + y \* stride)\[0\];  
134.                     **if** (p0 != 0)  
135.                     \{  
136.                         pTemp = p + (x - 1) \* 3 + (y - 1) \* stride;  
137.                         p2 = pTemp\[0\] + 4;  
138.                         pTemp = p + x \* 3 + (y - 1) \* stride;  
139.                         p3 = pTemp\[0\] + 3;  
140.                         pTemp = p + (x + 1) \* 3 + (y - 1) \* stride;  
141.                         p4 = pTemp\[0\] + 4;  
142.                         pTemp = p + (x - 1) \* 3 + y \* stride;  
143.                         p1 = pTemp\[0\] + 3;  
144.                         min = GetMin(p0, p1, p2, p3, p4);  
145.                         pTemp = p + x \* 3 + y \* stride;  
146.                         pTemp\[0\] = (**byte**)Math.Min(min, 255);  
147.                         pTemp\[1\] = (**byte**)Math.Min(min, 255);  
148.                         pTemp\[2\] = (**byte**)Math.Min(min, 255);  
149.                     \}  
150.                 \}  
151.                 **else**  
152.                 \{  
153.                     pTemp = p + x \* 3 + y \* stride;  
154.                     pTemp\[0\] = 0;  
155.                     pTemp\[1\] = 0;  
156.                     pTemp\[2\] = 0;  
157.                 \}  
158.   
159.             \}  
160.   
161.         \}  
162.         **for** (**int** y = h - 1; y > 0; y--)  
163.         \{  
164.             **for** (**int** x = w - 1; x > 0; x--)  
165.             \{  
166.                 **if** (x > 0 && x < w - 1 && y > 0 && y < h - 1)  
167.                 \{  
168.                     p0 = (p + x \* 3 + y \* stride)\[0\];  
169.                     **if** (p0 != 0)  
170.                     \{  
171.                         pTemp = p + (x + 1) \* 3 + y \* stride;  
172.                         p5 = pTemp\[0\] + 3;  
173.                         pTemp = p + (x + 1) \* 3 + (y + 1) \* stride;  
174.                         p6 = pTemp\[0\] + 4;  
175.                         pTemp = p + x \* 3 + (y + 1) \* stride;  
176.                         p7 = pTemp\[0\] + 3;  
177.                         pTemp = p + (x - 1) \* 3 + (y + 1) \* stride;  
178.                         p8 = pTemp\[0\] + 4;  
179.                         min = GetMin(p0, p5, p6, p7, p8);  
180.                         pTemp = p + x \* 3 + y \* stride;  
181.                         pTemp\[0\] = (**byte**)Math.Min(min, 255);  
182.                         pTemp\[1\] = (**byte**)Math.Min(min, 255);  
183.                         pTemp\[2\] = (**byte**)Math.Min(min, 255);  
184.                         mMax = (**int**)Math.Max(min, mMax);  
185.                     \}  
186.                 \}  
187.                 **else**  
188.                 \{  
189.                     pTemp = p + x \* 3 + y \* stride;  
190.                     pTemp\[0\] = 0;  
191.                     pTemp\[1\] = 0;  
192.                     pTemp\[2\] = 0;  
193.                 \}  
194.   
195.             \}  
196.   
197.         \}  
198.         mMin = 0;  
199.         **for** (**int** y = 0; y < h; y++)  
200.         \{  
201.             **for** (**int** x = 0; x < w; x++)  
202.             \{  
203.                 pTemp = p + x \* 3 + y \* stride;  
204.                 **if** (pTemp\[0\] != 0)  
205.                 \{  
206.                     **int** temp = pTemp\[0\];  
207.                     pTemp\[0\] = (**byte**)((temp - mMin) \* 255 / (mMax - mMin));  
208.                     pTemp\[1\] = (**byte**)((temp - mMin) \* 255 / (mMax - mMin));  
209.                     pTemp\[2\] = (**byte**)((temp - mMin) \* 255 / (mMax - mMin));  
210.                 \}  
211.             \}  
212.   
213.         \}  
214.         src.UnlockBits(srcData);  
215.         **return** src;  
216.     \}  
217. \}  
218. **private** **double** GetDistance(**int** x, **int** y, **int** dx, **int** dy, **int** mode)  
219. \{  
220.     **double** result = 0;  
221.     **switch** (mode)  
222.     \{  
223.         **case** 1:  
224.             result = EuclideanDistance(x, y, dx, dy);  
225.             **break**;  
226.         **case** 2:  
227.             result = CityblockDistance(x, y, dx, dy);  
228.             **break**;  
229.         **case** 3:  
230.             result = ChessboardDistance(x, y, dx, dy);  
231.             **break**;  
232.     \}  
233.     **return** result;  
234. \}  
235. **private** **double** EuclideanDistance(**int** x, **int** y, **int** dx, **int** dy)  
236. \{  
237.     **return** Math.Sqrt((x - dx) \* (x - dx) + (y - dy) \* (y - dy));  
238. \}  
239. **private** **double** CityblockDistance(**int** x, **int** y, **int** dx, **int** dy)  
240. \{  
241.     **return** Math.Abs(x - dx) + Math.Abs(y - dy);  
242. \}  
243. **private** **double** ChessboardDistance(**int** x, **int** y, **int** dx, **int** dy)  
244. \{  
245.     **return** Math.Max(Math.Abs(x - dx), Math.Abs(y - dy));  
246. \}  
247. **private** **double** GetMin(**double** a, **double** b, **double** c, **double** d, **double** e)  
248. \{  
249.     **return** Math.Min(Math.Min(Math.Min(a, b), Math.Min(c, d)), e);  
250. \}

最后，分享一个专业的图像处理网站（微像素），里面有很多源代码下载：

[http://www.zealpixel.com/portal.php][http_www.zealpixel.com_portal.php]

[Link 1]: http://blog.csdn.net/trent1985/article/details/18081761
[SouthEast]: /images/20220708/6d50827386b243ee8c5e42318eb2b44d.png
[SouthEast 1]: /images/20220708/6069e2f86a574c1db0d7e323c7fe49fa.png
[SouthEast 2]: /images/20220708/8d327dbdbb4a4023baace635d6346cc0.png
[SouthEast 3]: /images/20220708/35d660e866404d71963899779fb4ce1c.png
[SouthEast 4]: /images/20220708/a7d8fe8cc4a1448b8ae3be30006487c0.png
[SouthEast 5]: /images/20220708/20f51fe6d28b4661a0fe7996482f5e94.png
[Link 2]: http://lib.csdn.net/base/datastructure
[SouthEast 6]: /images/20220708/b9358f2ec38842fbb5167cdd461f13e8.png
[SouthEast 7]: /images/20220708/42bfb26fc2bd46f6ad992fc82a200306.png
[Link 3]: http://lib.csdn.net/base/softwaretest
[SouthEast 8]: /images/20220708/494c6874588342c780e3c362248fc5c2.png
[Center]: /images/20220708/5de2d94bef1a4999aa6539269d4bae9d.png
[SouthEast 9]: /images/20220708/f5cf7c1b97c04fca8e58c556fc0e9f27.png
[SouthEast 10]: /images/20220708/611621769710479087163ad2b9de540a.png
[SouthEast 11]: /images/20220708/8effdddd0c4f4329914bbe0de3af9b50.png
[SouthEast 12]: /images/20220708/8cf75b56e4d5406397a5914b528662a0.png
[http_download.csdn.net_detail_trent1985_6841125]: http://download.csdn.net/detail/trent1985/6841125
[view plain]: http://blog.csdn.net/trent1985/article/details/18081761#
[CODE]: https://code.csdn.net/assets/CODE_ico.png
[CODE_CODE]: https://code.csdn.net/snippets/153608
[ico_fork.svg]: https://code.csdn.net/assets/ico_fork.svg
[ico_fork.svg 1]: https://code.csdn.net/snippets/153608/fork
[http_www.zealpixel.com_portal.php]: http://www.zealpixel.com/portal.php