提取图像的骨架(Skeleton)算法

图像的骨架似乎没有严格的数学定义，可认为是图像细化（Thinning）的产物（中轴可以看作一种骨架，其有严格的数学定义）。目前已经有许多细化算法，这些算法得到的骨架可能略有差异。本文实现了Khalid Sheed 的 K3M算法。该算法属于迭代腐蚀边界的一类算法，该类算法的思想是，假定从二值图像中物体的边界处同时开始燃烧，物体就会被逐步细化，但在燃烧过程中要保证满足一定条件的点被保留或者被“烧掉”，以确定燃烧结束后，剩下的最后一像素宽的图像为图像的骨架。这些条件的确定没有统一的标准，各个算法采取了不同的方案。一般来讲，为了满足计算的速度要求和算法的准确，迭代中算法会对图像边界上某点的3*3邻域内进行检查，判断是否满足要求。

K3M算法在每次迭代中需要进行六次检查

Phase 0. 标记出图像的边界，

Phase 1. 如果该点的邻域中有3个点（非零，以下皆如此）相邻，删除该点

Phase 2. 如果该点的邻域中有3或4个点相邻，删除该点。

Phase 3. 如果该点的邻域中有3，4，5个点相邻，删除该点。

Phase 4. 如果该点的邻域中有3，4，5，6个点相邻，删除该点。

Phase 5. 如果该点的邻域中有3，4，5，6，7个点相邻，删除该点。

Phase 6. 剩余的边界点取消标记，如果Phase 5中没有点被修改，停止迭代，否则返回Phase 0。

具体的步骤可以阅读论文：http://matwbn./ksiazki/amc/amc20/amc2029.pdf。论文中算法实现的一个小技巧是，对邻域中的8个点的值看作二进制，即二进制编码，这样不同的值就对应邻域中不同的状态。迭代中通过计算即可判断该点是否满足条件，是否可以删除。具体细节请移步阅读论文。

算法的测试结果如下：

参考：

http://homepages.inf./rbf/HIPR2/skeleton.htm

http://home./~saeed/arts/2001%20CAIP.pdf

http://www.cs./~algorith/files/thinning.shtml

http://matwbn./ksiazki/amc/amc20/amc2029.pdf

[cpp] view plain copy

