黃童心， 王文珂， 張 慧， 宋征軒

（1. 清華大學(xué)軟件學(xué)院，北京 100084；2. 清華大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)系，北京 100084；3. 信息系統(tǒng)安全教育部重點(diǎn)實(shí)驗(yàn)室，北京 100084；4. 清華信息科學(xué)與技術(shù)國(guó)家實(shí)驗(yàn)室，北京 100084）

基于場(chǎng)分布的平面散亂點(diǎn)集B樣條曲線重建算法

黃童心1,3,4， 王文珂1,2,3,4， 張 慧1,3,4， 宋征軒1,3,4

（1. 清華大學(xué)軟件學(xué)院，北京 100084；2. 清華大學(xué)計(jì)算機(jī)科學(xué)與技術(shù)系，北京 100084；3. 信息系統(tǒng)安全教育部重點(diǎn)實(shí)驗(yàn)室，北京 100084；4. 清華信息科學(xué)與技術(shù)國(guó)家實(shí)驗(yàn)室，北京 100084）

平面散亂點(diǎn)集的曲線重建是逆向工程研究的核心問(wèn)題之一。該文在Goshtasby算法的基礎(chǔ)上，提出了一種基于場(chǎng)分布的平面散亂點(diǎn)集B樣條曲線重建算法。首先，通過(guò)估計(jì)場(chǎng)強(qiáng)基函數(shù)的邊界提高量子化效率，生成散亂點(diǎn)集場(chǎng)分布的數(shù)字圖像；然后，利用圖像細(xì)化結(jié)合改進(jìn)的BFS（Breadth-First-Search）算法來(lái)避免數(shù)字圖像中由于存在大量冗余分支像素而難以生成脊輪廓的問(wèn)題；最后，采用加權(quán)最小二乘法延長(zhǎng)重建曲線，改進(jìn)Goshtasby算法所得的開曲線在端點(diǎn)處收縮的缺點(diǎn)。實(shí)驗(yàn)表明，對(duì)于帶噪聲的平面稠密點(diǎn)集，該算法可有效地重建反映點(diǎn)集形狀和走向的B樣條曲線。

計(jì)算機(jī)應(yīng)用；B樣條曲線重建；場(chǎng)分布；散亂點(diǎn)集

由平面采樣點(diǎn)集重建曲線是逆向工程幾何造型領(lǐng)域中的一個(gè)重要問(wèn)題[1-3]。根據(jù)采樣數(shù)據(jù)點(diǎn)集的性質(zhì)可將曲線重建方法分為有序離散點(diǎn)集的曲線擬合和無(wú)序散亂點(diǎn)集的曲線擬合。對(duì)于前者而言，曲線重建效果往往受到數(shù)據(jù)點(diǎn)集參數(shù)化方法的影響，其代表有均勻參數(shù)化、弦長(zhǎng)參數(shù)化和向心參數(shù)化[4-7]等。在逆向工程領(lǐng)域獲取的測(cè)量數(shù)據(jù)一般是大規(guī)模的無(wú)序散亂點(diǎn)，有噪聲且采樣密度較為稠密，難以直接進(jìn)行曲線擬合，往往需要首先進(jìn)行定序[1-2]。特別地，當(dāng)數(shù)據(jù)點(diǎn)集反映的幾何拓?fù)浣Y(jié)構(gòu)較為復(fù)雜時(shí)，僅僅將其擬合成一條單獨(dú)的曲線是不合適的，應(yīng)按其結(jié)構(gòu)特征分別以各條連續(xù)的光滑曲線來(lái)表示。

目前無(wú)序散亂點(diǎn)集曲線重建的典型方法有如下幾類：Levin、Lee的移動(dòng)最小二乘法(MLS)[2,8]，Pottmann的圖像細(xì)化法[3]，F(xiàn)ang的彈力模型法[9]和Taubin的隱式單純模型法[10]，Amenta的基于Voronoi圖Crust算法的β骨架法[11]等。其中，移動(dòng)最小二乘法需要對(duì)數(shù)據(jù)點(diǎn)集進(jìn)行兩次局部最小二乘回歸，計(jì)算量較大；圖像細(xì)化法的準(zhǔn)確性受點(diǎn)集投影的網(wǎng)格分辨率的直接影響，重建效果較為粗糙；模型法依賴于給定優(yōu)化模型的目標(biāo)函數(shù)，往往需要進(jìn)行迭代求解，對(duì)于數(shù)據(jù)點(diǎn)集中的噪聲較為敏感；骨架法則不適用于稠密點(diǎn)集，而且抗噪性不夠理想。國(guó)內(nèi)也有學(xué)者提出了平面無(wú)序點(diǎn)集曲線重建的跟蹤算法[12]和場(chǎng)表示算法[13]，避免了上述方法的某些弊端，文獻(xiàn)[12]成功模擬連續(xù)曲線的跟蹤過(guò)程并引入質(zhì)量控制機(jī)制解決了簡(jiǎn)單光滑曲線的重建問(wèn)題；而文獻(xiàn)[13]則沿著構(gòu)造場(chǎng)函數(shù)的梯度方向運(yùn)動(dòng)以其極限位置為重建曲線。Lin的基于區(qū)間B樣條[14]的曲線重建方法則可以有效過(guò)濾噪聲，以及處理不均勻的采樣點(diǎn)集。

由于點(diǎn)集形狀千差萬(wàn)別，上述曲線重建方法的適用范圍各不相同，但是它們均難以實(shí)現(xiàn)對(duì)稠密、帶噪聲的散亂點(diǎn)集進(jìn)行有效區(qū)域劃分以重建出多條曲線。文獻(xiàn)[1]通過(guò)模擬散亂數(shù)據(jù)點(diǎn)集中各個(gè)數(shù)據(jù)點(diǎn)形成空間場(chǎng)分布的物理現(xiàn)象，有效解決了該問(wèn)題，但是還存在如下不足之處：量子化生成表示散亂點(diǎn)集場(chǎng)分布的數(shù)字圖像時(shí)運(yùn)算效率低下；數(shù)字圖像中含有大量冗余分支像素使得脊輪廓生成困難；以及重建曲線不能很好反映點(diǎn)集端點(diǎn)處的形狀特征。本文在文獻(xiàn)[1]算法的基礎(chǔ)上針對(duì)這些不足進(jìn)行了改進(jìn)，并給出實(shí)例證明了本文算法的有效性。

1 相關(guān)工作

