2D-3D配準(zhǔn)中的初始位姿估計(jì)方法

孫?濤，郭?可，劉傳耙，張?弢，宋軼民, 3，馬信龍

2D-3D配準(zhǔn)中的初始位姿估計(jì)方法

孫?濤1，郭?可1，劉傳耙1，張?弢2，宋軼民1, 3，馬信龍2

(1．天津大學(xué)機(jī)構(gòu)理論與裝備設(shè)計(jì)教育部重點(diǎn)實(shí)驗(yàn)室，天津 300354；2. 天津市天津醫(yī)院，天津 300211；3. 天津仁愛(ài)學(xué)院機(jī)械工程系，天津 301636)

手術(shù)導(dǎo)航；2D-3D配準(zhǔn)；初始位姿估計(jì)；正交雙平面投影

在骨科手術(shù)中，通常采用成像較快的熒光透視圖像為醫(yī)生提供術(shù)中手術(shù)器械與病患部位的相對(duì)位置信息[1-2]，但術(shù)中二維圖像難以準(zhǔn)確反映病患部位的空間位置和解剖信息，因此需要與術(shù)前三維CT數(shù)據(jù)結(jié)合，即通過(guò)2D-3D圖像配準(zhǔn)建立2D像素平面與3D模型空間的統(tǒng)一，實(shí)現(xiàn)術(shù)中、術(shù)前圖像在空間位置與解剖紋理上的對(duì)應(yīng)關(guān)系，從而引導(dǎo)手術(shù)的進(jìn)行，這是病灶精準(zhǔn)定位和施術(shù)成功的關(guān)鍵[3]．

綜上，本文提出一種新的初始位姿估計(jì)方法，通過(guò)正交融合投影原理，構(gòu)建空間與平面參數(shù)間的映射關(guān)系，結(jié)合區(qū)域生長(zhǎng)算法完成對(duì)多研究對(duì)象的5個(gè)三維空間參數(shù)的有效估計(jì)，以增大2D-3D配準(zhǔn)的初始變換范圍、降低配準(zhǔn)時(shí)間、提高配準(zhǔn)的精度和魯棒性，并進(jìn)行實(shí)驗(yàn)驗(yàn)證．

1?初始位姿估計(jì)方法

1.1?參數(shù)映射關(guān)系

大部分2D-3D醫(yī)學(xué)圖像配準(zhǔn)是將三維CT數(shù)據(jù)通過(guò)數(shù)字重建影像技術(shù)生成DRR(digitally reconstructed radiograph)圖像，比較DRR圖像與術(shù)中平面圖像的相似性，再不斷改變投影環(huán)境下三維物體的空間參數(shù)，使得在某個(gè)空間參數(shù)下投影生成的DRR圖像與平面圖像的相似性達(dá)到最大[17]，此時(shí)的空間參數(shù)就是相同相機(jī)坐標(biāo)投影環(huán)境下生成平面圖像的三維物體的位姿．因?yàn)槿S物體的空間位姿的改變會(huì)引起DRR圖像的改變，所以可建立空間與平面參數(shù)間的映射關(guān)系．

首先通過(guò)模擬術(shù)中C臂成像設(shè)備的成像過(guò)程，構(gòu)建相應(yīng)的投影環(huán)境來(lái)生成DRR圖像，其中正交融合投影過(guò)程的示意如圖1所示．所構(gòu)建的投影坐標(biāo)系原點(diǎn)為三維影像中心點(diǎn)，其坐標(biāo)軸與經(jīng)CT掃描后所成圖像自帶的坐標(biāo)軸平行．

圖1?正交融合投影示意

式(2)計(jì)算可得

其中

根據(jù)上述推導(dǎo)過(guò)程，可得空間與平面參數(shù)間的映射關(guān)系，但平面的參數(shù)信息包含在DRR圖像灰度信息之中，不能直接提取，因此需要在后續(xù)利用對(duì)應(yīng)目標(biāo)圖像與正交雙平面模板間的配準(zhǔn)來(lái)獲取平面參數(shù)信息．

1.2?初始位姿估計(jì)原理與方法

圖2?初始位姿估計(jì)原理

