Python morphology.medial_axis函数代码示例

    def convert_binary_image_to_sdf(self, binary_img, vis = False):
        binary_data = np.array(binary_img)
        skel, sdf_in = morph.medial_axis(binary_data, return_distance = True)
        useless_skel, sdf_out = morph.medial_axis(self.upper_bound_ - binary_data, return_distance = True)
        
        sdf = sdf_out - sdf_in

        # display the sdf and skeleton
        if vis:
            dist_on_skel = sdf * skel
            
            fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(8, 4))
            ax1.imshow(binary_data, cmap=plt.cm.gray, interpolation='nearest')
            ax1.axis('off')
            ax2.imshow(dist_on_skel, cmap=plt.cm.spectral, interpolation='nearest')
            ax2.contour(binary_data, [0.5], colors='w')
            ax2.axis('off')

            fig.subplots_adjust(hspace=0.01, wspace=0.01, top=1, bottom=0, left=0, right=1)
            plt.show()

            plt.imshow(sdf)
            plt.show()

            plt.imshow(skel)
            plt.show()

        return sdf, skel

    def run(self):

        print 'Gerando Caracteristicas....'
        images = self.__image_list(image_path('circinatum'), image_path('kelloggii'), image_path('negundo'))

        arquivo = open(data_path(),'w')
        count = 0
        for i in images:
            count+=1
            img, vetor, arclen, classe = self.__features(i)

            # Media Curvatura
            arquivo.write(str(np.mean(np.abs(vetor))) + ',')
            # Comprimento de Arco
            arquivo.write(str(arclen) + ',')            
            # Area
            arquivo.write(str(np.sum(img)) + ',')          
            # Numero Pixels Esqueleto
            arquivo.write(str(np.sum(morphology.medial_axis(img))) + ',')
            # Classe Folhas
            arquivo.write(classe)
            arquivo.write("\n")

        arquivo.close()

        print '100%'
        print 'Total Imagens: ' + str(count)

def label_nuclei(binary, min_size):
    '''Label, watershed and remove small objects'''

    distance = medial_axis(binary, return_distance=True)[1]

    distance_blured = gaussian_filter(distance, 5)

    local_maxi = peak_local_max(distance_blured, indices=False, labels=binary, min_distance = 30)

    markers = measure_label(local_maxi)

#    markers[~binary] = -1

#    labels_rw = segmentation.random_walker(binary, markers)

#    labels_rw[labels_rw == -1] = 0

#    labels_rw = segmentation.relabel_sequential(labels_rw)

    labels_ws = watershed(-distance, markers, mask=binary)

    labels_large = remove_small_objects(labels_ws,min_size)

    labels_clean_border = clear_border(labels_large)

    labels_from_one = relabel_sequential(labels_clean_border)

#    plt.imshow(ndimage.morphology.binary_dilation(markers))
#    plt.show()

    return labels_from_one[0]

def skeletonize_mitochondria(mch_channel):
    mch_collector = np.max(mch_channel, axis=0)  # TODO: check max projection v.s. sum
    skeleton_labels = np.zeros(mch_collector.shape, dtype=np.uint8)

    # thresh = np.max(mch_collector)/2.
    thresh = threshold_otsu(mch_collector)
    # use adaptative threshold? => otsu seems to be sufficient in this case

    skeleton_labels[mch_collector > thresh] = 1
    skeleton2 = skeletonize(skeleton_labels)
    skeleton, distance = medial_axis(skeleton_labels, return_distance=True)
    active_threshold = np.mean(mch_collector[skeleton_labels]) * 5

    # print active_threshold
    transform_filter = np.zeros(mch_collector.shape, dtype=np.uint8)
    transform_filter[np.logical_and(skeleton > 0, mch_collector > active_threshold)] = 1
    skeleton = transform_filter * distance

    skeleton_ma = np.ma.masked_array(skeleton, skeleton > 0)
    skeleton_convolve = ndi.convolve(skeleton_ma, np.ones((3, 3)), mode='constant', cval=0.0)
    divider_convolve = ndi.convolve(transform_filter, np.ones((3, 3)), mode='constant', cval=0.0)
    skeleton_convolve[divider_convolve > 0] = skeleton_convolve[divider_convolve > 0] / \
                                              divider_convolve[divider_convolve > 0]
    new_skeleton = np.zeros_like(skeleton)
    new_skeleton[skeleton2] = skeleton_convolve[skeleton2]
    skeleton = new_skeleton

    return skeleton_labels, mch_collector, skeleton, transform_filter

def execute_Skeleton(proxy,obj):

	from skimage.morphology import medial_axis

	threshold=0.1*obj.threshold

	try: 
		img2=obj.sourceObject.Proxy.img
		img=img2.copy()
	except: 
		sayexc()
		img=cv2.imread(__dir__+'/icons/freek.png')

	data = cv2.cvtColor(img, cv2.COLOR_BGR2GRAY)

	# Compute the medial axis (skeleton) and the distance transform
	skel, distance = medial_axis(data, return_distance=True)

	# Distance to the background for pixels of the skeleton
	dist_on_skel = distance * skel

	# entferne ganz duenne linien
	dist_on_skelw =(dist_on_skel >= threshold)* distance

	say("size of the image ...")
	say(dist_on_skelw.shape)
#	skel = np.array(dist_on_skelw,np.uint8) 
	skel = np.array(dist_on_skelw *255/np.max(dist_on_skelw),np.uint8) 
	obj.Proxy.img=cv2.cvtColor(skel*100, cv2.COLOR_GRAY2BGR)
	obj.Proxy.dist_on_skel=dist_on_skelw

def segment_cells(frame, mask=None):
    """
    Compute the initial segmentation based on ridge detection + watershed.
    This works reasonably well, but is not robust enough to use by itself.
    """
    
    blurred = filters.gaussian_filter(frame, 2)
    ridges = enhance_ridges(frame)
    
    # threshold ridge image
    thresh = filters.threshold_otsu(ridges)
    thresh_factor = 0.6
    prominent_ridges = ridges > thresh_factor*thresh
    prominent_ridges = morphology.remove_small_objects(prominent_ridges, min_size=256)
    prominent_ridges = morphology.binary_closing(prominent_ridges)
    prominent_ridges = morphology.binary_dilation(prominent_ridges)
    
    # skeletonize
    ridge_skeleton = morphology.medial_axis(prominent_ridges)
    ridge_skeleton = morphology.binary_dilation(ridge_skeleton)
    ridge_skeleton *= mask
    ridge_skeleton -= mask
    
    # label
    cell_label_im = measure.label(ridge_skeleton)
    
    # morphological closing to fill in the cracks
    for cell_num in range(1, cell_label_im.max()+1):
        cell_mask = cell_label_im==cell_num
        cell_mask = morphology.binary_closing(cell_mask, disk(3))
        cell_label_im[cell_mask] = cell_num
    
    return cell_label_im 

def skeletonize_mitochondria(mCh_channel):

    mch_collector = np.max(mCh_channel, axis=0)  # TODO: check how max affects v.s. sum
    labels = np.zeros(mch_collector.shape, dtype=np.uint8)

    # thresh = np.max(mch_collector)/2.
    thresh = threshold_otsu(mch_collector)
    # TODO: use adaptative threshold? => otsu seems to be sufficient in this case
    # http://scikit-image.org/docs/dev/auto_examples/xx_applications/plot_thresholding.html#sphx
    # -glr-auto-examples-xx-applications-plot-thresholding-py
    #  log-transform? => Nope, does not work
    # TODO: hessian/laplacian of gaussian blob detection?

    labels[mch_collector > thresh] = 1
    skeleton2 = skeletonize(labels)
    skeleton, distance = medial_axis(labels, return_distance=True)
    active_threshold = np.mean(mch_collector[labels]) * 5

    # print active_threshold
    transform_filter = np.zeros(mch_collector.shape, dtype=np.uint8)
    transform_filter[np.logical_and(skeleton > 0, mch_collector > active_threshold)] = 1
    skeleton = transform_filter * distance

    skeleton_ma = np.ma.masked_array(skeleton, skeleton > 0)
    skeleton_convolve = ndi.convolve(skeleton_ma, np.ones((3, 3)), mode='constant', cval=0.0)
    divider_convolve = ndi.convolve(transform_filter, np.ones((3, 3)), mode='constant', cval=0.0)
    skeleton_convolve[divider_convolve > 0] = skeleton_convolve[divider_convolve > 0] \
                                              / divider_convolve[divider_convolve > 0]
    new_skeleton = np.zeros_like(skeleton)
    new_skeleton[skeleton2] = skeleton_convolve[skeleton2]
    skeleton = new_skeleton

    return labels, mch_collector, skeleton, transform_filter

def _process_img_morph(img, threshold=.5, scale=1):
    if scale > 1:
        up_img = transform.pyramid_expand(img, upscale=scale, order=3)  # type: np.ndarray
        img = (255. * up_img).astype(img.dtype)
    img_min, img_max = img.min(), img.max()
    bin_img = (img >= img_min + (img_max - img_min) * threshold)
    skel, dist_map = morphology.medial_axis(bin_img, return_distance=True)
    return img, bin_img, skel, dist_map