假設(shè)給定平面散亂數(shù)據(jù)點(diǎn)集S = { pi= (xi, yi) | i = 1, …, N }，其中pi表示點(diǎn)集S中的第i個(gè)數(shù)據(jù)點(diǎn)，N表示S中數(shù)據(jù)點(diǎn)的個(gè)數(shù)(如圖1所示)，作者的目標(biāo)就是根據(jù)數(shù)據(jù)點(diǎn)集的特征擬合出一條或多條反映該點(diǎn)集形狀和走向的B樣條曲線。文獻(xiàn)[1]首先利用圖像梯度跟蹤散亂數(shù)據(jù)點(diǎn)集形成場(chǎng)分布的局部極大值點(diǎn)在平面的投影獲取點(diǎn)集的脊線(即曲線的初始逼近)，然后借助脊線將原始數(shù)據(jù)點(diǎn)集劃分成多個(gè)子集，以及對(duì)各個(gè)子集的數(shù)據(jù)點(diǎn)進(jìn)行排序，最后采用傳統(tǒng)的曲線擬合方法對(duì)已定序的數(shù)據(jù)點(diǎn)重建出相應(yīng)的曲線。

文獻(xiàn)[1]假設(shè)散亂數(shù)據(jù)點(diǎn)集 S中各個(gè)數(shù)據(jù)點(diǎn)pi都具有以其為中心的單調(diào)遞減線性場(chǎng)，并且該點(diǎn)在平面內(nèi)任意一點(diǎn)p(x, y)處的場(chǎng)強(qiáng)度由式(1)所示的逆多元二次（inverse multiquadric）徑向基函數(shù)所定義[15]

式(1)的幾何意義是數(shù)據(jù)點(diǎn) pi到平面內(nèi)任意一點(diǎn)p的反向距離，即兩點(diǎn)間歐氏距離的倒數(shù)，其中參數(shù)r是某個(gè)大于0的常量以避免被0除。不難看出，數(shù)據(jù)點(diǎn)pi在點(diǎn)p處的場(chǎng)強(qiáng)度會(huì)隨著p到pi距離越遠(yuǎn)而越小。而散亂數(shù)據(jù)點(diǎn)集S在平面內(nèi)任意一點(diǎn)p處的場(chǎng)強(qiáng)度則是點(diǎn)集中各個(gè)數(shù)據(jù)點(diǎn)pi在該點(diǎn)的場(chǎng)作用之和，由式(2)所示的場(chǎng)表面（field surface）函數(shù)表示

場(chǎng)表面函數(shù) g的局部極大值點(diǎn)構(gòu)成了其脊線，跟蹤場(chǎng)表面函數(shù)g的局部極大值點(diǎn)在xy平面的投影，就可獲取曲線的初始逼近。實(shí)驗(yàn)證明參數(shù)r還可作為獲取脊輪廓的平滑度調(diào)整因子。一般而言，隨著參數(shù)r的增大，對(duì)應(yīng)的場(chǎng)表面越趨于平滑，因此獲取脊輪廓也就越平滑。當(dāng)r取值為5時(shí)，可獲得一個(gè)較理想的平滑度。同時(shí)，注意到直接通過(guò)解析式方法求解函數(shù)g的局部極大值點(diǎn)是相當(dāng)困難的，因此文獻(xiàn)[1]給出了以下的重建算法基本步驟：

步驟1量子化生成表示散亂數(shù)據(jù)點(diǎn)集場(chǎng)分布的數(shù)字圖像（digital image），即將場(chǎng)表面g轉(zhuǎn)化成數(shù)字圖像。這里的圖像指的是一個(gè)二維數(shù)組，它的各個(gè)元素也稱為像素，是通過(guò)對(duì)顯式函數(shù)z = g(x, y)在xy平面上沿著x和y軸方向進(jìn)行等間隔劃分采樣得到的對(duì)應(yīng)像素強(qiáng)度值（如圖2所示）。數(shù)字圖像的分辨率則定義為像素的總數(shù)目。

步驟 2找到數(shù)字圖像中的次要脊（minor ridge）像素和主要脊（major ridge）像素，并且依據(jù)分支（branch）像素將這些次要脊像素劃分成各個(gè)次要脊線段（minor-ridge segment）。其中幾個(gè)重要概念的定義如下：

次要脊像素 其強(qiáng)度值高于它的兩個(gè)鄰接像素的強(qiáng)度值且這兩個(gè)鄰接像素在彼此相反的兩個(gè)方向上，即滿足式(3)所示的4個(gè)條件之一的像素便可被稱為次要脊，其中I(i, j)表示i行j列像素的強(qiáng)度值

主要脊像素 其強(qiáng)度值高于它在梯度方向（gradient direction，即強(qiáng)度值改變最大方向）上的鄰接像素的強(qiáng)度值。

次要脊線段 次要脊像素的連通集合，其中各個(gè)像素至多有兩個(gè)鄰接次要脊像素。

分支像素 兩條以上次要脊線段的公共端點(diǎn)，即它至少有3個(gè)鄰接次要脊像素。

步驟 3移除不含主要脊像素的次要脊線段，找出剩余的次要脊線段所形成的脊輪廓（ridge contour）。根據(jù)上述定義，主要脊是次要脊的子集，而且主要脊的像素點(diǎn)一般不會(huì)落在較小的、有噪聲的次要脊線段上（或稱為旁支），而是落在代表了原始數(shù)據(jù)點(diǎn)集脊線的輪廓上。脊輪廓就是由次要脊線段連接而成的集合，其中各個(gè)次要脊線段至少含有一個(gè)主要脊像素。

步驟4利用脊輪廓將原始數(shù)據(jù)點(diǎn)集劃分為多個(gè)子集，并通過(guò)跟蹤各個(gè)子集的數(shù)據(jù)點(diǎn)在對(duì)應(yīng)脊輪廓上的投影對(duì)它們進(jìn)行排序。

步驟5對(duì)于各個(gè)子集已排序的數(shù)據(jù)點(diǎn)，應(yīng)用現(xiàn)有的方法，如Piegl和Tiller在文獻(xiàn)[16]中介紹的B樣條擬合有序離散數(shù)據(jù)點(diǎn)集的方法，得到重建的B樣條曲線。

圖1 平面散亂數(shù)據(jù)點(diǎn)集

圖2 表示散亂數(shù)據(jù)點(diǎn)集場(chǎng)分布的數(shù)字圖像

文獻(xiàn)[1]提出的基于散亂數(shù)據(jù)點(diǎn)集場(chǎng)分布的曲線重建方法可基于數(shù)據(jù)點(diǎn)集的物理特征自動(dòng)完成對(duì)點(diǎn)集的區(qū)域劃分，因此適用于多條曲線的擬合，但是也存在如下問(wèn)題限制了其實(shí)用性：