然后，根據(jù)需要配準(zhǔn)的目標(biāo)投影生成正交雙平面模板，將目標(biāo)圖像和模板的灰度信息進(jìn)行平面配準(zhǔn)，計(jì)算得到正側(cè)位平面配準(zhǔn)參數(shù)．本文采用歸一化互相關(guān)[20]作為相似性測(cè)度來(lái)計(jì)算初始位姿估計(jì)以及后續(xù)2D-3D配準(zhǔn)過(guò)程的平面圖像相似性．其中，歸一化互相關(guān)系數(shù)的計(jì)算過(guò)程可以表示為

2?實(shí)?驗(yàn)

2.1?實(shí)驗(yàn)數(shù)據(jù)與測(cè)試環(huán)境

本文使用公開(kāi)的顱骨數(shù)據(jù)[19]和天津醫(yī)院提供的完整股骨、包含上下兩骨塊(近端骨和遠(yuǎn)端骨)的骨折股骨的CT數(shù)據(jù)來(lái)驗(yàn)證上文提出的初始位姿估計(jì)方法．為了方便獲得任意位姿下的平面影像，本實(shí)驗(yàn)采用DRR影像模擬術(shù)中熒光圖像進(jìn)行配準(zhǔn)．因?yàn)闊晒鈭D像成像時(shí)，成像設(shè)備發(fā)出的X射線透過(guò)病患部位后強(qiáng)度會(huì)減少，不同的組織結(jié)構(gòu)造成的X射線衰減程度有差異，這種差異在熒光圖像上顯示為各組織結(jié)構(gòu)明暗的不同．而DRR影像成像原理是模擬這種透視過(guò)程，故模擬生成的熒光影像與術(shù)中熒光影像具有一致性，這樣不僅減少了患者和醫(yī)生受到的輻射，也可以將圖像生成的真實(shí)位姿作為配準(zhǔn)結(jié)果的評(píng)估標(biāo)準(zhǔn)，具有可重復(fù)性．同時(shí)在實(shí)際應(yīng)用中，本文的數(shù)據(jù)處理算法對(duì)于術(shù)中熒光影像與DRR圖像的數(shù)據(jù)差異具有魯棒性．

本實(shí)驗(yàn)使用Visual Studio 2017作為軟件平臺(tái)，結(jié)合圖像開(kāi)發(fā)工具包ITK(insight segmentation and registration toolkit)編程實(shí)現(xiàn)初始位姿估計(jì)與2D-3D配準(zhǔn)過(guò)程，電腦配置為Intel?Xeon?Gold 5120 CPU @2.20GHz處理器，內(nèi)存為64GB．

2.2?評(píng)估方法

以投影生成仿術(shù)中平片時(shí)三維物體的空間位姿作為“金標(biāo)準(zhǔn)”，用生成的匹配模板進(jìn)行初始位姿估計(jì)，統(tǒng)計(jì)估計(jì)值與實(shí)際真值的誤差和不需要人工校正的比例．因?yàn)樵?D-3D圖像配準(zhǔn)中使用雙平面方法的配準(zhǔn)精度要優(yōu)于單平面[22]，且目前多采用Powell優(yōu)化方法[23]進(jìn)行配準(zhǔn)，所以本文用估計(jì)的參數(shù)值作為起始搜索位姿，采用雙平面Powell法進(jìn)行2D-3D配準(zhǔn)(以下簡(jiǎn)稱(chēng)結(jié)合法)，將得到的結(jié)果與不使用初始位姿估計(jì)的雙平面Powell法(以下簡(jiǎn)稱(chēng)單一法)進(jìn)行配準(zhǔn)的結(jié)果進(jìn)行統(tǒng)計(jì)，比較最終配準(zhǔn)結(jié)果與實(shí)際真值的誤差．

2.3?實(shí)驗(yàn)結(jié)果與分析

表1?初始位姿估計(jì)結(jié)果

Tab.1?Results of initial pose estimation

