一種新的等值線繪制方法研究

爨 瑩，朱航洲

（西安石油大學(xué) 計(jì)算機(jī)學(xué)院，陜西 西安 710065）

等值線圖是在石油勘探開發(fā)、采礦、地質(zhì)、地球物理、地球化學(xué)和氣象等工程和技術(shù)領(lǐng)域內(nèi)應(yīng)用極廣的一種圖形，是眾多領(lǐng)域成果表示的重要圖件之一。在等值線繪制過程中，首先要對(duì)散點(diǎn)數(shù)據(jù)進(jìn)行網(wǎng)格化，而常用的數(shù)據(jù)網(wǎng)格化方法（如距離加權(quán)法、方位法、趨勢(shì)面法、加權(quán)最小二乘法、疊加法）的精度和有效性都不高；Kriging估值技術(shù)雖然比傳統(tǒng)的估值技術(shù)具有更高的先進(jìn)性和有效性，但這種方法目前還沒有得到廣泛應(yīng)用，許多細(xì)節(jié)問題還處于研究階段。本文介紹利用Surfer軟件提供的可編程方法對(duì)散點(diǎn)數(shù)據(jù)進(jìn)行網(wǎng)格化，然后用矩形追蹤法對(duì)網(wǎng)格化后的數(shù)據(jù)進(jìn)行追蹤，并實(shí)現(xiàn)了放大、縮小、漫游、光滑等功能，繪制出的等值線具有更強(qiáng)的實(shí)用性。

1 繪制等值線圖原理

1.1 離散數(shù)據(jù)網(wǎng)格化

為了繪制等值線圖，首先需要把離散分布的數(shù)據(jù)點(diǎn)網(wǎng)格化。為此，要建立數(shù)字高程模型（DEM），其主要功能就是將離散數(shù)據(jù)網(wǎng)格化。其主要步驟為：首先，由原始數(shù)據(jù)點(diǎn)的橫、縱坐標(biāo)的最小值和最大值構(gòu)成矩形網(wǎng)格的邊界；然后，根據(jù)用戶需要將矩形再劃分為m×n個(gè)小矩形，即形成矩形網(wǎng)格；最后，用插值算法得到每個(gè)網(wǎng)格點(diǎn)的高程值。

DEM內(nèi)插有多種算法，常用的有距離加權(quán)法、樣條函數(shù)內(nèi)插、最小二乘配置法和Kriging方法等。本文使用Suffer軟件中提供的可編程方法對(duì)散點(diǎn)數(shù)據(jù)網(wǎng)格化，其具體算法如下。

（1）輸入?yún)?shù)；

（2）根據(jù)輸入?yún)?shù)從數(shù)據(jù)庫(kù)提取 X、Y坐標(biāo)及高程值Z，并保存到文本數(shù)據(jù)文件；

（3）在程序中引用 surfer并創(chuàng)建 surfer（srf）和 grid（grd）對(duì)象；

（4）用srf對(duì)象的GridData方法將文本文件轉(zhuǎn)換成surfer格式文件（.grd文件）；

（5）用grd對(duì)象的GetNode方法獲取 grd文件里生成的網(wǎng)格點(diǎn)高程值，寫入網(wǎng)格化后的數(shù)據(jù)文件。

1.2 等值點(diǎn)的計(jì)算

令任意一個(gè)網(wǎng)格邊兩端的數(shù)據(jù)值分別為 Z1、Z2，如圖 1所示，對(duì)于任意一條等值線 W，它在該邊上 Z（x，y）-W=0，即等值線所通過的點(diǎn)，可以有0～3個(gè)根。

當(dāng)?shù)戎稻€恰好等于網(wǎng)格節(jié)點(diǎn)的值時(shí)，這種網(wǎng)格點(diǎn)被稱作退化點(diǎn)，當(dāng)遇到退化點(diǎn)時(shí)，在網(wǎng)格點(diǎn)值上加上一個(gè)非常小的正數(shù)，以消除退化點(diǎn)[3]。

如圖1中的Z1、Z2值都不會(huì)等于等值線 W 的值，這時(shí)，網(wǎng)格邊上對(duì)于某一條等值線 W，存在根Z（x，y）-W=0的可能性有：

（1）Z1＞W(wǎng)，Z2＞W(wǎng) 或 Z1＜W，Z2＜W 時(shí)，沒有根，如圖 1中的（a）和（b），或存在兩個(gè)根，如圖 1 的（e）和（f）。

