从分割图像中去除白色边框

17

我正在尝试使用下面的代码通过Kmeans分割肺部CT图像:

def process_mask(mask):
    convex_mask = np.copy(mask)
    for i_layer in range(convex_mask.shape[0]):
        mask1  = np.ascontiguousarray(mask[i_layer])
        if np.sum(mask1)>0:
            mask2 = convex_hull_image(mask1)
            if np.sum(mask2)>2*np.sum(mask1):
                mask2 = mask1
        else:
            mask2 = mask1
        convex_mask[i_layer] = mask2
    struct = generate_binary_structure(3,1)
    dilatedMask = binary_dilation(convex_mask,structure=struct,iterations=10)

    return dilatedMask

def lumTrans(img):
    lungwin = np.array([-1200.,600.])
    newimg = (img-lungwin[0])/(lungwin[1]-lungwin[0])
    newimg[newimg<0]=0
    newimg[newimg>1]=1
    newimg = (newimg*255).astype('uint8')
    return newimg


def lungSeg(imgs_to_process,output,name):

    if os.path.exists(output+'/'+name+'_clean.npy') : return
    imgs_to_process = Image.open(imgs_to_process)
    
    img_to_save = imgs_to_process.copy()
    img_to_save = np.asarray(img_to_save).astype('uint8')

    imgs_to_process = lumTrans(imgs_to_process)    
    imgs_to_process = np.expand_dims(imgs_to_process, axis=0)
    x,y,z = imgs_to_process.shape 
  
    img_array = imgs_to_process.copy()  
    A1 = int(y/(512./100))
    A2 = int(y/(512./400))

    A3 = int(y/(512./475))
    A4 = int(y/(512./40))
    A5 = int(y/(512./470))
    for i in range(len(imgs_to_process)):
        img = imgs_to_process[i]
        print(img.shape)
        x,y = img.shape
        #Standardize the pixel values
        allmean = np.mean(img)
        allstd = np.std(img)
        img = img-allmean
        img = img/allstd
        # Find the average pixel value near the lungs
        # to renormalize washed out images
        middle = img[A1:A2,A1:A2] 
        mean = np.mean(middle)  
        max = np.max(img)
        min = np.min(img)
        
        kmeans = KMeans(n_clusters=2).fit(np.reshape(middle,[np.prod(middle.shape),1]))
        centers = sorted(kmeans.cluster_centers_.flatten())
        threshold = np.mean(centers)
        thresh_img = np.where(img<threshold,1.0,0.0)  # threshold the image
       
        eroded = morphology.erosion(thresh_img,np.ones([4,4]))
        dilation = morphology.dilation(eroded,np.ones([10,10]))
        
        labels = measure.label(dilation)
        label_vals = np.unique(labels)
        regions = measure.regionprops(labels)
        good_labels = []
        for prop in regions:
            B = prop.bbox
            if B[2]-B[0]<A3 and B[3]-B[1]<A3 and B[0]>A4 and B[2]<A5:
                good_labels.append(prop.label)
        mask = np.ndarray([x,y],dtype=np.int8)
        mask[:] = 0
       
        for N in good_labels:
            mask = mask + np.where(labels==N,1,0)
        mask = morphology.dilation(mask,np.ones([10,10])) # one last dilation
        imgs_to_process[i] = mask

    m1 = imgs_to_process
    
    convex_mask = m1
    dm1 = process_mask(m1)
    dilatedMask = dm1
    Mask = m1
    extramask = dilatedMask ^ Mask
    bone_thresh = 180
    pad_value = 0

    img_array[np.isnan(img_array)]=-2000
    sliceim = img_array
    sliceim = sliceim*dilatedMask+pad_value*(1-dilatedMask).astype('uint8')
    bones = sliceim*extramask>bone_thresh
    sliceim[bones] = pad_value


    x,y,z = sliceim.shape
    if not os.path.exists(output): 
        os.makedirs(output)
    
    img_to_save[sliceim.squeeze()==0] = 0
    
    im = Image.fromarray(img_to_save)

    im.save(output + name + '.png', 'PNG')

问题是分割后的肺部仍然包含像这样的白色边框:

分割后的肺部(输出):

segmented lung

未分割的肺部(输入):

unsegmented lung

完整代码可以在Google Colab笔记本中找到。code