（1） 在步驟1中，若數(shù)字圖像的采樣分辨率較高，或者原始數(shù)據(jù)點(diǎn)集中數(shù)據(jù)點(diǎn)的個(gè)數(shù)較多，則量子化生成數(shù)字圖像時(shí)運(yùn)算效率低下。

（2） 在步驟 2 、3中，數(shù)字圖像往往含有大量冗余分支像素，不利于生成形體較佳的脊輪廓，進(jìn)而影響后期散亂點(diǎn)的排序效果。

（3） 在步驟4中，對(duì)于開點(diǎn)集，由于點(diǎn)集端點(diǎn)附近的主要脊像素分布較少，根據(jù)步驟3，大部分端點(diǎn)處的脊像素都應(yīng)作為旁支而被刪除，導(dǎo)致獲取的脊輪廓在點(diǎn)集端點(diǎn)處收縮，因此重建出的開曲線不能很好地反映點(diǎn)集端點(diǎn)的形狀。

2 本文算法

針對(duì)上述文獻(xiàn)[1]算法的不足，本文作了以下改進(jìn)：首先，通過(guò)估計(jì)場(chǎng)強(qiáng)基函數(shù)的邊界，即各個(gè)數(shù)據(jù)點(diǎn)的場(chǎng)作用范圍，減少了大量冗余運(yùn)算從而提高了量子化生成數(shù)字圖像的效率；其次，在不改變圖像分辨率的前提下，采用圖像細(xì)化結(jié)合改進(jìn)的BFS算法解決了脊輪廓生成困難的問(wèn)題；最后，利用加權(quán)最小二乘法對(duì)開點(diǎn)集端點(diǎn)附近因脊輪廓縮短而無(wú)法有效投影到脊輪廓上的數(shù)據(jù)點(diǎn)進(jìn)行排序，進(jìn)而延長(zhǎng)了重建曲線的端點(diǎn)，改善了開曲線的重建效果。

2.1 估計(jì)場(chǎng)函數(shù)邊界提高量子化效率

量子化生成數(shù)字圖像時(shí)，一般首先將散亂數(shù)據(jù)點(diǎn)集的最小包圍盒沿著 x y平面坐標(biāo)軸方向按柵格尺寸λ進(jìn)行等間隔劃分，其中每個(gè)柵格就代表了數(shù)字圖像中的一個(gè)像素，柵格尺寸λ則決定了數(shù)字圖像的分辨率；然后取數(shù)字圖像中各個(gè)柵格的中心點(diǎn)為采樣點(diǎn)，分別代入式(2)定義的場(chǎng)表面函數(shù)g中計(jì)算出該柵格代表像素的場(chǎng)強(qiáng)度。為了確保原始數(shù)據(jù)點(diǎn)集中絕大部分?jǐn)?shù)據(jù)點(diǎn)都映射到數(shù)字圖像唯一的柵格單元中，柵格尺寸λ應(yīng)根據(jù)點(diǎn)集分布均勻性進(jìn)行設(shè)置，一般取點(diǎn)集分布密度ρ的k倍。點(diǎn)集分布密度ρ可由下列方法估計(jì)得到[17]：在原始數(shù)據(jù)點(diǎn)集中隨機(jī)取 m個(gè)數(shù)據(jù)點(diǎn)pi，計(jì)算離點(diǎn)pi最近點(diǎn)的距離di，則分布密度為，m一般取20～30。于是柵格尺寸

設(shè)散亂數(shù)據(jù)點(diǎn)集中最大、最小 x 、y坐標(biāo)分別為 xmax、ymax、xmin、ymin，同時(shí)為避免散亂點(diǎn)位于包圍盒外，給定容差ε放大包圍盒，則沿著x、y坐標(biāo)軸方向劃分的柵格數(shù)分別為其中符號(hào)表示向上取整。柵格總數(shù)a*b即為數(shù)字圖像的分辨率。為計(jì)算簡(jiǎn)便，也可以令分辨率的高和寬a=b=512。根據(jù)式(2)，要得到數(shù)字圖像中一個(gè)像素的場(chǎng)強(qiáng)度，需要計(jì)算數(shù)據(jù)點(diǎn)集中的所有數(shù)據(jù)點(diǎn)到該像素的反向距離，若點(diǎn)集中數(shù)據(jù)點(diǎn)的個(gè)數(shù)為N，則每個(gè)像素需要計(jì)算N次反向距離，于是總的計(jì)算次數(shù)為a*b*N。量子化運(yùn)算時(shí)間會(huì)隨著原始數(shù)據(jù)點(diǎn)集中數(shù)據(jù)點(diǎn)的個(gè)數(shù)或者數(shù)字圖像分辨率的增加而迅速增加，而點(diǎn)云數(shù)據(jù)通常規(guī)模較大，這樣龐大的計(jì)算量是不可接受的。

本文通過(guò)估計(jì)場(chǎng)強(qiáng)基函數(shù)的邊界減少大量冗余運(yùn)算，進(jìn)而提高了量子化生成數(shù)字圖像的效率，其重要依據(jù)在于散亂數(shù)據(jù)點(diǎn)集中各個(gè)數(shù)據(jù)點(diǎn)的場(chǎng)輻射范圍是有限的，可以近似認(rèn)為其作用范圍為距離該數(shù)據(jù)點(diǎn)半徑為R（以像素為單位）的圓形域內(nèi)，而圓形域外部的場(chǎng)強(qiáng)度幾乎衰減為0（如圖3所示）。只要圓形域半徑R足夠大，使得其外部的場(chǎng)強(qiáng)度I(R)小于一定的閾值δ，便可忽略圓形域外部區(qū)域場(chǎng)強(qiáng)度的計(jì)算，并且不會(huì)影響數(shù)字圖像的生成效果，實(shí)驗(yàn)證明δ可取0.01。因此，對(duì)于每個(gè)數(shù)據(jù)點(diǎn)，只需計(jì)算以其為中心的圓形域范圍內(nèi)πR2個(gè)像素的場(chǎng)強(qiáng)度，此時(shí)總的計(jì)算次數(shù)為πR2*N。當(dāng)R＜＜a或者R＜＜b時(shí)，可大幅度提高計(jì)算速度。由于R是通過(guò)估計(jì)場(chǎng)強(qiáng)基函數(shù)邊界得到的某個(gè)常量，于是量子化生成數(shù)字圖像的時(shí)間復(fù)雜度約為O(N)。

圖3 數(shù)據(jù)點(diǎn)的場(chǎng)輻射范圍

