結合廣義重心坐標與Voronoi剖分的函數分片逼近

肖艷陽， 涂錦燦， 陳中貴

（1. 廈門大學福建省智慧城市感知與計算重點實驗室，福建 廈門 361005；2. 廈門大學計算機科學系，福建 廈門 361005）

結合廣義重心坐標與Voronoi剖分的函數分片逼近

肖艷陽1,2， 涂錦燦1,2， 陳中貴1,2

（1. 廈門大學福建省智慧城市感知與計算重點實驗室，福建 廈門 361005；2. 廈門大學計算機科學系，福建 廈門 361005）

結合廣義重心坐標理論，提出了一個新方法，以解決在平面區域上的函數逼近問題。該方法通過構建基于廣義重心坐標的最優分片函數來逼近目標函數。采用Voronoi圖來劃分區域，并提出一個度量逼近誤差的能量函數。推導出該函數的導數后，采用一種高效的Voronoi節點更新方法來獲得區域的最優剖分，并通過最優剖分構建最優分片函數。由于該方法對不連續函數具有良好地逼近能力，因此將其應用在圖像逼近問題中。分別在解析函數和彩色圖像上對該方法進行實驗，均獲得了很好的逼近效果。

廣義重心坐標；函數逼近；Voronoi圖；圖像逼近

函數逼近問題是指任意給定一個目標函數，在該函數的定義域上，構造一個新的函數，使得新函數與目標函數盡可能一致。兩個函數之間不一致的程度可用逼近誤差來度量，誤差越小，兩個函數越接近。如果在目標函數的定義域上只構建一個函數來逼近它，誤差可能會較大。一般的想法總是在誤差最大的區域上細分成多個區域。但是這會造成過多的子域，如何減少剖分的數量也是一類值得研究的問題[1]。這種將定義域劃分成多個子域，再在每個子域上分別構建一個新函數來逼近目標函數的方法即分片的函數逼近，它是研究逼近問題的重要方法。

對于分片逼近來說，基于多項式的簡易性和靈活性特點[2]，Chen等[3]提出了一種在Voronoi剖分上構建最優分片多項式的方法。分片多項式需要同時考慮多項式的階次和對區域的剖分，一旦指定了階次和剖分的數量，分片多項式的逼近問題就轉化成如何劃分區域。因為得到剖分之后，每個子域上的多項式就可以通過求解線性最小二乘問題得到。因此，對區域良好地剖分是分片逼近問題的關鍵。

受此啟發，本文提出一種新的分片逼近方法，即結合廣義重心坐標和Voronoi剖分，通過構造最優分片函數來逼近二維定義域上的目標函數，其中，逼近誤差采用常用的 L2范數來衡量。因此本文方法只需要考慮如何獲得區域的最優剖分。劃分區域有許多方法，如三角剖分，本文采用Voronoi圖，它是由Voronoi節點集合唯一決定的結構。在形成的剖分中，每一個Voronoi多邊形都對應一個逼近函數，而由于節點集合唯一決定了Voronoi圖的形成，因此逼近誤差間接取決于節點集合的分布。把逼近誤差極小化對應的剖分稱為最優剖分，此時的點集分布為最優分布，該剖分對應的分片逼近函數稱為最優分片函數。

在此基礎上，本文提出了一個新的能量函數來衡量逼近誤差，并推導出相應的梯度。然后采用一種高效的節點優化方法來更新節點位置，使得能量函數值極小化，從而獲得最優剖分和最優分片函數。為了驗證方法的有效性，分別在連續的和不連續的解析函數上應用該方法，均能獲得很好地逼近效果。由于該方法對不連續函數有很強地逼近能力，而圖像可以看作是不連續的函數，所以將該方法應用到圖像逼近領域。實驗結果表明，該方法對彩色圖像也能產生較好的結果。

1 相關研究和背景

1.1 相關研究