数据集的示例在这里


第一个想法是尝试更改那些膨胀核。10像素似乎太多了。 - Matt Hall
@kwinkunks 我已将其更改为[5,5],但某些白色边框仍然存在。我不确定如何确认肺部的关键区域是否被剔除。 - Rawan
2
你有没有看过使用Canny之类的边缘检测算法时数据的变化?我认为你想要做的是找到灰色区域的边缘并保留内部,然后进行处理。 - Jay Obermark
2个回答

11
对于这个问题,我不建议使用Kmeans颜色量化技术,因为这种技术通常用于有各种颜色且希望将其分割成主要颜色块的情况。查看先前的答案以获取典型用例。由于你的CT扫描图像是灰度的,Kmeans的表现不会很好。以下是一种可能的解决方案,使用OpenCV进行简单的图像处理:
  1. 获取二进制图像。 加载输入图像,转换为灰度图像使用Otsu阈值法,并查找轮廓

  2. 创建一个空白掩模以提取所需对象。 我们可以使用np.zeros()创建与输入图像相同大小的空掩模。

  3. 使用轮廓面积和长宽比筛选轮廓。 通过确保轮廓在指定的面积阈值以及长宽比内搜索肺部对象。我们使用cv2.contourArea()cv2.arcLength()cv2.approxPolyDP()进行轮廓周长和轮廓形状近似。如果我们已经找到了肺部对象,我们利用cv2.drawContours()将我们的掩模填充为白色,以表示我们要提取的对象。

  4. 将掩模与原始图像进行逻辑与操作。 最后,我们将掩模转换为灰度图像,并使用cv2.bitwise_and()进行逻辑与操作,以获得我们的结果。


这是我们的图像处理流程,按步骤可视化:

灰度化 -> 奥茨阈值

检测到需要提取的对象会用绿色标出 -> 填充掩模

使用按位与运算得到我们的结果 -> 可选的带有白色背景的结果

代码

import cv2
import numpy as np

image = cv2.imread('1.png')
highlight = image.copy()
original = image.copy()

# Convert image to grayscale, Otsu's threshold, and find contours
gray = cv2.cvtColor(image, cv2.COLOR_BGR2GRAY)
thresh = cv2.threshold(gray, 0, 255, cv2.THRESH_BINARY_INV | cv2.THRESH_OTSU)[1]
contours = cv2.findContours(thresh.copy(), cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
contours = contours[0] if len(contours) == 2 else contours[1]

# Create black mask to extract desired objects
mask = np.zeros(image.shape, dtype=np.uint8)

# Search for objects by filtering using contour area and aspect ratio
for c in contours:
    # Contour area
    area = cv2.contourArea(c)
    # Contour perimeter
    peri = cv2.arcLength(c, True)
    # Contour approximation
    approx = cv2.approxPolyDP(c, 0.035 * peri, True)
    (x, y, w, h) = cv2.boundingRect(approx)
    aspect_ratio = w / float(h)
    # Draw filled contour onto mask if passes filter
    # These are arbitary values, may need to change depending on input image
    if aspect_ratio <= 1.2 or area < 5000:
        cv2.drawContours(highlight, [c], 0, (0,255,0), -1)
        cv2.drawContours(mask, [c], 0, (255,255,255), -1)

# Convert 3-channel mask to grayscale then bitwise-and with original image for result
mask = cv2.cvtColor(mask, cv2.COLOR_BGR2GRAY)
result = cv2.bitwise_and(original, original, mask=mask)

# Uncomment if you want background to be white instead of black
# result[mask==0] = (255,255,255)

# Display
cv2.imshow('gray', gray)
cv2.imshow('thresh', thresh)
cv2.imshow('highlight', highlight)
cv2.imshow('mask', mask)
cv2.imshow('result', result)

# Save images
# cv2.imwrite('gray.png', gray)
# cv2.imwrite('thresh.png', thresh)
# cv2.imwrite('highlight.png', highlight)
# cv2.imwrite('mask.png', mask)
# cv2.imwrite('result.png', result)
cv2.waitKey(0)

2
解决这个问题的更简单方法是使用形态学腐蚀。只需调整阈值即可。

网页内容由stack overflow 提供, 点击上面的
可以查看英文原文,
原文链接