技術(shù)領(lǐng)域

本發(fā)明涉及正電子發(fā)射斷層成像領(lǐng)域，尤其是涉及一種正電子發(fā)射斷層成像系統(tǒng)探測器晶體中心位置圖生成方法及晶體像素查找表生成方法。

背景技術(shù)：

正電子湮滅事例經(jīng)PET探測器像素陣列解碼之后產(chǎn)生事例解碼之后產(chǎn)生事例的位置分布，統(tǒng)計之后形成光子事件二維直方圖，由于探測器實現(xiàn)中編碼與解碼的非線性，生成的光子事件二維直方圖會形成不規(guī)則形變。以一個探測器block為例，晶體陣列為16*16的等距規(guī)則排列像素陣列，解碼之后會發(fā)生桶形、蝶形、旋轉(zhuǎn)、壓縮、擴(kuò)張等各種形式的形變。若不做任何處理，數(shù)據(jù)獲取時，在線處理無法判斷晶體像素的行列坐標(biāo)，直接導(dǎo)致圖像分辨率下降，嚴(yán)重的可以導(dǎo)致圖像信息錯誤。如果能確定之后的行列坐標(biāo)和解碼前的硬件晶體之間的對應(yīng)關(guān)系，就能判斷晶體像素的行列坐標(biāo)。目前業(yè)界廣泛使用的做法是使用晶體像素查找表(Crystal Lookup Table，CLT)來描述解碼前后的晶體對應(yīng)關(guān)系。

由于PET探測器的晶體漂移隨時間推移會不停惡化，使CLT每隔一段時間就要被重新繪制。現(xiàn)階段CLT靠手工繪制，但是一臺PET設(shè)備往往有幾萬個晶體，手工繪制過程耗時耗力，成為PET設(shè)備維護(hù)中的主要時間瓶頸。為了既保證矯正效果，又縮短維護(hù)時間，目前最常用的方法就是半自動繪制CLT，即計算機(jī)自動生成一個中間結(jié)果：往往不是分割好的CLT，而是晶體中心位置圖CCPM(Crystal Central Position Map)。自動生成的CCPM質(zhì)量越好，技術(shù)檢查和修正結(jié)果花費(fèi)的時間越少。

現(xiàn)有技術(shù)中，已經(jīng)公開的自動定位各晶體中心位置的算法，大致可以分為基于一維投影和直接在二維圖像域處理兩大類?；谝痪S投影的方法先把光子事件二維直方圖投影到兩個方向，最常見的是將圖像分別投影在X方向和Y方向分別投影，則二維檢測問題被簡化為一維尋峰問題。求解可以使用一維的極值檢測，一維求導(dǎo)定位法，甚至更加復(fù)雜的Gaussian fitting。但是投影算法過度依賴光子事件二維直方圖中晶體行列的規(guī)整性，當(dāng)晶體行列出現(xiàn)扭曲時，一維尋峰會失敗。二維圖像域處理包括特征檢測域聚類算法兩大類。特征檢測的實際性能往往取決于圖像質(zhì)量和預(yù)處理，容易受圖像噪聲、灰度、不均勻性等影響出現(xiàn)漏檢和錯檢。聚類算法缺少空間約束，更容易受偽影的影響，即使獨(dú)立的偽影也能使算法結(jié)果出錯，極大地影響探測結(jié)果，并導(dǎo)致后續(xù)的晶體映射錯位，嚴(yán)重的會導(dǎo)致整個晶體陣列的映射完全變形。綜上所述，現(xiàn)有眾多的晶體中心位置檢測算法，穩(wěn)定性差，出現(xiàn)誤檢和漏檢的可能性非常大。因此需要提出一種穩(wěn)定性較好的晶體中心位置檢測方法。

技術(shù)實現(xiàn)要素：

為了解決現(xiàn)有技術(shù)中晶體中心位置檢測方法錯檢率高，性能不穩(wěn)定的問題，本發(fā)明提供了一種晶體中心位置圖生成方法。

一種晶體中心位置圖生成方法，包括如下步驟：

