前往小程序,Get更优阅读体验!
立即前往
首页
学习
活动
专区
工具
TVP
发布
社区首页 >专栏 >医疗影像处理:去除医疗影像中背景的影响2D/3D【numpy-code】| CSDN博文精选

医疗影像处理:去除医疗影像中背景的影响2D/3D【numpy-code】| CSDN博文精选

作者头像
AI科技大本营
发布2019-11-19 22:20:26
8810
发布2019-11-19 22:20:26
举报
文章被收录于专栏:AI科技大本营的专栏

作者 | chestnut--

来源 | CSDN博客

在医疗影像中特别是CT影像,包含大量的背景,在进行器官分割时,首先去除背景对分割的效果有很好的提升。本博客使用Python处理医疗影像并去除背景的影像。

使用的样例数据来自Multimodal Brain Tumor Segmentation Challenge 2018(https://www.med.upenn.edu/sbia/brats2018/registration.html)

读取数据的基本信息

代码语言:javascript
复制
import nibabel as nib
import numpy as np
from nilearn.image.image import _crop_img_to as crop_img_to
file = 'flair.nii.gz'
image = nib.load(file)
print(image.shape)

(240, 240, 155)

代码语言:javascript
复制
print(np.max(image.get_data()),np.min(image.get_data()))

470, 0

注意:0为背景像素值

显示其中的一个slice

代码语言:javascript
复制
import matplotlib.pyplot as plt
data = image.get_data()[:,:,52]
plt.imshow(data)

得到裁剪的空间位置坐标

代码语言:javascript
复制
def get_slice_index(data, rtol=1e-8):

    infinity_norm = max(-data.min(), data.max())
    passes_threshold = np.logical_or(data < -rtol * infinity_norm,
                                     data > rtol * infinity_norm)  ##
    if data.ndim == 4:
        passes_threshold = np.any(passes_threshold, axis=-1)

    coords = np.array(np.where(passes_threshold))
    start = coords.min(axis=1)
    end = coords.max(axis=1) + 1

    # pad with one voxel to avoid resampling problems
    start = np.maximum(start - 1, 0)
    end = np.minimum(end + 1, data.shape[:3])

    slices = [slice(s, e) for s, e in zip(start, end)]
    return slices

使用True 和False标记背景区域,背景像素值为False。

代码语言:javascript
复制
def have_back(image):
    background_value=0
    tolerance=0.00001
    is_foreground = np.logical_or(image.get_data() < (background_value - tolerance),
                                  image.get_data()> (background_value + tolerance))
    foreground = np.zeros(is_foreground.shape, dtype=np.uint8)
    foreground[is_foreground] = 1
    return foreground

调用,得到背景图标记,通过背景图标记得到坐标位置

代码语言:javascript
复制
foreground = have_back(image)
crop = get_slice_index(foreground)
print(crop)

[slice(45, 192, None), slice(47, 210, None), slice(0, 137, None)]

分别代表X,Y,Z的裁剪坐标

裁剪

代码语言:javascript
复制
image_o = crop_img_to(image, crop, copy=True)

image_o shape (147, 163, 137)

显示裁剪之后的结果

代码语言:javascript
复制
import matplotlib.pyplot as plt
data_o = image_o.get_data()[:,:,52]
plt.imshow(data_o)

对比着看一下:

以上参考的是[1]中的代码。

[2]中相同的处理方式,同时能够处理2D和3D。

整体的代码如下,主要使用到了get_none_zero_region和crop_ND_volume_with_bounding_box这两个函数:

代码语言:javascript
复制
def get_none_zero_region(im, margin):
    """
    get the bounding box of the non-zero region of an ND volume
    """
    input_shape = im.shape
    if(type(margin) is int ):
        margin = [margin]*len(input_shape)
    assert(len(input_shape) == len(margin))
    indxes = np.nonzero(im)
    idx_min = []
    idx_max = []
    for i in range(len(input_shape)):
        idx_min.append(indxes[i].min())
        idx_max.append(indxes[i].max())

    for i in range(len(input_shape)):
        idx_min[i] = max(idx_min[i] - margin[i], 0)
        idx_max[i] = min(idx_max[i] + margin[i], input_shape[i] - 1)
    return idx_min, idx_max

def crop_ND_volume_with_bounding_box(volume, min_idx, max_idx):
    """
    crop/extract a subregion form an nd image.
    """
    dim = len(volume.shape)
    assert(dim >= 2 and dim <= 5)
    if(dim == 2):
        output = volume[np.ix_(range(min_idx[0], max_idx[0] + 1),
                               range(min_idx[1], max_idx[1] + 1))]
    elif(dim == 3):
        output = volume[np.ix_(range(min_idx[0], max_idx[0] + 1),
                               range(min_idx[1], max_idx[1] + 1),
                               range(min_idx[2], max_idx[2] + 1))]
    elif(dim == 4):
        output = volume[np.ix_(range(min_idx[0], max_idx[0] + 1),
                               range(min_idx[1], max_idx[1] + 1),
                               range(min_idx[2], max_idx[2] + 1),
                               range(min_idx[3], max_idx[3] + 1))]
    elif(dim == 5):
        output = volume[np.ix_(range(min_idx[0], max_idx[0] + 1),
                               range(min_idx[1], max_idx[1] + 1),
                               range(min_idx[2], max_idx[2] + 1),
                               range(min_idx[3], max_idx[3] + 1),
                               range(min_idx[4], max_idx[4] + 1))]
    else:
        raise ValueError("the dimension number shoud be 2 to 5")
    return output

使用:

代码语言:javascript
复制
margin = 5 
bbmin, bbmax = get_none_zero_region(arr, margin)
bbmin, bbmax

([34, 17], [171, 191])

arr为2D或者3D的矩阵。

代码语言:javascript
复制
volume = crop_ND_volume_with_bounding_box(arr, bbmin, bbmax)

参考

3DUnetCNN

https://github.com/ellisdg/3DUnetCNN

3DUnet-Tensorflow-Brats18

https://github.com/tkuanlun350/3DUnet-Tensorflow-Brats18/blob/master/utils.py

本文参与 腾讯云自媒体同步曝光计划,分享自微信公众号。
原始发表:2019-11-18,如有侵权请联系 cloudcommunity@tencent.com 删除

本文分享自 AI科技大本营 微信公众号,前往查看

如有侵权,请联系 cloudcommunity@tencent.com 删除。

本文参与 腾讯云自媒体同步曝光计划  ,欢迎热爱写作的你一起参与!

评论
登录后参与评论
0 条评论
热度
最新
推荐阅读
领券
问题归档专栏文章快讯文章归档关键词归档开发者手册归档开发者手册 Section 归档