1. #include <opencv2/imgproc/imgproc.hpp>  
2. #include <opencv2/objdetect/objdetect.hpp>  
3. #include <opencv2/highgui/highgui.hpp>  
4. #include<vector>  
5. #include<iostream>  
6. #include<algorithm>  
7.   
8. using std::vector;  
9.   
10. vector<int> GetFlags(int a[],int length)  
11. {  
12.     vector<int> vec;  
13.     int neighbour[]={1,2,4,8,16,32,64,128,1,2,4,8,16,32,64};  
14.     for(int i=0;i<length;i++)  
15.     {  
16.         for(int j=0;j<8;j++)  
17.         {  
18.             int sum=0;  
19.             for(int k=j;k<=j+a[i];k++)  
20.                 sum+=neighbour[k];  
21.             vec.push_back(sum);  
22.             std::cout<<sum<<" ";  
23.         }  
24.     }  
25.     std::cout<<std::endl;  
26.     return vec;  
27. }  
28. void skeleton(cv::Mat &Input) //Input-binary image  
29. {  
30.     int a0[]={1,2,3,4,5,6};  
31.     int a1[]={2};  
32.     int a2[]={2,3};  
33.     int a3[]={2,3,4};  
34.     int a4[]={2,3,4,5};  
35.     int a5[]={2,3,4,5,6};  
36.     vector<int> A0=GetFlags(a0,6);  
37.   
38.     vector<int> A1=GetFlags(a1,1);  
39.   
40.     vector<int> A2=GetFlags(a2,2);  
41.     vector<int> A3=GetFlags(a3,3);  
42.     vector<int> A4=GetFlags(a4,4);  
43.     vector<int> A5=GetFlags(a5,5);  
44.     vector<cv::Point2i> border;  
45.     bool modify=true;  
46.     int neighbour[3][3]={  
47.         {128,1,2},  
48.         {64,0,4},  
49.         {32,16,8}  
50.     };  
51.     int row=Input.rows;  
52.     int col=Input.cols;  
53.     while(modify)  
54.     {  
55.         modify=false;  
56.         // flag the border Pharse 0  
57.         for(int m=1;m<row-1;++m)  
58.         {  
59.             for(int n=1;n<col-1;++n)  
60.             {  
61.                 int weight=0;  
62.                 for(int j=-1;j<=1;++j)  
63.                 {  
64.                     for(int k=-1;k<=1;k++)  
65.                     {  
66.                         weight+=neighbour[j+1][k+1]*Input.at<uchar>(m+j,n+k);  
67.                     }  
68.                 }  
69.                 if(std::find(A0.begin(),A0.end(),weight)!=A0.end())  
70.                     border.push_back(cv::Point2i(m,n));  
71.             }  
72.         }  
73.         //Pharse 1  
74.         vector<cv::Point2i>::iterator first=border.begin();  
75.         while(first!=border.end())  
76.         {  
77.             int weight=0;  
78.             for(int j=-1;j<=1;++j)  
79.             {  
80.                 for(int k=-1;k<=1;k++)  
81.                 {  
82.                     weight+=neighbour[j+1][k+1]*Input.at<uchar>((*first).x+j,(*first).y+k);  
83.                 }  
84.             }  
85.             if(std::find(A1.begin(),A1.end(),weight)!=A1.end())  
86.             {  
87.                 Input.at<uchar>((*first).x,(*first).y)=0;  
88.                 first=border.erase(first);  
89.             }  
90.             else  
91.             ++first;  
92.         }  
93.         //Pharse2  
94.         first=border.begin();  
95.         while(first!=border.end())  
96.         {  
97.             int weight=0;  
98.             for(int j=-1;j<=1;++j)  
99.             {  
100.                 for(int k=-1;k<=1;k++)  
101.                 {  
102.                     weight+=neighbour[j+1][k+1]*Input.at<uchar>((*first).x+j,(*first).y+k);  
103.                 }  
104.             }  
105.             if(std::find(A2.begin(),A2.end(),weight)!=A2.end())  
106.             {  
107.                 Input.at<uchar>((*first).x,(*first).y)=0;  
108.                 first=border.erase(first);  
109.             }  
110.             else  
111.             ++first;  
112.         }  
113.         //Pharse3  
114.         first=border.begin();  
115.         while(first!=border.end())  
116.         {  
117.             int weight=0;  
118.             for(int j=-1;j<=1;++j)  
119.             {  
120.                 for(int k=-1;k<=1;k++)  
121.                 {  
122.                     weight+=neighbour[j+1][k+1]*Input.at<uchar>((*first).x+j,(*first).y+k);  
123.                 }  
124.             }  
125.             if(std::find(A3.begin(),A3.end(),weight)!=A3.end())  
126.             {  
127.                 Input.at<uchar>((*first).x,(*first).y)=0;  
128.                 first=border.erase(first);  
129.             }  
130.             else  
131.             ++first;  
132.         }  
133.         //Pharse4  
134.         first=border.begin();  
135.         while(first!=border.end())  
136.         {  
137.             int weight=0;  
138.             for(int j=-1;j<=1;++j)  
139.             {  
140.                 for(int k=-1;k<=1;k++)  
141.                 {  
142.                     weight+=neighbour[j+1][k+1]*Input.at<uchar>((*first).x+j,(*first).y+k);  
143.                 }  
144.             }  
145.             if(std::find(A4.begin(),A4.end(),weight)!=A4.end())  
146.             {  
147.                 Input.at<uchar>((*first).x,(*first).y)=0;  
148.                 first=border.erase(first);  
149.             }  
150.             else  
151.             ++first;  
152.         }  
153.         //Pharse5  
154.         first=border.begin();  
155.         while(first!=border.end())  
156.         {  
157.             int weight=0;  
158.             for(int j=-1;j<=1;++j)  
159.             {  
160.                 for(int k=-1;k<=1;k++)  
161.                 {  
162.                     weight+=neighbour[j+1][k+1]*Input.at<uchar>((*first).x+j,(*first).y+k);  
163.                 }  
164.             }  
165.             if(std::find(A5.begin(),A5.end(),weight)!=A5.end())  
166.             {  
167.                 Input.at<uchar>((*first).x,(*first).y)=0;  
168.                 first=border.erase(first);  
169.                 modify=true;  
170.             }  
171.             else  
172.             ++first;  
173.         }  
174.         //Pharse6  
175.         border.clear();  
176.     }  
177.     for(int m=1;m<row-1;++m)  
178.     {  
179.         for(int n=1;n<col-1;++n)  
180.         {  
181.             int weight=0;  
182.             for(int j=-1;j<=1;++j)  
183.             {  
184.                 for(int k=-1;k<=1;k++)  
185.                 {  
186.                     weight+=neighbour[j+1][k+1]*Input.at<uchar>(m+j,n+k);  
187.                 }  
188.             }  
189.             if(std::find(A0.begin(),A0.end(),weight)!=A0.end())  
190.                     Input.at<uchar>(m,n)=0;;  
191.         }  
192.     }  
193.   
194. }  
195. int main()  
196. {  
197.     cv::Mat raw=cv::imread("test.bmp");  
198.     cv::Mat image(raw.rows,raw.cols,CV_8UC1);  
199.     cv::cvtColor(raw,image,CV_RGB2GRAY);  
200.     cv::Mat binaryImage(image.rows,image.cols,CV_8UC1);  
201.     cv::threshold(image,binaryImage,150,1,CV_THRESH_BINARY_INV);  
202.   
203.     cv::imshow("input",image);  
204.     skeleton(binaryImage);  
205.         for(int p=0;p<binaryImage.rows;p++)  
206.         {  
207.             for(int q=0;q<binaryImage.cols;q++)  
208.             {  
209.                 if(binaryImage.at<uchar>(p,q)==1)  
210.                     binaryImage.at<uchar>(p,q)=0;  
211.                 else  
212.                     binaryImage.at<uchar>(p,q)=255;  
213.             }  
214.         }  
215.     cv::imshow("output",binaryImage);  
216.   
217.     cv::waitKey(0);  
218.     return 0;  
219. }