式(1)所示的逆多元二次徑向基函數(shù)可轉(zhuǎn)化成一元函數(shù)[R2+r2]-1/2的形式，其中表示歐幾里德(Euclidean)范數(shù)。由于該一元函數(shù)表示的場(chǎng)強(qiáng)值I(R) = [R2+r2]-1/2是關(guān)于半徑R的單調(diào)遞減函數(shù)，當(dāng)參數(shù)r取經(jīng)驗(yàn)值5時(shí)，令I(lǐng)(R)小于閾值 δ = 0.01，可計(jì)算得到 R 的最小值為100，則散亂數(shù)據(jù)點(diǎn)的反向距離場(chǎng)強(qiáng)值可表示為如式(4)所示分段的形式

如果數(shù)字圖像分辨率為512*512，則對(duì)于逆多元二次徑向基函數(shù)，通過(guò)估計(jì)場(chǎng)強(qiáng)基函數(shù)邊界可以將生成數(shù)字圖像的運(yùn)算效率提高一個(gè)量級(jí)。徑向基函數(shù)實(shí)際上是利用一個(gè)一元函數(shù)作用在歐幾里德距離上，然后做平移，因此具有物理上各向同性的性質(zhì)，適用于模擬數(shù)據(jù)點(diǎn)的場(chǎng)分布。于是也可以選擇其他衰減速度更快的徑向基函數(shù)，如高斯基函數(shù)、緊支柱正定基函數(shù)[18]等來(lái)模擬數(shù)據(jù)點(diǎn)形成的場(chǎng)，從而更加高效地實(shí)現(xiàn)數(shù)字圖像的量子化。以高斯基函數(shù)為例，其定義如式(5)所示當(dāng)式中參數(shù)σ仍取經(jīng)驗(yàn)值5，且該式表示的場(chǎng)強(qiáng)值小于閾值δ = 0.01時(shí)，可計(jì)算得到R的最小值為11。仍然取分辨率為512*512，則采用高斯徑向基函數(shù)時(shí)，運(yùn)算時(shí)間約是原始數(shù)字圖像量子化方法的πR2/ (512*512) ≈ 1/690，大幅提高了生成數(shù)字圖像的效率。

2.2 圖像細(xì)化結(jié)合改進(jìn)BFS算法去除旁支

圖1所示的散亂點(diǎn)集經(jīng)量子化后，生成的數(shù)字圖像如圖4所示：圖中左側(cè)表示的是根據(jù)式(3)定義找到的數(shù)字圖像中的次要脊像素，它們實(shí)際上就是使得場(chǎng)表面函數(shù)g具有局部極大值的那些像素點(diǎn)；右側(cè)表示的是指定局部像素區(qū)域內(nèi)各種類型的像素排列情況的放大示意，其中方形標(biāo)記表示主要脊像素、半徑較大的圓形標(biāo)記表示次要脊像素、半徑較小的圓形標(biāo)記表示分支像素。

圖4表明，初始的數(shù)字圖像不僅含有大量冗余分支像素，而且其次要脊像素排列情況具有明顯的不規(guī)則性，還帶有一些噪聲脊像素，難以生成用于散亂數(shù)據(jù)點(diǎn)集區(qū)域劃分以及各子集數(shù)據(jù)點(diǎn)排序的脊輪廓。本文采用圖像細(xì)化結(jié)合改進(jìn)的BFS算法解決這一問(wèn)題：首先，利用圖像細(xì)化算法對(duì)數(shù)字圖像進(jìn)行初步細(xì)化，抽取出只有一個(gè)像素寬度的圖像骨架，以此大幅減少冗余分支像素；然后，通過(guò)廣度優(yōu)先遍歷（BFS）數(shù)字圖像中的所有次要脊像素，將數(shù)字圖像轉(zhuǎn)化為廣度優(yōu)先生成森林；最后，對(duì)生成森林中各棵生成樹進(jìn)行剪支運(yùn)算，去除不含主要脊像素的次要脊線段（即旁支）以獲取所需的脊輪廓。

圖4 數(shù)字圖像中各種類型的像素

（1） 利用圖像細(xì)化初步細(xì)化

一般說(shuō)來(lái)，柵格尺寸λ的設(shè)定對(duì)于生成數(shù)字圖像的質(zhì)量和精度有直接影響，尺寸越小數(shù)字圖像越精細(xì)，反之則越粗糙。在實(shí)際應(yīng)用中，為保證重建曲線的精度，有必要適當(dāng)減小柵格尺寸以增大數(shù)字圖像的分辨率，然而實(shí)驗(yàn)證明較高的圖像分辨率往往會(huì)形成大量的冗余分支像素（如圖4所示）。這些冗余分支像素也滿足至少含有 3個(gè)鄰接次要脊像素的基本條件，但是不同于作為兩條以上次要脊線段公共端點(diǎn)的正常分支像素，它們是由于許多次要脊像素并列排列，從而導(dǎo)致超出一個(gè)像素的線條寬度所造成的，因此并不能作為將數(shù)字圖像中次要脊像素劃分成各個(gè)次要脊線段的依據(jù)。作者考慮采用圖像細(xì)化算法對(duì)數(shù)字圖像進(jìn)行初步細(xì)化，提取出數(shù)字圖像的線形骨架以剔除冗余分支像素。

圖像細(xì)化算法種類較多，代表性的有Hilditch細(xì)化算法[19]、Pavlidis細(xì)化算法[20]、Rosenfeld細(xì)化算法[21]、索引表細(xì)化算法[22]等。細(xì)化算法一般遵循以下準(zhǔn)則：圖像骨架必須保持原圖像的連通性和拓?fù)浣Y(jié)構(gòu)；圖像骨架盡可能是原圖像的中心線；細(xì)化的結(jié)果應(yīng)是一個(gè)像素寬度的線條圖像，即非特殊點(diǎn)（交點(diǎn)和端點(diǎn)）在8連通意義下只有兩個(gè)鄰接像素等。根據(jù)數(shù)字圖像中次要脊像素排列特征的不同可以選擇其中任意一種算法進(jìn)行初步細(xì)化，本文采用了文獻(xiàn)[22]介紹的索引表細(xì)化算法。

