工業(yè) CT 系統(tǒng)幾何參數(shù)校正

陳志祥 須穎 石錦洋

摘要：工業(yè)計算機斷層成像技術(shù)（computer tomography，CT）系統(tǒng)幾何參數(shù)不準確會影響重建圖像的質(zhì)量和尺寸測量的精度。為校正工業(yè) CT 幾何位姿失配參數(shù)，提出基于標準球板模體的幾何校正方法。首先，利用球板在兩個不同位置的投影圖求解出射線源到探測器的距離以及射線源到運動平臺的距離。再在360°范圍內(nèi)間隔一定角度均勻采集球板投影信息，擬合球心投影軌跡，解析求解探測器偏角η、φ和探測器中心點坐標。最后根據(jù)重建出的體素模型的球心距測量結(jié)果，利用牛頓迭代法，以工業(yè) CT 尺寸測量最大允許誤差為約束條件迭代求解探測器偏角θ。實驗結(jié)果顯示：該文幾何校正方法可以求解出所有工業(yè) CT 幾何位姿失配參數(shù)。參數(shù)校正后，提升圖像質(zhì)量，消除重建圖像偽影，且工業(yè) CT 球心距測量最大誤差由44.8μm 降低到16.3μm，表明該幾何校正方法可以滿足實際工業(yè) CT 系統(tǒng)校正的要求。

關(guān)鍵詞：工業(yè) CT;幾何校正;校正模體;圖像重建

中圖分類號： TP391.4文獻標志碼： A文章編號：1674–5124（2021）12–0023–06

Geometric parameter calibration of industrial CT system

CHEN Zhixiang，XU Ying，SHI Jinyang

（School of Mechanical and Electrical Engineering， Guangdong University of Technology， Guangzhou 510006， China）

Abstract： Incorrect geometric parameters of the industrial CT （computer tomography） system will affect the qualityof thereconstructedimageandtheaccuracyof dimensionalmeasurement. Inordertocalibrate geometric mismatching parameters of industrial CT， a geometric calibration method based on the standard spherical plate phantom was proposed. Firstly， the distance from the ray source to the detector and the distance from the ray source to the rotating platform were calculated by using the projection of the spherical plate at two different positions. Then， the projection information of the spherical plate was uniformly collected at a certain angle interval within 360° range， fitted the sphere center projection trajectory， and the detector deflection angle η ， φ and the detector center point coordinates were calculated analytically. Finally， based on the spherical- center distance measurement results of the voxel model， the Newton iteration method was used to solved the detectordeflectionangle θ bytakingtheMaximumPermissibleErrorof theindustrialCTdimensional measurementastheconstraintcondition. Theexperimentalresultsshowedthatthegeometriccorrection method could solved all the geometric position mismatch parameters of industrial CT. After the parameters were calibrated， the image quality was improved， and the reconstructed image artifacts were eliminated. The maximum error of industrial CTspherical distance measurement was reduced from 44.8μm to 16.3μm，indicated that the geometric calibration method could met the requirements of actual industrial CT system calibration.

Keywords： industrial computer tomography; geometry calibration; calibration phantom; image reconstruction

0引言

工業(yè)計算機斷層成像技術(shù)（ computer tomography，CT）是一門集光學(xué)、機械、控制和計算機科學(xué)于一體的數(shù)字化 X 射線成像技術(shù)，其利用物體對 X 射線吸收率的不同來獲取物體的多角度投影圖，并根據(jù)這些投影圖重建出物體的三維圖像[1]。目前這項技術(shù)正在向尺寸測量領(lǐng)域發(fā)展，但因為工業(yè) CT 通常采用 FDK 算法對投影數(shù)據(jù)進行重建，而 FDK（ Feldkamp-Davis-Kress）算法要求工業(yè) CT系統(tǒng)幾何位姿滿足理想條件[2]，實際情況下，工業(yè) CT 系統(tǒng)難于滿足理想條件，系統(tǒng)幾何參數(shù)會存在偏差。在圖像重建過程中，如果不進行校正，會影響工業(yè) CT 尺寸測量的精度。

