domenica 9 settembre 2012

Topological skeleton


For my final work I needed to "skeletonize" an image. 

If you don't know what I'm talking about, you can see the wikiedia page:

en.wikipedia.org/wiki/Topological_skeleton
and
en.wikipedia.org/wiki/Straight_skeleton

I looked for a SIMPLE code, just to know i fthis solution could work on my images and... I didn't find it...


Everytime i found a code, it was ultra-optimized-at least 250 lines of code. I can understand that a people prefer to post just an optimized code, but if you need just something working out-of-the-box is tiringto not understand what the code is doing.

anyway, after tryed twice a code without results, I implemented mine. I know it is not optimized, but it is SIMPLE and IT WORKS.

It's based on the article:
T. Y. Zhang and C. Y. Suen. 1984. A fast parallel algorithm for thinning digital patterns. Commun. ACM 27, 3 (March 1984), 236-239.

and YES, it's old and not perfect, but I repeat, it works simple and my implementation is very similar to the one used in the article, so you can understand what the code does.

if you need an explanation on how the code works, you can look at this page:
http://jeff.cs.mcgill.ca/~godfried/teaching/projects97/azar/skeleton.html
even if the algorithm is a bit different, the idea behind the conditions to switch a pixel from black to white are the same.

this is the c++ code, I use OpenCV 2.4.something. this is my 90-lines-of-code NOT OPTIMIZED code. Just the first sub-iteration is commented, but the code is very simple.

//****************************************************************************//
// Copyright (C) 2012, Michele (MikyMouse) Martinelli                      //
//****************************************************************************//
// This file is free software; you can redistribute it and/or modify it    //
// under the terms of the GNU Lesser General Public License as published by   //
// the Free Software Foundation; either version 2.1 of the License, or (at    //
// your option) any later version.                                            //
//****************************************************************************//
int skeletonize(Mat img){
    int iter;
    for (iter=0;first_subiteration(img)+second_subiteration(img)>0;iter++);
   
    //this return the number of iteration made
    return iter;
}

int first_subiteration(Mat imgMat){
    int changed=0;
   
    //for each non-zero pixel p1
    for(int i=2; i<imgMat.rows-2; i++)
            for(int j=2; j<imgMat.cols-2; j++){
                unsigned char p1 = imgMat.at<uchar>(i,j)>0?1:0;
                if(p1==0)
                    continue;
                //compute pixel surrounding p1
                unsigned char p2 = imgMat.at<uchar>(i-1,j)>0?1:0;
                unsigned char p3 = imgMat.at<uchar>(i-1,j+1)>0?1:0;
                unsigned char p4 = imgMat.at<uchar>(i,j+1)>0?1:0;
                unsigned char p5 = imgMat.at<uchar>(i+1,j+1)>0?1:0;
                unsigned char p6 = imgMat.at<uchar>(i+1,j)>0?1:0;
                unsigned char p7 = imgMat.at<uchar>(i+1,j-1)>0?1:0;
                unsigned char p8 = imgMat.at<uchar>(i,j-1)>0?1:0;
                unsigned char p9 = imgMat.at<uchar>(i-1,j-1)>0?1:0;
               
                //compute B(p1)
                int Bp1=p2+p3+p4+p5+p6+p7+p8+p9;
               
                //compute A(p1)
                unsigned char points[10]={p2,p3,p4,p5,p6,p7,p8,p9,p2};
                unsigned char Ap1=0;
                for(int ind=0;ind<9;ind++)
                    if(points[ind]==0 && points[ind+1]==1)
                        Ap1++;
               
                //conditions as reported in the article
                if ( Bp1>=2 && Bp1<=6 )
                    if (Ap1==1)
                        if (p2*p4*p6==0)
                            if(p4*p6*p8==0){
                                imgMat.at<uchar>(i,j)=0;
                                changed++;   
                            }
            }
        return changed;
}

int second_subiteration(Mat imgMat){
    int changed=0;
   
    for(int i=2; i<imgMat.rows-2; i++)
            for(int j=2; j<imgMat.cols-2; j++){
                unsigned char p1 = imgMat.at<uchar>(i,j)>0?1:0;
                if(p1==0)
                    continue;
                unsigned char p2 = imgMat.at<uchar>(i-1,j)>0?1:0;
                unsigned char p3 = imgMat.at<uchar>(i-1,j+1)>0?1:0;
                unsigned char p4 = imgMat.at<uchar>(i,j+1)>0?1:0;
                unsigned char p5 = imgMat.at<uchar>(i+1,j+1)>0?1:0;
                unsigned char p6 = imgMat.at<uchar>(i+1,j)>0?1:0;
                unsigned char p7 = imgMat.at<uchar>(i+1,j-1)>0?1:0;
                unsigned char p8 = imgMat.at<uchar>(i,j-1)>0?1:0;
                unsigned char p9 = imgMat.at<uchar>(i-1,j-1)>0?1:0;
               
                int Bp1=p2+p3+p4+p5+p6+p7+p8+p9;
                unsigned char points[10]={p2,p3,p4,p5,p6,p7,p8,p9,p2};
                unsigned char Ap1=0;
                for(int ind=0;ind<9;ind++)
                    if(points[ind]==0 && points[ind+1]==1)
                        Ap1++;
               
                if ( Bp1>=2 && Bp1<=6 )
                   
                    if (Ap1==1)if (p2*p4*p8==0)
                            if(p2*p6*p8==0){
                                imgMat.at<uchar>(i,j)=0;
                                changed++;
                                }
            }
   
    return changed;       
}

Nessun commento:

Posta un commento