獲取光子事件二維直方圖I，將所述光子事件二維直方圖I歸一化至[0，1]之間，獲得圖像I₀；對所述圖像I₀進(jìn)行迭代處理得到圖像I_i；使用自動閾值對所述圖像I_i進(jìn)行二值化處理，獲得二值圖像集S_B；提取所述二值圖像集S_B中的每幅二值圖像的連通域及其對應(yīng)的中心位置；判斷所述連通域彼此是否有重疊；

保留有重疊的連通域中迭代次數(shù)最大的連通域；生成晶體中心位置圖，所述晶體中心位置圖包括有重疊的連通域中迭代次數(shù)最大的連通域?qū)?yīng)的中心位置以及沒有重疊的連通域?qū)?yīng)的中心位置。

優(yōu)選地，所述判斷所述連通域彼此是否有重疊的步驟具體是：

將一個連通域內(nèi)的任意位置與另一個連通域內(nèi)的任意位置進(jìn)行比較，如果相同即表示這兩個連通域彼此有重疊。

優(yōu)選地，所述有重疊的連通域中僅保留迭代次數(shù)最大的連通域的步驟具體為：比較所述彼此有重疊的連通域的對應(yīng)圖像I_i的迭代次數(shù)i，剔除對應(yīng)迭代次數(shù)較小的連通域。

優(yōu)選地，所述有重疊的連通域中僅保留迭代次數(shù)最大的連通域的步驟具體為：比較所述彼此有重疊的連通域的面積，剔除面積較大的連通域。

優(yōu)選地，所述步驟還包括設(shè)置面積閾值，剔除連通域面積大于所述面積閾值的連通域。

優(yōu)選地，所述將光子事件二維直方圖I歸一化至[0，1]之間采用線性收放的方式進(jìn)行，即max(I)是所述光子事件二維直方圖中的最大值，min(I)是所述光子事件二維直方圖中的最小值。

優(yōu)選地，所述光子事件二維直方圖最大值max(I)按照比例ε獲得，即取所述光子事件二維直方圖的灰度直方圖面積的p處的灰度為最大值max(I)，p＝1-ε，圖中超出最大值的max(I)的值強(qiáng)制轉(zhuǎn)換成max(I)。

優(yōu)選地，所述比例其中n是所述光子事件二維直方圖中實際包含的晶體數(shù)量，N_p是光子事件二維直方圖的像素數(shù)量。

優(yōu)選地，所述自動閾值的選取基于光子事件二維直方圖的灰度直方圖，所述閾值將所述光子事件二維直方圖的灰度直方圖分成兩部分，所述兩部分的質(zhì)心到所述自動閾值的距離相等。

本發(fā)明還提供了一種晶體像素查找表生成方法，采用上述生成的晶體中心位置圖，根據(jù)所述晶體中心位置圖劃分各相鄰晶體條的分界線，生成晶體像素查找表。

與現(xiàn)有技術(shù)相比，本發(fā)明提供的晶體中心位置圖生成方法，采用連續(xù)非極大值衰退的方法進(jìn)行晶體中心位置定位，本方法本質(zhì)上是對圖像中的非波峰部分進(jìn)行連續(xù)衰減，從而突出波峰部分，對于一維或者二維尋峰算法而言，極大地提高了尋峰的準(zhǔn)確度。同時，本發(fā)明還提供了晶體像素查找表生成方法。

附圖說明

圖1為本發(fā)明晶體中心位置圖生成方法的流程示意圖；

圖2為本發(fā)明的連續(xù)非極大值衰退的示意圖；

圖3為本發(fā)明生成的晶體像素查找表。

具體實施方式

在下面的描述中闡述了很多具體細(xì)節(jié)以便于充分理解本發(fā)明。但是本發(fā)明能夠以很多不同于在此描述的其它方式來實施，本領(lǐng)域技術(shù)人員可以在不違背本發(fā)明內(nèi)涵的情況下做類似推廣，因此本發(fā)明不受下面公開的具體實施的限制。其次，本發(fā)明利用示意圖進(jìn)行詳細(xì)描述，在詳述本發(fā)明實施例時，為便于說明，所述示意圖只是實例，其在此不應(yīng)限制本發(fā)明保護(hù)的范圍。

如背景技術(shù)所述，現(xiàn)有技術(shù)中的晶體中心位置檢測方法，不管是一維投影法，還是直接在二維圖像域處理法，其穩(wěn)定性和精確性都有待提高。因此本發(fā)明提出一種穩(wěn)定性較好的晶體中心位置檢測方法。

如圖1所示，一種晶體中心位置圖生成方法，包括如下步驟：

步驟S10，獲取光子事件二維直方圖I，將所述光子事件二維直方圖I歸一化至0到1之間，獲得圖像I₀。

對入射光子事件位置進(jìn)行二維分布統(tǒng)計，獲取光子事件二維直方圖I，將所述光子事件二維直方圖進(jìn)行歸一化處理至[0，1]之間。在一個實施例中，所述將光子事件二維直方圖歸一化至[0，1]之間采用線性收放的方式進(jìn)行，即max(I)是所述光子事件二維直方圖中的最大值，min(I)是所述光子事件二維直方圖中的最小值。為了避免極端值對歸一化的影響，對圖像面積的ε處做截斷處理，所述比例ε取決于晶體陣列中晶體的數(shù)量n和光子事件二維直方圖的像素數(shù)N_p。在本發(fā)明中，采用的是16*16的晶體陣列，則每個陣列有256個晶體，而光子事件二維直方圖的分辨率(即像素數(shù)N_p)是256*256。因此，ε＝256/(256*256)＝0.4％，那么圖像直方圖面積的99.6％處的值為max(I)，圖像中大于max(I)的值都強(qiáng)制轉(zhuǎn)成max(I)，即全圖最大值就是max(I)。

在一個實施例中，在對圖像進(jìn)行處理歸一化處理之前，對所述光子事件位置直方圖進(jìn)行均衡化矯正。本技術(shù)方案采用CLAHE算法對圖像進(jìn)行均衡化矯正，但是其他矯正的方法均可使用，例如用原圖減去或除以圖像中緩慢變化的背景內(nèi)容。圖像中緩慢變化的背景內(nèi)容可以通過在圖像域中使用大尺度平滑，在頻域中壓制高頻部分，或使用曲面擬合等方法獲得。

步驟S20，對所述圖像I₀進(jìn)行迭代處理得到圖像I_i。

所述迭代處理即為連續(xù)非極大值衰退。所述圖像I₀經(jīng)過多次迭代處理，圖像中的非波峰部分進(jìn)了連續(xù)的衰減，突顯出波峰部分，突顯的波峰部分更易于被探測。

迭代之后的圖像I_i＝I_i-1×I₀，其中i∈{1…n}，i為迭代的次數(shù)，n為最大迭代次數(shù)，即終止條件。在一個實施例中，對圖像I₀進(jìn)行50次迭代處理，如圖2所示，分別為迭代1次、5次、10次、15次、20次、25次、30次、35次、40次的中間圖像。從圖中可以看出，經(jīng)歷迭代的次數(shù)越多，非波峰部分信號越弱，波峰部分越明顯。

步驟S30，使用自動閾值對所述圖像I_i進(jìn)行二值化處理，獲得二值圖像B_i。

自動閾值算法有許多種，比如著名的otsu分割，最大熵分割，isodata等等。在一個實施例中采用isodata自動閾值法，但是其他自動閾值方法均可使用。具體實現(xiàn)方式如下：

先提取圖像I_i中的非零最小值為閾值t的，開始迭代：

a)計算圖像I_i中灰度大于t的所有像素的平均值m_sup；

b)計算圖像I_i中灰度小于t的所有像素的平均值m_inf；

c)檢驗t是否等于m_sup和m_inf的平均值，如是，則結(jié)束，如否，則增加t的值(t＝t+1)，從回步驟a)。

迭代結(jié)束后t的值，正好是圖像I_i的灰度直方圖中大于t部分的質(zhì)心和小于t部分的質(zhì)心的平均值，即所述兩部分的質(zhì)心到所述自動閾值的距離相等。此時的t即為最后的閾值。以t對圖像I_i進(jìn)行二值化，大于t的位置賦為1，小于t的位置賦為0，得到圖像B_i。

步驟S40，對所述二值圖像集S_B中的每幅二值圖像B_i提取連通域及其對應(yīng)的中心位置。

在一個實施例中，通過如下方式提取連通域及其對應(yīng)的中心位置：建立一個空的堆棧T，遍歷B_i上的像素位置(x，y)，當(dāng)其值為1時，進(jìn)行如下步驟：

1)將當(dāng)前位置設(shè)為(x，y)；2)將(x，y)壓入堆棧T；3)遍歷B_i上(x，y)的四鄰域((x+1，y)，(x-1，y)，(x，y+1)，(x，y-1))，如果其中的某個位置值為1，針對此位置，重復(fù)步驟1)-3)；4)將T中所有位置取出，并將B_i上這些位置對應(yīng)的值改為0。取出的所有位置即為一個連通域，位置的個數(shù)即為此連通域面積，所有位置的平均值既為此連通域的中心。

步驟S50，判斷所述連通域彼此是否有重疊。

將一個連通域內(nèi)的任意位置與另一個連通域內(nèi)的任意位置進(jìn)行比較，如果相同即表示這兩個連通域彼此有重疊。

在一個實施例中，如果一個連通域的中心位置與另一個連通域的任意位置相同，即表示這兩個連通域彼此有重疊。

在一個實施例中，將新提取的連通域內(nèi)的所有位置與已確定保留的所有連通域的中心位置進(jìn)行一一比較，如果相同，即表示這新提取的連通域與已確定保留的某一個連通域有重疊。

步驟S60，保留有重疊的連通域中迭代次數(shù)最大的連通域。

由于最終的晶體中心位置圖是由連通域集組成，重疊的連通域?qū)嵸|(zhì)上是對同一個晶體重復(fù)記錄，因此需要有選擇性地過濾有重疊的連通域。連通域來源于二值圖像B_i。來源于相同二值圖像的連通域彼此之間是沒有重疊的，但不同二值圖像所產(chǎn)生的連通域會出現(xiàn)重疊。故相同迭代次數(shù)的連通域彼此沒有重疊，彼此有重疊的連通域迭代次數(shù)必不相同?？梢灾苯油ㄟ^比較重疊連通域的迭代次數(shù)，保留迭代次數(shù)多的連通域來過濾連通域中的重疊。也可以通過比較連通域面積，保留面積較小的連通域達(dá)到同樣的目的，因為迭代次數(shù)多的連通域的面積必然小于迭代次數(shù)少的連通域的面積。

實際應(yīng)用中，為了減少重復(fù)運(yùn)算，可以將步驟S40，S50和S60同時進(jìn)行(即在提取連通域時就決定此連通域是保留還是丟棄)。同時可以通過設(shè)置連通域面積閾值c進(jìn)一步加快運(yùn)行速度。具體的，需要如下步驟：

1)建立確定保留的連通域集合S_c以及其對應(yīng)的中心位置的集合S_p，初始使S_c和S_p為空集；

2)二值圖像B_i＝n得到的連通域彼此之間不會有重疊。且B_i＝n是所有二值圖像中迭代次數(shù)最大的，因此如果與其他二值圖像得到的連通域有重疊，都應(yīng)該保留二值圖像B_i＝n得到的連通域。因此將所有二值圖像B_i＝n得到的連通域加入S_c中，其中心位置加入S_p；

3)按從i＝n-1，i＝n-2，…直到i＝1的順序從二值圖像B_i中提取連通域。提取每個連通域時，如果其面積小于面積閾值c，且此連通域不包括任何S_p已有的中心位置，則將此連通域加入S_c中，其中心位置加入S_p，否則丟棄。得到的S_c即為篩選后的連通域集，S_p即為晶體中心位置集。

步驟S70，生成晶體中心位置圖，所述晶體中心位置圖包括有重疊的連通域中迭代次數(shù)最大的連通域?qū)?yīng)的中心位置以及沒有重疊的連通域?qū)?yīng)的中心位置。在一個實施例中，將晶體中心位置集S_p中所有位置在二維圖像上標(biāo)出，即生成晶體中心位置圖。

