[發(fā)明專利]一種基于ICP算法和B樣條的錐形束CT和CT圖像配準(zhǔn)方法

有效

• 申請?zhí)枺?/td>: 201710148067.3
• 申請日: 2017-03-13
• 公開（公告）號: CN106934821B
• 公開（公告）日: 2020-06-23
• 發(fā)明（設(shè)計）人: 李海;劉帥
• 申請（專利權(quán)）人: 中國科學(xué)院合肥物質(zhì)科學(xué)研究院
• 主分類號: G06T7/30
• 分類號: G06T7/30;G06T7/10
• 代理公司: 北京科迪生專利代理有限責(zé)任公司 11251
• 代理人: 楊學(xué)明;顧煒

1.一種基于ICP算法和B樣條的錐形束CT和CT圖像配準(zhǔn)方法，其特征在于，所述方法包括：

步驟1：獲取待配準(zhǔn)圖像，即CT圖像I₁和參考圖像，即錐形束CT圖像I₂；

步驟2：對所述待配準(zhǔn)圖像和參考圖像進(jìn)行圖像閾值分割，將所述待配準(zhǔn)圖像中的目標(biāo)物分割出來并提取其坐標(biāo)點，然后對坐標(biāo)點進(jìn)行采樣生成點集數(shù)據(jù)P＝{p₁,p₂,…,p_t}，其中p₁,p₂,…,p_t表示目標(biāo)物的坐標(biāo)點，t表示點集P中坐標(biāo)點個數(shù)，將所述參考圖像中的目標(biāo)物分割出來并提取其坐標(biāo)點，然后對坐標(biāo)點進(jìn)行采樣生成點集數(shù)據(jù)Q＝{q₁,q₂,…,q_t},其中q₁,q₂,…,q_t表示目標(biāo)物的坐標(biāo)點，點集Q中坐標(biāo)點個數(shù)同點集P中個數(shù)相同；

步驟3：對所述待配準(zhǔn)圖像點集數(shù)據(jù)P和參考圖像點集數(shù)據(jù)Q進(jìn)行基于點集的配準(zhǔn)即ICP仿射配準(zhǔn)，其步驟主要為：首先對所述仿射變換陣M_affine進(jìn)行初始化；利用仿射變換陣M_affine對所述待配準(zhǔn)點集P進(jìn)行變換得到點集P’即P’＝P*M_affine；然后使用公式計算所述參考圖像點集Q和P’之間的均方誤差，其中p_i和q_i分別表示第i個坐標(biāo)點；使用LBFGS算法不斷優(yōu)化M_affine對均方誤差求最小化，得到最優(yōu)解M_affine；

步驟4：根據(jù)所述仿射變換陣M_affine對所述待配準(zhǔn)圖像進(jìn)行仿射配準(zhǔn)得到ICP配準(zhǔn)圖像I_icp；

步驟5：將獲取的ICP配準(zhǔn)圖像I_icp作為基于體素的配準(zhǔn)初值，即對所述仿射配準(zhǔn)結(jié)果I_icp和所述參考圖像I₂進(jìn)行B樣條配準(zhǔn)得到最終結(jié)果I_final；

所述步驟2中的對所述待配準(zhǔn)圖像和參考圖像進(jìn)行圖像分割的步驟，包括：

對所述待配準(zhǔn)圖像和參考圖像的像素灰度進(jìn)行歸一化為[0,255]；

對所述歸一化處理后的待配準(zhǔn)圖像和參考圖像進(jìn)行閾值分割；

對所述閾值分割后的待配準(zhǔn)圖像和參考圖像進(jìn)行像素坐標(biāo)點提取和采樣并生成點集數(shù)據(jù)；

所述的對所述待配準(zhǔn)圖像和參考圖像的像素灰度進(jìn)行歸一化，采用如下公式：

其中g(shù)和g’分別表示歸一化后灰度值和原始圖像CT值，A和B分別用公式A＝wc-ww/2，B＝wc+ww/2，其中ww表示圖像上CT值的窗口范圍，wc表示圖像上CT值窗口中心；

所述對所述歸一化處理后的待配準(zhǔn)圖像和參考圖像進(jìn)行閾值分割的目標(biāo)物為骨骼區(qū)域，其中骨骼的CT值范圍為[100,1000],使用歸一化公式計算其灰度值；

所述待配準(zhǔn)圖像和參考圖像的目標(biāo)進(jìn)行閾值分割后提取其坐標(biāo)點，然后對坐標(biāo)點進(jìn)行采樣生成點集數(shù)據(jù)的方法為：首先提取出256灰度級中骨骼范圍的坐標(biāo)點，然后對提取出的坐標(biāo)點進(jìn)行1/50等間距采樣得到點集圖像；

所述步驟5中，B樣條配準(zhǔn)采用多分辨率B樣條彈性配準(zhǔn)方法，多分辨率層數(shù)分為1到3層，步驟為：使用公式：

對輸入圖像進(jìn)行高斯平滑處理，其中I₂為錐形束CT圖像像素的坐標(biāo)位置矩陣，I_icp為ICP配準(zhǔn)結(jié)果圖像像素的坐標(biāo)位置矩陣，σ為正態(tài)分布的標(biāo)準(zhǔn)偏差，根據(jù)圖像大小進(jìn)行降采樣實施多分辨率配準(zhǔn)，先將圖像在低分辨進(jìn)行B樣條粗配準(zhǔn)，首先初始化控制點，然后對兩幅圖像進(jìn)行相似度測量，通過對相似度值的判斷修改控制點，將修改后的控制點代入B樣條變換函數(shù)對待配準(zhǔn)圖像進(jìn)行變換,使用LBFGS優(yōu)化方法不斷修改控制點對相似度值進(jìn)行最小化求解，從而得到最終配準(zhǔn)結(jié)果I_final；

所述B樣條變換函數(shù)為：采用公式求取像素點在(x，y，z)位置處的形變場，其中：

其中表示向下取整，φ_i,j,k表示x，y，x軸控制點間距分別為δ_i,δ_j,δ_k，大小為n_x×n_y×n_z的網(wǎng)格中序號為ijk的控制點的值,其中i，j，k為了求像素點在(x,y,z)處的位移值T(x,y,z)僅與相鄰的4³個控制點有關(guān)；B_l,B_m,B_n均為三次B樣條函數(shù)基函數(shù)，B_l表示B樣條函數(shù)中第l個基函數(shù)，其計算公式為B₀(u)＝(1-u)³/6，B₁(u)＝(3u³-6u²+4)/6，B₂(u)＝(-3u³+3u²+3u+1)/6，B₃(u)＝u³/6,其中0≤u≤1， B_m和B_n的計算同B_l，起到權(quán)函數(shù)的作用，根據(jù)控制點到(x,y,z)的距離大小來加權(quán)每個控制點對T(x,y,z)的貢獻(xiàn)；l,m,n分別為公式中的求和序號，取值為0,1,2,3；

所述多分辨率層數(shù)分為1到3層，每層XYZ軸最小分辨率參數(shù)分別設(shè)為128×128×10,256×256×20,512×512×40；每層分辨率的B樣條網(wǎng)格間距分別設(shè)為：32×32×16,16×16×8,8×8×4。

2.根據(jù)權(quán)利要求1所述的基于ICP算法和B樣條的錐形束CT和CT圖像配準(zhǔn)方法，其特征在于：對兩幅圖像進(jìn)行相似度測量時，采用灰度均方差作為相似度測量。