之前有解釋預處理部分的函數,不過覺得還不夠詳細,同時文字解釋還不夠直觀,所以現在想一步步運行下,打印輸出
首先讀取原始數據,包括相應的注釋(即結節標簽)【注意】注釋文件中的標簽是按x,y,z的順序給的,但是origin以及spacing都是按照z,y,x的順序,所以要逆序處理一下([:,::-1])
raw_data,origin,spacing,isflip = load_itk_image("/home/dataset/LUNA16/subset3/1.3.6.1.4.1.14519.5.2.1.6279.6001.126631670596873065041988320084.mhd")
annos = np.array(pandas.read_csv("/home/dataset/LUNA16/CSVFILES/annotations.csv"))
this_annos = np.copy(annos[annos[:,0]==("1.3.6.1.4.1.14519.5.2.1.6279.6001.126631670596873065041988320084")])
raw_label = (this_annos[:,1:-1][:,::-1]-origin)/spacing
然后輸出下原始數據的shape,以及標簽的內容
print(raw_data.shape) print(raw_label)
輸出
(272, 512, 512)
[[83.84986792 232.12420469969365 348.18196764820215] [229.169256536 308.7302064585643 158.2345410262112]]
然后可視化其中一張切片
plt.imshow(raw_data[229])
效果
現在開始進入預處理流程
之前讀取過原始數據,此處不再讀取,除原始數據及標簽外,還需要讀取LUNA16為每個CT提供的掩碼,用於剔除肺部以外區域
Mask,origin,spacing,isflip = load_itk_image("/home/dataset/LUNA16/seg-lungs-LUNA16/1.3.6.1.4.1.14519.5.2.1.6279.6001.126631670596873065041988320084.mhd")
掩碼是與CT圖像同樣大小的三維圖像,區別在於掩碼只有0,1兩種值,這張CT的大小是(272,512,512),掩碼也是(272,512,512)
值得注意的是此處的掩碼卻並不是0,1二值,而是0,3,4三種值,0依然代表肺部以外區域,3代表左肺,4代表右肺,為便於將肺部一起處理,還需要將左右兩個肺合並一下
m1 = Mask==3 #LUNA16的掩碼有兩種值,3和4 m2 = Mask==4 Mask = m1+m2
然后這里做了一個我認為沒那么重要但是給處理帶來很多麻煩的環節,那就是對掩碼求取邊界(肺部邊界),只保留邊界內的數據
不妨想象一個正方形,里面有一個小圓,小圓就是掩碼,那么此處做的就是求小圓的最小外接矩形,將矩形外的部分砍掉
xx,yy,zz= np.where(Mask) box = np.array([[np.min(xx),np.max(xx)],[np.min(yy),np.max(yy)],[np.min(zz),np.max(zz)]])
打印輸出看一下box
array([[ 21, 264], [130, 407], [ 62, 441]])
做到這,預處理已經差不多了,再加幾步
box = box*np.expand_dims(spacing,1)/np.expand_dims(resolution,1) #對邊界即掩碼的最小外部長方體應用新分辨率 box = np.floor(box).astype('int')
打印輸出
array([[ 26, 330], [ 96, 302], [ 46, 327]])
對這個邊界向外擴展一點,為了處理邊緣的像素
margin = 5 extendbox = np.vstack([np.max([[0,0,0],box[:,0]-margin],0),np.min([newshape,box[:,1]+2*margin],axis=0).T]).T
打印輸出
array([[ 21, 340], [ 91, 312], [ 41, 337]])
然后對掩碼進行一點處理
convex_mask = m1
dm1 = process_mask(m1) #對掩碼采取膨脹操作,去除肺部黑洞
dm2 = process_mask(m2)
dilatedMask = dm1+dm2
Mask = m1+m2
extramask = dilatedMask ^ Mask #異或操作,求出相比於原始掩碼膨脹后多出來的區域
這里mask的大小沒有變過,仍然是(272,512,512)
bone_thresh = 210
pad_value = 170
sliceim = lumTrans(sliceim) #對原始數據閾值化,並歸一化
sliceim = sliceim*dilatedMask+pad_value*(1-dilatedMask).astype('uint8') #170對應歸一化話后的水,掩碼(膨脹過后)外的區域補充為水
bones = (sliceim*extramask)>bone_thresh #210對應歸一化后的骨頭,凡是大於骨頭的區域都填充為水
sliceim[bones] = pad_value
此時CT數據即sliceim的大小也沒有變過,現在要變化一下,進行重采樣
sliceim1,_ = resample(sliceim,spacing,resolution,order=1)
查看下此時的大小
(340, 380, 380)
最后還記得之前求的box嗎,我們只需要box內的數據即可
sliceim2 = sliceim1[extendbox[0,0]:extendbox[0,1], #將extendbox內數據取出作為最后結果
extendbox[1,0]:extendbox[1,1],
extendbox[2,0]:extendbox[2,1]]
處理完數據,還需要處理標簽
之前已經將世界坐標轉換為體素坐標,現在要對其應用新的分辨率(這里取[1,1,1])
raw_label = raw_label*spacing/resolution
輸出
array([[104.8123349, 172.279793861, 258.41647013999994], [286.46157067, 229.13584731999998, 117.43977386999997]], dtype=object)
最后的最后,減去box的下界
得到
[[83.84986792 232.12420469969365 348.18196764820215] [229.169256536 308.7302064585643 158.2345410262112]]
上面處理過的數據和標簽與完整預處理后的clean.npy和label.npy是一樣的,證明這個分解的過程沒什么紕漏
完結,撒花