本文采用的索引表細(xì)化算法是在應(yīng)用上述細(xì)化準(zhǔn)則的基礎(chǔ)上，依據(jù)次要脊像素p的8個(gè)鄰接像素是否也為次要脊像素，構(gòu)造一張包含256種情況的剔除表，水平、豎直方向交替掃描數(shù)字圖像中各行各列的次要脊像素，通過(guò)查表即可判斷該次要脊像素是否應(yīng)該被剔除。以次要脊像素為中心像素的3*3的像素模板為例，若該次要脊像素為內(nèi)部點(diǎn)、直線端點(diǎn)，或除去后增加連通分量則不能被刪除。圖5給出了次要脊像素p八鄰域情況的幾個(gè)示例：(a)中心像素點(diǎn)為一個(gè)內(nèi)部點(diǎn)，不能被刪除，否則骨架被掏空；同樣道理，(b)中心像素點(diǎn)也不能被刪除；(c)中心像素點(diǎn)不是骨架，可以被刪除；(d)中心像素點(diǎn)不能被刪除，否則破壞圖像的連通性；(e)中心像素點(diǎn)也不是骨架，可以被刪除；(f)中心像素點(diǎn)為直線的端點(diǎn)，不能被刪除，否則整條直線都將被刪除；(g)中心像素點(diǎn)為孤立點(diǎn)，一般作為噪聲可以被刪除。利用索引表細(xì)化算法初步細(xì)化后，獲得了理想的結(jié)果，無(wú)論是冗余分支像素，還是數(shù)字圖像中的噪聲都大大減少了（如圖 6 所示）。不妨設(shè)數(shù)字圖像中的像素總數(shù)為n，圖像中包含的次要脊像素個(gè)數(shù)為 m，則該細(xì)化算法時(shí)間復(fù)雜度約為O(nlogm)。

圖5 次要脊像素八鄰域情況

圖6 初步細(xì)化后的數(shù)字圖像

（2） 利用BFS生成森林去除旁支

通過(guò)上面的步驟，大量冗余分支像素已被剔除，接下來(lái)可采用迭代的方式去除數(shù)字圖像中的次要脊旁支。旁支又可分為普通分支（branches）和環(huán)（loops）兩種類型。普通分支的一個(gè)端點(diǎn)是分支像素，另一個(gè)端點(diǎn)則是只有一個(gè)鄰接像素的次要脊像素，因此它是開的。與此相反，環(huán)則是封閉的，它只有分支像素作為其端點(diǎn)，如圖6右側(cè)的局部放大示意所示。為了有效識(shí)別這兩種類型的旁支，并且統(tǒng)一處理次要脊像素排列的各種復(fù)雜情況，本文采用了改進(jìn)的BFS算法。

改進(jìn)的 BFS算法基本思想是通過(guò)廣度優(yōu)先遍歷數(shù)字圖像中的所有次要脊像素，將數(shù)字圖像轉(zhuǎn)化為多棵以次要脊像素為結(jié)點(diǎn)、相鄰次要脊像素的鄰接關(guān)系為邊的廣度優(yōu)先生成樹所構(gòu)成的廣度優(yōu)先生成森林。具體轉(zhuǎn)化方式如下：逐行逐列掃描數(shù)字圖像中的像素，若遇到第一個(gè)未被訪問(wèn)的次要脊像素 s ，就將其作為初始結(jié)點(diǎn)（生成樹的根）插入到搜索隊(duì)列中，然后依次移出隊(duì)首結(jié)點(diǎn)u，并將隊(duì)首結(jié)點(diǎn)u所有未訪問(wèn)過(guò)的鄰接次要脊像素v加入搜索隊(duì)列，同時(shí)將v作為新結(jié)點(diǎn)、鄰接關(guān)系（u, v）作為新邊插入到本次循環(huán)創(chuàng)建的生成樹中，直到隊(duì)列為空循環(huán)搜索過(guò)程才結(jié)束。一方面，由于生成森林具有不含環(huán)（或稱回路）的特性，因此環(huán)類型的旁支在計(jì)算得到生成森林后也被轉(zhuǎn)化為普通分支類型，有利于對(duì)這兩種類型的旁支進(jìn)行統(tǒng)一處理；另一方面，由于非連通的生成森林是由多棵連通的生成樹所構(gòu)成，每棵生成樹就代表了數(shù)字圖像中次要脊像素形成的一個(gè)連通分量，因此在將數(shù)字圖像轉(zhuǎn)化為生成森林后，實(shí)際上就完成了對(duì)數(shù)字圖像中的次要脊像素的分塊工作，有助于后期根據(jù)點(diǎn)集的拓?fù)浣Y(jié)構(gòu)特征生成多條脊輪廓，以及順利解決點(diǎn)集的區(qū)域劃分問(wèn)題。將數(shù)字圖像轉(zhuǎn)化為生成森林的算法描述如下：

輸入細(xì)化后的數(shù)字圖像I

輸出以次要脊像素為結(jié)點(diǎn)、相鄰次要脊像素的鄰接關(guān)系為邊的廣度優(yōu)先生成森林Gπ

GENERATE-BFS-FOREST(I )

for each pixel p of image I // 逐行逐列掃描數(shù)字圖像I中的像素p

do if p is minor ridge // 如果p為次要脊像素，初始化p的訪問(wèn)標(biāo)記為WHITE，即未訪問(wèn)

then color[p]←WHITE

create a BFS forest Gπ// 創(chuàng)建廣度優(yōu)先生成森林Gπ

for each pixel s of image I // 逐行逐列掃描數(shù)字圖像I中的像素s

do if s is minor ridge and color[s]=WHITE // 如果s是次要脊像素且未被訪問(wèn)

then create a BFS tree T // 創(chuàng)建一棵廣度優(yōu)先生成樹T

color[s]←BLACK // 設(shè)置 s 的訪問(wèn)標(biāo)記為BLACK，表示已被訪問(wèn)

insert pixel s as a vertex into the tree T // 將像素s作為結(jié)點(diǎn)插入生成樹T中

root[T ]←s // 以初始結(jié)點(diǎn)s作為生成樹T的根

Q←Ф // 初始化隊(duì)列Q為空

ENQUEUE(Q,s) // 將像素s插入Q隊(duì)尾

while Q≠Ф // 若隊(duì)列Q非空，則執(zhí)行while循環(huán)；否則，不執(zhí)行

do u←DEQUEUE(Q) // 移出Q隊(duì)首u(yù)

for each v in u’s eight adjoint pixels // 逐個(gè)查找u的8個(gè)鄰接像素v

do if v is minor ridge and color[v]=WHITE // 如果v是次要脊像素且未被訪問(wèn)

then color[v]←BLACK // 設(shè)置v的訪問(wèn)標(biāo)記為BLACK

insert pixel v as a vertex into the tree T // 將像素v作為結(jié)點(diǎn)插入生成樹T中

insert edge(u, v) into the tree T // 將(u, v)作為邊添加進(jìn)生成樹T中