（2）Z1＜W，Z2＞W(wǎng) 或 Z1＞W(wǎng)，Z2＜W 時(shí)，存在一個(gè)根，如圖 1中的（c），或存在兩個(gè)根，如圖 1中的（d），其中 A根是重根，如圖 1中的（g）。

對(duì)于立方函數(shù)，利用數(shù)值法求根。這時(shí)可以將網(wǎng)格邊長(zhǎng)nx和ny再分成等間距的細(xì)分段，如圖1中的（d）所示，然后求出每個(gè)細(xì)分段兩端的函數(shù)值Zk和Zk+1。當(dāng)Zk或Zk+1等于等值線值W時(shí)，同樣需要加上一個(gè)非常小的正數(shù)，以消除退化點(diǎn)。這樣，當(dāng)（Zk-W）（Zk+1-W）＜0時(shí)，說明該細(xì)分段上存在一個(gè)根。該細(xì)分段的長(zhǎng)度可以小到以機(jī)器步長(zhǎng)來計(jì)算，最小可以到描繪儀機(jī)器步長(zhǎng)的2倍，這樣可以保證在曲面函數(shù)上追蹤所得到的等值線連續(xù)光滑而不會(huì)發(fā)生方向上的突變。在找到根值的細(xì)分段上，再以Zk和Zk+1作線性內(nèi)插求得所在根及等值點(diǎn)的位置。

數(shù)值法找根時(shí)有時(shí)會(huì)失根，尤其對(duì)于重根情況，如圖 1（d）中的 A根，由于它兩端的 Zk和 Zk+1都小于 W，因此（Zk-W）（Zk+1-W）＞0，A 根就會(huì)失去，但這對(duì)于繪制等值線并無多大妨礙，至多漏繪一些極小的封閉等值線或等值線應(yīng)該封閉而不完全封閉。

1.3 等值線的追蹤

把網(wǎng)格縱邊和橫邊的細(xì)分段所組成的小網(wǎng)格稱為單元，它們的4個(gè)節(jié)點(diǎn)稱為單元節(jié)點(diǎn)，以區(qū)別于原始的網(wǎng)格與網(wǎng)格節(jié)點(diǎn)。這樣，對(duì)于每一個(gè)單元在其單元節(jié)點(diǎn)上求出函數(shù)值后（若遇到有退化點(diǎn)時(shí)，同樣需要通過加一個(gè)很小的正數(shù)來處理），有圖2所示的8種連接等值線的可能（其中+號(hào)表示節(jié)點(diǎn)值大于等值線值，-號(hào)表示節(jié)點(diǎn)值小于等值線值）。

（1）單元節(jié)點(diǎn)上函數(shù)值同時(shí)大于或同時(shí)小于等值線值時(shí)，該單元內(nèi)沒有等值線通過，如圖2（a）所示。

（2）單元節(jié)點(diǎn)上函數(shù)值大于和小于等值線時(shí)，單元邊上要么存在兩個(gè)根，要么存在4個(gè)根。存在兩個(gè)根時(shí)，兩根在對(duì)邊上，如圖 2（b）、圖 2（c）所示，兩根在鄰邊上，如圖 2（d）、圖 2（e）、圖 2（g）、圖 2（h）所示；存在 4個(gè)根時(shí)，4個(gè)根分別在單元的4個(gè)邊上。對(duì)于兩個(gè)根的情況，只要找到根的位置就可以相連；對(duì)于4個(gè)根的情況，則需規(guī)定底邊根與左邊根相連，右邊根與上邊根相連或底邊根與右邊根相連，右邊根與上邊根相連，如圖2（f）所示，這樣就能保證等值線間不會(huì)出現(xiàn)相交現(xiàn)象[4]。

圖3表示在一個(gè)矩形網(wǎng)格（i，j）內(nèi)，將縱邊 ny等分成6個(gè)細(xì)分段，將橫邊nx等分成5個(gè)細(xì)分段的連法示意圖。在等分后的每個(gè)單元節(jié)點(diǎn)上算出函數(shù)值Z。Z＞W(wǎng)（W 為等值線值）時(shí)，以“+”表示；Z＜W 時(shí)，以“-”表示。這樣在每個(gè)“+”、“-”節(jié)點(diǎn)之間必定存在一個(gè)等值線根，可用線性內(nèi)插求出它的位置，連成的等值線如圖3中的折線所示。由于網(wǎng)格邊上的函數(shù)是立方函數(shù)，它在網(wǎng)格邊上最多有 3個(gè)根，如圖 3中，底邊上有 3個(gè)根 A1、A2、A3。在理論上只要對(duì)全部網(wǎng)格內(nèi)的單元節(jié)點(diǎn)掃描一遍，把求得的根值用圖2中的規(guī)則相連即可畫出等值線圖。

但是這樣計(jì)算會(huì)花費(fèi)計(jì)算機(jī)大量時(shí)間，而且大部分單元內(nèi)是不會(huì)有等值線通過的。因此，需要設(shè)計(jì)一種總等值線的方法，使得只在等值線會(huì)通過的單元節(jié)點(diǎn)上才求出函數(shù)值，而大部分不通過等值線的單元不必求節(jié)點(diǎn)函數(shù)值，以節(jié)省機(jī)器運(yùn)算時(shí)間、加快速度。這里還需指出，圖3是一個(gè)復(fù)雜化的例子，實(shí)際上，對(duì)于特定的等值線，一個(gè)網(wǎng)格內(nèi)的等值線可能是很少的，或者只有一條或少數(shù)幾條，還有很多網(wǎng)格根本就沒有。此外，當(dāng)分段加密時(shí)，等值線將逐漸變成光滑曲線，這是因?yàn)閿M合網(wǎng)格點(diǎn)的函數(shù)是一個(gè)連續(xù)光滑的函數(shù)。因此，利用上述方法，就不需要對(duì)等值線再作光滑處理。

將等值線分為從邊界出發(fā)到邊界結(jié)束的等值線以及在內(nèi)部封閉的等值線兩種類型，如圖4所示。無論哪種類型的等值線，都必須已知有順序相鄰的 A1、A2兩個(gè)等值點(diǎn)。 設(shè)它們的坐標(biāo)值分別為（a1x，a1y）及（a2x，a2y），分別以 nx和 ny為單位長(zhǎng)度，其坐標(biāo)的整數(shù)值為（a1j，a1i）及（a2j，a2i）。任意設(shè)兩個(gè)相鄰的單元如圖 5中的單元Ⅰ、Ⅱ所示，A2落在單元Ⅰ和Ⅱ的鄰邊上，A1落在單元Ⅰ的其他三邊的任意一邊上。用（i，j）表示單元Ⅰ的序號(hào)，用（i+1，j）、（i，j+1）、（i-1，j）、（i，j-1）表示 單 元Ⅱ 的 序號(hào)，分別相當(dāng)于圖 5 中的（a）、（b）、（c）、（d）4 種情況，則下一等值點(diǎn)追蹤方向依次由下述4個(gè)判別式確定。

（1）若 a1i＜a2i，即圖 5（a）的情況，下一點(diǎn)在單元（i+1，j）的另外三邊上；

（2）若 a1j＜a2j，即圖 5（b）的情況，下一點(diǎn)在單 元 （i，j+1）的另外三邊上；

（3）若 a2j＜a2x，即圖 5（c）的情況，下一點(diǎn)在單元（i-1，j）的另外三邊上；

（4）若 a2i＜a2y，即圖 5（d）的 情 況 ，下 一 點(diǎn) 在 單 元 （i，j-1）的另外三邊上。

