用于分割CT体模CT图像的直径测量的简单算法实例。

现在列出了三个程序，演示了应用于计算机断层扫描检查的基本形态图像处理。这些程序以MATLAB，Mathematica和Python计算机语言提供。请参阅其他存储库。程序读取DICOM格式文件以获得元信息和图像数据。分割图像以将感兴趣区域与背景分开。基于元信息中包括的像素尺寸，计算等效圆盘直径。然后使用包含在体像素内的像素衰减值计算等效水盘直径。

这些程序用于计算CT体模的等效圆盘直径。CT体模具有已知直径。因此，可以将结果与已知的真实值进行比较。在这种情况下，用于形态学算法的Mathematica参数基于体模的测量。MATLAB程序和Python程序都依赖于Otsu方法来确定形态分割的阈值。

这些程序旨在作为医学物理学家可以用来编写用于从CTDIvol计算SSDE的工具的示例。SSDE工具可以演示如何计算尺寸特定剂量估算是更好的辐射暴露量度，以优化辐射暴露和图像质量之间的权衡。大多数CT设施使用较旧的CT辐射测量。通过演示SSDE（可以自动化）的计算，可以说服CT设施采用更好的方法来优化辐射暴露。SSDE为预测图像噪声和辐射暴露提供了更好的参数。

如果用户计划使用这些程序，则需要进行一些修改。必须修改程序以反映发现轴向CT图像的文件夹（通常在光盘上提供）。还可能必须修改程序以反映DICOM图像头中包括的元信息中的不同值。不同CT制造商的DICOM标题中的元信息不同。DICOM标头元信息应包括有关各个像素尺寸的信息。Mathematica程序确实包含根据等效磁盘直径从CTDIvol（SSDE转换因子）计算SSDE的公式。这应该为用户提供足够的信息来修改程序供自己使用。

预备工作

这个Python（版本3.6，使用Spyder开发）程序演示了医学图像的处理，用于基本分割和直径测量。直径测量可用于计算尺寸特定剂量估计（SSDE）。SSDE是用于优化辐射剂量和图像诊断质量之间折衷的有价值的参数。
。
该算法是一种非常基本的方法，但对于此应用程序而言是有效且高效的。该算法在用于图像处理的Python程序的示例中详述。在这种情况下，该算法应用于CT体模的CT图像。体模具有已知的直径。这有助于确认算法的准确性。

介绍了使用Mathematica，MATLAB和Python的类似算法。人们对编程语言有自己的偏好。通过提出简单的方法（由于它在教科书和免费提供的代码示例中提供，因此不应具有可专利性），有些人可能会发现这些示例很有用。

这些程序可能是医学物理学家感兴趣的。该程序可用于向提供CT成像的设施演示使用SSDE优化CT检查的辐射暴露的优势。

代码

# -*- coding: utf-8 -*-
"""
Created on Wed Dec 26 13:28:16 2018
@author: John
"""

import os
import sys
import math
import numpy as np
import matplotlib.pyplot as plt 
import pathlib
import pydicom
from itertools import chain
from skimage.filters import threshold_otsu
from skimage.measure import label, regionprops
from scipy import ndimage as ndi
from skimage import morphology
from skimage.filters import threshold_otsu
from skimage.segmentation import clear_border
from skimage.measure import label, regionprops

ds=pydicom.dcmread('c:\\Users\\Public\\PhtmTest.dcm')
image=ds.pixel_array
pspx,pspy=ds.PixelSpacing
rescale_intercept=int(ds.RescaleIntercept)
threshold=threshold_otsu(image)
bw=clear_border(morphology.closing(image>threshold,morphology.square(3)))
plt.figure(dpi=300)
plt.title('image')
plt.imshow(image,cmap=plt.cm.bone)
plt.show()
plt.figure(dpi=300)
plt.title('bw')
plt.imshow(bw)
plt.show()
labelimage=label(bw,connectivity =1)
props=regionprops(labelimage)
area=[ele.area for ele in props ]
largest_blob_indx=np.argmax(area)
largest_blob_label=props[largest_blob_indx].label
body_pixels=np.zeros(image.shape,dtype=np.uint8)
body_pixels[labelimage == largest_blob_label] = 1
plt.figure(dpi=300)
plt.title('body pixels')
plt.imshow(body_pixels)
plt.show()
filled_body_pixels=ndi.binary_fill_holes(body_pixels)
plt.figure(dpi=300)
plt.title('filled body pixels')
plt.imshow(filled_body_pixels)
plt.show()
pixel_count=np.sum(filled_body_pixels)
pixel_radius=math.sqrt(pixel_count/math.pi)
eq_diam=.2*pspx*pixel_radius
pixel_attenuations= image[filled_body_pixels==1]
atest2=[(((i+rescale_intercept)/1000. )+1.0) for i in np.nditer(pixel_attenuations)]
tot_pix_area=sum(atest2)
pa_radius=math.sqrt(tot_pix_area/math.pi)
pa_diam=.2*pa_radius*pspx