def skeleton(seg):
    skel, dist = skmorph.medial_axis(seg, return_distance=True)
    node, edge, leaf = (spim.label(g, np.ones((3, 3), bool))[0] for g in skel2graph(skel))

    trim_edge = (edge != 0) & ~(skmorph.binary_dilation(node != 0, np.ones((3, 3), bool)) != 0)
    trim_edge = spim.label(trim_edge, np.ones((3, 3), bool))[0]

    leaf_edge_vals = skmorph.binary_dilation(leaf != 0, np.ones((3, 3), bool)) != 0
    leaf_edge_vals = np.unique(trim_edge[leaf_edge_vals])
    leaf_edge_vals = leaf_edge_vals[leaf_edge_vals > 0]
    leaf_edge = leaf != 0

    trim_edge = ndshm.fromndarray(trim_edge)
    leaf_edge = ndshm.fromndarray(leaf_edge)
    Parallel()(delayed(set_msk)(leaf_edge, trim_edge, l) for l in leaf_edge_vals)
    trim_edge = np.copy(trim_edge)
    leaf_edge = np.copy(leaf_edge)

    leaf_edge[(skmorph.binary_dilation(leaf_edge, np.ones((3, 3), bool)) != 0) & (edge != 0)] = True
    leaf_edge = spim.label(leaf_edge, np.ones((3, 3), bool))[0]

    leaf_edge_node = skmorph.binary_dilation(leaf_edge != 0, np.ones((3, 3), bool)) != 0
    leaf_edge_node = ((node != 0) & leaf_edge_node) | leaf_edge
    leaf_edge_node = spim.label(leaf_edge_node, np.ones((3, 3), bool))[0]

    cand_node = leaf_edge_node * (node != 0)
    cand_node = cand_node.nonzero()
    cand_node = np.transpose((leaf_edge_node[cand_node],) + cand_node + (2 * dist[cand_node],))

    cand_leaf = leaf_edge_node * (leaf != 0)
    cand_leaf = cand_leaf.nonzero()
    cand_leaf = np.transpose((leaf_edge_node[cand_leaf],) + cand_leaf)

    if len(cand_node) > 0 and len(cand_leaf) > 0:
        cand_leaf = ndshm.fromndarray(cand_leaf)
        cand_node = ndshm.fromndarray(cand_node)
        pruned = Parallel()(delayed(prune_leaves)(cand_leaf, cand_node, j) for j in np.unique(cand_node[:, 0]))
        cand_leaf = np.copy(cand_leaf)
        cand_node = np.copy(cand_node)

        pruned_ind = []
        for p in pruned:
            pruned_ind.extend(p)
        pruned_ind = tuple(np.transpose(pruned_ind))

        pruned = ~skel

        pruned = ndshm.fromndarray(pruned)
        leaf_edge = ndshm.fromndarray(leaf_edge)
        Parallel()(delayed(set_msk)(pruned, leaf_edge, l) for l in np.unique(leaf_edge[pruned_ind]))
        pruned = np.copy(pruned)
        leaf_edge = np.copy(leaf_edge)

        pruned = ~pruned
    else:
        pruned = skel

    return pruned

def get_text_image(text, font, point_size):
    text = get_text(text, font, point_size)
    io.imsave(IMAGE_FILENAME, text)
    thresh = skfilter.threshold_otsu(text)
    binary = text > thresh
    skel, distance = morphology.medial_axis(binary, return_distance=True)
    distance = distance.astype(np.uint16)
    skel = skel.astype(np.uint16)
    return skel*distance