分片函數逼近。Powell[2]詳細介紹了函數逼近的概念和基本方法。分片逼近是該領域的一種重要方法，其中，分片多項式的方法已經有較多的研究成果[1-8]。大多數相關研究是針對單變量函數的逼近問題，而多變量的多項式逼近問題則更加復雜。如果多項式的基函數是正交的，Fedele和Ferrise[4]的方法可以用來處理最小二乘逼近問題。單變量的正交多項式可以使用 Gram-Schmidt過程來構建，而多變量的情況則面臨諸多數學和算法上的困難。然而在一些特殊區域中已經可以構造二變量的正交多項式，如Dunkl[5]的正六邊形區域和Farouki等[6]的三角形區域。在任意區域中構建多變量正交多項式的問題仍然值得研究。Lecot和 Lévy[7]則使用冪基多項式構建二變量的逼近函數，區域的剖分是根據底層函數值進行聚類得到的。Chen等[3]的方法中也使用了冪基多項式，但是采用了更為簡單的Voronoi圖作為區域的剖分結構，通過獲取最優剖分得到最優逼近多項式。Nivoliers和Lévy[8]的方法在平面區域上即可看做Chen等[3]方法的特殊情況，即多項式取常數值。本文與Chen等[3]的方法類似，但是可使用幾何意義更鮮明的重心坐標來構造分片函數，并獲得了很好的效果，而且數值計算更加穩定。

圖像逼近。該領域已經有很多研究成果和高效的算法被提出來[3,7-14]。許多方法是通過自動提取圖像特征的幾何信息實現的，如通過邊緣檢測算法提取圖像的特征曲線，再根據特征線構建矢量曲線和像素簇。比較常見的用幾何信息來表達圖像的方法，如 數據 相關 的三 角化 (data-dependent triangulation, DDT)，通過構建圖像的三角剖分，再線性插值出像素的顏色值。Dyn等[9]結合三角剖分的頂點位置分布和對應的值進行三角剖分的構建，同時討論了不同的定義準則的優、缺點。Li和Adams[10]在此基礎上提出了一種用來表達圖像的三角網格生成框架，根據不同的規則，該方法不斷地向網格中誤差最大的三角片上添加頂點，從而逐漸細分網格。DDT方法的另一個方向是通過簡化復雜網格并在達到最大誤差時結束，如 Yin等[11]的方法中使用了網格簡化。Su和Willis[12]則提出了像素級別的DDT方法，可將低分辨率的圖像放大到更高分辨率。但是大多數方法都不能顯式地表達出圖像的不連續特點，為了捕捉到圖像的特征線，DDT方法只能在特征線兩邊放置大量的點。Tu和Adams[13]的方法則根據圖像不連續特征，提出了一種新的網格模型生成方法。首先提取圖像的特征線，根據特征線定義約束的Delaunay邊并產生初始網格，然后再不斷地向網格誤差大的片上添加頂點。Kreylos和Hamann[14]則應用模擬退火的方法優化網格頂點位置來獲取最優三角剖分，繼而得出擬合圖像。Lecot和Lévy[7]將Cohen-Steiner等[15]的框架擴展到圖像處理中，從而提出了一種新方法，可以將光柵圖像矢量化。該方法將圖像劃分成多個不規則的子區域，每個區域的邊界用矢量曲線來表達，內部由事先構造的三角像素組成，再根據能量函數的梯度更新子區域。Yin等[11]結合一種新的二次誤差度量(quadric error metrics)，提出了將光柵圖像參數化的新方法，并將其應用在圖像重建中。Chen等[3]的方法則通過獲取圖像的最優剖分，再從該剖分中構建最優多項式的方法重建圖像。類似地，將本文方法應用在圖像逼近中則通過構建基于廣義重心坐標的最優分片函數來擬合原圖像，能獲得很好的效果。李海洋等[16]應用 K-均值聚類的方法產生了類似于文獻[7]中的圖像分割效果，它們的剖分都是不規則的，本文方法應用在圖像上的剖分在一定程度上也可看做一種分割結果。

1.2 背景知識

Voronoi圖是計算幾何中重要的基礎概念，在計算機圖形學、圖像處理、地理信息系統和生命科學領域都有廣泛地應用[17]。假設在二維緊致區域Ω?R2中，有N個位置互異的節點組成的集合則在歐式距離度量下，節點xi對應的Voronoi區域為：

本文方法基于廣義重心坐標，其包括均值坐標(mean value coordinates, MVC)[18]、調和坐標(harmonic coordinates, HC)[19]、格林坐標(Green coordinates, GC)[20]和Wachspress坐標[21]。MVC將多邊形內的任意點表達為多邊形各個頂點的線性組合，因此在多邊形內部具有良好的光滑性，在頂點處則是C0連續。假設在一個任意的平面多邊形中，各個頂點表示為｛p1,p2,…,pmi｝，mi表示該多邊形的頂點個數，則在多邊形內的任意一點可以表示成：