針對這個問題，Noo[3]利用其設(shè)計的兩點模體，通過獲取模體多角度的投影數(shù)據(jù)，擬合橢圓運動軌跡，根據(jù)橢圓的相關(guān)參數(shù)求解系統(tǒng)幾何參數(shù)，在忽略探測器偏角θ對其他參數(shù)的影響的情況下，精確求解其余參數(shù)，但實際工業(yè) CT 系統(tǒng)中偏角θ是一定存在的，且也會影響圖像的重建質(zhì)量。Jacobson[4]設(shè)計了一種金屬絲校準模型，利用兩條不共面的直線在不同角度的投影平面內(nèi)的交點位置的不同，根據(jù)這些交點的軌跡來計算出系統(tǒng)幾何參數(shù)，但該方法需要準確知道校正模體中每一根金屬絲準確的空間位置，在實際過程中，很難精確定位空間中金屬絲的位置，且對空間中金屬絲線位置測量的難度明顯高于球板球心位置測量的難度。Zhao[5]利用含有九個標準球的正方形模體，采集模體單個角度的投影數(shù)據(jù)，根據(jù)標準球的投影位置關(guān)系迭代求出系統(tǒng)幾何參數(shù)，雖能通過迭代求解出系統(tǒng)所有的失配參數(shù)，但該方法僅根據(jù)球板一個位置的投影關(guān)系求解，其求解精度難以保證。目前，工業(yè) CT 系統(tǒng)參數(shù)校正的方法[6-14]都是利用特定校正模體，采集模體不同角度的投影數(shù)據(jù)，根據(jù)投影數(shù)據(jù)建立起與工業(yè)CT 系統(tǒng)幾何位姿參數(shù)的映射關(guān)系，再求解出相應(yīng)的參數(shù)。

本文采用高精度球板模體，采集球板模體在不同角度，不同位置下的二維投影圖，求出不同情況下球板各球體球心在探測器上的投影坐標，擬合運動軌跡，再根據(jù)投影關(guān)系，利用解析法精確求解除θ外的失配參數(shù)。之后，用帶參數(shù)的 FDK 重建算法軟件重建出的球板三維圖像，將其導(dǎo)入可視化軟件 VG Studio MAX 3.0中并提取邊緣，測量球體球心距，分析工業(yè) CT 測量誤差，以工業(yè) CT 尺寸測量最大允許誤差為約束條件，利用牛頓迭代法求解系統(tǒng)最后一個失配參數(shù)。相較于其他校正方法，本方法可求解出工業(yè) CT 系統(tǒng)所有的失配參數(shù)，且用工業(yè) CT 尺寸測量誤差來對比校正前后的效果，更為直觀。

1工業(yè) CT 系統(tǒng)和幾何參數(shù)

工業(yè) CT 系統(tǒng)主要由射線源、旋轉(zhuǎn)平臺和探測器構(gòu)成，根據(jù) FDK 算法的要求，理想的工業(yè) CT 系統(tǒng)幾何結(jié)構(gòu)應(yīng)滿足：射線源焦點和探測器中心點的連線經(jīng)過旋轉(zhuǎn)軸，且垂直于探測器平面;旋轉(zhuǎn)軸與探測器中心列平行。工業(yè) CT 系統(tǒng)理想位姿關(guān)系如圖1所示。其中f 為射線源到旋轉(zhuǎn)軸的距離;d 為射線源到探測器中心的距離;（u0 ， v0）為射線源焦點在探測器上的垂直投影點。

為更好地研究工業(yè) CT 系統(tǒng)的幾何位姿關(guān)系，建立三維坐標系O? XYZ，如圖1所示。X 軸沿射線方向;Y 軸沿水平方向垂直于射線方向且與探測器中心行平行;Z 軸由右手坐標系原則確定。同樣在探測器所在平面上建立二維坐標系O?uv，u 軸沿探測器中心行方向;v 軸沿探測器中心列方向。工業(yè) CT 系統(tǒng)在安裝調(diào)試過程中，系統(tǒng)位姿關(guān)系一定會存在偏差，工業(yè) CT 的校正即校正這些參數(shù)。經(jīng)國內(nèi)外眾多研究人員對工業(yè) CT 校正的研究，需要校正的參數(shù)可概括為：探測器繞 X、Y、Z 軸的轉(zhuǎn)動角度η、θ、φ;射線源焦點在探測器上的垂直投影點（u0 ， v0）和在 u、v 方向上的偏移Δu、Δv;源到旋轉(zhuǎn)中心的距離以及源到探測器中心的距離。其中和探測器有關(guān)的5個參數(shù)如圖2所示。

2參數(shù)校正

2.1校正模體

本文使用的校正模體是球板標準器，如圖3所示，由12個高精度的紅寶石球體和碳纖維平板組成，部分球體之間的球心距已知。

2.2計算d、f

在保持管電壓、管電流和曝光時間等掃描參數(shù)不變的情況下，僅改變探測器到運動平臺的距離，采集兩組不同位置時的球板二維投影圖。擬合球板各球體在二維投影圖上的圓心位置，再計算圓心距 L，L 為像素數(shù)量與像素大小的乘積，即

式中：L——被測量，mm;

n——像素數(shù)量;

M——放大比;

D——像元尺寸，mm。

