使用pydicom实现Dicom文件读取与CT图像窗宽窗位调整

使用pydicom实现Dicom文件读取与CT图像窗宽窗位调整1.前言为了能够在Labelme上对Dicom图像进行编辑,这里对python环境下Dicom文件的读取进行了研究。在Dicom图像中CT的窗宽窗位是一个很重要的概念,但是找了半天在pydicom中没有相关设置函数,这里跟DCMTK还不一样。但是可以根据两个tag得到CT图像的CT值,那就是(0028|1052):rescaleintercept和(0028|1053):rescales…

大家好,又见面了,我是你们的朋友全栈君。

1. 前言

为了能够在Labelme上对Dicom图像进行编辑,这里对python环境下Dicom文件的读取进行了研究。在Dicom图像中CT的窗宽窗位是一个很重要的概念,但是找了半天在pydicom中没有相关设置函数,这里跟DCMTK还不一样。但是可以根据两个tag得到CT图像的CT值,那就是(0028|1052):rescale intercept和(0028|1053):rescale slope。则按照下面的算子得到CT图像,进而就可以调整窗宽窗位了

Hu = pixel * slope + intercept

至于那个部位的窗宽窗位是多少各位看官就可以自行百度了。

2. 代码实现

# -*- coding=utf-8 -*-
import matplotlib.pyplot as plt
import pydicom
import pydicom.uid
import sys
import PIL.Image as Image
from PyQt5 import QtGui
import os

step1:读取Dicom图像数据与得到CT值图像(CT图)

have_numpy = True

try:
    import numpy
except ImportError:
    have_numpy = False
    raise

sys_is_little_endian = (sys.byteorder == 'little')

NumpySupportedTransferSyntaxes = [
    pydicom.uid.ExplicitVRLittleEndian,
    pydicom.uid.ImplicitVRLittleEndian,
    pydicom.uid.DeflatedExplicitVRLittleEndian,
    pydicom.uid.ExplicitVRBigEndian,
]

# 支持的传输语法
def supports_transfer_syntax(dicom_dataset):
    """ Returns ------- bool True if this pixel data handler might support this transfer syntax. False to prevent any attempt to try to use this handler to decode the given transfer syntax """
    return (dicom_dataset.file_meta.TransferSyntaxUID in
            NumpySupportedTransferSyntaxes)


def needs_to_convert_to_RGB(dicom_dataset):
    return False


def should_change_PhotometricInterpretation_to_RGB(dicom_dataset):
    return False