def make_skeleton(image, mindist):
    """Return a skeletonization of the image, filtered and normalized."""
    skel, distance = morphology.medial_axis(image, return_distance=True)

    dist_on_skel = distance * skel

    dist_on_skel = filter_skeleton_distances(dist_on_skel)
    normalized = normalize_skeleton(dist_on_skel)

    return normalized

 def test_01_01_rectangle(self):
     '''Test skeletonize on a rectangle'''
     image = np.zeros((9, 15), bool)
     image[1:-1, 1:-1] = True
     #
     # The result should be four diagonals from the
     # corners, meeting in a horizontal line
     #
     expected = np.array([[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
                          [0,1,0,0,0,0,0,0,0,0,0,0,0,1,0],
                          [0,0,1,0,0,0,0,0,0,0,0,0,1,0,0],
                          [0,0,0,1,0,0,0,0,0,0,0,1,0,0,0],
                          [0,0,0,0,1,1,1,1,1,1,1,0,0,0,0],
                          [0,0,0,1,0,0,0,0,0,0,0,1,0,0,0],
                          [0,0,1,0,0,0,0,0,0,0,0,0,1,0,0],
                          [0,1,0,0,0,0,0,0,0,0,0,0,0,1,0],
                          [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]], bool)
     result = medial_axis(image)
     assert np.all(result == expected)
     result, distance = medial_axis(image, return_distance=True)
     assert distance.max() == 4

def loadImage(filename):
    global im
    global distanceImage
    global lapImage
    global skeleton
    global resultImage
    global im1
    global skimage

    im = np.asarray(Image.open(filename))
    distanceImage, lapImage, resultImage = medial(im)
    skeleton = lapImage < 0
    im1 = ax.imshow(im, cmap = "gray", interpolation="nearest")
    skimage = medial_axis(im, return_distance=False)
    lapImage = lapImage - np.max(lapImage)

 def test_01_02_hole(self):
     '''Test skeletonize on a rectangle with a hole in the middle'''
     image = np.zeros((9, 15), bool)
     image[1:-1, 1:-1] = True
     image[4, 4:-4] = False
     expected = np.array([[0,0,0,0,0,0,0,0,0,0,0,0,0,0,0],
                          [0,1,0,0,0,0,0,0,0,0,0,0,0,1,0],
                          [0,0,1,1,1,1,1,1,1,1,1,1,1,0,0],
                          [0,0,1,0,0,0,0,0,0,0,0,0,1,0,0],
                          [0,0,1,0,0,0,0,0,0,0,0,0,1,0,0],
                          [0,0,1,0,0,0,0,0,0,0,0,0,1,0,0],
                          [0,0,1,1,1,1,1,1,1,1,1,1,1,0,0],
                          [0,1,0,0,0,0,0,0,0,0,0,0,0,1,0],
                          [0,0,0,0,0,0,0,0,0,0,0,0,0,0,0]],bool)
     result = medial_axis(image)
     assert np.all(result == expected)

    def __transform(self):
        img_gray = io.imread(self.image_path, True)
        #Aplicar otsu para binarizar a imagem
        img_otsu = filter.threshold_otsu(img_gray) 
        self.__img = img_gray < img_otsu

        # Procura contornos da imagem binarizada
        self.__contours = measure.find_contours(self.__img, 0.5)
        
        arclen=0.0
        for n, contour in enumerate(self.__contours):
            arclenTemp=0.0
            for indice, valor in enumerate(contour):
                if indice > 0:
                    d1 = math.fabs(round(valor[0]) - round(contour[indice-1,0]))
                    d2 = math.fabs(round(valor[1]) - round(contour[indice-1,1]))
                    if d1+d2>1.0:
                        arclenTemp+=math.sqrt(2)
                    elif d1+d2 == 1:
                        arclenTemp+=1
            if arclenTemp > arclen:
                arclen = arclenTemp
                bestn = n
        
        #Transforma a lista contours[bestn] em uma matriz[n,2]
        aux = np.asarray(self.__contours[bestn])

        #---------------------------  Curvatura --------------
        vetor = [] #vetor que irá receber as curvaturas k(t)
        
        for i in range(len(aux)-2):    #Percorrer ate -2 para não pegar elementos inexistentes
            #---------------------------  Curvatura --------------
            #Inverter as posições em relação a fórmula pois o x esta no lugar do y
            b1 =  ( (aux[i-2,1]+aux[i+2,1]) + (2*(aux[i-1,1] + aux[i+1,1])) - (6*aux[i,1]) ) / 12
            b2 = ( (aux[i-2,0]+aux[i+2,0]) + (2*(aux[i-1,0] + aux[i+1,0])) - (6*aux[i,0]) ) / 12
            c1 =  ( (aux[i+2,1]-aux[i-2,1]) + (4*(aux[i+1,1] - aux[i-1,1])) ) / 12
            c2 = ( (aux[i+2,0]-aux[i-2,0]) + (4*(aux[i+1,0] - aux[i-1,0])) ) / 12

            k =  (2*(c1*b1 - c2*b2)) / ((c1**2 + c2**2)**(3/2))

            vetor.append(k) #append: insere objeto no final da lista

        self.__media_curvatura = np.mean(np.abs(vetor))
        self.__comprimento_arco = arclen
        self.__area = np.sum(self.__img)
        self.__esqueleto_pixel = np.sum(morphology.medial_axis(self.__img))
        self.__bestn = bestn