其中，Ck≥0,k=1,…,mi，Ck（q）代表點q在該多邊形中相對頂點pk的重心坐標值（權值）。

圖 1是均值重心坐標的示例圖。在該多邊形中，內部的任意點相對多邊形的每個頂點都有一個權值C。權值越大，代表頂點對其影響越大。以圖1為例，該多邊形有6個頂點，其中左上的子圖代表多邊形內部的點對于頂點1的權值分布，中上的子圖代表多邊形內部的點對于頂點 2的權值分布，以此類推。

Wachspress坐標[21]最先用于有限元計算，在平面中，它的計算可采用Meyer等[22]的方法，但是計算結果可能會產生極點，在非凸的多邊形中會產生負值坐標。因此Wachspress坐標適用于單個凸多邊形。而MVC由于在整個平面上具有良好地閉形式定義，因此適用于不自交的任意多邊形以及多邊形的嵌套[18]。與MVC不同，HC的定義并不是閉形式，但是在多邊形內部它也是非負的，并且在整個平面沒有極點，因此HC能夠很好地保持以往重心坐標的性質[19]。Lipman等[20]結合格林第三定理提出了用于變形的GC，它在計算中考慮了多邊形邊界的法向量，因此變形過程可以保持細節特征。

圖1 多邊形內的點相對各個頂點的權值分布

2 基于廣義重心坐標的分片函數逼近

2.1 能量函數的定義

設函數f（x）是定義在二維緊致區域Ω中的函數。在該區域上有Voronoi剖分取定函數空間Φ，基函數取為廣義重心坐標基。根據廣義重心坐標理論，在任意一個 Voronoi多邊形區域內都可以構建一個定義在Φ中的函數其中λk,k=1,…,mi表示多邊形的頂點pk在該函數中的取值，其余符號的意義同重心坐標部分。λk具有良好的幾何意義，即在該平面多邊形中，頂點pk的高度值，k=1,…,mi。

于是在每個子區域Ωi中，f（x）都可以用一個新函數Hi（x）來逼近。建立一個衡量逼近誤差的函數，即每個子域上誤差的累加和，表達為：

將該函數稱為本文的能量函數。將目標轉化成求解該函數的極小值，以獲得逼近f（x）的最優分片函數。

對能量函數來說，最優的概念包含兩部分：對區域Ω的剖分最優和在每個子域Ωi上對應的函數最優。當剖分確定時，子域Ωi上對應的最優函數滿足下式：

在子域Ωi上的逼近誤差：

分別對λk,k=1,…,mi求偏導，并令其為0，得：

整理得：

因此Ωi上的最優函數僅依賴于剖分。而唯一決定Voronoi剖分的就是節點的位置，所以可以將節點集合看做唯一變量，重寫能量函數為：

下面探討如何找到點集X的最優位置分布。現有的許多優化方法，如梯度下降法、共軛梯度法等，都是基于梯度信息。因此，先推導出能量函數的梯度。

2.2 E（X）梯度的計算

用 Ji表示與Ωi相鄰的Voronoi多邊形對應的節點的序號集合。為了推導E（X）對節點xi的導數，只需要考慮在E（X）中與xi相關的項。因此，能量函數對 xi求偏導：

注意到xi出現在被積函數和積分區域當中，需考慮積分區域中變量的不同情況。

首先，應用Leibniz規則[23]來簡化這個公式。假設Dt是一個在時間t上光滑連續的二維區域，g（x,t）,x∈Dt是定義在Dt上的一個函數。將一個在區域邊界?Dt上的點的速度矢量表示為v=?x?t，用n代表該邊界上指向外的單位法向量。一般的Leibniz規則[23]是：

其中，ds是在封閉的邊界曲線tD?上的弧長單元。因此，在式(4)中應用Leibniz規則，得到：

其中，Ωij=?Ωi∩?Ωj，是Ωi與Ωj的公共邊。

式(5)右邊的第一項，根據包絡定理[24]可以得到：

則等式兩邊對xi求偏導，得到：

因此，得到：