# 加载Dicom图像数据
def get_pixeldata(dicom_dataset):
    """If NumPy is available, return an ndarray of the Pixel Data. Raises ------ TypeError If there is no Pixel Data or not a supported data type. ImportError If NumPy isn't found NotImplementedError if the transfer syntax is not supported AttributeError if the decoded amount of data does not match the expected amount Returns ------- numpy.ndarray The contents of the Pixel Data element (7FE0,0010) as an ndarray. """
    if (dicom_dataset.file_meta.TransferSyntaxUID not in
            NumpySupportedTransferSyntaxes):
        raise NotImplementedError("Pixel Data is compressed in a "
                                  "format pydicom does not yet handle. "
                                  "Cannot return array. Pydicom might "
                                  "be able to convert the pixel data "
                                  "using GDCM if it is installed.")

    # 设置窗宽窗位
    #dicom_dataset.

    if not have_numpy:
        msg = ("The Numpy package is required to use pixel_array, and "
               "numpy could not be imported.")
        raise ImportError(msg)
    if 'PixelData' not in dicom_dataset:
        raise TypeError("No pixel data found in this dataset.")

    # Make NumPy format code, e.g. "uint16", "int32" etc
    # from two pieces of info:
    # dicom_dataset.PixelRepresentation -- 0 for unsigned, 1 for signed;
    # dicom_dataset.BitsAllocated -- 8, 16, or 32
    if dicom_dataset.BitsAllocated == 1:
        # single bits are used for representation of binary data
        format_str = 'uint8'
    elif dicom_dataset.PixelRepresentation == 0:
        format_str = 'uint{}'.format(dicom_dataset.BitsAllocated)
    elif dicom_dataset.PixelRepresentation == 1:
        format_str = 'int{}'.format(dicom_dataset.BitsAllocated)
    else:
        format_str = 'bad_pixel_representation'
    try:
        numpy_dtype = numpy.dtype(format_str)
    except TypeError:
        msg = ("Data type not understood by NumPy: "
               "format='{}', PixelRepresentation={}, "
               "BitsAllocated={}".format(
                   format_str,
                   dicom_dataset.PixelRepresentation,
                   dicom_dataset.BitsAllocated))
        raise TypeError(msg)

    if dicom_dataset.is_little_endian != sys_is_little_endian:
        numpy_dtype = numpy_dtype.newbyteorder('S')

    pixel_bytearray = dicom_dataset.PixelData

    if dicom_dataset.BitsAllocated == 1:
        # if single bits are used for binary representation, a uint8 array
        # has to be converted to a binary-valued array (that is 8 times bigger)
        try:
            pixel_array = numpy.unpackbits(
                numpy.frombuffer(pixel_bytearray, dtype='uint8'))
        except NotImplementedError:
            # PyPy2 does not implement numpy.unpackbits
            raise NotImplementedError(
                'Cannot handle BitsAllocated == 1 on this platform')
    else:
        pixel_array = numpy.frombuffer(pixel_bytearray, dtype=numpy_dtype)
    length_of_pixel_array = pixel_array.nbytes
    expected_length = dicom_dataset.Rows * dicom_dataset.Columns
    if ('NumberOfFrames' in dicom_dataset and
            dicom_dataset.NumberOfFrames > 1):
        expected_length *= dicom_dataset.NumberOfFrames
    if ('SamplesPerPixel' in dicom_dataset and
            dicom_dataset.SamplesPerPixel > 1):
        expected_length *= dicom_dataset.SamplesPerPixel
    if dicom_dataset.BitsAllocated > 8:
        expected_length *= (dicom_dataset.BitsAllocated // 8)
    padded_length = expected_length
    if expected_length & 1:
        padded_length += 1
    if length_of_pixel_array != padded_length:
        raise AttributeError(
            "Amount of pixel data %d does not "
            "match the expected data %d" %
            (length_of_pixel_array, padded_length))
    if expected_length != padded_length:
        pixel_array = pixel_array[:expected_length]
    if should_change_PhotometricInterpretation_to_RGB(dicom_dataset):
        dicom_dataset.PhotometricInterpretation = "RGB"
    if dicom_dataset.Modality.lower().find('ct') >= 0:  # CT图像需要得到其CT值图像
        pixel_array = pixel_array * dicom_dataset.RescaleSlope + dicom_dataset.RescaleIntercept  # 获得图像的CT值
    pixel_array = pixel_array.reshape(dicom_dataset.Rows, dicom_dataset.Columns*dicom_dataset.SamplesPerPixel)
    return pixel_array, dicom_dataset.Rows, dicom_dataset.Columns

step2:对于CT图像设置窗宽窗位

# 调整CT图像的窗宽窗位
def setDicomWinWidthWinCenter(img_data, winwidth, wincenter, rows, cols):
    img_temp = img_data
    img_temp.flags.writeable = True
    min = (2 * wincenter - winwidth) / 2.0 + 0.5
    max = (2 * wincenter + winwidth) / 2.0 + 0.5
    dFactor = 255.0 / (max - min)

    for i in numpy.arange(rows):
        for j in numpy.arange(cols):
            img_temp[i, j] = int((img_temp[i, j]-min)*dFactor)

    min_index = img_temp < 0
    img_temp[min_index] = 0
    max_index = img_temp > 255
    img_temp[max_index] = 255

    return img_temp

step3:获取Dicom中的tag信息

第一种方式:

# 加载Dicom图片中的Tag信息
def loadFileInformation(filename):
    information = {}
    ds = pydicom.read_file(filename)
    information['PatientID'] = ds.PatientID
    information['PatientName'] = ds.PatientName
    information['PatientBirthDate'] = ds.PatientBirthDate
    information['PatientSex'] = ds.PatientSex
    information['StudyID'] = ds.StudyID
    information['StudyDate'] = ds.StudyDate
    information['StudyTime'] = ds.StudyTime
    information['InstitutionName'] = ds.InstitutionName
    information['Manufacturer'] = ds.Manufacturer
    print(dir(ds))
    print(type(information))
    return information

第二种方式

dcm = pydicom.dcmread(fileanme)  # 加载Dicom数据

print(dcm[0x0008, 0x0060])
>>(0008, 0060) Modality                            CS: 'MR'
print(dcm[0x0008, 0x0060].VR)
>>CS
print(dcm[0x0008, 0x0060].value)
>>MR

step4:Dicom图像数据转换为PIL.Image

dcm = pydicom.dcmread(fileanme)  # 加载Dicom数据
dcm_img = Image.fromarray(img_data)  # 将Numpy转换为PIL.Image
dcm_img = dcm_img.convert('L')

# 保存为jpg文件,用作后面的生成label用
dcm_img.save('temp.jpg')
# 显示图像
dcm_img.show()

3. 结果展示

调整了窗宽窗位的脑部CT图像:
这里写图片描述

4. 参考资料

  1. Pydicom User Guide
  2. 【医学影像】窗宽窗位与其处理方法
版权声明:本文内容由互联网用户自发贡献,该文观点仅代表作者本人。本站仅提供信息存储空间服务,不拥有所有权,不承担相关法律责任。如发现本站有涉嫌侵权/违法违规的内容, 请联系我们举报,一经查实,本站将立刻删除。

发布者:全栈程序员-站长,转载请注明出处:https://javaforall.net/148747.html原文链接:https://javaforall.net

(1)
全栈程序员-站长的头像全栈程序员-站长


相关推荐

  • CTK插件框架学习4-创建跨平台插件工程「建议收藏」

    CTK插件框架学习4-创建跨平台插件工程「建议收藏」在上一篇博文中已经实现了一个简单的插件和测试程序的编写,但是插件跟应用是分开独立的工程。实际应用开发中需要把相关的库和头文件打包到一个工程中,如下图所示,这样比较方便调试开发,也为创建跨平台工程提供了便利。此节我们将创建一个初步完整的工程,工程文件中包含应用程序以及要使用的各个插件,同时将各个平台编译后的ctk插件库文件也整合到一起。目前支持如下三个平台:系统CPU编译器说明…

    2022年5月20日
    40
  • docker部署jenkins安装使用教程_docker搭建python开发环境

    docker部署jenkins安装使用教程_docker搭建python开发环境前言使用docker安装jenkins环境,jenkins构建的workspace目录默认是在容器里面构建的,如果我们想执行python3的代码,需进容器内部安装python3的环境。进jenki

    2022年7月31日
    3
  • php中PHPMailer发送带附件的电子邮件方法

    php中PHPMailer发送带附件的电子邮件方法

    2021年9月21日
    56
  • 物理讨论题复习

    物理讨论题复习请简要回答避雷针的工作原理?避雷针由于曲率半径小,电荷面密度大,从而产生尖端放电现象,导致自身与带电云层形成回路。导致自身电荷放出从而不会被雷击中,当带电云层密度过大,避雷针通过接地把电引下大地“分子电流假说“是谁提出的?请解释“分子电流”。安培。在原子、分子等物质微粒内部,存在着一种环形电流—-分子电流。分子电流使每个物质都成为微小的磁体,他的两侧相当于两个磁极请解释”磁偶极子“。磁偶极子是类比电偶极子而建立的物理模型。具有等值异号的两个点磁荷构成的系统称为磁偶极子。磁偶极子的物理

    2022年10月29日
    0
  • node.js常用npm命令

    node.js常用npm命令本文主要介绍npm的常用命令,如果用过淘宝镜像cnpm同样适用。特别注意,此处的指令多为node.js的依赖包,所以node.js是必不可少。一、安装node.js的依赖包Tips:每次都要打开cmd,进行指令操作,后续就不再提醒了。npminstallname>如:npminstallgulp默认安装express的最新版本如:npminstallgulp

    2022年7月16日
    12
  • Python 二进制,十进制,十六进制转换「建议收藏」

    Python 二进制,十进制,十六进制转换「建议收藏」十六进制到十进制使用int()函数,第一个参数是字符串’0Xff’,第二个参数是说明,这个字符串是几进制的数。 转化的结果是一个十进制数。>>>int(‘0xf’,16) 15二进制到十进制>>>int(‘10100111110′,2)   1342八进制到十进制>>>int(’17’,8)  15其实可以

    2022年5月19日
    41

发表回复

您的邮箱地址不会被公开。 必填项已用 * 标注

关注全栈程序员社区公众号