ENQUEUE(Q, v) // 將像素v插入Q隊(duì)尾

add tree T to the forest Gπ// 將本次循環(huán)得到的生成樹T添加到生成森林Gπ中

return Gπ// 返回生成森林Gπ

當(dāng)數(shù)字圖像中所有的次要脊像素被訪問(wèn)并作為結(jié)點(diǎn)插入生成森林的各棵生成樹后，建生成森林的過(guò)程結(jié)束。由于設(shè)置了訪問(wèn)標(biāo)記，因此保證了每個(gè)次要脊像素至多只有一次入隊(duì)和出隊(duì)。入隊(duì)和出隊(duì)操作需要O(1)的時(shí)間，則隊(duì)列操作所需的全部時(shí)間為 O(m)。而只有當(dāng)每個(gè)次要脊像素出隊(duì)時(shí)，才會(huì)查找其8個(gè)鄰接像素，故查找鄰接像素花費(fèi)時(shí)間也為 O(m)。此外，由于初始化及建生成森林均需遍歷一次數(shù)字圖像，相應(yīng)的時(shí)間復(fù)雜度為O(n)。于是，建生成森林總的運(yùn)行時(shí)間為O(m+n)，通常m＜＜n。

建完生成森林，接下來(lái)需要對(duì)生成森林的各棵生成樹進(jìn)行剪支運(yùn)算。從葉子結(jié)點(diǎn)開始向樹根方向搜索，記錄搜索路徑path上經(jīng)過(guò)的各個(gè)次要脊像素并設(shè)置訪問(wèn)標(biāo)記，如果在path路徑上遇到第一個(gè)分支像素，或者通過(guò)path路徑無(wú)法找到新的未訪問(wèn)的結(jié)點(diǎn)，則此次搜索結(jié)束。然后根據(jù)此次搜索的path路徑上是否含有主要脊像素，判斷是否應(yīng)將該路徑上的像素作為旁支而從生成樹中予以摘除。對(duì)于每棵生成樹，剪支算法描述如下：

輸入 生成森林中各棵生成樹 T=(V，E)，其中 V是T中結(jié)點(diǎn)的集合，E是T中邊的集合

輸出 無(wú)

REMOVE-BRANCHES(T )

for each vertex p∈V [T ] // 逐個(gè)遍歷生成樹T中的結(jié)點(diǎn)p，初始化p的訪問(wèn)標(biāo)記為WHITE

do color[p]←WHITE

for each vertex u∈V [T ] // 逐個(gè)遍歷生成樹T中的結(jié)點(diǎn)u

do if u is a leaf and color[u]=WHITE // 如果u是葉子結(jié)點(diǎn)(即度為1，此處度表示鄰接邊數(shù))且未被訪問(wèn)

then make a path // 開始記錄搜索路徑path

repeat

color[u]←BLACK // 設(shè)置 u的訪問(wèn)標(biāo)記為BLACK

add u to the path // 將u添加到搜索路徑path中

for each edge(u,v)∈E[T ] // 逐條查找生成樹T的邊(u,v)

do if color[v]=WHITE and v is not branch // 如果邊的另一個(gè)端點(diǎn)v未被訪問(wèn)且不是分支像素

then u←v // 更新u為v，跳出for循環(huán)

break

until color[u]=BLACK // 若找到新結(jié)點(diǎn)u，則繼續(xù)執(zhí)行repeat循環(huán)；否則，不執(zhí)行

if path is not empty and doesn’t contain major ridge // 如果搜索路徑path非空且不存在主要脊像素

then remove vertexes of path from the tree T // 將path路徑上的像素從生成樹T中刪除

return

剪支運(yùn)算中，若令生成樹中的結(jié)點(diǎn)數(shù)為V，邊數(shù)為E，則初始化和尋找葉結(jié)點(diǎn)的運(yùn)行時(shí)間為O(V )，搜索路徑的運(yùn)行時(shí)間為O(E)，總的時(shí)間復(fù)雜度為O(V+E)。由于結(jié)點(diǎn)數(shù)V也表示生成樹代表的數(shù)字圖像中某個(gè)連通分量所包含次要脊像素的個(gè)數(shù)，因此V一般遠(yuǎn)小于數(shù)字圖像中的像素個(gè)數(shù)n，并且對(duì)樹形結(jié)構(gòu)而言，邊數(shù)E = V-1，則O(V+E ) = O(2V-1) = O(V )。由于非連通的生成森林是由多棵連通的生成樹所構(gòu)成，因此生成森林中的結(jié)點(diǎn)個(gè)數(shù)就為數(shù)字圖像中次要脊像素的個(gè)數(shù)m，則所有生成樹的剪支運(yùn)算時(shí)間復(fù)雜度應(yīng)為 O(m)，可見剪支運(yùn)算的開銷是較低的。綜上所述，BFS算法去除旁支總的時(shí)間復(fù)雜度為O(m+n)+ O(m)=O(m+n)。本文利用圖像細(xì)化結(jié)合改進(jìn)的 BFS算法去除旁支后得到了形體較好的脊輪廓（如圖7所示）。

圖7 去除旁支后生成的脊輪廓

對(duì)于幾何拓?fù)浣Y(jié)構(gòu)較為復(fù)雜的多連通數(shù)據(jù)點(diǎn)集，通過(guò)上述改進(jìn)的BFS算法可獲得多棵代表該點(diǎn)集拓?fù)浣Y(jié)構(gòu)的生成樹，經(jīng)剪支運(yùn)算后生成多條對(duì)應(yīng)的脊輪廓。這些脊輪廓被用于對(duì)散亂數(shù)據(jù)點(diǎn)集進(jìn)行區(qū)域劃分以及各個(gè)子集的數(shù)據(jù)點(diǎn)排序：首先，將原始數(shù)據(jù)點(diǎn)集中的每個(gè)數(shù)據(jù)點(diǎn)pi賦給離它最近的一條脊輪廓，以此將點(diǎn)集劃分為各個(gè)子集；然后，確定各子集的數(shù)據(jù)點(diǎn)到代表該子集的脊輪廓中的最近像素（或稱之為投影）；最后，跟蹤輪廓上的脊像素，以投影像素被訪問(wèn)的順序?qū)ιy點(diǎn)進(jìn)行排序。

2.3 加權(quán)最小二乘法延長(zhǎng)重建曲線

將平面無(wú)序散亂點(diǎn)集轉(zhuǎn)化為各個(gè)有序的離散點(diǎn)子集后，就可采用文獻(xiàn)[1]介紹的參數(shù)化方法和Piegl和Tiller在文獻(xiàn)[16]中介紹的B樣條最小二乘逼近方法將每個(gè)離散點(diǎn)子集擬合成滿足一定擬合誤差限制的B樣條曲線。