顱骨實(shí)驗(yàn)中，單一法與結(jié)合法都有較好的配準(zhǔn)結(jié)果，因此統(tǒng)計(jì)每次迭代后所有平移參數(shù)和所有旋轉(zhuǎn)參數(shù)配準(zhǔn)誤差絕對(duì)數(shù)的平均值(以下分別簡(jiǎn)稱(chēng)為平移配準(zhǔn)誤差和旋轉(zhuǎn)配準(zhǔn)誤差)并進(jìn)行繪制，如圖3所示；在完整股骨、骨折股骨近端骨和遠(yuǎn)端骨的實(shí)驗(yàn)中，有10%的實(shí)驗(yàn)配準(zhǔn)結(jié)果嚴(yán)重偏離“金標(biāo)準(zhǔn)”，且全來(lái)自于單一法，因此本文僅對(duì)剩余實(shí)驗(yàn)結(jié)果進(jìn)行統(tǒng)計(jì)，繪制每次迭代后的平移配準(zhǔn)誤差和旋轉(zhuǎn)配準(zhǔn)誤差圖，如圖4～圖6所示．由于大部分實(shí)驗(yàn)經(jīng)過(guò)6次迭代后，配準(zhǔn)結(jié)果趨于相對(duì)穩(wěn)定，因此只統(tǒng)計(jì)各實(shí)驗(yàn)前6次迭代過(guò)程中的配準(zhǔn)誤差．

圖3?顱骨實(shí)驗(yàn)中每次迭代的配準(zhǔn)誤差

圖4?股骨實(shí)驗(yàn)中每次迭代的配準(zhǔn)誤差

圖5?近端骨實(shí)驗(yàn)中每次迭代的配準(zhǔn)誤差

圖6?遠(yuǎn)端骨實(shí)驗(yàn)中每次迭代后的配準(zhǔn)誤差

各研究對(duì)象配準(zhǔn)統(tǒng)計(jì)結(jié)果表明，使用單一法進(jìn)行2D-3D配準(zhǔn)時(shí)，平均迭代次數(shù)為7.5次，平均配準(zhǔn)時(shí)間為8089s；而基于結(jié)合法的2D-3D配準(zhǔn)平均迭代次數(shù)為3次，平均配準(zhǔn)時(shí)間(包括初始位姿估計(jì)時(shí)間)為4330s．與單一法相比，結(jié)合法在平均迭代次數(shù)上降低60%，平均配準(zhǔn)時(shí)間上縮短46.5%，并且具有精度高、魯棒性好的優(yōu)點(diǎn)，最終2D-3D配準(zhǔn)結(jié)果中各參數(shù)的配準(zhǔn)誤差的平均值±標(biāo)準(zhǔn)差統(tǒng)計(jì)結(jié)果如表2?所示．

表2?2D-3D參數(shù)配準(zhǔn)誤差

Tab.2?2D/3D parameter registration errors

圖7?各實(shí)驗(yàn)的參數(shù)配準(zhǔn)誤差精度提升率

3?結(jié)?論

本文提出了一種正交雙平面模板匹配初始位姿估計(jì)方法，該方法利用正交融合投影原理構(gòu)建了正交投影下三維空間參數(shù)與平面圖像坐標(biāo)參數(shù)間的映射關(guān)系，將已知位姿投影下的DRR圖像作為匹配模板與待配準(zhǔn)的圖像進(jìn)行比較，通過(guò)平面坐標(biāo)參數(shù)的變化來(lái)對(duì)空間參數(shù)值進(jìn)行估計(jì)，并且基于區(qū)域生長(zhǎng)算法實(shí)現(xiàn)多個(gè)研究對(duì)象的配準(zhǔn)．為了驗(yàn)證本文提出的初始位姿估計(jì)方法的有效性，采用了顱骨、完整股骨、骨折股骨的近端骨和遠(yuǎn)端骨作為研究對(duì)象，各目標(biāo)通過(guò)20組實(shí)驗(yàn)進(jìn)行初始位姿估計(jì)與后續(xù)2D-3D配準(zhǔn)實(shí)驗(yàn)對(duì)比．得到全文結(jié)論如下．

(2) 在2D-3D配準(zhǔn)實(shí)驗(yàn)中，相較于單一法，結(jié)合法在整體配準(zhǔn)時(shí)間上平均縮短46.5%，在各實(shí)驗(yàn)中單個(gè)平移參數(shù)精度最高提升46.5%，單個(gè)旋轉(zhuǎn)參數(shù)精度最高提升91.1%，同時(shí)實(shí)驗(yàn)結(jié)果也兼具了良好的魯?棒性.

(3) 通過(guò)各組實(shí)驗(yàn)，驗(yàn)證了在較大的初始位姿偏差下，本文提出的初始位姿估計(jì)方法對(duì)不同研究對(duì)象都具有很好的適應(yīng)性，這也將有助于為骨科手術(shù)中的多模影像導(dǎo)航系統(tǒng)提供最佳或接近最佳的搜索參數(shù)選擇．