根據(jù)式（1），利用上述采集的兩組投影圖數(shù)據(jù)，利用球板各球體間球心距，即可建立多個關(guān)于 d 和 f 的二元一次方程，在精確知道運動平臺沿射線源方向上的運動距離的情況下，可求解出多個 d、f值，最后取 d、f 的平均值作為最終校正結(jié)果。

2.3計算η、u0、v0

采集 N（N 為偶數(shù)）張1號球和12號球每隔一定角度的投影圖，計算出1號球和12號球在不同角度投影圖下的球心位置（ui ， vi）（i =1， 2，3， ··· ， N）。并根據(jù)這些球心坐標擬合運動軌跡，可以得到兩個橢圓（橢圓1和橢圓2），擬合出的橢圓可用橢圓方程式表示如下：

式中：（ ， ）——橢圓中心點坐標;

a、b、c——橢圓方程參數(shù)。

由式（3）分別計算出1號球和12號球球心連線在每隔180°投影圖下的中心點位置（1 ， 1）、（2 ， 2），如圖4所示。

其中j = i+ N/2;i =1， 2， ··· ， N/2。

探測器面內(nèi)扭轉(zhuǎn)角η即可由式（4）計算得出。探測器平面內(nèi)射線源焦點的垂直投影點（u0 ， v0）可由式（5）、（6）計算出。

2.4計算φ

Sun[15]等利用四球模體校正系統(tǒng)幾何參數(shù)，提出平行于轉(zhuǎn)軸的兩球投影圓心距之間滿足以下關(guān)系：

式中：——平行于轉(zhuǎn)軸的兩球投影圓心距，

——實際圓心距，mm。

根據(jù)這一原理，將 1 號球到 10 號球的球心距和 3 號球到 12 球的球心距 SD3-12 代入式 （7），可計算出。

2.5 計算

根據(jù)上述計算出的探測器失配參數(shù)對在理想坐標系下的 1 號球球心坐標進行坐標變換得，則：

其中：

經(jīng)過坐標變換后得到：

射線源焦點坐標S =（0， 0， 0）經(jīng)坐標變換后，射線源焦點在理想坐標系下的坐標應(yīng)為S＼=（xs ，ys ，zs）;將 S＼與p 相連交探測器于點p = v ＼），p ＼即為1號球球心經(jīng) X 射線投射后在探測器上的投影點。其中：

同理可計算出其他球體球心經(jīng) X 射線投射后在探測器上的投影點坐標。

此時使用帶參數(shù)的 FDK 重建算法軟件，利用上面所計算出的6個系統(tǒng)失配參數(shù)，重建出此時球板模體的三維圖像，再利用 VG Studio MAX 3.0圖像可視化軟件測量各球體間的球心距，找出此時球心距誤差最大的兩球體。以1號球和12號球為例，此時兩球球心在探測器上的投影點間距離為：

式（12）中，除θ外的其他參數(shù)均已知，故S D1?12= f （θ）是一個關(guān)于θ的函數(shù)，解析求解過于復(fù)雜，可用牛頓迭代法求解，以工業(yè) CT 尺寸測量最大允許誤差為迭代約束條件，迭代求解出θ。

3實驗驗證

3.1實驗設(shè)備及參數(shù)

為驗證上述方法的可行性，在實際工業(yè) CT 系統(tǒng)中進行了實驗，實驗所使用的設(shè)備為天津三英精密儀器有限公司生產(chǎn)的 nanoVoxel-4000，如圖5所示。射線源為開管發(fā)射式，最大管電壓為220 kV，最大管電流為1.0 mA;探測器成像面積為427 mm ×427 mm，像素矩陣為3072×3072;像素細節(jié)分辨力為500 nm，空間分辨率為2?m。

3.2校正過程

首先驗證 d 和f，在 d 固定的情況下，改變f其他參數(shù)如表1所示，采集兩組投影圖。表1 初始參數(shù)

通過計算可得，射線源到探測器的距離誤差為–5.2029 mm，射線源到旋轉(zhuǎn)軸的距離誤差為–4.3907 mm。

在調(diào)整工業(yè) CT 系統(tǒng) d 和f后，再次掃描校正模體，保持 d=1500 mm，f=500 mm，在其他參數(shù)不變的情況下360°范圍內(nèi)采集球板模體的投影圖，然后根據(jù)本文方法求解圓心，擬合橢圓軌跡如圖6所示，橢圓參數(shù)計算出工業(yè) CT 其余幾何參數(shù)，最后利用迭代法經(jīng)4次迭代后計算出工業(yè) CT 需校正的最后一個參數(shù)。

至此，所有參數(shù)校正完畢，最終校正結(jié)果如表2所示。

3.3校正前后效果對比

實際重建過程中，利用帶參數(shù)的 FDK 重建算法軟件重建球板模體，校正前后的重建圖像對比如圖7所示，可以看出，校正后球板模體的重建圖像質(zhì)量有了明顯提升，消除了圖像偽影，提升了圖像清晰度。