由于分布在開點(diǎn)集端點(diǎn)附近的主要脊像素?cái)?shù)目較少，因此根據(jù)Goshtasby曲線重建方法的步驟3，大部分端點(diǎn)處的脊像素都應(yīng)作為旁支而被刪除，從而導(dǎo)致生成的脊輪廓在點(diǎn)集端點(diǎn)處收縮，因此重建出的開曲線無(wú)法很好地反映點(diǎn)集端點(diǎn)處的形狀（如圖8所示）。

圖8 開點(diǎn)集重建曲線端點(diǎn)處收縮

對(duì)開點(diǎn)集進(jìn)行曲線重建時(shí)，縮短的脊輪廓會(huì)使得點(diǎn)集端點(diǎn)附近沿著脊輪廓的延長(zhǎng)線方向上的散亂點(diǎn)都投影到脊輪廓的端點(diǎn)（像素）上，因此，僅通過(guò)跟蹤脊輪廓上的投影無(wú)法對(duì)點(diǎn)集端點(diǎn)附近的數(shù)據(jù)點(diǎn)進(jìn)行有效排序，進(jìn)而造成重建出的開曲線在端點(diǎn)處縮短。由于在對(duì)散亂數(shù)據(jù)點(diǎn)集進(jìn)行區(qū)域劃分時(shí)，已經(jīng)計(jì)算過(guò)散亂數(shù)據(jù)點(diǎn)集的各個(gè)數(shù)據(jù)點(diǎn)到相應(yīng)脊輪廓的投影，可找出投影到脊輪廓端點(diǎn)的數(shù)據(jù)點(diǎn)形成的端點(diǎn)鄰域點(diǎn)集，并對(duì)其進(jìn)行局部擬合得到一條局部回歸線，通過(guò)對(duì)端點(diǎn)鄰域點(diǎn)集中的數(shù)據(jù)點(diǎn)到該局部回歸線的投影對(duì)其進(jìn)行排序，進(jìn)而延長(zhǎng)重建曲線的端點(diǎn)。記脊輪廓端點(diǎn)為p*，投影到p*的某個(gè)鄰域范圍內(nèi)的數(shù)據(jù)點(diǎn)pj(xj, yj)構(gòu)成的點(diǎn)集為A，同時(shí)不妨設(shè)A中數(shù)據(jù)點(diǎn)的個(gè)數(shù)為k，顯然A?S且k＜＜N。考慮先采用加權(quán)最小二乘法（WLS）對(duì)端點(diǎn)鄰域點(diǎn)集A中的所有數(shù)據(jù)點(diǎn)擬合出一條關(guān)于點(diǎn)p*的局部回歸線L*：y = ax + b，即極小化式(6)所示的二次函數(shù)

圖9 延長(zhǎng)端點(diǎn)后的重建曲線

3 實(shí)例與分析

作者在微機(jī)平臺(tái)（P4 3.06G CPU，512M RAM）上用C++語(yǔ)言實(shí)現(xiàn)了本文提出的算法，下面給出應(yīng)用該算法重建曲線的幾個(gè)實(shí)例。圖 10所示的是開點(diǎn)集的重建曲線（包括(a)波浪線和(b)螺旋線），可以看出重建曲線較好地反映了點(diǎn)集的形狀，保持了尖點(diǎn)和端點(diǎn)的特征，有效避免了開曲線在點(diǎn)集端點(diǎn)處收縮退化的現(xiàn)象。圖 1 1所示的是封閉點(diǎn)集的重建曲線，對(duì)于封閉點(diǎn)集，利用 BFS生成樹去除旁支后得到的脊輪廓并非形成封閉的環(huán)（生成樹的基本性質(zhì)），但是由于脊輪廓的兩個(gè)端點(diǎn)十分接近，甚至很可能是彼此鄰接的，因此，可以根據(jù)其兩個(gè)端點(diǎn)之間距離是否小于給定閾值來(lái)判斷點(diǎn)集是否封閉，進(jìn)而重建出所需的封閉曲線。圖 1 2所示的是多連通點(diǎn)集的重建曲線，利用多連通點(diǎn)集的場(chǎng)分布也具有多連通性的特點(diǎn)，可實(shí)現(xiàn)對(duì)點(diǎn)集的區(qū)域劃分，然后針對(duì)各個(gè)連通子集數(shù)據(jù)點(diǎn)的開閉特征重建出相應(yīng)的曲線。以上實(shí)驗(yàn)表明，對(duì)于稠密且?guī)г肼暤钠矫嫔y點(diǎn)集，本文算法可以有效重建出反映點(diǎn)集特征和形狀的B樣條曲線。

圖10 開點(diǎn)集的重建曲線

圖11 封閉點(diǎn)集的重建曲線

圖12 多連通點(diǎn)集的重建曲線

表1分別給出了幾個(gè)曲線重建實(shí)例的擬合誤差和運(yùn)行時(shí)間。其中，Error_Ave和Error_Max分別表示B樣條重建曲線的平均誤差和最大誤差，其定義如式(7)所示

其中 di是各子集中數(shù)據(jù)點(diǎn)qi到對(duì)應(yīng)的B樣條重建曲線上參數(shù)點(diǎn)C）的投影距離。

對(duì)各種類型的點(diǎn)集，初始化生成數(shù)字圖像的條件一致，即其分辨率都是根據(jù)點(diǎn)集分布密度進(jìn)行自適應(yīng)選擇。對(duì)于本文提出的改進(jìn)步驟，步驟1是用于提高效率的，步驟 2和步驟 3相比Goshtasby的算法會(huì)消耗一定的時(shí)間。從表1可以看出，步驟2和步驟3所消耗的時(shí)間大概占總重建時(shí)間的 15%～30%，時(shí)間消耗在可以接受的范圍內(nèi)，但曲線重建效果得到了顯著改善。另外，值得注意的是，本文提出的曲線重建算法的主要目的是追蹤散亂點(diǎn)集的脊線，所得到的擬合曲線也只是一次擬合的結(jié)果，并沒有進(jìn)行迭代計(jì)算以減小誤差以及對(duì)曲線進(jìn)行光順處理。如果加上這些附加操作，則本文算法提出的兩個(gè)步驟消耗的相對(duì)時(shí)間將會(huì)更少。

表1 本文各個(gè)重建實(shí)例擬合誤差和運(yùn)行時(shí)間的比較

4 總 結(jié)