式(6)即為能量函數關于Voronoi節點的導數公式。計算得到函數的梯度之后，采用一種新穎的基于梯度信息的優化方法來優化點集X的分布。

2.3 優化算法

Nivoliers和Lévy[8]介紹了一種新的基于梯度的優化算法，Chen等[3]則從自適應的角度對迭代步長作了改進。這里引述文獻[3]中的算法過程。從初始位置開始，節點的位置根據下式進行迭代更新，

i迭代時的步長。需要特別說明的是，在這個優化式子中，由于下降方向是通過將各自梯度分量除模得到的，所以不再是梯度的方向，該優化過程也不是梯度下降法。每一次迭代的步長由下式控制：

其中，maxJ是指定的最大迭代次數，是節點xi的初始步長。其中k是縮放參數，Ωi是節點 xi對應的Voronoi多邊形區域。

這種方法比經典的優化算法，如傳統的梯度下降法、共軛梯度法和L-BFGS法，更加高效，能夠優化經典算法無法優化的情況。基于此，提出本文算法如下：

算法1. 基于重心坐標的分片函數逼近

輸入. 函數f（x），初始點集X，最大迭代次數Jmax

輸出. 點集X的最優Voronoi剖分和對應的最優函數集合

(1) 計算點集X的Voronoi剖分

(2) 計算點集的初始步長

(3)DO

計算當前點集X的Voronoi剖分和對應的最優分片函數

計算能量函數的梯度

根據優化方法更新點的位置，仍記為X

WHILE迭代次數小于Jmax

(4) 計算結果點集X的Voronoi剖分和對應的最優分片函數

由于本文中能量函數的連續性不高，特別是對于不連續的圖像來說，更達不到 C2連續。梯度下降法、共軛梯度法和L-BFGS方法等傳統的優化算法，無法徹底地優化能量函數。在實驗中，傳統的優化算法運行了幾次迭代之后就停止了，Voronoi剖分并沒有得到較好地優化。本文采用的優化算法將變量的梯度信息除模，為其指定了優化方向，但該方向已不再是梯度下降的方向。因此給定合適的步長之后，該算法可以強制更新Voronoi節點的位置，直到達到指定的迭代次數。

對于算法中的初始化，本文采用簡單的隨機算法。在產生隨機節點時，拋棄重合的點，以保證剖分個數為指定的數量。

將該算法應用在圖像逼近問題時，把圖像作為底層函數，將像素顏色值正規化。先產生初始剖分，根據算法流程，在達到最大迭代次數時，得到最優剖分。據此，在每個分片上計算最優逼近函數。根據這個函數，對位于該分片內的像素重新計算顏色值。新的顏色值可能會越界（大于255或小于0），因此將其截斷在0～255范圍內，即大于255的值設為255，小于0的值設為0。再把計算得到的顏色值填充到該像素中。重復該過程即可得到擬合圖像。擬合圖像與原圖像的相似程度可以用峰值信噪比(peak signal to noise ratio, PSNR)來衡量。理論上，PSNR值越高，兩幅圖相似程度越大。

3 實驗結果

本文算法使用了C++編程語言。所有的實驗都是在Intel I5 3.1 GHz CPU和4.0 GB內存的PC機進行的。其中，使用計算幾何(computational geometry algorithms library, CGAL)算法庫[25]來計算 Voronoi剖分。

采用文獻[26]的程序來計算重心坐標。另外由于本文算法需要計算很多的積分，為了提高效率，對于一般的解析函數積分，使用 Dunavant積分規則[26]來簡化計算。具體的做法是，對每個 Voronoi多邊形，將其劃分成三角形剖分，每個三角形由多邊形中心與多邊形每條邊的兩個頂點連接而成。然后在每個三角形中利用Dunavant規則[27]計算積分。而對于圖像上的積分，不再使用普通的正方形像素，而是將每個多邊形劃分成三角像素的集合，每個三角像素的面積約為原像素面積的2倍。對于梯度的計算，邊界積分使用經典的兩點高斯積分規則（連續函數）和一點高斯積分規則（圖像）。

圖2和圖3是本文算法對一般解析函數的剖分圖。圖2所示為算法對連續函數：

的剖分結果。該函數在x+y=0處是C0連續的，從圖中可以看到，Voronoi剖分的分布具有規律性，雖然函數是連續的，但在x+y=0處函數值變化較大，剖分很敏銳地捕捉到了這一特點。其他地方基本上沿著函數的梯度方向分布。

