剖面平面圖反算數(shù)據(jù)的實現(xiàn)

畢 武 袁小龍 段新力 黃顯義 彭仲秋 向詩強 張 恒（①新疆地礦局物化探大隊昌吉831100②烏魯木齊金維圖文信息科技有限公司烏魯木齊830091）

剖面平面圖反算數(shù)據(jù)的實現(xiàn)

畢 武 袁小龍 段新力 黃顯義 彭仲秋 向詩強 張 恒
（①新疆地礦局物化探大隊昌吉831100②烏魯木齊金維圖文信息科技有限公司烏魯木齊830091）

為了獲取矢量圖件上剖面平面圖的數(shù)據(jù)，可以將剖面平面圖線文件分解成基線和剖面線，根據(jù)數(shù)據(jù)比例尺、成圖比例尺等參數(shù)及線段相交原理，將原始數(shù)據(jù)反算出來。本文主要介紹了線段相交原理在剖面平面圖反算原始數(shù)據(jù)中的應(yīng)用，并編程實現(xiàn)了該方法，給出了實例。

剖面平面圖矢量化數(shù)字化線段相交。

1 引言

在航磁數(shù)據(jù)庫建庫工作中，接觸的原始資料很大一部分是記錄在厘米紙上的剖面曲線，早期的工作是利用CAD軟件和數(shù)字化儀，分別輸入剖面基線和剖面曲線，然后進(jìn)行計算處理得到剖面數(shù)據(jù)坐標(biāo)和值，其數(shù)字化精度常因不同的人操作而出現(xiàn)人為誤差，所以一般是采用兩個人分別數(shù)字化后，再取平均值的方式處理，這樣既費時又效率不高。特別是對現(xiàn)在的大比例尺磁測剖面平面圖，剖面線交錯密集，如果仍用以上方法數(shù)字化，將是一項很復(fù)雜的工作。本文作者采用MAPGIS進(jìn)行剖面圖矢量化，利用Qt Designer這一跨平臺圖形工具包和Fortran編程語言，反算出數(shù)據(jù)和坐標(biāo)，實現(xiàn)了剖面圖數(shù)字化程序的開發(fā)，集成在GeoIPAS V2.8系統(tǒng)中。程序操作簡單、方便、實用，生成的數(shù)據(jù)結(jié)果與原測數(shù)據(jù)吻合度好。

2 算法模型

《地理信息系統(tǒng)算法基礎(chǔ)》[1]第24頁，算法二。

定義A、B、C、D為二維空間的點，則有向線段AB和CD的參數(shù)方程為：

AB=A+r(B-A)，r∈[0，1]，CD=C+s(D-C)，s∈[0，1]

如果AB和CD相交，則：

A+r(B-A)=C+s(D-C)=＞Ax+r(Bx-Ax)=Cx+s(Dx-Cx)

Ay+r(By-Ay)=Cy+s(Dy-Cy)，r，s∈[0，1]

解方程，求r和s：

r=((Ay-Cy)(Dx-Cx)-(Ax-Cx)(Dy-Cy))∕((Bx-Ax) (Dy-Cy)-(By-Ay)(Dx-Cx))

s=((Ay-Cy)(Bx-Ax)-(Ax-Cx)(By-Ay))∕((Bx-Ax) (Dy-Cy)-(By-Ay)(Dx-Cx))

設(shè)P為直線AB和CD的交點，則：

P=A+r(B-A)=＞Px=Ax+r(Bx-Ax)，Py=Ay+r(By-Ay)

如果(0≤r≤1)and(0≤s≤1)，則有向線段AB和CD的交點存在，否則交點不存在。

如果(Bx-Ax)(Dy-Cy)-(By-Ay)(Dx-Cx)為0，則AB和CD平行。

如果(Ay-Cy)(Dx-Cx)-(Ax-Cx)(Dy-Cx)為0，則AB和CD共線。

如果直線AB和CD相交，而交點不位于線段AB和CD之間，則交點位置可以通過如下條件判斷：

如果r＞1，則P位于有向線段AB的延長線上；

如果r＜0，則P位于有向線段BA的延長線上；

如果s＞1，則P位于有向線段CD的延長線上；

如果s＜0，則P位于有向線段DC的延長線上。

3 功能設(shè)置

程序主要功能，見圖1：

⑴按輸入線順序?qū)?yīng)：基線和剖面線的明碼文件是一一對應(yīng)的，即可以按輸入順序來進(jìn)行轉(zhuǎn)換；

⑵按屬性順序?qū)?yīng)：根據(jù)兩個屬性文件的編號來進(jìn)行對應(yīng)，這種情況適用于當(dāng)在MAPGIS編制圖時，基線和剖面線沒有按照一定順序輸入時，就需要使用屬性文件來進(jìn)行對應(yīng)；

⑶圖幅比例尺：剖面圖成圖比例；

⑷數(shù)據(jù)比例尺：剖面線縱向的顯示比例；

⑸方位角選擇：在計算的過程中，可以根據(jù)剖面線的角度來進(jìn)行值的運算；

⑹正負(fù)值鏡像：如果計算結(jié)果與實際正負(fù)值相反，可選擇鏡像選項改正。

圖1 剖面圖數(shù)字化界面

4 程序?qū)崿F(xiàn)

一個完整的計算過程為：剖面圖數(shù)字化→轉(zhuǎn)為明碼文件→輸入剖面基線明碼文件→輸入剖面線明碼文件→輸入基線屬性文件→輸入剖面線屬性文件→輸入保存結(jié)果文件名→輸入剖面圖參數(shù)→計算輸出數(shù)字化剖面數(shù)據(jù)文本結(jié)果。

