基于CT圖像的骨骼缺損修復方法研究

劉暢 周競宇

摘要：在3D打印領域中，針對找到受損骨骼待修補部分的問題上，提出一種利用模板進行配準找出需要修補部分的方法。這種模板可以是相似度高的模型，也可以是病人的鏡像骨。在完成方法過程中，利用到了最小盒算法和ICP配準，前者是對模型的粗配準，后者是對模型的精配準;在最后進行剔除錯誤部分和平滑，生成打印模型。實驗表明，該方法能完整精確打印出待修補處三維模型，與缺損模型平滑連接。

關鍵詞：3D打印;骨骼修復;三維重建;剛性配準;點云

中圖分類號：TP301.6 文獻標識碼：A

文章編號：1009-3044（2019）31-0231-03

1概述

隨著3D打印技術的飛速發(fā)展，醫(yī)學領域和計算機領域的快速結合，也帶來了更多優(yōu)秀人才去加入研究骨骼缺損修復這一領域。多方面領域人才的結合更快的推動了骨骼修復的發(fā)展，生物醫(yī)療研究人員致力于生產(chǎn)優(yōu)質(zhì)的仿生材料去貼合人體構造，而在計算機領域則是構造更完美的模型去進行打印。因此，如何設計更適應骨骼缺陷的三維模型成了如今研究的重點。

在構建契合破損骨骼的三維模型時，關鍵在于我們要得到精確的形狀，并使其能夠與正常骨平滑連接。在本文的研究中，我們應需要找到合適的待匹配模板，而這種模板是自于模型庫或病人的對側骨。所有的后續(xù)研究都是在這樣的基礎上進行的。

2總體流程

該修復方法首先重建受損骨骼模型嘲由已知的CT序列得到;接著對模型預處理，選取合適的模板;再通過對應的點云模型求取變換矩陣組，最后進行模型配準和計算缺損部分，具體流程如圖1所示。重建過程在win10系統(tǒng)環(huán)境下完成，采用

2.1骨骼模型重建

根據(jù)原始CT圖像的橫斷面、冠狀面及矢狀面逐層分割，能更好地提高圖像的完整性和準確性。CT數(shù)值反應的是組織的密度，我們可以參照正常人體組織CT值表，通過調(diào)窗去除大部分非骨骼組織，實現(xiàn)粗分割;然后再利用分割算法[3]進行精分割，最后對分割結果平滑處理。

2.1.1骨骼粗分割

通過查表可知，骨的平掃CT值在150～1000，因此我們對原始數(shù)據(jù)進行范圍規(guī)定化，在此過程中，若CT值小于150，則灰度值置為0，;若CT值高于1000，則灰度值置為255;若CT值大于150而小于1000，則灰度值按照（1）得出，在（1）中g（i，j，k）為灰度值，f（i，j，k）為CT值。

g（i，j，k）=（f（i，j，k）-150/（1000-150）*255 （1）

2.1.2骨骼精分割

由于目標具有相同的灰度特征，我們使用Otsu全局閾值分割方法方法是根據(jù)圖像上灰度值的分布將圖像分成背景和前景兩部分看待，前景就是我們要按照閾值分割出來的部分。背景和前景的分界值就是我們要求出的閾值。遍歷不同的閾值，計算不同閾值下對應的背景和前景之間的類間方差，當類間方差取得極大值時，此時對應的閾值就是OTSU算法所求的閾值。

2.1.3平滑去噪

在進行骨骼精分割以后，圖像中也存在少許噪點。此外，骨骼邊緣比較粗糙，會直接影響后續(xù)點云模型的精度。所以我們還需要進行平滑處理。

步驟1形態(tài)學開運算去噪：可選擇小半徑的圓形結構元，對精分割結果進行開運算?？稍诒３止趋啦蛔兊那闆r下去除小噪點，并斷開不同骨骼間的偽連接。

步驟2高斯濾波平滑邊緣：選擇小半徑高斯濾波器對去噪結果進行平滑，然后以閾值T=0.5分割平滑結果。

2.2求取缺損部分模型

求取缺損模型首先應該是在適配模型庫中找到與缺損模型最相近的完整模型，由于我們需要得到最契合的缺損模型，因此本文提出的方法是在獲取模型點云的基礎上進行配準的。利用MarchingCubes算法獲取體模表面三角片，并計算三角片重心作為點云集。接著進行點云降采樣，這樣做的目的是減少點云數(shù)量，提高配準速率。