圖3所示則是算法對不連續函數：

的剖分結果。分析結果可知，Voronoi剖分在函數的斷層處可捕捉到函數的不連續特征。圖中由Voronoi區域邊界形成的近似正方形的內外區域分別對應函數的不同定義域的表達。

圖2 本文算法逼近連續函數的結果

圖3 本文算法逼近不連續函數的結果

在圖像逼近中，普通的灰度圖像可視為一個不連續的函數，由于本文方法具有很強地逼近不連續函數的能力（圖3），所以將本文算法直接應用在灰度圖像上是可行的，實驗結果也證明了這一點。再做一些擴展，對于一幅RGB的彩色圖像，將紅、綠和藍色通道分離出來，視為3個不連續的函數，r（x），g（x）和b（x）。為了在圖像中同時逼近這 3個函數，把能量函數修改為：

其中，R*（x），G*（x）和B*（x）分別是對函數 r（x），g（x）和b（x）的最優函數。這個能量函數的導數推導過程與2.2節類似，算法1仍然可以執行。圖4是本文算法的流程，圖5～8是算法的結果，所有列出的PSNR值均為擬合圖像與原圖像的PSNR。

圖4采用簡單的隨機算法初始化500個Voronoi節點，如圖4(b)所示。通過不斷地迭代更新節點位置，使得能量函數值逐漸減小，Voronoi剖分趨于最優。實驗結果表明，絕大部分的Voronoi區域的邊界都能夠與圖像的特征線重合，因此能夠很好地擬合出原圖像。圖 4(b)～(d)的底層圖像都是利用剖分對應的逼近函數重新擬合出來的。

本文方法可看做是文獻[3]類型方法的擴充。一般的Voronoi多邊形的頂點個數在5～7之間，以自由度來看，本文方法與文獻[3]中的二次逼近（6個自由度）相近。圖5是一個比較結果，從PSNR的角度看，本文方法與二次逼近方法得到的效果相近，而比線性逼近效果更好。而且，本文所使用的基函數具有更直觀的幾何意義，在數值計算方面也更加穩定。這主要體現在以下兩方面：①本文算法和 Chen等[3]的算法，在求解最優逼近函數的過程中，都需要解一個線性方程組。Chen等[3]利用冪基多項式構建逼近函數，在求解該過程中的線性方程組時，有可能出現退化情況（即該方程組的秩小于方程的數量），如二次退化成線性的情況。而本文方法并不會出現這種退化情況。②Chen等[3]的算法中，多項式的階次越高，求解最優多項式的線性方程組的系數矩陣越龐大，計算的復雜程度越大。而本文算法在該過程中的系數矩陣只跟Voronoi多邊形的頂點個數相關。

本文方法還對算法內部的參數，即剖分數量和不同類型的重心坐標做了比較。圖6分別是500、1 000和2 000個剖分執行算法后生成的擬合圖，從圖6(a)～(c)的PSNR值變化可知，逼近誤差與剖分數量呈反比。但是過多的分片不利于算法效率，需平衡選擇。圖7所示為在相同的初始剖分下，分別使用均值坐標和Wachspress坐標執行算法的結果，對應的部分最優剖分分別顯示在圖像的下半部分?？梢钥吹?，使用不同的坐標產生的最優剖分雖然不同，但是均能產生良好地逼近效果。在實驗中，由于多邊形都是凸的，不需要考慮坐標值為負的情況，因此除了計算方式不一樣外，算法的其他部分不需要做任何修改。算法的更多結果見圖8。

圖4 本文算法的流程圖

圖5 本文方法與文獻[3]方法的比較

圖6 不同剖分數量對結果的影響

圖7 相同的初始化

圖8 本文算法實驗結果

表1是算法的運行時間統計。圖像分辨率、剖分數量和計算梯度的精度都會影響算法中積分的計算，從而影響算法的效率。最大迭代次數也是影響算法運行時間的重要參數。將本文算法歸結為以下四部分：構建三角像素、求解最優逼近函數、計算梯度和更新節點。分別對其進行時間統計，得到本文算法主要耗時在于構建三角像素和求解最優逼近函數。

表1 運行時間統計

