区域生长 python版 dicom图像

#!\Anaconda3\envs\py35 python

# -*- coding: utf-8 -*-
#!@Time   : 2018/6/1 15:39
#!@Author : python
#!@File   : .py



from operator import eq
import numpy as np
import pydicom as dicom
from skimage.measure import label
from skimage.measure import regionprops
from skimage.morphology import thin

import matplotlib.pyplot as plt
import cv2


def regiongrowth（orig_image seed_matrix threshold）:
    “““区域生长算法

    # Arguments：
        orig_image：待分割图像
        seed_matrix：种子矩阵
        threshold：生长条件，即待分割图中的像素值与目标像素值之间的可承受的差值
    # Returns
        seg：最终分割后的mask
        num_region：分割区域的个数

    # Example
        regiongrow（orig_image seed_matrix 50）

    “““
    orig_image = np.array（orig_image dtype=np.float32）                         # 将图像数值格式转为float

    markers = thin（seed_matrix）                                                 # 将种子所在区域缩小为一个点，后面限定区域时使用
    coordinate_nonzero = markers.nonzero（）                                      # 当种子为矩阵时，获取种子点的坐标值
    coordinate_nonzero = np.transpose（coordinate_nonzero）                       # 获取种子点的坐标值的转置
    seed_value = []
    for coordinate_i in coordinate_nonzero:
        seed_value.append（orig_image[coordinate_i[0] coordinate_i[1]]）         # 获取种子点的灰度值

    size_x size_y = orig_image.shape
    seg_image = np.array（np.zeros（（size_x size_y）） dtype=bool）                # 初始化一个bool格式的矩阵储存分割结果

    for i each_seed_value in enumerate（seed_value）:
        all_satisfied = abs（orig_image - each_seed_value） <= threshold          # 所有满足单个种子点阈值条件的像素点
        seg_image = seg_image | all_satisfied                                   # 所有满足多有种子点阈值条件的像素点

    label_seg_image = label（seg_image）                                          # 将所有连通区域打label
    area_region = regionprops（label_seg_image）                                  # 统计被标记的区域的面积分布，返回值为显示区域总数

    seg = np.array（np.zeros（（size_x size_y）） dtype=np.float32）                # 初始化一个float格式的矩阵储存分割结果
    for masker in coordinate_nonzero:                                           # 提取maskers所在的连通区域
        for each_area_region in area_region:
            for coord in each_area_region.coords:
                if sum（eq（coord masker）） == 2:                                 # 若连通区域中包含masker，则保留此连通区域
                    for idx in each_area_region.coords: