《數(shù)字圖像處理與圖像通信》朱秀昌

1第7章圖像重建

7.1計(jì)算機(jī)斷層掃描技術(shù)

7.2投影定理

7.3傅立葉投影重建

7.4卷積逆投影重建

7.5代數(shù)重建

7.6三維圖像重建的體繪制

7.7三維圖像重建的面繪制2圖像重建：由一系列沿直線投影圖來(lái)重建二維圖像，由一系列二維圖像重建三維物體。成像方式：透射斷層成像發(fā)射斷層成像反射斷層成像射線種類：X射線成像、核磁共振成像、正電子發(fā)射成像、

超聲成像、微波成像、激光共焦成像、……3射線投影成像的基本原理：人體組織對(duì)X射線吸收和散射，造成衰減，人體內(nèi)的不同結(jié)構(gòu)，比如脂肪、胰、骨骼對(duì)X射線吸收能力有所不同。入射線圖7.1組織對(duì)射線的吸收散射線散射線4投射斷層成像：射線穿過(guò)物體，在檢測(cè)器上得到的遭受衰減的值＝＝射線的投影，根據(jù)投影可以了解物體對(duì)射線的吸收程度。發(fā)射斷層成像：發(fā)射源在物體內(nèi)部，將具有放射性的離子（放射元素）注入物體內(nèi)部，在物體外部檢測(cè)其經(jīng)過(guò)物體吸收之后放射量。反射斷層成像：將入射信號(hào)（通常是單色平面波）入射到物體上，通過(guò)檢測(cè)經(jīng)物體散射（反射）后的信號(hào)強(qiáng)度來(lái)重建。5透射投影成像，圖7.2表示等強(qiáng)度的射線透過(guò)不同密度分布時(shí)的情況，每塊上的數(shù)字表示每塊的密度或衰減，總的衰減是疊加的，一條射線束通過(guò)均勻密度物質(zhì)的厚塊，另一射線通過(guò)不等密度的厚塊組合，但檢測(cè)器的記錄相同，因此，投影重建時(shí)需要一系列投影才能重建二維圖像。入射線6222入射線6141入射線少透射高密度體多透射入射線低密度體圖7.2等強(qiáng)度射線穿透不同組織的情況6

發(fā)射投影成像如，正電子發(fā)射成像（PET：PositronEmissionTomography）采用在衰減時(shí)放出正電子的放射性離子，放出的正電子很快與負(fù)電子相撞湮滅而產(chǎn)生一對(duì)相背運(yùn)動(dòng)的光子。相對(duì)放置的兩個(gè)檢測(cè)器接收到這兩個(gè)光子就可以確定一條射線，檢測(cè)器圍繞物體呈環(huán)形分布，相對(duì)的兩個(gè)檢測(cè)器構(gòu)成一組檢測(cè)器對(duì)，檢測(cè)由一對(duì)正負(fù)電子產(chǎn)生的光子。正電子負(fù)電子光子光子PET成像系統(tǒng)示意圖檢測(cè)器檢測(cè)器77.1

計(jì)算機(jī)斷層掃描技術(shù)計(jì)算機(jī)斷層掃描技術(shù)又稱為計(jì)算機(jī)層析或CT（ComputedTomography)利用數(shù)字圖像處理技術(shù)來(lái)獲取三維圖像。CT機(jī)通常包括X射線管、X射線檢測(cè)器、掃描機(jī)架、病人床、用來(lái)重建圖像結(jié)構(gòu)的工作站。圖7.4CT掃描成像的示意圖8醫(yī)學(xué)影像領(lǐng)域：ComputedTomography(CT)：獲1979年諾貝爾獎(jiǎng)（NobelPrice）布爾赫、珀塞爾，獲1952年諾貝爾獎(jiǎng)，發(fā)現(xiàn)了核磁共振現(xiàn)象勞特布爾（美）、P·曼斯菲爾德（英）獲2003年年諾貝爾獎(jiǎng)，核磁共振的研究（英）G.N.Hounsfield

（美）AllanM.Cormack9CT實(shí)例：掃描系統(tǒng)的X射線源和檢測(cè)器，始終保持嚴(yán)格的相對(duì)靜止；射線管發(fā)出的是直線形波束，掃描機(jī)構(gòu)圍繞人體作旋轉(zhuǎn)加平移運(yùn)動(dòng)；以檢測(cè)器的位置為自變量，就構(gòu)成如圖7.5(b)的電流—位置函數(shù)曲線。圖7.5CT一次平移掃描所獲得的輸出信號(hào)10第一次直線平移掃描完畢后，掃描系統(tǒng)旋轉(zhuǎn)一個(gè)小角度，再作第二次直線式平移掃描，獲得另一組投影數(shù)據(jù)；重復(fù)以上過(guò)程，便得到很多組投影數(shù)據(jù)；對(duì)這些數(shù)據(jù)進(jìn)行處理形成三維圖像。圖7.6頭顱CT掃描成像示意圖117.2

投影定理一個(gè)N維函數(shù)在第N

-1維上的映射稱為函數(shù)

f在第N

-1維的投影。二維：函數(shù)f(x,y)在x軸上（沿y方向）的投影函數(shù)f(x,y)在y軸上（沿x方向）的投影設(shè)f(x,y)

的傅立葉變換為F(u,v)，可得：（7.1）（7.2）（7.3）12把式（7.3）代入到式（7.1）可得：可知gy(x,y)是F(u,0)的傅氏反變換，或gy(x,)

的傅氏變換G(u)與F(u,0)相同。結(jié)論，函數(shù)f(x,y)在x軸上投影gy(x,y)的傅立葉變換等于f(x,y)的傅立葉變換F(u,v)在(u,v)平面上沿u軸平面上的切片F(xiàn)(u,0)。（7.4）13沿y軸的投影圖示沿y軸的的投影示意圖f(x,y)(a)二維函數(shù)f(x,y)在x軸上投影yxgy(x)

(b)f(x,y)傅立葉變換F(u,v)在u軸上切片F(xiàn)(u,v)vuF(u,0)14假設(shè)函數(shù)f(x,y)投影到一條經(jīng)過(guò)旋轉(zhuǎn)的直線上t1，t是一條與t1平行經(jīng)過(guò)原點(diǎn)的直線，與t垂直經(jīng)過(guò)原點(diǎn)的直線為s，該直線s與x軸的夾角為θ，直線t1離開(kāi)原點(diǎn)的距離為s1，如圖7.9所示。以s和t可用θ為極坐標(biāo)：函數(shù)f(x,y)沿著t1方向s投影為：圖7.9坐標(biāo)旋轉(zhuǎn)關(guān)系(7.6)15將投影式(7.6)只對(duì)s1作一維傅立葉變換，將指數(shù)項(xiàng)作變換，得uv(u,v)θ0r（7.8）（7.9）16投影定理（切片定理）：

f(x,y)在一條與x軸夾角為θ，離開(kāi)原點(diǎn)距離為r的直線上的投影的傅立葉變換＝＝二維傅立葉變換在與u軸成θ方向上的切片三維圖像重建基礎(chǔ)：若投影變換G(r,θ)中對(duì)所有的r和θ值都已知，則圖像的二維傅立葉變換也可以完全確定，進(jìn)行二維傅立葉反變換，就可以得到f(x,y)。圖7.10投影定理示意圖f(x,y)yxθ

vuF(u,v)F(r,θ)tsθ

17三維投影定理：令f(x,y,z)表示一個(gè)三維物體，它的三維傅立葉變換為如ω=0

其中g(shù)z(x,y)正是f(x,y,z)在(x,y)平面上的投影，即表明f(x,y,z)在(x,y)平面上投影的傅立葉變換＝＝

f(x,y,z)的三維傅立葉變換F(u,v,ω)在ω=0平面上的切面F(u,v,0)。與(x,y)平面成夾角為θ的平面上投影的傅立葉變換＝＝三維傅立葉變換F(u,v,ω)在與(u,v)平面成θ角的切面F(u,v,θ)。（7.13）（7.12）（7.11）187.3

傅立葉投影重建傅立葉投影重建的基礎(chǔ)：傅立葉投影定理。根據(jù)投影定理，就可以得到F(u,v)分別在相應(yīng)角度位置上的切片；當(dāng)切片趨向無(wú)窮多，就可獲得在(u,v)平面上的所有F(u,v)值；由F(u,v)進(jìn)行傅立葉反變換就可以重建圖像f(x,y)。將f(x,y)沿s方向的投影表達(dá)式及其一維傅立葉變換式改寫(xiě)為：令u=Rcosθ,v=Rsinθ，根據(jù)投影定理，用極坐標(biāo)(R,θ)來(lái)表示:結(jié)論：如果知道所有R和θ的投影變換值G(R,θ)，則變換域的二維函數(shù)將全部確定，取傅立葉反變換就可以得到圖像函數(shù)。（7.13）（7.14）（7.15）19利用傅立葉變換的共軛對(duì)稱性，積分限由0～2π換成0～π，R→|R|，積分限由0～∞換成－∞～＋∞，上式：

記傅立葉投影重建圖像為：（7.18）（7.17）（7.16）20以上是理想的情況－－可以獲得無(wú)窮多個(gè)投影，對(duì)連續(xù)圖像的傅立葉重建。如果只有有限個(gè)角度的投影g(ρ,θn)，θn表示nθ，G(R,θ)可用在一系列采樣點(diǎn)(m△s，θn)

上對(duì)g(?)求和得到，△s為沿著射線方向采樣點(diǎn)的間距，采樣點(diǎn)數(shù)為M，式(7.14)可寫(xiě)成：

令R=k△R，k為整數(shù)，△R為頻率域上采樣間距，采樣點(diǎn)數(shù)為M，取則有

根據(jù)極坐標(biāo)上點(diǎn)(k△R,n△θ)的值G(k△R,θn)插值出在直角坐標(biāo)上點(diǎn)(k△u,n△v)的值F(k△u,n△v)，從而反傅立葉變換得到f(k△x,n△y)。（7.20）（7.19）21注意：必須得到所有投影數(shù)據(jù)后再能重建圖像，不能根據(jù)所獲得的部分投影數(shù)據(jù)重建圖像，重建圖像需要進(jìn)行傅立葉反變換。傅立葉變換法步驟：(1)對(duì)N個(gè)不同θ方向上投影進(jìn)行一維傅立葉變換。

(2)在傅立葉變換空間從極坐標(biāo)向直角坐標(biāo)插值。

(3)利用式(7.17)或離散形式的傅立葉頻譜進(jìn)行反變換得到重建圖像。圖7.11傅立葉空間的直角和極坐標(biāo)網(wǎng)格227.4

卷積逆投影重建卷積逆投影重建法：以投影切片定理為基礎(chǔ)；傅立葉變換重建法：計(jì)算量比較小，但要二維插值，在射電天文學(xué)研究中得到應(yīng)用廣泛；卷積逆投影法：能快速實(shí)現(xiàn)，噪聲較小時(shí)可重建出準(zhǔn)確清晰的圖像，在X射線CT成像中應(yīng)用廣泛。23在計(jì)算投影的一維傅立葉變換F(R,θ)時(shí)，R為頻域極軸變量。投影數(shù)據(jù)g(ρ,θ)總是被有限截?cái)唷．?dāng)ρ的取樣間隔為d時(shí)，在頻率R的變化范圍將是-d/2～d/2，投影切片定理可近似成:

記因?yàn)?，上式又可?xiě)成（7.21）（7.22）（7.23）24由(7.17)可知：

卷積逆投影法重建圖像為：由式（7.24）可知，右邊正是投影數(shù)據(jù)g(ρ,θ)與脈沖響應(yīng)h(ρ)所表示的濾波器的卷積，h(·)為卷積函數(shù)。求f’(x,y,θ)

則是在θ角方向上卷積了的投影，因此從式（7.22）求f(x,y)可被認(rèn)為是求逆投影過(guò)程，即卷積逆投影重建法。（7.24）（7.25）25式（7.22）所表示的h(ρ)正是頻率響應(yīng)為|R|的濾波器，又稱為重建濾波器。卷積逆投影重建的關(guān)鍵－－設(shè)計(jì)重建濾波器。重建濾波器的設(shè)計(jì)對(duì)式(7.22)求積分：在離散情況下，重建濾波器被有限截?cái)啵O(shè)在ρ的變化范圍內(nèi)，對(duì)其取M點(diǎn)，采樣間隔為d，即ρ=md，。26用極坐標(biāo)表示，用離散值代替連續(xù)積分，則式(7.23)可以寫(xiě)成：用極坐標(biāo)和直角坐標(biāo)之間的關(guān)系，求出與極坐標(biāo)點(diǎn)相對(duì)應(yīng)的直角坐標(biāo)系上的點(diǎn)，從而得到用離散值表示的，最后，根據(jù)求逆投影式(6.22)得到重建的圖像：(7.27)和(7.28)是一組便于計(jì)算機(jī)快速運(yùn)算的表達(dá)式。（7.27）（7.28）27卷積逆投影重建法的重建濾波器。1/2dvu-1/2d|R|H(R)圖7.12重建濾波器的頻率響應(yīng)287.5

代數(shù)重建投影重建的傅立葉變換法和濾波器逆投影都在變換域（頻率域）內(nèi)處理；在連續(xù)域內(nèi)進(jìn)行解析處理，為便于計(jì)算機(jī)實(shí)現(xiàn)，引入離散化和有限近似。另一類方法－－代數(shù)重建，逐次逼近的迭代算法，也稱級(jí)數(shù)展開(kāi)法。代數(shù)重建法：需要在重建的目標(biāo)上加一柵格，將目標(biāo)劃分為許多等大小的體積單元；計(jì)算每個(gè)體元的衰減系數(shù)，寫(xiě)成矩陣形式如下：

Y=AX（7.30）注意：在許多位置，射線束只是部分地通過(guò)一些體元。29代數(shù)法重建：切片＋柵格＝小體積元（小立方體）圖7.13掃描重建柵格30迭代方法來(lái)求解式(7.30)：根據(jù)物體已有的先驗(yàn)知識(shí)，先對(duì)未知圖像的各像素都賦予一個(gè)初始估值，利用這些假設(shè)數(shù)據(jù)去計(jì)算各射線穿過(guò)對(duì)象時(shí)可能得到的投影值,將算得的值和實(shí)測(cè)投影值進(jìn)行比較，按照差異獲得一個(gè)修正值，根據(jù)這些修正值，修正各對(duì)應(yīng)射線穿過(guò)的諸像素值，如此反復(fù)迭代，直到計(jì)算值和實(shí)測(cè)值接近所要求的精度為止。31具體實(shí)施步驟如下：(1)對(duì)于未知圖像各像素均給予一個(gè)假定的初始值，從而得到一幅初始計(jì)算圖像，如設(shè)各像素的初始值均為0。(2)根據(jù)假設(shè)圖像，求對(duì)應(yīng)第i條射線穿過(guò)時(shí)，各體單元應(yīng)得到的各個(gè)相應(yīng)投影值Z1*，Z2*，…，Zn*。(3)將計(jì)算值Z1*，Z2*，…，Zn*和對(duì)應(yīng)的實(shí)測(cè)值Z1，Z2，…，Zn進(jìn)行比較，然后取對(duì)應(yīng)差值ΔZi=Zi－Zi*作為修正值。(4)用每條射線的修正值來(lái)修正和該射線相交的諸像素值。(5)用修正后的像素值重復(fù)(1)～(4)步，直到計(jì)算值和實(shí)測(cè)值之差，即修正值小到所允許的精度為止。32四單元矩陣實(shí)例：第一次運(yùn)算：水平照射后，將射線和放入圖像單元，第二次運(yùn)算：垂直照射后，再將射線和加到上面那個(gè)圖像單元中去，021322442437591533第三次運(yùn)算：取向右上方的對(duì)角線照射，獲得的射線和加到前圖像單元中，第四次運(yùn)算：取向左上方的對(duì)角線照射，并將此射線和加到前面所獲得取的圖像單元中去，最后一次運(yùn)算：從每個(gè)圖像單元中減去最小數(shù)6，然后用3去除每個(gè)單元值，便可解出該問(wèn)題開(kāi)始的那個(gè)矩陣。310812033231612915347.6

三維圖像重建的體繪制前面：傅立葉投影重建、卷積逆投影重建、代數(shù)法重建－－進(jìn)行降維處理，重建的基礎(chǔ)是投影切片定理，三維物體的重建是從大量的二維投影圖像得到，二維圖像的重建是從大量的一維投影圖像得到，這類重建方法在CT、核磁共振中應(yīng)用廣泛，重建圖像質(zhì)量較高。本節(jié)：三維圖像重建－－“直接體繪制法”－－“三維圖像顯示技術(shù)”，“三維空間數(shù)據(jù)可視化”實(shí)際應(yīng)用中，要求在獲得大量二維切片圖像之后，希望能觀察到三維圖像，即三維圖像顯示，獲得更加直觀、逼真、形象的視覺(jué)效果。關(guān)鍵：已知物體的各層切片圖像，通過(guò)三維重建在屏幕上顯示物體。35對(duì)于分布在三維空間的大量數(shù)據(jù)來(lái)說(shuō)，比較簡(jiǎn)單的一種是這些數(shù)據(jù)均勻地分布在三維網(wǎng)格點(diǎn)上，即在x，y，z三個(gè)方向上，網(wǎng)格點(diǎn)之間的距離均相等，如圖所示。圖7.15規(guī)則數(shù)據(jù)場(chǎng)示意圖及掃描得到的頭部切片圖36在三維空間的規(guī)則數(shù)據(jù)場(chǎng)中，無(wú)需給出各數(shù)據(jù)點(diǎn)空間位置，只要給出三維網(wǎng)格某一角點(diǎn)的空間位置和某一數(shù)據(jù)點(diǎn)的序號(hào)（第k層的i，j值），即可根據(jù)網(wǎng)格間距所對(duì)應(yīng)的距離求出該點(diǎn)的空間位置。每個(gè)網(wǎng)格是結(jié)構(gòu)化數(shù)據(jù)的一個(gè)元素，通常叫做“體元”(Voxel)，數(shù)據(jù)場(chǎng)的函數(shù)值f(i,j,k)

對(duì)應(yīng)于三維空間(i,j,k)

的位置。圖7.15(b)是掃描得到的一系列頭部切片圖像的其中三幅切片。直接體繪制：將已采集到計(jì)算機(jī)上的三維離散數(shù)據(jù)場(chǎng)重新采樣（先重構(gòu)再采樣），按照一定的規(guī)則轉(zhuǎn)換為圖像顯示緩存中的二維離散信號(hào)。注：將三維離散（圖像）數(shù)據(jù)場(chǎng)集合稱為物體空間，將顯示圖像的屏幕稱為圖像空間，如圖7.16所示。37體繪制方法中最常見(jiàn)的有兩種：（1）圖像空間掃描的體繪制法：如圖，從圖像空間到物體空間的方法。圖7.16重構(gòu)和重采樣圖7.17光線投射體繪制的重采樣38（2）物體空間掃描的體繪制法：如圖，從物體空間到圖像空間的方法。圖7.18空間采樣點(diǎn)對(duì)屏幕的貢獻(xiàn)示意圖39光線投射體繪制方法：是一個(gè)三維離散數(shù)據(jù)場(chǎng)的重采樣和圖像合成的過(guò)程(1)對(duì)物體空間的三維離散數(shù)據(jù)場(chǎng)進(jìn)行預(yù)處理；(2)從圖像空間——顯示屏幕上的每一個(gè)像素點(diǎn)根據(jù)設(shè)定的觀察方向發(fā)出一條投射光線；(3)為了增強(qiáng)三維逼真效果，突出顯示不同組織的邊界面，可以采樣表面明暗計(jì)算；(4)計(jì)算每條射線對(duì)屏幕像素點(diǎn)的貢獻(xiàn)。用該方法對(duì)40幅頭部切片重建結(jié)果如圖所示。圖7.19頭部切片數(shù)據(jù)光線投射重建后顯示407.7

三維圖像重建的面繪制面繪制技術(shù)：由已知物體的三維網(wǎng)格離散數(shù)據(jù)，在三維空間數(shù)據(jù)場(chǎng)中構(gòu)造出中間幾何圖元，然后再由傳統(tǒng)的計(jì)算機(jī)圖形學(xué)技術(shù)實(shí)現(xiàn)面繪制。面繪制方法構(gòu)造出的可視化圖形當(dāng)中間幾何圖元較小時(shí)，可以得到光滑的表面和清晰的圖像，在構(gòu)造出中間圖元之后，可以利用現(xiàn)有的圖形硬件實(shí)現(xiàn)繪制功能，速度比體繪制要快，因而得到了廣泛的應(yīng)用。面繪制的典型算法：立方體步進(jìn)法（MC法：MarchingCubes)，四面體步進(jìn)法（MT法：MarchingTetrahedra)，剖分立方體法(DividingCubes)等。41對(duì)于三維數(shù)據(jù)場(chǎng)，確定的是物體表面輪廓信息，－－類似于二維圖像的邊緣輪廓線信息。（1）確定包含表面輪廓面的體元：離散的三維空間數(shù)據(jù)場(chǎng)中的一個(gè)體元如圖7.15(a)所示，8個(gè)數(shù)據(jù)點(diǎn)位于該體元的8個(gè)角點(diǎn)位置上。一個(gè)體元的8個(gè)角點(diǎn)位置用A，B，…，F(xiàn)表示，其8個(gè)角點(diǎn)對(duì)應(yīng)的狀態(tài)為1，0，0，0，0，0，0，0，則等值面必定與AB、AD、AE相交。每個(gè)體元有8個(gè)角點(diǎn)，每個(gè)角點(diǎn)有0，1兩種狀態(tài)，每個(gè)體元共有個(gè)不同的狀態(tài)。圖7.22角點(diǎn)的狀態(tài)及其體元與等值面的相交42我們可以利用不同種對(duì)稱性將256