由于最優逼近函數求解出來之后，梯度的計算只涉及到線積分，而且相關參數也已經求解得到，因此計算量較小，耗時較少。重心坐標的計算貫穿于整個算法過程，在實驗中，對于512×512的灰度圖像，在多邊形頂點個數均為6的情況下，本文算法對所有像素執行一次重心坐標的計算耗時為0.185 s。然而有些像素可能需要多次計算重心坐標，如某個像素同時出現在求解最優函數過程中和計算梯度過程中。另外需要注意的是，對于RGB圖像來說，不同通道是分別計算的，重心坐標的計算量為灰度圖像的3倍，因此算法運行時間更久。總體來說本文算法計算量較大，相比于Chen等[3]的算法，本文算法的運行速度也較慢。

4 結 束 語

本文提出了一種基于重心坐標的分片函數逼近方法。在廣義重心坐標和Voronoi剖分的基礎上，首先給出了衡量逼近質量的能量函數，該能量函數僅依賴于Voronoi節點的位置。 推導出該函數的導數后，采用文獻[3]中改進的基于梯度的優化方法來優化節點位置，從而獲取最優剖分和最優分片函數。實驗結果表明，該方法對解析函數非常有效，并且能夠應用在圖像逼近問題上。關于與其他方法的比較，在相同自由度下，本文結果與其相近甚至更好，而且文獻[3]中使用的基函數，其幾何意義不明顯，因此使用重心坐標更加合理。

不足之處在于，本文方法的逼近能力有限，當多項式階次取更高時，Chen等[3]方法的逼近能力遠遠高于本文方法。而且，本文算法的運行時間相對更久（表1），對于某些應用可能是不可接受的。該算法與許多因素相關，主要耗時在于構建三角像素和梯度的計算。以后需考慮采用并行算法或圖形處理器(graphic processing unit, GPU)執行來加速本算法。研究本文方法在不同基函數下的實現和內在聯系，或者將其拓展到高維空間中，獲取高維的最優剖分，或者應用在視頻相關領域中，都是本文方法的后續工作。

[1] Prentice J S C. Range and domain partitioning in piecewise polynomial approximation [J]. Studies in Mathematical Sciences, 2011, 2(2): 67-77.

[2] Powell M J D. Approximation theory and methods [M]. Cambridge: Cambridge University Press, 1981: 25-30.

[3] Chen Zhonggui, Xiao Yanyang, Cao Juan.Approximation by piecewise polynomials on Voronoi tessellation [J]. Graphical Models, 2014, 76(5): 522-531.

[4] Fedele G, Ferrise A. Explicit solution of the finite time L2-norm polynomial approximation problem [J]. Applied Mathenatics and Computation, 2011, 217(21): 8354-8359.

[5] Dunkl C F. Orthogonal polynomials on the hexagon [J]. SIAM J. Applied Mathenatics, 1987, 47(2): 343-351.

[6] Farouki R T, Goodman T N T, Sauer T. Construction of orthogonal bases for polynomials in bernstein form on triangular and simplex domains [J]. Computer Aided Geometric Design, 2003, 20(4): 209-230.

[7] Lecot G, Lévy B. Ardeco: Automatic region detection and conversion [C]//Proceedings of the 17th Eurographics Conference on Rendering Techniques. Nicosia, Cyprus: Eurographics Association, 2006: 349-360.

[8] Nivoliers V, Lévy B. Approximating functions on a mesh with restricted Voronoi diagrams [J]. Computer Graphics Forum, 2013, 32(5): 83-92.

[9] Dyn N, Levin D, Rippa S. Data dependent triangulations for piecewise linear interpolation [J]. IMA Journal of Numerical Analysis, 1990, 10(1): 137-154.

[10] Li Ping, Adams M D. A tuned mesh-generation strategy for image representation using data-dependent triangulation [J]. IEEE Transactions on Image Processing, 2013, 22(5): 2004-2018.

[11] Yin Xuetao, Femiani J, Wongka P, et al. A new QEM for parametrization of raster images [J]. Computer Graphics Forum, 2011, 30(8): 2440-2451.

[12] Su Dan, Willis P. Image interpolation by pixel-level data-dependent triangulation [J]. Computer Graphics Forum, 2004, 23(3): 189-201.

