【发布时间】:2021-01-17 13:22:16
【问题描述】:
我有一个dicom 图像,但图像被填充。我有代码从图像中删除填充,以便只留下扫描,但我必须使用 ImageJ 打开图像并手动查找图像开始和结束位置的 x 和 y 轴的最小值和最大值。扫描的灰度值范围为-3000 to 2000。填充区域的值为0。有没有办法在不手动的情况下找到这些最小值和最大值?
原图:
所需图像:
【问题讨论】:
标签: python opencv padding medical-imaging
我有一个dicom 图像,但图像被填充。我有代码从图像中删除填充,以便只留下扫描,但我必须使用 ImageJ 打开图像并手动查找图像开始和结束位置的 x 和 y 轴的最小值和最大值。扫描的灰度值范围为-3000 to 2000。填充区域的值为0。有没有办法在不手动的情况下找到这些最小值和最大值?
原图:
所需图像:
【问题讨论】:
标签: python opencv padding medical-imaging
下面是一个使用 SimpleITK 裁剪背景的 Python 脚本。
基本思想是创建一个由不是背景值的像素组成的蒙版图像。然后它使用 SimpleITK 的 LabelShapeStatisticsImageFilter 来查找该掩码图像中非零像素的边界框。
import SimpleITK as sitk
img = sitk.ReadImage("padded-image.png")
# Grey background in this example
bg_value = 161
# Create a mask image that is just non-background pixels
fg_mask = (img != bg_value)
# Compute shape statistics on the mask
lsif = sitk.LabelShapeStatisticsImageFilter()
lsif.Execute(fg_mask)
# Get the bounds of the mask.
# Bounds are given as [Xstart, Ystart, Xwidth, Ywidth]
bounds = lsif.GetBoundingBox(1)
print(bounds)
Xmin_crop = bounds[0]
Ymin_crop = bounds[1]
Xmax_crop = img.GetWidth() - (bounds[0]+bounds[2])
Ymax_crop = img.GetHeight() - (bounds[1]+bounds[3])
# Crop parameters are how much to crop off each side
cropped_img = sitk.Crop(img, [Xmin_crop, Ymin_crop], [Xmax_crop, Ymax_crop])
sitk.Show(cropped_img)
sitk.WriteImage(cropped_img, "cropped-image.png")
因为我使用了您的 8 位 PNG 图像,所以背景值设置为 161。如果您使用原始的 16 位 DICOM CT,您将使用背景值 0。SimpleITK 可以读取 DICOM,以及其他图像格式的数量。
有关 LabelShapeStatisticsImageFilter 类的更多信息,这里是文档:https://simpleitk.org/doxygen/latest/html/classitk_1_1simple_1_1LabelShapeStatisticsImageFilter.html#details
【讨论】:
这是 Python/OpenCV 中使用颜色阈值和轮廓来查找边界框的另一种方法。
输入:
import cv2
import numpy as np
# read image
img = cv2.imread('scan.png')
# threshold on gray color (161,161,161)
lower = (161,161,161)
upper = (161,161,161)
thresh = cv2.inRange(img, lower, upper)
# invert threshold image so border is black and center box is white
thresh = 255 - thresh
# get external contours (presumably just one)
contours = cv2.findContours(thresh, cv2.RETR_EXTERNAL, cv2.CHAIN_APPROX_SIMPLE)
contours = contours[0] if len(contours) == 2 else contours[1]
cntr = contours[0]
x,y,w,h = cv2.boundingRect(cntr)
# crop to bounding rectangle
crop = img[y:y+h, x:x+w]
# save cropped image
cv2.imwrite('scan_thresh.png',thresh)
cv2.imwrite('scan_crop.png',crop)
cv2.imshow("THRESH", thresh)
cv2.imshow("CROP", crop)
cv2.waitKey(0)
cv2.destroyAllWindows()
裁剪结果:
【讨论】:
无需借助复杂图像分析的 SITK 或 CV 之类复杂(并且导入大/慢)的东西 - 您只需使用 numpy 即可轻松完成。
恕我直言,这将更快更可靠:
# if a is your image:
same_cols = np.all(a == a[0, :], axis=0)
same_cols_index = np.where(same_cols==False)[0]
C0,C1 = same_cols_index[0], same_cols_index[-1] + 1
same_rows = np.all(a == a[:, 0], axis=1)
same_rows_index = np.where(same_rows==False)[0]
R0,R1 = same_rows_index[0], same_rows_index[-1] + 1
print('rows', R0, R1)
print('cols', C0, C1)
a_snipped = a[R0:R1, C0:C1]
这里的逻辑是
例子
# make a sample image
a = np.zeros((512,512), dtype=np.int32)
r0, r1 = 53, 421
c0, c1 = 43, 470
rnd = np.random.randint(-3000, 2000, (r1-r0, c1-c0))
a[r0:r1, c0:c1] = rnd
plt.imshow(a, cmap='gray', vmin=-50, vmax=50)
same_cols = np.all(a == a[0, :], axis=0)
same_cols_index = np.where(same_cols==False)[0]
C0,C1 = same_cols_index[0], same_cols_index[-1] + 1
same_rows = np.all(a == a[:, 0], axis=1)
same_rows_index = np.where(same_rows==False)[0]
R0,R1 = same_rows_index[0], same_rows_index[-1] + 1
print('rows', R0, R1)
print('cols', C0, C1)
a_snipped = a[R0:R1, C0:C1]
plt.imshow(a_snipped, cmap='gray', vmin=-3000, vmax=2000)
第 53 421 行 43 470 列
【讨论】: