人腦的NURBS曲面三維建模研究

張 斌，謝松云，李寧飛

（西北工業大學 電子信息學院，陜西 西安 710129）

EEG與MEG都能夠提供具有良好時間分辨率 （通常為毫秒級）的數據，但空間分辨率有限。與之相對應，功能磁共振（function Magnetic Resonance Imaging， fMRI）提供了良好的空間分辨率卻相對較差的時間分辨率。研究人員們開始從將腦電與功能磁共振技術進行結合的角度尋求靈感和算法[1]，其目的就是結合兩者的優勢，并最終得到更為合理精確的源定位結果。人腦建模是腦磁圖（Magnetoencephalogram，MEG）與腦電圖（Electroencephalogram，EEG）計算中正逆問題的先決條件。同時它也是基于EEG/fMRI的多模態醫學數據融合中約束方法首先要解決的問題。實際應用中廣泛應用的往往是球頭模型，通過將人腦部簡化為球使得正逆問題得到簡化。最為常用的為三層球頭模型，該模型包含頭皮、顱骨以及大腦三層，其它的類似球模型包括五層球模型以及重疊球模型等。然而，在使用偶極子方法進行求逆時，由于球模型與真實頭部的差異，在相關的數值運算以及激活源定位時會引入巨大的誤差[1]。這一缺陷也使得球模型的應用受到了限制。為了提高新模型的準確程度，研究者開始考慮在構建模型時結合人腦MRI數據。

隨后產生了幾種建立真實頭模型的數學方法，如基于體元的有限元方法（Finite Element Method，FEM）以及基于面元的邊界元方法（Boundary Element Method， BEM）[2]，這兩種方法都是通過使用三角元來構建模型。為了保證模型的準確性，使用的三角元必須很小，因此需要大量的節點坐標，隨之而來的就是龐大的運算量，致使真實頭模型無法得到廣泛應用[3]。本文提出了一種基于非均勻有理B樣條（Non-Uniform Rational B-Splines，NURBS）的構建真實頭模型的新方法。NURBS曲面是一種在工業上早已得到廣泛應用的方法，其它領域也有許多基于NURBS的建模方法研究。插值與擬合是最為基礎的兩種方法[4]。這里使用對MRI圖像進行插值的方法重建真實頭模型。該方法在保證結果良好的情況下建模速度更快且所使用的點更少。

1 總體流程設計

整個建模流程如圖1所示，首先對MRI數據進行圖像分割、邊緣提取等預處理，接著將圖像邊緣數據信息進行整合并規范化，最后對預處理得到的數據插值，反求出真實頭模型。

圖1 插值法構建頭模型流程圖Fig.1 Flow chart of head modeling by interpolation

2 圖像預處理

2.1 圖像分割

圖像預處理部分的第一步是對MRI圖像進行分割[5]。使用的數據為Compressed NifTI（.nii.gz）格式，是磁共振常用的數據格式。掃描圖像為217×181×181矢狀圖，即從如圖2所示的矢狀面角度看，每層圖像大小為217×181，共掃描了181層。

由于實際圖像不可能在特定結構區域有穩定的灰度值，因此需要應用一定的形態學方法對圖像加以處理，去除可能的噪點，得到邊界后才能進行下一步的初始化。

圖2 圖像分割Fig.2 Image segmentation

2.2 邊緣提取

該模塊的第二部是邊界提取，這里使用變形方法完成頭部MRI圖像邊緣提取。該方法首先對原始MRI圖像進行閾值處理并二值化如圖3（a）所示。二值化后求得該圖像的重心和半徑，利用其重心和半徑構造一個位于頭內部的初始球面，再通過不斷變形迭代更新曲面，最終得到所要求的邊緣結果。接下來再進行形態學處理，使用圖像閉合方法，即對圖像進行先膨脹后腐蝕，再進行填孔運算，最終得到圖3（b）所示的大腦邊界圖像

圖3 邊界提取Fig.3 Boundary extractioin

3 基于NURBS的曲面重構

3.1 數據初始化

經過圖像的預處理模塊，我們已經得到了MRI圖像的頭部邊界，要以這些邊界點作為已知數據點進行曲面插值。要進行曲面插值，首先要求這些數據點必須是有序的，所以在對MRI圖像進行分割得到邊緣點后，必須按照一定規則，通過邊界跟蹤的方法按順序排列這些邊緣點，從而為插值做準備。此外，由于邊緣提取方法所限，不可避免的會出現一些孤立的噪聲點，可能導致跟蹤方法失敗，因此在邊緣跟蹤開始前，需要先對提取的邊緣數據去噪，即將孤立點以及與邊緣相連的一些異常點刪除。

邊界跟蹤即從圖像的一個邊緣點出發，按照一定的規則搜索下一邊緣點直到滿足特定終止條件，由此得到目標邊界。具體步驟為：首先，確定起始搜索點。本文取圖像左上角第一個有值點作為起始搜索點S，將其位置存入邊界點序列并作為第一個邊界點。將所有點搜索狀態置0，存儲在一個與圖像等大小的零矩陣中，且該點不是已搜索得到的邊界點。

然后，通過對提取的腦邊緣進行八鄰域搜索，將各掃描層大腦邊界點的坐標數據按順序以矩陣形式存儲在一個元胞數組內。由于總的數據點數目依然很大，如果直接將這些點用于構建模型仍然會給計算機很大的運算負載，與此同時，生成的模型也缺乏光順性。因此需要在進行曲面重構前對數據進行預處理。

由于MRI數據中各層腦切面大小不同，因此所存儲的邊界數據點也不相同，首先必須將數據矩陣規范為相同大小，并在這一過程中降低矩陣中元素的數目。本文所用MRI數據中共有146層矢狀面包含大腦數據，故邊界矩陣中包含146行數據點。各行對應數據由數十列到數百列不等。如果進行曲面重構時需要一個大小為m×n的矩陣，那么對于數據長度超過n列的行來說就可以進行采樣，而對于少于n列數據點的則按順序插入一些點。通過這個方法就可以生成一個m×n的邊界點矩陣，表示為 Qij（i=0，1，…，m；j=0，1，…，n）。 由于最終建模時所需要的n小于多數切面對應的邊界數據點，因此在對數據進行預處理后的總數據量有大幅減少。這里n的選取可以根據各層所需要的重構曲線質量自由選擇。

3.2 曲面重構

接下來的問題就是如何根據邊界數據使用NURBS重構大腦曲面。由于數據點已經被規范化為一個矩陣，因此這個問題可以通過一般的曲面插值得到解決。

一個（p，q）次的B樣條曲面表示形式如下[6]：

其中 Pi，j為控制點或定義點，Ni，q（u）與 Nj，q（v）為 p 與 q 次下的B樣條基函數，u和v為節點，對應節點向量分別為：

曲面插值問題可以表示為

其中

這里使用雙三次樣條曲面插值，即p=3，q=3。插值過程包含以下幾步：

1）數據點參數化，采用積累弦長參數化方法

首先令d為總弦長

那么有u0=0 un=1

該方法可用于大多數情況，且對應產生的插值曲線結果光順性較好，符合通常的要求，因此是應用最為廣泛的參數化方法。它在數據點分布不均勻或控制多邊形弦長不等長的情況下，能夠避免如上均勻參數化方法的問題，給出符合數據點或弦長分布的參數化結果，所以它在某種程度上也被粗略的認為是弧長參數化。這里需要強調的是，曲線插值最后生成的結果光順性好壞并不僅僅取決于所選的參數化方法，同時還由所選的插值方法所決定。此外，這里采用弦長參數化方法也并不等同于說最終的插值曲線參數是積累弦長。

2）計算節點向量U

節點通常可以進行等距離分割

然而通常并不推薦使用這種方法，這是因為當結合積累弦長參數化方法使用時，有可能會產生一個奇異方程組。因此使用如下的平均方法

通過使用這一方法可以使得節點的分布反映出參數化后點uk的分布。

3）利用公式（4）對U向使用NURBS曲線插值，求得一組“中間”控制點。

4）根據公式（5）對V向使用NURBS曲線插值，求得控制點 Pi，j。

5）使用 U，V 以及 Pi，j計算反求 NURBS 曲面。

4 重構結果

整個建模過程在MATLAB環境下完成。最終重構真實頭模型中的大腦部分如圖4所示。

圖4 大腦重構結果Fig.4 Reconstructed result

計算機仿真結果如圖4所示，最終的生成結果良好。除此之外，在數據預處理時，數據點的個數由超過70 000點降至7 000余點。而如果使用三角元構建模型，數字將大幅提升。上述建模過程在一臺處理器為Intel Core（TM）2 Duo Processor，內存為2G的普通臺式機上僅耗時十秒左右。

5 結 論

該建模系統使用模塊化思想結合NURBS新方法構建真實腦模型，具有高速、少數據建模的特點。對于該方法有效性的更進一步評估可以通過在腦電的正逆問題計算中進行量化分析。此外，未來開展的研究中可以使用更為復雜的NURBS方法。NURBS作為真實頭模型構建中一種很有潛力的方法，是值得投入更多精力的。

[1]Lalancette M，Quraan M，Cheyne D.Evaluation of multiplesphere head models for MEG source localization[J].Physics in Medicine and Biology，2011，56（17）:5622-5636.

[2]Akahn Z，Acar C E，Gencer N G.Development of realistic head models for electromagnetic source imaging of the human brain engineering[J].Proceedings of the 23rd Annual International Conference of the IEEE Medicine and Biology Society，2001，1:899-902.

[3]He J J，Shen H，Hu D W.A survey of the EEG/fMRI fusion analysis:head models，methods and applications[J].Computer Engineering&Science，2007，29（12）:74-81.

[4]Piegl L，Tiller W.The NURBS Book[M].New York:Springer，1997.

[5]Smith S M.Fast robust automated brain extraction[J].Human Brain Mapping，2002，17（3）:143-155.

[6]Piegl L，Tiller W.Reducing control points in surface interpolation[J].IEEE Computer Graphics and Applications，2000，20（5）：70-74.

[7]李焱，時芝勇，海曉濤.基于改進蟻群算法的配電網重構[J].陜西電力，2010（9）：22-25.LI Yan,SHI Zhi-yong,HAI Xiao-tao.Distribution Network Reconfiguration Based on Improved Ant Colony Algorithm[J].Shaanxi Electric Power，2010（9）：22-25.