［1］ Russakoff D B，Rohlfing T，Mori K，et al. Fast generation of digitally reconstructed radiographs using attenuation fields with application to 2D-3D image registration[J]. IEEE Transactions on Medical Imaging，2005，24(11)：1441-1454.

［2］ Khamene A，Bloch P，Wein W，et al. Automatic registration of portal images and volumetric CT for patient positioning in radiation therapy[J]. Medical Image Analysis，2006，10(1)：96-112.

［3］ Markelj P，Tomazevic D，Likar B，et al. A review of 3D/2D registration methods for image-guided interventions[J]. Medical Image Analysis，2012，16(3)：642-661.

［4］ Frederik M，Vandermeulen D，Paul S，et al. Comparative evaluation of multiresolution optimization strategies for multimodality image registration by maximization of mutual information[J]. Medical Image Analysis，1999，3(4)：373-386.

［5］ Kaiser M，John M，Heimann T，et al. 2D/3D registration of TEE probe from two non-orthogonal c-arm directions[C]//International Conference on Medical Imaging Computing and Computer-Assisted Intervention. Boston，USA，2014：283-290.

［6］ Muhit A，Pickering M R，Ward T，et al. A comparison of the 3D kinematic measurements obtained by single-plane 2D-3D image registration and RSA[C]//32nd Annual International Conference of the IEEE Engineering in Medicine and Biology Society. Buenos Aires，Argentina，2010：6288-6291.

［7］ Miao S，Wang Z J，Liao R. A convolutional neural network approach for 2D/3D medical image registration[J]. Optoelectronics and Advanced Materials-Rapid Communications，2015，4(3)：636-648.

［8］ Pickering M R，Muhit A，Scarvell J M，et al. A new multi-modal similarity measure for fast gradient-based 2D-3D image registration[C]//Annual International Conference of the IEEE Engineering in Medicine and Biology-Society. Minneapolis，USA，2009：5821-5824.

［9］ Akter M，Lambert A J，Pickering M R，et al. Robust initialisation for single-plane 3D CT to 2D fluoroscopy image registration[J]. Computer Methods in Biomechanics and Biomedical Engineering：Imaging & Visualization，2015，3(3)：147-171.

［10］ Livyatan H，Yaniv Z，Joskowicz L. Gradient-based 2D/3D rigid registration of fluoroscopic X-ray to CT[J]. IEEE Transactions on Medical Imaging，2003，22(11)：1395-1406.

［11］ Van der Bom M J，Bartels L W，Gounis M J，et al. Robust initialization of 2D-3D image registration using the projection-slice theorem and phase correlation[J]. Medical Physics，2010，37(4)：1884-1892.

［12］ Andreas V S，Tom C，Graeme P. Increasing the automation of a 2D-3D registration system[J]. IEEE Transactions on Medical Imaging，2013，32(2)：387-399.

［13］ Akter M，Lambert A J，Pickering M R，et al. Improved robustness of 3D CT to 2D fluoroscopy image registration using log polar transforms[C]//International Conference on Electrical. Rajshahi，Bangladesh，2016：1-5.

［14］ Miao S，Wang Z J，Liao R. A CNN regression approach for real-time 2D/3D registration[J]. IEEE Transactions on Medical Imaging，2016，35(5)：1352-1363.

［15］ Piat S，Miao S，Liao R，et al. Dilated Fully Convolutional Network for Multi-Agent 2D/3D Medical Image Registration：US，2019050999A1[P]. 2019-02-14.

［16］ Litjens G，Kooi T，Bejnordi B E，et al. A survey on deep learning in medical image analysis[J]. Medical Image Analysis，2017，7(5)：60-88.

［17］ 王?雷. 影像導(dǎo)航手術(shù)中2D/3D圖像配準(zhǔn)[D]. 長(zhǎng)春：中國(guó)科學(xué)院，2015.

Wang Lei. 2D/3D Image Registration in Image-Guided Surgery[D]. Changchun：Chinese Academy of Sciences，2015(in Chinese).

［18］ 孫?濤，閆?巍，孫振輝，等. 外固定支架模型重建軟件開(kāi)發(fā)與安裝參數(shù)識(shí)別研究[J]. 天津大學(xué)學(xué)報(bào)(自然科學(xué)與工程技術(shù)版)，2020，53(6)：593-600.

