To transform a closed digital curve, remove one pixel, remember it and submit the opened curve to the following python code.
For example, the minimal closed curve:
For example, the minimal closed curve:
gives:
closed curve
[[0 0 0 0 0]
[0 0 1 0 0]
[0 1 0 1 0]
[0 1 0 1 0]
[0 0 1 0 0]
[0 0 0 0 0]]
index list [array([1, 2]), array([2, 1]), array([3, 1]), array([4, 2]), array([3, 3]), array([2, 3]), array([1, 2])] The fourth pixel: [4 2]
Using functions, the code is a little bit cleaner, it relies only on numpy and scipy.ndimage for mathematical morphology (hit or miss), since only endpoints are needed here.If branched points are required, I will continue to use mahotas implementation of hit or miss operator because it handles the "don't care" pixels written with a "2" in the structuring element. Thus some of my scripts may not be run on windows, because mahotas has some installation problems due to freeimage on that plateform. I don't know if the scipy implementation of the hit or miss operator handle the "don't care" pixels.
# -*- coding: utf-8 -*-
"""
Created on Mon Jun 20 09:08:38 2011
@author: Jean-Patrick
"""
import time
import numpy as np
from scipy import ndimage as nd
def opencurveToPixelsList(im):
def scipyEndPoints(imcurve):
se1=np.array([[0,0,0],[0,1,0],[0,1,0]])
se2=np.array([[0,0,0],[0,1,0],[0,0,1]])
EndPoints=np.zeros(imcurve.shape,dtype=np.uint)
#perform hit & miss
for i in (0,3):
hit1=nd.morphology.binary_hit_or_miss(imcurve,se1)
se1=np.rot90(se1)
EndPoints=EndPoints+hit1
hit2=nd.morphology.binary_hit_or_miss(imcurve,se2)
se2=np.rot90(se1)
EndPoints=EndPoints+hit2
return EndPoints
#total number of pixels in the curve
pixel=np.array([-1,-1])
pixN=np.sum(im==1)
pixelsList=[]
#get end point by H&Miss op by ndimage
#print scipyEndPoints(im)
#get end points by hit&miss operator provided by mahotas
#ep=bep.get_endpoints(im)
ep=scipyEndPoints(im)
lab,n=nd.label(ep)
#where the first endpoint is
first_indices=np.where(lab==1)
pixel[0]=first_indices[0][0]
pixel[1]=first_indices[1][0]
pixelsList.append(np.copy(pixel))
#print "first point",first_indices," vector",walkingPixel
firstEndPoint=np.uint8(lab==1)
#keep an image of the second end point
lastEndPoint=np.uint8(lab==2)
last_Indices=np.where(lastEndPoint==1)
#print "last point",last_Indices
######################################################
## walk on curve with a 3x3 neighborhood
######################################################
#Init
current_curve=np.copy(im)
current_point=np.copy(firstEndPoint)
###start to walk on curve
#remove the second last point
#current_curve=current_curve-lastEndPoint##too heavy?, try indices
current_curve[last_Indices[0][0],last_Indices[1][0]]=0
c_point_ind=np.where(current_point==1)
#print c_point_ind
li=c_point_ind[0][0]
col=c_point_ind[1][0]
for i in range (0,pixN-2):
#3x3 neighborhood arround the endpoint
neighbor=current_curve[li-1:li+2,col-1:col+2]
neighbor[1,1]=0
#################
#Can only handle a curve
#such np.where(neihgbor==1) must gives only one pixel
#############
nextPointIndices=np.where(neighbor==1)##vectO'M=nextPointIndices
##remove the first point from the curve
current_curve[li,col]=0
pixN=pixN-1
#print current_curve
##compute nextPoint indices in the original base vectOM=OO'+O'M
li=(li-1)+nextPointIndices[0][0]
col=(col-1)+nextPointIndices[1][0]
pixel[0]=li
pixel[1]=col
pixelsList.append(np.copy(pixel))
#don't forget the last pixel
pixel[0]=last_Indices[0][0]
pixel[1]=last_Indices[1][0]
pixelsList.append(np.copy(pixel))
return pixelsList
############################
def closedcurveToPixelsList(im):
'''
given a closed curve (no test),
not touching the image border(no test) as:
[0,0,0,0,0],
[0,0,1,0,0],
[0,1,0,1,0],
[0,1,0,1,0],
[0,0,1,0,0],
[0,0,0,0,0]
The function returns an array list of ordered
points along the curve:
[1,2],[2,3],[3,3],[4,2],[3,1],[2,1],[1,2]
'''
firstpoint=np.copy(np.where(im==1))
print np.where(im==1)[1][0]
print "first point",firstpoint
openedcurve=np.copy(im)
openedcurve[firstpoint[0][0],firstpoint[1][0]]=0
print "now it's open"
print openedcurve
openlist=opencurveToPixelsList(openedcurve)
#add the first point of the closed curved at the
#begiining and the end of the list:
pixel=np.array([-1,-1])
pixel[0]=firstpoint[0][0]
pixel[1]=firstpoint[1][0]
openlist.append(np.copy(pixel))
openlist.insert(0,np.copy(pixel))
return openlist
#tests
im=np.array([[0,0,0,0,0,0,0],
[0,1,0,0,1,0,0],
[0,0,1,0,0,1,0],
[0,0,0,1,1,0,0],
[0,0,0,0,0,0,0]])
liste=opencurveToPixelsList(im)
closedcurve=np.array([
[0,0,0,0,0],
[0,0,1,0,0],
[0,1,0,1,0],
[0,1,0,1,0],
[0,0,1,0,0],
[0,0,0,0,0]])
liste2=closedcurveToPixelsList(closedcurve)
print liste
print liste[0]
print "closed curve"
print closedcurve
print liste2
print liste2[3]