本發(fā)明還提供了晶體像素查找表生成方法，采用上述生成的晶體中心位置圖，根據(jù)所述晶體中心位置圖劃分各相鄰晶體條的分界線，生成晶體像素查找表。

根據(jù)晶體中心位置圖劃分各相鄰晶體條的分界線，利用得到的分界線完成對晶體條區(qū)域的劃分填充，將填充數(shù)據(jù)回朔到圖像存儲中，完成晶體像素查找表的生成。

對于實際的晶體越是到中間能量越大，即像素值越大；而越到周邊，其能量越小。因此，分界線所在位置像素值是全局最小，以及分界線本身光滑連續(xù)的這兩個條件，這兩個條件互相制約，形成了傳統(tǒng)的內(nèi)能和外能。本發(fā)明使用更精準(zhǔn)的約束條件：由于晶體實際上是以行列形式排列的，因此一條連續(xù)的分界線應(yīng)該是相鄰兩行或兩列晶體之間的分界線。由于在建立網(wǎng)格模板時使用樣條把晶體中心點(diǎn)的行列都串了起來，這等于為分界線確立了范圍。而且，每條分界線從模型上來說就是一條最優(yōu)化路徑。

本實施例使用動態(tài)規(guī)劃確定各晶體條的分界線，對于一條相鄰兩行之間的分界線可以分為w段，其中w為晶體中心圖像的寬度。動態(tài)規(guī)劃中所需要的能量函數(shù)為E_n(x,x′)＝a*d(x,x′)+b*I(x)。其中，E_n(x,x′)為相鄰兩段之間的能量遞增,x和x'分別是相鄰兩段的位置，d為取距離差函數(shù)。由于我們要求分界線盡可能平滑，因此相鄰兩段之間的距離差被限制在2個像素以下。I(x)為0-1均衡化后的圖像像素值，a和b都是能量函數(shù)的權(quán)重，其數(shù)值大小根據(jù)實際需要可以調(diào)節(jié)，在此具體實施例中取a＝0.1，b＝1.5，由于圖像中心信息比較準(zhǔn)確，我們采用從中間向左右兩邊的順序?qū)ふ易疃搪窂健?/p>

具體地包括：(1)在正確的晶體中心位置圖上以兩條串起相鄰兩行晶體中心點(diǎn)的樣條為邊界線并將上述邊界線分為w段，令邊界線所圍的圖像中心位置為階段0，記階段0以左為階段1至階段e，記階段0以右為階段e+1至階段w-1，記E_i(x)為階段i的能量，1≤i≤w-1，令E₀(x)＝1.5I(x)。(2)獲取階段1至階段e所有階段累積能量最小值，并記錄其回朔位置x₁’，其中每一階段的累積能量E_i(x)＝E_i-1(x′)+E_n(x,x′)。同樣階段e+1至階段w-1所有階段累積能量最小值，并記錄其回朔位置x₂’。(3)根據(jù)起始位置進(jìn)行左右路徑連接的所有組合，計算總累積能量。(4)總累積能量最小的路徑即為分界線。通過上述操作即可完成晶體分割，生成如圖3所示的最終的晶體像素查找表。

本發(fā)明涉及的晶體中心位置圖生成方法，采用連續(xù)非極大值衰退的方法進(jìn)行晶體中心位置定位，本方法本質(zhì)上是對圖像中的非波峰部分進(jìn)行連續(xù)衰減，從而突出波峰部分，對于一維或者二維尋峰算法而言，極大地提高了尋峰的準(zhǔn)確度。同時，本發(fā)明還提供了晶體像素查找表生成方法。

雖然本發(fā)明披露如上，但本發(fā)明并非限定于此。任何本領(lǐng)域技術(shù)人員，在不脫離本發(fā)明的精神和范圍內(nèi)，均可作各種更動與修改，因此本發(fā)明的保護(hù)范圍應(yīng)當(dāng)以權(quán)利要求所限定的范圍為準(zhǔn)。