以圖4的網(wǎng)格數(shù)組為例，首先從網(wǎng)格數(shù)組四邊尋找第一類等值點(diǎn)，當(dāng)在底邊、左邊、右邊（或上邊）的邊界上找到一個(gè)等值點(diǎn)時(shí)，令其為A2。因?yàn)樽粉櫼粭l等值線必須先要知道兩個(gè)等值點(diǎn)，所以必須再假設(shè)一個(gè)等值點(diǎn)A1。當(dāng) A2在底邊邊界時(shí)，令 A1、A2的位置關(guān)系滿足：a1i=-1，a1j=a2j；當(dāng) A2在左邊邊界時(shí)，令 a1i=a2i，a1j=-1；當(dāng)A2在右邊（或上邊）邊界時(shí)，令 a1i=a2i，a1j=a2j。 然后利用 4個(gè)判別式順次進(jìn)行判斷，即可追蹤到下一等值點(diǎn)A3。一旦追蹤到 A3，再把 A2充到 A1中，A3充到 A2中，就可繼續(xù)追蹤，一直追蹤到圖形邊界為止，這樣就完成了第一類等值線的繪制。對(duì)于第二類等值線，從內(nèi)部網(wǎng)格的橫邊或縱邊上先尋找它的第一個(gè)等值點(diǎn)，令其為A2。對(duì)于橫邊上的A2，再設(shè)一個(gè)等值點(diǎn)A1，使其和A2的位置關(guān)系滿足a1i=-1，a1j=a2j；對(duì)于位于縱邊上的 A2，令其a1i=a2i，a1j=-1，然后應(yīng)用4個(gè)判別式順次進(jìn)行判斷，即可追蹤到下一個(gè)等值點(diǎn)A3。按照第一類等值點(diǎn)的追蹤方法繼續(xù)追蹤下去，直到重新追蹤到首先找到的那個(gè)A2點(diǎn)的位置為止，這時(shí)便完成了一條封閉等值線的繪制。必須指出，無論屬于哪種類型的等值線，每當(dāng)追蹤到等值點(diǎn)落在網(wǎng)格邊上（而不是網(wǎng)格內(nèi)部的單元邊上）時(shí)，必須將其位置記在一個(gè)表內(nèi)。對(duì)于封閉性等值線，還需把第一個(gè)等值點(diǎn)的位置記住，以便每次開始找到第一個(gè)A2點(diǎn)時(shí)，與表中的位置作比較，如果已經(jīng)在表中，則不再重繪。對(duì)于封閉等值線，每追蹤到一個(gè)位于網(wǎng)格上的點(diǎn)時(shí)，還需同記錄下來的第一個(gè)等值點(diǎn)位置進(jìn)行比較，以判斷是否回到原來的出發(fā)位置。

按上述方法追蹤等值線，可以避免計(jì)算許多沒有等值線通過的單元的節(jié)點(diǎn)函數(shù)值，但也有可能漏掉一些只在一個(gè)網(wǎng)格之內(nèi)的小封閉等值線，如圖4中用虛線表示的等值線。但這些等值線在整個(gè)等值線圖形中并不起重要作用，忽略不繪并不影響圖形質(zhì)量，卻可以節(jié)省許多機(jī)器時(shí)間。

2 等值線繪制

關(guān)于等值線的繪制，主要有數(shù)據(jù)處理和圖形生成兩部分。在數(shù)據(jù)處理部分，主要是將離散數(shù)據(jù)網(wǎng)格化；而在圖形生成部分，系統(tǒng)將采用矩形網(wǎng)格追蹤法來生成等值線，然后就要對(duì)所生成的等值線進(jìn)行光滑處理和標(biāo)注。

離散數(shù)據(jù)網(wǎng)格化主要調(diào)用Surfer提供的方法進(jìn)行，具體程序如下。

網(wǎng)格追蹤程序主要封裝在CContourTrace類中。

利用本文方法所生成的等值線如圖6所示。

本文是提出了一種新的方法來實(shí)現(xiàn)等值線的繪制，它首先用成熟軟件Surfer將原始數(shù)據(jù)進(jìn)行網(wǎng)格化，再利用矩形網(wǎng)格追蹤法追蹤網(wǎng)格，并按照實(shí)際需求實(shí)現(xiàn)了對(duì)等值線的縮放、平移、光滑和標(biāo)注等功能。利用本方法所繪制的等值線主要在石油勘探開發(fā)、采礦等領(lǐng)域得到應(yīng)用。使用本方法所繪制的等值圖較其他方法效果更好、精確度更高、速度更快，是一個(gè)實(shí)用有效的等值線繪制方法。

[1]王永會(huì)，宋曉宇，欒方軍.基于網(wǎng)格的等值線圖快速生成算法[J].計(jì)算機(jī)工程與應(yīng)用，2001（17）：124-125.

[2]陳學(xué)工，劉凱敏.一種基于格網(wǎng)法快速生成等值線的算法[J].電腦與信息技術(shù)，2007，15（3）：4-6.

[3]宋麗娟，龔曉峰，鐘猛.基于網(wǎng)格法的等值線繪制方法[J].現(xiàn)代電子技術(shù)，2005，28（14）：65-67.

[4]韋美雁，杜丹蕾.基于規(guī)則網(wǎng)格的等值線的生成研究[J].湖南科技學(xué)院學(xué)報(bào)，2007，28（4）：39-41.

[5]于嘉，吳旭.一種改進(jìn)的矩形網(wǎng)格等值線追蹤算法[J].河南師范大學(xué)學(xué)報(bào)，2008，36（6）：34-36.

[6]張利，王俊彪，張賢杰.基于矩形網(wǎng)格追蹤法的曲面主曲率等值線生成算法 [J].計(jì)算機(jī)應(yīng)用研究，2009，26（8）：3179-3181.