2.2.1點云粗配準

在得到匹配模板和缺損模型的點云后，利用最小包圍盒算法進行點云的粗配準。最小包圍盒也叫外接最小矩陣，是一種求解離散點集最優(yōu)包圍空間的算法，基本思想是用體積稍大且特性簡單的幾何體（稱為包圍盒）來近似地代替復雜的幾何對象。

首先，獲取已知點云的三個主方向，獲取質(zhì)心，計算協(xié)方差，獲得協(xié)方差矩陣，求取協(xié)方差矩陣的特征值和特長向量，特征向量即為主方向。接著，利用已獲得的主方向和質(zhì)心，將輸入點云轉(zhuǎn)換至原點，且主方向與坐標系方向重回，建立變換到原點的點云的包圍盒。最后，給輸入點云設置主方向和包圍盒，通過輸人點云到原點點云變換的逆變換實現(xiàn)。

通俗來說，最小包圍盒算法是求解點集XOY平面內(nèi)投影點的最小面積包圍盒，即把空間化問題轉(zhuǎn)化為平面問題來求解。

2.2.2點云精配準

由蘭州交通大學楊軍等人提出的改進的ICP算法是本文用于點云精配準的主要方法，在該方法中主要在初始值的敏感性以及整體與部分模型的配準精度上做出了改進。選擇改進的ICP算法適用于精配準是由于在初始配準后還未達到想要的效果，因此選用KD-Tree搜索大幅度地提高了搜索效率，制定自適應最有閾值則是在經(jīng)典的ICP算法上提高了配準精度。

假如我們手上的數(shù)據(jù)有原始點云數(shù)據(jù)P和經(jīng)過平移加旋轉(zhuǎn)后的點云Q，那么ICP算法就是求出最優(yōu)的從P到Q的坐標變換矩陣。從（2）中我們可以看到我們從P到Q的過程，其中旋轉(zhuǎn)為R.平移為T。

2.3計算缺損部分

為解決無法通過點云直接得到缺損部分的問題，我們以點云模型為基準建立初始采樣網(wǎng)格，后利用變換矩陣對網(wǎng)格進行逆運算，再對原始體數(shù)據(jù)進行采樣，即可得到配準后的體模型。

把適配模型網(wǎng)格稱為模板，缺損模型網(wǎng)格稱為實例，那么現(xiàn)在用模板減去實例就可以得到粗結果。當然我們應該在粗結果上進行進行錯誤部分的提出，其中包括刪除扁平、細長、顆粒狀部分，最后我們需要利用三維高斯核和球形中值濾波器對剩余部分進行平滑處理。

3實驗

在本實驗中選用羊的脛骨作為實驗對象，數(shù)據(jù)集一共83張516KB的CT圖像，CT數(shù)據(jù)如圖2所示。讀取整個序列，并根據(jù)dicom參數(shù)（ResaleSlope、RescaleIntercept、WindowCenter、WindowWidth、ShceLocationl進行調(diào)窗與調(diào)序，獲得數(shù)值為0～255的灰度圖像;接著利用閾值分割提取所有骨骼，利用連通分量篩選目標骨骼;再根據(jù)dicom參數(shù)（PixelSpacing、SliceTh-ickness），重建序列體數(shù)據(jù)，過程示例如圖3所示。

在選取到適匹配模板后，那就可以生成缺損模型和匹配到的完整模型的點云，如圖4所示，經(jīng)過點云的粗配準即最小包圍盒算法以及精配準ICP算法后的點云如圖5所示，在經(jīng)過這一步操作后我們可以得到變換矩陣，這一矩陣的得到是為了下一步做鋪墊。

在得到缺損模型前，需要用采樣網(wǎng)格對點云進行包圍，在圖6中三維坐標下豎直采樣網(wǎng)格代表缺損模型，傾斜采樣網(wǎng)格代表匹配到的完整模型，在已知變換矩陣的條件下對采樣網(wǎng)格進行逆運算，從而得到配準后模型。

那么在這樣的情況下，就可以進行兩模型的直接相減得到粗略的缺損模型，如圖7所示，在經(jīng)過剔除錯誤部分和平滑后，就可以生成打印模型。MarchingCubes面繪制算法獲取三角片，并以光順優(yōu)化獲得更好的三角片，最后以sTL的格式導出，結果的缺損模型如圖8所示。