基于gpu的快速三維ct迭代重建系統(tǒng)的制作方法
【專利摘要】本發(fā)明公開基于GPU的快速三維CT迭代重建系統(tǒng)，包括數(shù)據(jù)輸入模塊，預(yù)處理模塊，正/反投影模塊，變量更新模塊，迭代終止模塊，結(jié)果輸出模塊；數(shù)據(jù)輸入模塊主要包括輸入投影數(shù)據(jù)；預(yù)處理模塊包括對投影數(shù)據(jù)作相關(guān)的預(yù)處理，并將處理后的數(shù)據(jù)以及與重建相關(guān)的參數(shù)傳入到GPU；正/反投影模塊進行包括正投影步驟、記錄信息步驟和修正步驟、反投影步驟操作，分別計算得到正投影系統(tǒng)矩陣和反投影系統(tǒng)矩；變量更新模塊根據(jù)反投影結(jié)果對當(dāng)前迭代過程中的變量值進行更新；迭代終止模塊包括判斷當(dāng)前迭代是否滿足迭代的終止條件；結(jié)果輸出模塊將迭代結(jié)果輸出。本發(fā)明只需計算一次射線與體素的相交情況，減少了計算系統(tǒng)矩陣所需的計算量，加快迭代重建的速度。
【專利說明】
基于GPU的快速三維CT迭代重建系統(tǒng)
技術(shù)領(lǐng)域
[0001] 本發(fā)明屬于X-射線斷層成像(CT)技術(shù)領(lǐng)域，特別涉及基于GPU的快速三維CT迭代 重建系統(tǒng)。
【背景技術(shù)】
[0002] 計算機斷層成像技術(shù)(Computed Tomography，CT)已經(jīng)廣泛應(yīng)用于人體組織成像、 工業(yè)無損檢測等領(lǐng)域。CT重建算法可以分為迭代重建算法和解析重建算法兩大類。解析重 建算法的優(yōu)點是算法簡單、重建速度快，缺點是對數(shù)據(jù)的完備性要求較高，并且投影數(shù)據(jù)中 的噪聲對重建圖像的質(zhì)量影響非常大。相比解析重建算法，迭代重建算法的優(yōu)點是重建圖 像質(zhì)量高，且適用于各種形式的采集數(shù)據(jù)，即使是在有限角度投影數(shù)據(jù)的情況下，也能重建 出較好的圖像，典型的迭代算法有SART、I SRA、WLS和EM-ML等。但是由于迭代重建算法計算 量大、收斂速度慢、重建時間長，目前已成為制約其廣泛應(yīng)用的最大瓶頸。
[0003] 迭代重建算法的迭代過程可分解為四個步驟:正投影、修正、反投影和體素更新。 在迭代算法中，與正投影、反投影采用同一驅(qū)動方式相比，正投影使用射線驅(qū)動方式、反投 影使用體素驅(qū)動方式的組合能有效去除環(huán)狀偽影。因此，在迭代算法中，正投影采用基于射 線的驅(qū)動方式計算系統(tǒng)矩陣；而反投影采用基于體素驅(qū)動方法計算系統(tǒng)矩陣。兩種驅(qū)動方 式的區(qū)別在于前者只考慮射線穿過體素的密度，后者還需考慮周圍的體素對當(dāng)前射線穿過 的體素的影響。在代碼實現(xiàn)過程中，通常將正投影和修正步驟合并，一起進行處理；反投影 和體素更新步驟合并，同時處理。
[0004]在迭代重建算法其中，正投影和反投影的計算量非常大，占迭代重建算法總的計 算量的比例很高，因此提高迭代重建算法的運行速度，關(guān)鍵在于提高正、反投影的計算效 率。其中，求解射線與體素的相交情況占用了正投影和反投影的絕大部分的計算量。對于傳 統(tǒng)的方法，正投影和反投影是分別進行的正投影需要確定射線與每個體素的相交情況;在 進行反投影的過程，該操作還得重新進行一次。
[0005] 目前，GPU(圖形處理器，英語:Graphics Processing Unit，縮寫:GPU)技術(shù)的出現(xiàn) 減少了傳統(tǒng)迭代算法的運行時間，因為射線驅(qū)動方式的正投影和使用體素驅(qū)動的反投影是 支持GPU并行加速的。但重建效率依然較低，難以滿足重建效率要求高的應(yīng)用。因此，提出基 于GPU的快速三維CT迭代重建方法和系統(tǒng)，在理論和實際應(yīng)用中都具有較高的價值。