利用 VG Studio MAX 3.0圖像可視化軟件分別對校正前后重建的球板模體三維圖像進行尺寸測量，兩次測量結(jié)果如圖8所示?？芍Ｕ蠊I(yè) CT 尺寸測量的最大誤差由44.8?m 降低到16.3?m，效果明顯，本文工業(yè) CT 系統(tǒng)參數(shù)校正的方法降低了工業(yè) CT 尺寸測量的誤差。

4結(jié)束語

本文基于高精度球板模體，提出了一種工業(yè) CT 幾何參數(shù)校正的方法，利用球板在不同角度下的投影圖像，獲取球心投影運動軌跡，擬合橢圓，再根據(jù)橢圓參數(shù)和工業(yè) CT 尺寸測量的最大誤差，結(jié)合解析法與迭代法，可求解出工業(yè) CT 系統(tǒng)所有幾何失配參數(shù)。實驗表明，校正后重建圖像的質(zhì)量有了明顯的提升，減少了圖像偽影，校正后工業(yè) CT 尺寸測量最大誤差由44.8?m 降低到16.3?m，能滿足實際工業(yè) CT 系統(tǒng)參數(shù)的校正要求。但是，本文幾何校正方法需要高精度的校正模板，需要準確知道球板模體各球心距，也受限于運動平臺回轉(zhuǎn)精度和沿射線源方向的運動精度，否則會影響校正參數(shù)的準確性。

參考文獻

[1]張朝宗.工業(yè) CT 技術(shù)和原理[M].北京：科學(xué)出版社， 2009.[2] FELDKAMP L A， DAVIS L C， KRESS J W. Practical cone-beam algorithm[J]. Journal of the Optical Society of America A， 1984， 1（6）：612-619.

[3] NOO F， CLACKDOYLE R， MENNESSIER C， et al. Analyticmethodbasedonidentificationofellipseparametersfor scanner calibration in cone-beam tomography[J]. Physics in Medicine and Biology， 2000， 45（11）：3489-3508.

[4] JACOBSON M W， KETCHA M D， CAPOSTAGNO S， et al.A line fiducial method for geometric calibration of cone-beam CTsystemswithdiversescantrajectories[J]. Physicsin Medicine and Biology， 2018， 63（2）：025030.

[5] ZHAO J T， HU X D， ZOU J，etal. Geometric parametersestimationandcalibrationincone-beammicro-CT[J]. Sensors， 2015， 15（9）：22811-22825.

[6]周凌宏， 李翰威， 徐圓， 等.錐束 CT 圓軌道掃描的幾何校正[J].光學(xué)精密工程， 2014， 22（10）：2847-2854.

[7] SMEKALLV，KACHELRIEM，STEPINAE，etal.Geometricmisalignmentandcalibrationincone-beam tomography[J]. Medical Physics， 2004， 31（12）：3242.

[8] CHO Y， MOSELEY D J，SIEWERDSEN J H， et al. Accuratetechniqueforcompletegeometriccalibrationof cone-beam computedtomographysystems[J]. MedicalPhysics， 2005，32（4）：968-983.

[9] XIANG Q， WANG J， CAI Y F. A geometric calibrationmethodforconebeamCTsystem[C]//EighthInternational ConferenceonDigitalImageProcessing （ICDIP 2016）. International Society for Optics and Photonics， 2016.

[10] PATEL V， CHITYALA R N， HOFFMANN K R， et al. Self- calibrationofacone-beammicro-CTsystem[J]. MedicalPhysics， 2008， 36（1）：48-58.

[11]張才鑫， 黃魁東， 陳華.基于標定量具的快速錐束 CT 幾何校準方法[J].儀器儀表學(xué)報， 2019， 40（1）：175-182.

[12] ZHOU K， HUANG Y X， MENG X X， et al. A new method for cone-beamcomputedtomographygeometricparameters estimation[J]. JournalofComputerAssistedTomography， 2016， 40（4）：639-648.

[13] LI G， LUO S H， YOU C Y， et al. A novel calibration method incorporatingnon-linearoptimizationandball-bearing markers for cone-beam CT with a parameterized trajectory[J]. Medical Physics， 2019， 46（1）：152-164.

[14] LI H C， HONG W， LIU Y， et al. A novel method of micro- tomographygeometricanglecalibrationwithrandom phantom[J]. Journal of X-ray science and technology， 2017，25（4）：641-652.

[15] SUN Y， HOU Y， ZHAO F Y， et al. A calibration method for misalignedscannercogeometryinne-beamcomputed tomography[J]. NDT and E International， 2006， 39（6）：499-513.

（編輯：莫婕）