剖面圖參數(shù)包括：成圖比例尺、數(shù)據(jù)比例尺、剖面方位角的選擇。剖面方位角考慮了航測和地面測量的多種情況，如果剖面線為統(tǒng)一的方位角，就選擇取“第一條線方位角”；如果各剖面線方位角稍有變化，還可以選擇“取所有線的平均值”；對于一個工區(qū)如果設(shè)計了多個方位角或者環(huán)形方位角的情況，還可以選擇按“每條線取方位角”。

程序計算如下：

!基線點和剖面點不相等，計算基線的方位角：

do 500 i=1,nline-1

if(iselectline.eq.2)then

fwj=atan((y1(i,1)-y1(i,n1(i)))∕(x1(i,1)-x1(i,n1(i))))

else

endif

!循環(huán)計算，先判斷過剖面線上每個點的垂直基線的線與基線是否相交，如果相交，再計算交點坐標(biāo)：

do 600 j=1,n2(i)

!直線CD(xc,yc,xd,yd)是以垂直基線且過剖面線上的點的直線：

xc=x2(i,j)-zmax*cos(fwj+pi∕2.0)

yc=y2(i,j)-zmax*sin(fwj+pi∕2.0)

xd=x2(i,j)+zmax*cos(fwj+pi∕2.0)

yd=y2(i,j)+zmax*sin(fwj+pi∕2.0)

do 700 j1=2,n1(i)

xa=x1(i,j1-1)

ya=y1(i,j1-1)

xb=x1(i,j1)

yb=y1(i,j1)

st=(xb-xa)*(yc-yd)-(xc-xd)*(yb-ya)

t=((xc-xa)*(yc-yd)-(xc-xd)*(yc-ya))∕st

s=((xb-xa)*(yc-ya)-(xc-xa)*(yb-ya))∕st

if(s.ge.0.and.s.le.1.and.t.ge.0.and.t.le.1)then

isEmpty(i,j)=1

x3(i,j)=xa+t*(xb-xa)

y3(i,j)=ya+t*(yb-ya)

dx=x2(i,j)-x3(i,j)

dy=y2(i,j)-y3(i,j)

dxy=sqrt(dx*dx+dy*dy)

if(dxy.gt.0.0001)then

if(dy.gt.0.0)then

z3(i,j)=dxy

else

z3(i,j)=-1.0*dxy

endif

else

z3(i,j)=0.0

endif

goto 600

else

!基線比剖面線短

endif

700continue

600continue

500 continue

5 例子

有一航磁數(shù)據(jù)，經(jīng)MAPGIS矢量化后如圖2原始剖面平面圖，經(jīng)過剖面平面圖數(shù)字化處理后，用GeoIPAS系統(tǒng)的剖面平面圖成圖，結(jié)果如圖3數(shù)字化結(jié)果成剖面平面圖，將原始的基線和剖面線與結(jié)果的剖面平面圖相互疊加，可以看到兩線有部分未完全重合，但相對于整個剖面線值來說誤差已經(jīng)很小，分析原因是：計算基線斜率的精度由有效位數(shù)引起，在剖面線與基線的距離變大時其累積的誤差也相應(yīng)變大。

圖2 原始剖面平面圖

圖3 數(shù)字化結(jié)果成剖面平面圖

6 結(jié)束語

做剖面平面圖數(shù)字化時應(yīng)注意以下幾點：

⑴剖面平面圖矢量時，基線和剖面線方向要一致，順序要對應(yīng)。如有線號，可輸入其線號作為屬性對應(yīng)基線和剖面線。

⑵數(shù)字化的精度主要取決于矢量圖精度，所以在做矢量化時要提高底圖的分辨率。

⑶對多幅圖拼接的剖面平面圖，如果在接圖部分比較復(fù)雜，要注意基線和剖面線拼接對應(yīng)關(guān)系。

⑷對剖面平面圖原圖矢量后要在MAPGIS中要做誤差校正，應(yīng)利用所有的圖框校正點。

⑸對線文件編輯完成后，需要進(jìn)行壓縮保存文件，然后再做明碼文件轉(zhuǎn)換。

⑹剖面平面圖數(shù)字化結(jié)果可成圖與原剖面平面圖對比檢查，保證各基線與剖面線的對應(yīng)正確，并可檢查矢量化精度。

剖面平面圖數(shù)字化是一個簡單方便將剖平圖轉(zhuǎn)換為數(shù)據(jù)的工具，尤其是用在對航磁剖平圖老資料的數(shù)字化建庫中，對提高對老資料的再利用上有很大方便性。

［1］張宏，溫永寧，劉愛利.地理信息系統(tǒng)算法基礎(chǔ).北京：科學(xué)出版社，2006.

［2］劉天佑.地球物理勘探概論.北京：地質(zhì)出版社，2007.

［3］張勝業(yè)，潘玉玲.應(yīng)用地球物理學(xué)原理.武漢：中國地質(zhì)大學(xué)出版社，2004.

［4］吳信才.MapGIS地理信息系統(tǒng)[M].北京：電子工業(yè)出版社，2004.

［5］吳立新，史文中.地理信息系統(tǒng)原理與算法.北京：科學(xué)出版社，2003.

［6］林曉彤.Fortran 90編程基礎(chǔ).青島：中國海洋大學(xué)出版社，2006.

［7］彭國倫.Fortran 95程序設(shè)計.北京：中國電力出版社，2002.

［8］（美）Papheal Pender，蘇劍，等.標(biāo)準(zhǔn)C++編程電子工業(yè)出版社，2002.

收稿：2014-12-29

DOI∶10.16206∕j.cnki.65-1136∕tg.2015.01.015