基于cuda架構(gòu)的gpu并行加速雙能譜ct重建方法
【技術(shù)領(lǐng)域】
[0001] 本發(fā)明涉及X射線CT技術(shù)領(lǐng)域雙能CT解析重建算法，尤其涉及一種基于CUDA架 構(gòu)的GPU并行加速雙能譜CT解析重建方法。
【背景技術(shù)】
[0002] 對于CT加速領(lǐng)域，采用GPU加速已經(jīng)成為重要的方法，而具體的加速方法是提高 CT產(chǎn)品的處理速度的關(guān)鍵環(huán)節(jié)，既是本領(lǐng)域技術(shù)人員研宄的核心技術(shù)，也直接決定產(chǎn)品的 市場競爭力。
[0003] 與單能CT重建不同，雙能CT重建并不是前者的簡單重復(fù)執(zhí)行，而具有其特殊性。 雙能CT解析重建算法分為后處理重建和預(yù)處理重建兩種算法。其中，后處理重建和單能 CT解析重建方法一致，是對投影數(shù)據(jù)直接重建，獲得被檢物的高低能線衰減系數(shù)信息，重建 的形狀精度較好，但由于不能從根本上消除射術(shù)硬化效應(yīng)，重建的被檢物材料精度難以保 證；而預(yù)處理重建與后處理算法相比，預(yù)處理重建算法要對雙能投影數(shù)據(jù)進(jìn)行分解，預(yù)處理 重建算法分為"基效應(yīng)分解模型"和"基材料分解模型"兩種模型，兩種模型是統(tǒng)一的，對分 解后的基材料（或基效應(yīng)）投影數(shù)據(jù)進(jìn)行重建，再根據(jù)基材料的重建圖像計算等效原子序 數(shù)和電子密度圖像。這種預(yù)處理算法理論上可以得到不受能譜影響的被檢物的等效原子序 數(shù)和電子密度信息，重建的材料精度較高。雙能CT重建想要同時獲得上述四種重建圖像信 息，時間因素是算法能否實際應(yīng)用的關(guān)鍵。
[0004] 如現(xiàn)有技術(shù)中，公開號為：101283913,名稱：CT圖像重建的GPU加速方法，該發(fā)明 專利申請公開了一種基于GPU硬件加速的CT圖像重建方法，該發(fā)明軟件部分通常包括基 于GPU的CT數(shù)據(jù)預(yù)處理模塊、基于GPU的CT數(shù)據(jù)濾波模塊、基于GPU的CT正投影模塊和 基于GPU的CT圖像重建反投影模塊，該發(fā)明利用GPU硬件實現(xiàn)了對CT圖像重建算法的加 速，重建部分和數(shù)據(jù)處理部分均在GPU上完成，通過對大數(shù)據(jù)處理采用分塊處理方法，對每 塊待重建體重建時只需要用到部分投影數(shù)據(jù)，減少了數(shù)據(jù)傳輸，從而可提高整個體重建速 度。
[0005] 又如現(xiàn)有技術(shù)中，公開號為：101596113,名稱：一種CT并行重建系統(tǒng)及成像方法， 該發(fā)明專利申請公開了一種CT并行重建系統(tǒng)及其成像方法，該發(fā)明所述系統(tǒng)包括前端采 樣器、中心節(jié)點、以及與所述中心節(jié)點連接的多個子節(jié)點，每個所述子節(jié)點均安裝有圖形 處理器。該發(fā)明所述成像方法包括如下步驟：1)將重建圖像劃分為多個BLOCK區(qū)域；2)中 心節(jié)點接收前端采樣器所采集的各角度下的原始投影數(shù)據(jù)；3)中心節(jié)點為各節(jié)點分配計 算任務(wù)，每個分配了計算任務(wù)的子節(jié)點計算一個或多個BLOCK區(qū)域的重建值；4)各子節(jié)點 完成所分配的重建計算；5)中心節(jié)點組合出完整的重建圖像。該發(fā)明采取劃分重建圖像子 區(qū)域的方法充分挖掘CT掃描及重建的并行特性，在數(shù)據(jù)采集的同時進(jìn)行GPU的并行重建， 重建時間由現(xiàn)有技術(shù)的分鐘級降低為毫秒級，達(dá)到準(zhǔn)實時的被測對象的斷層圖像重建顯 不〇
[0006] 由此可知，現(xiàn)有技術(shù)基本上都是利用GPU并行實現(xiàn)單一的CT圖像重建過程，以達(dá) 到加速圖像重建的目的，由于CT重建本身計算量大，耗時高，算法速度一直是CT相關(guān)算法 能否在實際中使用的關(guān)鍵影響因素，對于雙能CT解析重建算法中預(yù)處理重建算法的特殊 性，現(xiàn)有技術(shù)中并沒有給出有針對性的完善的并行加速解決方法，基于上述現(xiàn)有技術(shù)中的 不足，這也成為雙能CT解析重建算法發(fā)展過程中亟待解決的問題。

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