【發(fā)明內(nèi)容】

[0006] 本發(fā)明提出了基于GPU的快速三維CT迭代重建系統(tǒng)，主要是針對傳統(tǒng)三維CT迭代 重建所需計算時間長的問題，提高重建效率。
[0007] 本發(fā)明通過以下技術(shù)方案實現(xiàn)：
[0008] 基于GPU的快速三維CT迭代重建系統(tǒng)，包括數(shù)據(jù)輸入模塊，預(yù)處理模塊，正/反投影 模塊，變量更新模塊，迭代終止模塊，結(jié)果輸出模塊；
[0009] 所述的數(shù)據(jù)輸入模塊包括輸入投影數(shù)據(jù)、模體位置和尺寸、旋轉(zhuǎn)中心的位置、射線 源到旋轉(zhuǎn)中心的距離、探測器尺寸、探測器到旋轉(zhuǎn)中心的距離和投影角度；
[0010]所述的預(yù)處理模塊包括對投影數(shù)據(jù)作相關(guān)的預(yù)處理，并將處理后的數(shù)據(jù)以及與重 建相關(guān)的參數(shù)傳入到GPU;
[0011] 所述的正/反投影模塊進行包括正投影步驟、記錄信息步驟和修正步驟、反投影步 驟操作，該模塊在GPU上實現(xiàn)正投影步驟和反投影步驟，分別計算得到正投影系統(tǒng)矩陣和反 投影系統(tǒng)矩陣，具體包括:首先，選出一部分待計算正投影的射線，通過多線程分別對這些 選中的射線進行正投影，記錄每條射線正投影的結(jié)果，以及射線與體素相交的信息;然后， 在選中射線的正投影步驟和信息記錄步驟完成之后，根據(jù)記錄的信息進行修正步驟操作和 反投影步驟操作，其中，修正過程只與投影數(shù)據(jù)和每條射線的正投影有關(guān)，根據(jù)之前正投影 步驟記錄的信息可以得到當(dāng)前迭代過程中每條射線的投影值，從而可以快速地計算每條射 線的修正后投影值，在反投影過程中，根據(jù)之前正投影步驟記錄的信息可以快速地得到每 射線穿過體素的位置和長度，從而可以根據(jù)每條射線與體素的相交情況，快速地將每條射 線的修正后投影值分配到相交位置處相鄰的8個體素上;最后，繼續(xù)從剩余的射線中選取一 部分射線重復(fù)之前的操作，直到所有的射線都被計算過而終止循環(huán)；
[0012] 所述的變量更新模塊根據(jù)反投影結(jié)果對當(dāng)前迭代過程中的變量值進行更新；
[0013] 所述的迭代終止模塊包括判斷當(dāng)前迭代是否滿足迭代的終止條件，如果滿足迭代 的終止條件，終止迭代，輸出迭代結(jié)果到結(jié)果輸出模塊;否則，繼續(xù)之前的正/反投影操作和 變量更新操作；
[0014] 所述的結(jié)果輸出模塊包括存儲和顯示迭代重建后的數(shù)據(jù)，將迭代結(jié)果輸出。
[0015] 進一步，所述的正/反投影模塊反投影步驟中反投影操作通過多線程操作進行加 速，每個線程都需要分配一塊數(shù)據(jù)來存儲各體素的值。
[0016] 本發(fā)明與現(xiàn)有技術(shù)相比，優(yōu)點如下：
[0017] 1、正/反投影過程中，只需計算一次射線與體素的相交情況，減少了計算系統(tǒng)矩陣 所需的計算量，加快迭代重建的速度。
[0018] 2、正/反投影模塊將正投影步驟和反投影步驟合并，將減少確定射線與體素相交 情況的次數(shù)，加快迭代重建的速度。
【附圖說明】：
[0019] 下面結(jié)合附圖對本發(fā)明的【具體實施方式】作進一步詳細(xì)說明。
[0020] 圖1是本發(fā)明基于GPU的快速三維CT迭代重建方法的具體實施流程圖；
[0021 ]圖2是本發(fā)明正/反投影的具體實施示意圖；
[0022]圖3是本發(fā)明正投影和記錄信息的具體實施示意圖；
[0023]圖4是本發(fā)明修正和反投影的具體實施示意圖；
[0024]圖5是本發(fā)明基于GPU的快速三維CT迭代重建系統(tǒng)框圖。
【具體實施方式】
[0025]為使本發(fā)明實施的目的、特點和優(yōu)點能夠更為明顯易懂，下面結(jié)合附圖對本發(fā)明 的具體實施例作詳細(xì)的說明。
[0026]本發(fā)明提供了基于GPU的快速三維CT迭代重建方法和系統(tǒng)，如圖5所示，基于GPU的 快速三維CT迭代重建系統(tǒng)100包括數(shù)據(jù)輸入模塊I，預(yù)處理模塊2,正/反投影模塊3,變量更 新模塊4,迭代終止模塊5,結(jié)果輸出模塊6;如圖1~圖5所示，為該系統(tǒng)的具體實施流程圖： [0027]步驟SlOl為數(shù)據(jù)輸入，需要通過CT設(shè)備采集投影數(shù)據(jù)并輸入到數(shù)據(jù)輸入模塊1，以 及獲取該投影數(shù)據(jù)時的模體位置和尺寸、旋轉(zhuǎn)中心的位置、射線源到旋轉(zhuǎn)中心的距離、探測 器尺寸、探測器到旋轉(zhuǎn)中心的距離和投影角度；
[0028]步驟S102為預(yù)處理，預(yù)處理模塊2對投影數(shù)據(jù)作相關(guān)的預(yù)處理，例如，為了降低投 影數(shù)據(jù)中的噪聲，可以對投影數(shù)據(jù)進行非線性濾波、基于貝葉斯統(tǒng)計理論的濾波;為了克服 重建后圖像中的條形偽影、環(huán)狀偽影、金屬偽影和杯形偽影，可以采用基于字典學(xué)習(xí)、形態(tài) 分量分析等新興的圖像處理方法對投影數(shù)據(jù)進行處理，并將處理后的數(shù)據(jù)以及與重建相關(guān) 的參數(shù)傳入到GPU，所述處理后的數(shù)據(jù)以及與重建相關(guān)的參數(shù)具體指模體位置和尺寸、旋轉(zhuǎn) 中心的位置、射線源到旋轉(zhuǎn)中心的距離、探測器尺寸、探測器到旋轉(zhuǎn)中心的距離和投影角 度，由于這些處理后的數(shù)據(jù)和參數(shù)是在CPU上，這些數(shù)據(jù)還需傳入到GPU上才能使用GPU進行 加速；
[0029] 步驟S103為正/反投影，正/反投影模塊3進行包括正投影步驟S301、記錄信息步驟 S302和修正步驟S401、反投影步驟S402操作（如圖2所示），該模塊在GPU上實現(xiàn)正投影步驟 S301和反投影步驟S301，分別計算得到正投影系統(tǒng)矩陣和反投影系統(tǒng)矩陣，如圖2、圖3、圖4 所示，該正/反投影模塊3內(nèi)進行了正投影步驟S301和反投影步驟S401，分別得到正投影系 統(tǒng)矩陣和反投影系統(tǒng)矩陣，具體包括：
[0030]首先，選出一部分待計算正投影的射線，通過多線程分別對這些選中的射線進行 正投影(如圖3的S301)，記錄每條射線正投影的結(jié)果，以及射線與體素相交的信息(如圖3的 S302)，
[0031] 然后，在選中射線的正投影步驟S301和信息記錄步驟S302完成之后，根據(jù)記錄的 信息進行修正步驟S401和反投影步驟S402操作，
[0032]其中，修正過程只與投影數(shù)據(jù)和每條射線的正投影有關(guān)，根據(jù)之前正投影步驟記 錄的信息可以得到當(dāng)前迭代過程中每條射線的投影值，從而可以快速地計算每條射線的修 正后投影值，在反投影過程中，根據(jù)之前正投影步驟記錄的信息可以快速地得到每射線穿 過體素的位置和長度，從而可以根據(jù)每條射線與體素的相交情況，快速地將每條射線的修 正后投影值分配到相交位置處相鄰的8個體素上，這時的反投影操作避免了重新確定射線 和體素的相交情況，減少了算法的計算量，反投影操作可以通過多線程操作進行加速，但是 同時進行基于體素的操作容易產(chǎn)生"寫沖突"，每個線程都需要分配一塊數(shù)據(jù)來存儲各體素 的值，因此反投影操作不能啟用過多的線程，
[0033]最后，繼續(xù)從剩余的射線中選取一部分射線重復(fù)之前的操作，直到所有的射線都 被計算過而終止循環(huán)；
[0034]所述的正/反投影模塊3中正/反投影過程的修正步驟S401由所使用的具體迭代重 建算法決定，為了明顯易懂，以EM-ML算法為例進行詳細(xì)說明，以下是EM-ML算法的迭代公 式：
[0035]
[0036] 其中，模體內(nèi)有N個體素，在第k次迭代中，每個體素位置處的劑量值為利用M條 射線照射模體，投影數(shù)據(jù)為y，在公式中
b正投影過程
1修正因子， 對的值進行修正，將修正后的值進行反投曼
i后利用反投影值更新當(dāng)前 的xf值，在正投影過程中，使用基于射線的驅(qū)動方式計算系統(tǒng)矩陣a^，而反投影過程采用基 于體素的驅(qū)動方式計算系統(tǒng)矩陣b^，兩種驅(qū)動方式的差別在于:射線j穿過體素i的線長度 為1，對于基于射線的驅(qū)動方式，值為1和體素i的密度值的乘積;對于基于體素的驅(qū)動方 式，bij值為1和1中間位置處對應(yīng)密度值的乘積，通過對中間位置處周圍8個體素的密度值進 行插值得到該密度值，在正/反投影過程中使用非匹配的系統(tǒng)矩陣是為了減少重建圖像中 的環(huán)狀偽影，先計算正投〗
記錄每條射線的正投影值和的值及其位置;然后， 計算每條射線的修正因子，在前面兩個過程中，可以用多線同時計算多條射線的正投影值 及修正因子(如圖3所示）；最后，利用插值算法，將修正因子與的乘積賦給體素i及周圍其 它7個體素，在反投影過程中同樣可以使用多線程進行加速；
[0037]步驟S104為變量更新，變量更新模塊4按具體迭代算法所對應(yīng)的迭代公式更新當(dāng) 前的#值；
[0038]步驟S105為迭代終止，迭代終止模塊5判斷當(dāng)前迭代是否滿足迭代的終止條件，如 果滿足迭代的終止條件，終止迭代，輸出迭代結(jié)果到結(jié)果輸出模塊，否則，繼續(xù)之前的正/反 投影操作和變量更新操作；
[0039] 步驟S106為結(jié)果輸出，結(jié)果輸出模塊6將迭代結(jié)果輸出，包括存儲和顯示迭代重建 后的數(shù)據(jù)。
[0040] 本發(fā)明書中未作詳細(xì)描述的內(nèi)容屬于本領(lǐng)域?qū)I(yè)技術(shù)人員公知的現(xiàn)有技術(shù)。
[0041]以上所述僅是本發(fā)明的優(yōu)選實施方式，應(yīng)當(dāng)指出，對于本技術(shù)領(lǐng)域的普通技術(shù)人 員來說，在不脫離本發(fā)明原理的前提下，還可以做出若干的改進和潤飾，這些改進和潤飾也 應(yīng)視為本發(fā)明的保護范圍。
【主權(quán)項】
1. 基于GPU的快速三維CT迭代重建系統(tǒng)，其特征在于:包括數(shù)據(jù)輸入模塊，預(yù)處理模塊， 正/反投影模塊，變量更新模塊，迭代終止模塊，結(jié)果輸出模塊； 所述的數(shù)據(jù)輸入模塊包括輸入投影數(shù)據(jù)、模體位置和尺寸、旋轉(zhuǎn)中心的位置、射線源到 旋轉(zhuǎn)中心的距離、探測器尺寸、探測器到旋轉(zhuǎn)中心的距離和投影角度； 所述的預(yù)處理模塊包括對投影數(shù)據(jù)作相關(guān)的預(yù)處理，并將處理后的數(shù)據(jù)以及與重建相 關(guān)的參數(shù)傳入到GPU; 所述的正/反投影模塊進行包括正投影步驟、記錄信息步驟和修正步驟、反投影步驟操 作，該模塊在GHJ上實現(xiàn)正投影步驟和反投影步驟，分別計算得到正投影系統(tǒng)矩陣和反投影 系統(tǒng)矩陣，具體包括:首先，選出一部分待計算正投影的射線，通過多線程分別對這些選中 的射線進行正投影，記錄每條射線正投影的結(jié)果，以及射線與體素相交的信息;然后，在選 中射線的正投影步驟和信息記錄步驟完成之后，根據(jù)記錄的信息進行修正步驟操作和反投 影步驟操作，其中，修正過程只與投影數(shù)據(jù)和每條射線的正投影有關(guān)，根據(jù)之前正投影步驟 記錄的信息可以得到當(dāng)前迭代過程中每條射線的投影值，從而可以快速地計算每條射線的 修正后投影值，在反投影過程中，根據(jù)之前正投影步驟記錄的信息可以快速地得到每射線 穿過體素的位置和長度，從而可以根據(jù)每條射線與體素的相交情況，快速地將每條射線的 修正后投影值分配到相交位置處相鄰的8個體素上;最后，繼續(xù)從剩余的射線中選取一部分 射線重復(fù)之前的操作，直到所有的射線都被計算過而終止循環(huán)； 所述的變量更新模塊根據(jù)反投影結(jié)果對當(dāng)前迭代過程中的變量值進行更新； 所述的迭代終止模塊包括判斷當(dāng)前迭代是否滿足迭代的終止條件，如果滿足迭代的終 止條件，終止迭代，輸出迭代結(jié)果到結(jié)果輸出模塊;否則，繼續(xù)之前的正/反投影操作和變量 更新操作； 所述的結(jié)果輸出模塊包括存儲和顯示迭代重建后的數(shù)據(jù)，將迭代結(jié)果輸出。2. 如權(quán)利要求1所述的基于GPU的快速三維CT迭代重建系統(tǒng)，其特征在于:所述的正/反 投影模塊反投影步驟中反投影操作通過多線程操作進行加速，每個線程都需要分配一塊數(shù) 據(jù)來存儲各體素的值。
【文檔編號】G06T15/10GK105844690SQ201610216968
【公開日】2016年8月10日
【申請日】2016年3月29日
【發(fā)明人】桂志國, 張鵬程, 韓建寧, 尚禹, 楊民, 舒華忠
【申請人】中北大學(xué)