[13] Tu Xi, Adams M D. Improved mesh models of images through the explicit representation of discontinuities [J]. IEEE Canadian Journal of Electrical and Computer Engineering, 2013, 36(2): 78-86.

[14] Kreylos O, Hamann B. On simulated annealing and the construction of linear spline approximations for scattered data [J]. IEEE Transactions on Visualization and Computer Graphics, 2001, 7(1): 17-31.

[15] Cohen-Steiner D, Alliez P, Desbrun M. Variational shape approximation [J]. ACM Transactions on Graphics (TOG), 2004, 23(3): 905-914.

[16] 李海洋, 文永革, 何紅洲, 等. 基于隨機權重粒子群和K-均值聚類的圖像分割[J]. 圖學學報, 2014, 35(5): 755-761.

[16] Okabe A, Boots B, Sugihara K. Spatial tessellations: concepts and applications of Voronoi diagrams [M]. New York: Wiley, 1992: 1-3.

[17] Hormann K, Floater M. Mean value coordinates for arbitrary planar polygons [J]. ACM Transactions on Graphics (TOG), 2006, 25(4): 1424-1441.

[18] DeRose T, Meyer M. Harmonic coordinates [M]. Pixar Animation Studios, Emeryville, CA, Pixar Technical Momo, 2006: 06-12.

[19] Lipman Y, Levin D, Cohen-Or D. Green coordinates [J]. ACM Transactions on Graphics (TOG), 2008, 27(3): 1-10.

[20] Wachspress E. A rational finite element basis [M]. New York: Academic Press, 1975: 49-51.

[21] Meyer M, Lee H, Barr A, et al. Generalized barycentric coordinates on irregular polygons [J]. Journal of Graphics Tools, 2002, 7(1): 13-22.

[22] Flanders H. Differentiation under the integral sign [J]. The American Mathematical Monthly, 1973, 80(6): 615-627.

[23] Silberberg E, Suen W. 經濟學的結構——數量分析方法[M]. 3版. 張研, 譯. 北京: 清華大學出版社, 2004: 125-140.

[24] Fabri A. CGAL-the computational geometry algorithm library [C]//Proceedings of 10th International Meshing Roundtable, Newport Beach, California, 2001: 137-142.

[25] Gleicher M. Generalized barycentric coordinates for mesh deformation [EB/OL]. (2011-05-11) [2014-03-10]. http://pages.cs.wisc.edu/～csverma/CS777/bary.html.

[26] Dunavant D A. High degree efficient symmetrical Gaussian quadrature rules for triangle [J]. International Journal for Numerical Methods in Engineering, 1985, 21(6): 1129-1148.

Approximation by Piecewise Function Based on Generalized Barycentric Coordinates and Voronoi Tessellation

Xiao Yanyang1,2, Tu Jincan1,2, Chen Zhonggui1,2

(1. Fujian Key Laboratory of Sensing and Computing for Smart City, Xiamen University, Xiamen Fujian 361005, China; 2. Department of Computer Science, Xiamen University, Xiamen Fujian 361005, China)

Under the generalized barycentric coordinates theory, we propose a new method to solve the problem of approximating a given function on the planar domain. To accomplishing this, an optimal piecewise function which based on the generalized barycentric coordinates is constructed. We use the Voronoi tessellation to create a partition of the domain, then an energy function that measures the approximation error is built. After deriving the gradient of the energy function, an efficient optimization method is adopted to update the tessellation. The optimal piecewise function will be constructed from the optimal tessellation. Due to its good ability of approximating discontinuous functions, our method can be applied to image approximation field. In order to demonstrate its efficacy, some experiments on analytic functions and color images are designed, which have produced good results.

generalized barycentric coordinates; function approximation; Voronoi tessellation; image approximation

TP 391

A

2095-302X(2015)03-0367-09

2014-10-08；定稿日期：2014-10-24

國家自然科學基金資助項目(61100107, 61472332)；中央高?；究蒲袠I務費專項基金(廈門大學基礎創新科研基金)資助項目(20720140520)

肖艷陽(1991-)，男，江西贛州人，碩士研究生。主要研究方向為計算機圖形學。E-mail：jxndxyy@gmail.com

陳中貴(1982-)，男，浙江溫州人，副教授，博士。主要研究方向為計算機圖形學。E-mail：chenzhonggui@xmu.edu.cn