[0007] 為了彌補現(xiàn)有技術(shù)的缺陷，本發(fā)明提供一種基于CUDA架構(gòu)的GPU并行加速雙能譜 CT解析重建方法，結(jié)合雙能CT重建的特點，提供一種基于CUDA架構(gòu)利用GPU的并行能力， 進(jìn)行雙能CT快速解析重建同時獲得四幅重建圖像的方法，S卩：后處理重建的高低能線衰減 系數(shù)圖像和預(yù)處理重建的等效原子序數(shù)和電子密度圖像。
[0008] 為了實現(xiàn)上述目的，本發(fā)明所述采取的技術(shù)解決方案是：本發(fā)明所述一種基于 CUDA架構(gòu)的GPU并行加速雙能譜CT解析重建方法，所述方法包括以下基本步驟：
[0009] S1 :高低能采樣數(shù)據(jù)快速配準(zhǔn)；
[0010] S2 :并行計算原始投影數(shù)據(jù)及投影分解和投影加權(quán)；
[0011] S3 :雙能投影或者分解投影同步快速卷積濾波；
[0012] S4 :同時并行反投影重建幾類圖像；
[0013] 進(jìn)一步的，也可以包括：
[0014] S5 :由分解重建圖像計算等效原子序數(shù)圖像和電子密度圖像。
[0015] 根據(jù)不同應(yīng)用，上述步驟也可以調(diào)整為S1對原始投影數(shù)據(jù)進(jìn)行快速配準(zhǔn)；S2雙能 投影分解和投影加權(quán)。
[0016] 進(jìn)一步的，其中所述步驟S1的具體特征為：
[0017] S11 :所述高低能采樣數(shù)據(jù)由雙能譜CT系統(tǒng)采樣獲得，其雙能譜為通過快速切換 射線源產(chǎn)生的真實雙能譜或者通過使用雙層夾心探測器獲得的偽雙能譜；
[0018] S12 :所述雙能譜CT系統(tǒng)的掃描模式為圓周掃描或者螺旋掃描；雙能均采用順時 針掃描或者逆時針掃描；
[0019] S13 :所述雙能譜CT系統(tǒng)的探測器為單排探測器或者多排探測器或者面陣探測 器；
[0020] S14 :所述的CT系統(tǒng)相關(guān)參數(shù)在GPU設(shè)備端定義為常量參數(shù)，并在初始化時賦值， 以便加快GPU內(nèi)核函數(shù)對其訪問的速度；
[0021] S15:所述的配準(zhǔn)是為保證雙能CT預(yù)處理投影分解重建算法的分解精度，先對雙 能采樣數(shù)據(jù)進(jìn)行視角配準(zhǔn)，其特征在于：
[0022] S16:所述的雙能快速配準(zhǔn)使用計算最大互相關(guān)系數(shù)方法，為加速將計算雙能二維 采樣數(shù)據(jù)的互相關(guān)系數(shù)簡化為對采樣數(shù)據(jù)進(jìn)行一維投影，再計算雙能一維投影的互相關(guān)系 數(shù)；
[0023] S17 :采用GPU并行方式對高低能采樣進(jìn)行局部探測器單元數(shù)據(jù)二值化，并對二值 化后的數(shù)據(jù)進(jìn)行行方向的一維投影；
[0024] S18 :通過雙能一維投影的卷積實現(xiàn)計算互相關(guān)系數(shù)加速，使用CUFFT對雙能采樣 的一維投影變換到頻域進(jìn)行，獲得互相關(guān)系數(shù)數(shù)組，尋找數(shù)組中最大互相關(guān)系數(shù)所在位置 即為配準(zhǔn)需要調(diào)整的行數(shù)。
[0025] 進(jìn)一步地，其中所述步驟S2的具體特征在于：
[0026] S21 :計算原始雙能投影為加快讀取數(shù)據(jù)的速度，將采樣數(shù)據(jù)和空采明場數(shù)據(jù)都存 儲在紋理存儲器中；
[0027] S22:基于查找表實現(xiàn)雙能投影分解，初始化時將預(yù)先離線生成的"投影分解查找 表"存儲到顯存的二維紋理存儲器以便提高讀取速度，該查找表紋理采用非歸一化的浮點 型拾取坐標(biāo)，將高低能投影值作為橫縱坐標(biāo)查找對應(yīng)的分解投影值，并利用紋理存儲器硬 件實現(xiàn)的線性模式濾波功能對讀取的浮點型返回值進(jìn)行插值，實現(xiàn)快速高精度雙能投影分 解；
[0028] S23:在同一內(nèi)核函數(shù)中按照像素級并行先計算高低能原始投影數(shù)據(jù)，再利用高低 能原始投影數(shù)據(jù)到查找表紋理中查找分解投影數(shù)據(jù)，并對原始投影和分解投影的同時進(jìn)行 加權(quán)操作。根據(jù)具體應(yīng)用場景不同，加權(quán)值可以離線計算后，初始化時導(dǎo)入顯存紋理中，在 加權(quán)操作時讀取，或者在該內(nèi)核函數(shù)中實現(xiàn)同步計算該像素點的加權(quán)值。在同一內(nèi)核函數(shù) 中同時實現(xiàn)上述計算的優(yōu)點在于節(jié)省存取中間結(jié)果的時間，加快處理速度。
[0029] 進(jìn)一步地，其中所述步驟S3的具體特征在于：
[0030] S31 :濾波過程中使用CUDA架構(gòu)提供的基于GPU的快速傅里葉變換函數(shù)庫CUFFT;
[0031] S32 :初始化時將濾波核用FFT變換到頻域，并存儲到一維紋理存儲器中，以便在 頻域運算時加快讀取速度；
[0032] S33 :同步進(jìn)行雙能和雙分解后投影的FFT變換，利用復(fù)數(shù)的實部和虛部同時進(jìn)行 高低能投影的FFT變換或者分解后的雙投影的FFT變換；
[0033] S34 :利用FFT同時處理一批一維離散傅里葉變換，并將變換后的復(fù)數(shù)值存儲到復(fù) 數(shù)紋理以便計算頻域點積時加快讀取速度；
[0034] S35:使用同一內(nèi)核函數(shù)同步對雙能和雙分解后頻域投影進(jìn)行與濾波核的點積計 算，將濾波后的結(jié)果數(shù)據(jù)綁定到紋理存儲器，以便反投影重建時加快讀取速度。
[0035] 進(jìn)一步地，其中所述步驟S4的具體特征在于：
[0036] S41 :采用像素驅(qū)動的反投影重建，不存在對存儲重建圖像的全局顯存的"寫沖 突"；
[0037] S42 :同一內(nèi)核函數(shù)中按照像素級并行進(jìn)行幾類圖像反投影重建，同一線程完成四 種類型重建圖像的同一像素位置的反投影計算，即同時完成高低能線衰減系數(shù)圖像和兩種 分解圖像中各自的同一像素位置的重建；
[0038] S43 :重建圖像存儲在全局存儲器中，內(nèi)核函數(shù)內(nèi)部使用寄存器變量臨時存儲計算 用的變量，計算過程中不訪問全局存儲器，每個內(nèi)核函數(shù)完成X個均勻分布視角的反投影 計算，計算完X個視角后再給全局存儲器的相應(yīng)像素位置賦值，以便節(jié)省反投影重建的訪 存時間，根據(jù)具體應(yīng)用環(huán)境不同，該X=N/M。N為掃描一圈采集的全部視角；M為用來并行 處理視角的線程數(shù)量，且能夠被N整除，M等于1時，對重建的像素點不存在"寫沖突"，可以 不使用原子加操作進(jìn)行寫入全局存儲器，當(dāng)M大于1時，