Sun Tao，Yan Wei，Sun Zhenhui，et al. Software development and pose acquisition of the Taylor spatial frame[J]. Journal of Tianjin University(Science and Technology)，2020，53(6)：593-600(in Chinese).

［19］ Akter M，Lambert A J，Pickering M R，et al. A 2D-3D image registration algorithm using log-polar transforms for knee kinematic analysis[C]//International Conference of the IEEE Engineering in Medicine and Biology Society. Fremantle，Australia，2012：1-5.

［20］ Wu J，Kim M，Peters J，et al. Evaluation of similarity measures for use in the intensity-based rigid 2D-3D registration for patient positioning in radiotherapy[J]. Medical Physics，2009，36(12)：5391-5403.

［21］ Zhu Z，Li G. An automatic 2D-3D image matching method for reproducing spatial knee joint positions using single or dual fluoroscopic images[J]. Computer Methods in Biomechanics and Biomedical Engineering，2012，15(11)：1245-1256.

［22］ Zhu Z，Massimini D F，Wang G，et al. The accuracy and repeatability of an automatic 2D-3D fluoroscopic image-model registration technique for determining shoulder joint kinematics[J]. Medical Engineering & Physics，2012，34(9)：1303-1309.

［23］ Ding H Y，Bian Z F. A sub-pixel registration approach based on powell algorithm[J]. Acta Photonica Sinica，2009，38(12)：3322-3327.

Initial Pose Estimation Method in 2D/3D Registration

Sun Tao1，Guo Ke1，Liu Chuanba1，Zhang Tao2，Song Yimin1, 3，Ma Xinlong2

(1. Key Laboratory of Mechanism Theory and Equipment Design of Ministry of Education，Tianjin University，Tianjin 300354，China；2. Tianjin Hospital，Tianjin 300211，China；3. Department of Mechanical Engineering，Tianjin Renai College，Tianjin 301636，China)

In this paper，we propose a new initial pose estimation approach based on parameter mapping and image segmentation，focusing on the problems of small initial parameter transformation range and single registration object in the 2D/3D medical image registration process. This approach expands the range of initial motion deviation of 3D object spatial parameters and solves the issues regarding few estimable parameters and multibone image registration difficulty of the surgical navigation system during the registration process. First，the region growth algorithm is used to segment the multibone image and extract the target of interest to complete the target-data acquisition to be registered. Second，a projection imaging model is constructed based on the principle of orthogonal fusion projection to decompose the spatial rigid body transformation parameters into the positive lateral plane and establish the mapping relationship between spatial and plane parameters. Then，the orthogonal biplane template generated using the projection is matched with the corresponding target image，and the registration parameters of the obtained positive and lateral plane are converted into spatial parameters，to achieve the effective estimation of the initial pose of the 3D object. The method was verified based on computed tomography data of the skull，intact femur，and fractured femur. Then，the method was compared with the traditional 2D/3D image registration method. The experimental results showed that the proposed method could make the initial transformation range of the five spatial parameters reach ?±30mm or ±20°，and different research objects within this range show a good initial pose estimation effect. Compared with the traditional image registration method，the overall registration time of the 2D/3D iterative optimization based on the above method is 46.5% shorter on average，the registration accuracy of a single translation parameter is up to 69.2%，and the registration accuracy of a single rotation parameter is up to 91.1%. Simultaneously，the final registration results have adequate robustness.

surgical navigation；2D/3D registration；initial pose estimation；orthogonal biplane projection

the National Key Research and Development Program of China(No.2018YFB1307800)，the National Natural Science Foundation of China(No.51775367)，Tianjin Science and Technology Plan Project(No.18PTLCSY00080，No.20201193，No.18YFSDZC00010)．

10.11784/tdxbz202010007

TP391

A

0493-2137(2022)02-0143-08

2020-10-07；

2021-02-25.

孫?濤（1983—??），男，博士，教授，stao@tju.edu.cn.Email：m_bigm@tju.edu.cn

宋軼民，ymsong@tju.edu.cn.

國(guó)家重點(diǎn)研發(fā)計(jì)劃資助項(xiàng)目(2018YFB1307800)；國(guó)家自然科學(xué)基金資助項(xiàng)目(51775367)；天津市科技計(jì)劃資助項(xiàng)目(18PTLCSY00080，20201193，18YFSDZC00010)．

(責(zé)任編輯：王曉燕)