本文在文獻(xiàn)[1]算法的基礎(chǔ)上，提出了一種更加實(shí)用的平面散亂點(diǎn)集的B樣條曲線重建算法。相比文獻(xiàn)[1]的算法，本文算法通過(guò)估計(jì)場(chǎng)強(qiáng)基函數(shù)的邊界，提高了量子化生成數(shù)字圖像的運(yùn)算效率；通過(guò)引入圖像細(xì)化以及 BFS生成樹去除旁支，解決了冗余分支像素導(dǎo)致脊輪廓生成困難的問(wèn)題；通過(guò)對(duì)開點(diǎn)集端點(diǎn)附近沿著脊輪廓延長(zhǎng)線方向上的散亂點(diǎn)排序，改進(jìn)了文獻(xiàn)[1]在重建開曲線時(shí)端點(diǎn)處縮短的缺陷。如何將本文的重建方法推廣到三維的情形，進(jìn)而解決三維空間曲線的擬合問(wèn)題，將是今后進(jìn)一步研究的內(nèi)容。

[1] Goshtasby A A. Grouping and parameterizing irregularly spaced points for curve fitting [J]. ACM Transaction on Graphics, 2000, 19(3): 185-203.

[2] Lee I K. Curve reconstruction from unorganized points [J]. Computer Aided Geometric Design, 2000, 17(2): 161-177.

[3] Pottmann H, Randrup T. Rotational and helical surface approximation for reverse engineering [J]. Computing, 1998, 60(4): 307-322.

[4] Hastie T, Stuetzle W. Principal curves [J]. Statistical Association, 1989, 84(13): 502-516.

[5] Lee E T Y. Choosing nodes in parametric curve interpolation [J]. Computer Aided Design, 1989, 21(6): 363-370.

[6] Hoschek J, Lasser D. Fundamentals of computer aided geometric design [M]. Wellesley, MA: A.K. Peters, 1993. 118-212.

[7] Farin G. Curves and surfaces for computer aided geometric design: a practical guide (5th edition) [M]. NewYork: Academic Press, 2002. 161-166.

[8] Levin D. The approximation power of moving least-squares [J]. Mathematics of Computation, 1998, 67(224): 1517-1531.

[9] Fang L, Gossard D C. Multidimensional curve fitting to unorganized data points by nonlinear minimization [J]. Computer Aided Design, 1995, 27(1): 48-58.

[10] Taubin G, Rondfard R. Implicit simplicial models for adaptive curve reconstruction [J]. IEEE Transactions on Pattern Analysis and Machine Intelligence, 1996, 18(3): 321-325.

[11] Amenta N, Bern M, Eppstein D. The crust and the β-skeleton: combinatorial curve reconstruction [J]. Graphical Models and Image Processing, 1998, 60(2): 125-135.

[12] 鐘 綱, 楊勛年, 汪國(guó)昭. 基于場(chǎng)表示的平面無(wú)序點(diǎn)集曲線重建算法[J]. 計(jì)算輔助設(shè)計(jì)與圖形學(xué)學(xué)報(bào), 2002, 14(11): 1074-1079.

[13] 鐘 綱, 楊勛年, 汪國(guó)昭. 平面無(wú)序點(diǎn)集曲線重建的跟蹤 算 法 [ J]. 軟件學(xué)報(bào), 2002, 13(11): 2188-2193.

[14] Lin H W, Chen W, Wang G J. Curve reconstruction based on an interval B-spline curve [J]. The Visual Computer, 2005, 21(6): 418-427.

[15] Hardy R L. Theory and applications of the multiquadrics-biharmonic method [J]. Applied Mathematics and Computation, 1990, 19(8/9): 163-208.

[16] Piegl L, Tiller W. The NURBS book (2nd ed.)[M]. New York: Springer-Verlag, 1997. 405-413.

[17] 柯映林. 反求工程CAD建模理論、方法和系統(tǒng)[M].北京: 機(jī)械工業(yè)出版社, 2005. 16-17.

[18] 吳宗敏. 徑向基函數(shù)、散亂數(shù)據(jù)擬合與無(wú)網(wǎng)格偏微分方程數(shù)值解[J]. 工程數(shù)學(xué)學(xué)報(bào), 2002, 19(2): 1-12.

[19] Hilditch C J. Comparison of thinning algorithms on a parallel processor [J]. Image Vision Compute, 1983, 1(3): 115-132.

[20] Pavdis T A. Thinning algorithm for discrete binary image [J]. Computer Vision, Graphics and Image Processing, 1980, 13(2): 142-157.

[21] Rosenfeld A. A characterization of parallel thinning algorithms [J]. Information and Control, 1975, 29(3): 286-291.

[22] 張宏林. Visual C++數(shù)字圖像模式識(shí)別技術(shù)及工程實(shí)踐[M]. 北京: 人民郵電出版社, 2003. 235-251.

B-Spline Curve Reconstruction from Planar Unorganized Points Based on Field Distribution

HUANG Tong-xin1,3,4, WANG Wen-ke1,2,3,4, ZHANG Hui1,3,4, SONG Zheng-xuan1,3,4
( 1. School of Software, Tsinghua University, Beijing 100084, China; 2. Department of Computer Science and Technology, Tsinghua University, Beijing 100084, China; 3. Key Laboratory for Information System Security, Ministry of Education, Beijing 100084, China; 4. Tsinghua National Laboratory for Information Science and Technology, Beijing 100084, China )

Curve reconstruction from planar unorganized points is one of the most important problems in reverse engineering. A practical B-spline curve fitting algorithm based on Goshtasby’s approach is presented. The digital image representing field distribution by estimating the bound of field strength basis function is generated at first, and then an algorithm of image thinning associated with improved BFS is proposed to overcome difficulties of obtaining the ridge contour under the situation of redundant branch pixels, and finally the weighted least squares method is used to extend the reconstructed curve for overcoming the deficiency that reconstructed curve may be shortened. Experiments show that the algorithm is valid and practical for B-spline curve reconstruction, especially when the given points are dense and noisy.

computer application; B-spline curve reconstruction; field distribution; unorganized points

TP 391

A

1003-0158(2010)02-0073-11

2008-08-16

國(guó)家“973”計(jì)劃資助項(xiàng)目（2004CB719404）；國(guó)家自然科學(xué)基金重點(diǎn)資助項(xiàng)目（60635020）；國(guó)家“863”重點(diǎn)資助項(xiàng)目（2007AA040401）

黃童心（1984-），男，湖北咸寧人，碩士研究生，主要研究方向?yàn)橛?jì)算機(jī)輔助幾何設(shè)計(jì)。