python x y 坐标 z 代表值_luna16标签数据里的xyz,以及CT的dicom.ImagePositionPatient里的三个值分别代表哪些轴的初始点...

在来看看dicom.ImagePositionPatient中三个值代表的是哪个的那些轴注意看第一段代码中的x对应Origin[0],y对应Origin[1],z对应Origin[2].而Origin的获取是以下代码:Ps:这里面的一些道道给了我很大的困扰,这里在我脑袋清晰的时候将理解记录下来防止以后再次弄混.

先来看张图看看是怎么照CT:

所以先不看Z轴,只看X轴和Y轴这个面,看出X是列,而Y是行.这跟我们平时在图上画点很类似.图参考:https://www.cnblogs.com/h2zZhou/p/9072967.html

而dicom.ImagePositionPatient是直接从设备生成的dicom里读取的故而dicom.ImagePositionPatient里的三个值分别代表X轴,Y轴,Z轴的初始点.一定注意X是列,而Y是行.

x, y, z = int((nodules[idx, 0]-Origin[0])/SP[0]), int((nodules[idx, 1]-Origin[1])/SP[1]),

int((nodules[idx, 2]-Origin[2])/SP[2])

# 注意 y代表纵轴,x代表横轴

data[max(0, y - radius):min(data.shape[0], y + radius),

max(0, x - radius - pad):max(0, x - radius)] = 3000  # 竖线

data[max(0, y - radius):min(data.shape[0], y + radius),

min(data.shape[1], x + radius):min(data.shape[1], x + radius + pad)] = 3000  # 竖线

data[max(0, y - radius - pad):max(0, y - radius),

max(0, x - radius):min(data.shape[1], x + radius)] = 3000  # 横线

data[min(data.shape[0], y + radius):min(data.shape[0], y + radius + pad),

max(0, x - radius):min(data.shape[1], x + radius)] = 3000  # 横线

你可注意到x对应nodules[idx, 0],y对应nodules[idx,1],z对应nodules[idx,2].而后后面这段标出标签的方式就是赋一个很大值.但是有注释y代表纵轴,x代表横轴.在看上面那段代码第一行y相关的参数作为data赋值时的第一维度更加印证了上面的注释.

比如我们知道二维矩阵A,赋值取到点(x1,y1)代表A[x1,y1].而画图时想要标记出点(x1,y1)则需要如下代码:

plt.ion()

plt.scatter(y1, x1, marker='o', color='red', s=10, label='First')

plt.imshow(A, cmap=plt.cm.bone)

如果还模棱两可,再来看看第一段代码中的x对应Origin[0],y对应Origin[1],z对应Origin[2].而Origin的获取是以下代码:

filename='E:\\JLS\\dcm_data\\luna\\subset1\\1.3.6.1.4.1.14519.5.2.1.6279.6001.173106154739244262091404659845.mhd'

itkimage = sitk.ReadImage(filename)#读取.mhd文件

Origin=itkimage.GetOrigin()

而itkimage.GetOrigin()与这个CT序列的第一张dicom的ImagePositionPatient是一样的.这样就可以看出luna16的标签数据与dicom.ImagePositionPatient代表的轴是一样的,也是coordX代表列,coordX代表行.

PS:2019.11.6 博客的思路结论写了三遍,大体过程:是这样的,不反了反了是那样的,还是这样的.最后这次为了确保正确我还再次写了按照标签分别画出结果用luna16数据和目前在做项目的数据.具体代码如下:

# coding=UTF-8

import cv2

import numpy as np

from PIL import Image

from scipy import misc

import os

import sys

import cv2

from skimage import measure, morphology, segmentation

import matplotlib.pyplot as plt

import SimpleITK as sitk

import pandas

MIN_BOUND = -1000.0

MAX_BOUND = 400.0

def normalize(image):

image = (image - MIN_BOUND) / (MAX_BOUND - MIN_BOUND)

image[image > 1] = 1.

image[image < 0] = 0.

return image

def load_itk_image(filename):

with open(filename) as f:

contents = f.readlines()

line = [k for k in contents if k.startswith('TransformMatrix')][0]

transformM = np.array(line.split(' = ')[1].split(' ')).astype('float')

transformM = np.round(transformM)

if np.any(transformM != np.array([1, 0, 0, 0, 1, 0, 0, 0, 1])):

isflip = True

else:

isflip = False

itkimage = sitk.ReadImage(filename)

numpyImage = sitk.GetArrayFromImage(itkimage)

numpyOrigin = np.array(list(reversed(itkimage.GetOrigin())))

numpySpacing = np.array(list(reversed(itkimage.GetSpacing())))

return numpyImage, numpyOrigin, numpySpacing, isflip

F_path='Tannotations.csv'

L_path='annotations.csv'

F='/FORNEWTEST/change/NewForTest'

L='dcm_data/luna/subsets'

annosF = np.array(pandas.read_csv(F_path))

annosL = np.array(pandas.read_csv(L_path))

nameF = annosF[0][0]

nameL = annosL[55][0]

coordXF,coordYF,coordZF=annosF[0][1],annosF[0][2],annosF[0][3]

coordXL,coordYL,coordZL=annosL[55][1],annosL[55][2],annosL[55][3]

print(nameF,nameL)

for root, dirs, files in os.walk(F):

for file in files:

if (file[:-4] == nameF and file[-4:] == '.mhd'):

sliceim, origin, spacing, isflip = load_itk_image(

os.path.join(root, file))

ZF=np.absolute(coordZF-origin[0])/spacing[0]

YF=np.absolute(coordYF-origin[1])/spacing[1]

XF = np.absolute(coordXF - origin[2]) / spacing[2]

imgF=sliceim[ZF]

plt.ion()

plt.scatter(XF, YF, marker='o', color='red', s=10, label='First')

plt.imshow(imgF, cmap=plt.cm.bone)

plt.show()

for root, dirs, files in os.walk(L):

for file in files:

if (file[:-4] == nameL and file[-4:] == '.mhd'):

sliceim, origin, spacing, isflip = load_itk_image(os.path.join(root, file))

ZL=int(np.absolute(coordZL-origin[0])/spacing[0])

YL=np.absolute(coordYL-origin[1])/spacing[1]

XL = np.absolute(coordXL - origin[2]) / spacing[2]

imgL=sliceim[ZL]

imgC=imgL[YL-20:YL+20,XL-20:XL+20]

plt.ion()

plt.scatter(XL, YL, marker='o', color='red', s=2, label='First')

plt.imshow(imgL, cmap=plt.cm.bone)

plt.show()

plt.imshow(imgC, cmap=plt.cm.bone)

plt.show()

现在思路终于十分清晰了.真是要眼见为实!

  • 0
    点赞
  • 1
    收藏
    觉得还不错? 一键收藏
  • 0
    评论
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值