裴光菊,唐德利,魏德照
(1.國家測繪地理信息局重慶測繪院,重慶 40001)
基于ArcGIS、Visual C++的DEM數據接邊檢驗方法
裴光菊1,唐德利1,魏德照1
(1.國家測繪地理信息局重慶測繪院,重慶 40001)

為提高DEM數據接邊檢驗的高效性、準確性、可靠性,比較各種檢驗方法后,本文結合ArcGIS組件,利用Visual C++開發平臺開發了一套DEM檢驗軟件,并在實際檢驗項目中驗證了該軟件的使用效果。檢驗結果表明,該軟件在DEM接邊檢驗工作中高效、準確、可靠。
數字高程模型(DEM);數據接邊;質量檢驗
數字高程模型(digital elevation model,DEM)是用一組有序數值陣列形式表示地面高程的一種實體地面模型,是數字地形模型(digital terrain model,DTM)的一個分支。DEM的建立方法和存儲格式有很多種,如 ARC ASCII Grid(*.asc)、ARC Info GRID(*.grd)、ESRI FLOAT BIL(*.bil *.hdr *.blw)、DEM GB(*.dem)、GRID GB(*.grd)、DEM VZ(*.dem)等,其數據內容的表現形式有兩種:規則矩形格網和不規則三角網。由于DEM描述的是地面高程信息,在測繪、水文、氣象、地貌、地質、土壤、工程建設、通訊、氣象、軍事等國民經濟、國防建設以及人文和自然科學領域都有著廣泛的應用[1-6]。因此,DEM產品質量檢驗的智能化和數據拼接質量的檢驗十分重要。
在地理信息數據中,DEM按國家標準經緯度梯形分幅編號存放,其數據范圍由該梯形在高斯平面的最小外接矩形經過一定外擴來確定,相鄰圖幅間存在重疊區域,重疊區域的同名格網點高程要保持一致,這就要求對圖幅數據進行接邊處理。DEM數據的質量要素中接邊質量要求尤為嚴格,最終成果均要求同名格網的高程值必須相同。生產過程數據接邊的同名格網點高程值較差不能超過兩倍高程中誤差,符合要求時,以平均值確定格網點高程值,如不符合要求,則要查明原因并處理后再作接邊。DEM最終成果不允許相鄰圖幅同名點高程值存在差異[2]。
DEM接邊質量檢驗就是查看相鄰圖幅DEM數據重疊部分同名格網點的高程值是否能達到規定要求[7]。
在ArcGIS上打開相同投影帶中相鄰的兩幅DEM數據,利用該軟件的DEM視圖功能,直接比較同名格網點灰度值,獲取同名格網點高程差值,也可以利用ArcGIS的柵格計算功能直接以同名格網點灰度差值顯示重疊部分視圖,以獲取同名格網點高程差值。若開發個人軟件實現接邊檢驗,則嵌入ArcGIS的COM組件同時加載相鄰兩幅DEM數據,獲取重疊區域,并在重疊區域內循環讀取同名格網點高程數據對,并計算較差,完成循環后計算接邊中誤差,最后輸出檢驗報表。
跨投影帶圖幅數據接邊檢驗不能依靠手工完成,只能開發專業軟件來實現。由于不同投影帶下的相鄰圖幅在拼接前需要將其中一幅數據進行換帶處理,也就是需要按照很復雜的數學規則進行計算,使兩幅數據處于同一投影帶下,手工方法是不可能實現的,必須借助計算機程序進行批量處理。數據跨帶接邊檢驗與相同帶接邊檢驗不同,必須先把相鄰DEM數據投影到同一投影帶,然后再獲取同名格網點高程數據對,并計算較差,完成循環后計算接邊中誤差,最后輸出檢驗報表。
現有的DEM數據接邊質量檢驗方法主要有:
1)在ArcGIS下,采用單點查詢方式檢查同名格網高程。
2)在ArcGIS下,利用其柵格計算功能手工逐個文件檢查同名格網高程。具體方法是首先加載相鄰兩幅DEM數據,加載后圖幅鄰接處如圖1。然后在Spatial Anlyst工具欄下運行Raster Calculator工具,如圖2。在彈出的對話框中選擇一個數據文件→點擊所需運算符“<>”→選擇另一數據文件→點擊Evaluate按鈕。完成后,在圖層瀏覽窗口會自動產生一個新的柵格圖層Calculation,并在主窗口顯示該圖層,效果如圖3。圖中亮顯的星云狀部分即為灰度值為1的格網點,表示這些格網點高程值不同。

圖1 相鄰DEM數據加載后鄰接處

圖2 DEM柵格疊加工具

圖3 DEM柵格疊加計算結果
3)在ArcGIS下,利用其卷簾功能,觀察同名格網點像素灰度值。
4)在ArcGIS下利用VB進行內嵌式編程檢查。
5)利用ArcGIS組件,在Visual C++/Visual Basic等開發平臺下,開發獨立運行檢查軟件[3]。
前3種方法主要是基于手工方式,只能用于相同投影帶數據檢驗,不能用于跨帶數據接邊檢驗,且不能直接形成檢驗結果統計報表。第4種方法可以實現自動,但其在運行界面、運行速度及用戶使用舒適性等方面均有一定的局限性。實際數據生產中,DEM數據的接邊不一致現象非常普遍。然而,一個圖幅范圍的DEM數據柵格有上百萬個,數據量十分龐大,接邊檢驗是一項十分繁重的任務。這種情況下,采用手工檢驗方法是不切實際的。因此,本文探討了基于ArcEngine 9.3新增的柵格數據接口,在Visual C++6.0開發平臺下實現DEM數據的自動接邊檢驗,并運用到實際生產過程中,以提高DEM接邊檢驗的工作效率。
4.1 軟件功能要求
軟件設計以方便實用、高效準確、界面友好、易于掌握為原則。這就要求軟件界面能設置用戶參數(如坐標系統、高程系統),自動搜索指定路徑下的所有DEM數據文件,自動判斷鄰接圖幅DEM數據文件名及其存在性,判斷本圖幅與鄰接圖幅是否在同一投影帶下,進行數據讀入并獲取同名格網點的高程值,統計同名格網點高程較差及其中誤差,輸出檢驗報表。由于DEM數據量大,往往需要較長時間的數據處理,所以還需要設計適時刷新且不會卡死的進度顯示條。
4.2 軟件界面設計
由于軟件要實現批量自動處理和結果輸出,所以,除主窗口菜單、工具條外,還需要選擇數據文件路徑界面、軟件運行進度界面、運行結束提示界面和運行結果顯示界面。
4.3 設計框圖

圖4 設計框圖
4.4 核心功能實現
4.4.1 導入所需ArcGIS組件
#import文件夾所在路徑+"/com/esriGeometry.olb" raw_interfaces_only raw_native_types no_namespace named_guids;其他所需的組件同樣導入。com組件導入后,ArcGIS對象及方法才能被正常使用。
4.4.2 自動搜索DEM文件
自動搜索指定文件夾下的按比例尺分幅存放的DEM數據(其文件名為標準圖幅編號),并將結果存放在CStringList列表對象中。實現模塊為GetFileByFolder,傳遞參數為文件擴展名CString ContainFileEx,返回值為CStringList對象。關鍵代碼[5]:
ContainFileEx.MakeUpper(); CFileFind finder;
CString fileName,ls;
CStringList strList,strList0; strList.AddTail(BeginFolder);
POSITION ps=strList.GetHeadPosition();
while(ps)
{
CString strName=strList.GetNext(ps);
BOOL bWorking = finder.FindFile(strName+"\*.*");
while(bWorking)
{
bWorking = finder.FindNextFile( );
if (finder.IsDirectory()&&(!finder.IsDots()))//文件夾
……
}
}
4.4.3 計算并搜索鄰接DEM文件
循環讀取放在CStringList列表對象中的DEM文件,并獲取全部鄰接圖幅DEM文件名稱,驗證其存在。實現模塊為GetRoundNum,傳遞參數為當前DEM圖幅文件名CString strThisNum。由于不同比例尺圖幅的鄰接關系不一樣,所以編程過程中需要判斷DEM的分幅比例尺,按照國家標準分幅進行判斷和計算。
4.4.4 加載DEM文件數據并讀取同名格網點高程
加載DEM數據需要用到ArcGIS的IRaster組件接口IRasterPtr、IrasterLayer的IrasterLayerPtr接口[3]。代碼為:
IRasterLayerPtr pRasterLayer(CLSID_RasterLayer);
pRasterLayer->CreateFromFilePath((CComBSTR)strPath);其中,strPath是包含DEM數據的文件路徑。由于DEM采用高斯投影平面坐標系統,涉及到投影分帶的問題,如果鄰接圖幅在不同投影帶,則需要對加載數據進行投影換帶,以保證接邊檢驗的鄰接數據在相同投影帶下。需要引入ArcGIS的ISpatialReferenceFactory3、ISpatialReference組建接口:ISpatialReferenceFactory3Ptr、ISpatialReferencePtr。代碼為:
ipSRFactory->CreateESRISpatialReferenceFromPRJFile(str PrjP,&ipPSR);
pRasterLayer ->putref_SpatialReference(ipPSR);
其中,strPrjP為高斯投影坐標系統的定義文件(包含路徑)。注意,這里要先定義pRasterLayer的原坐標系統,再加載DEM數據,之后重定義pRasterLayer的坐標系統為目標系統,實現跨帶轉換。讀取格網點高程代碼為:
IRaster2Ptr ipRa2=ipRaster;// ipRaster為讀取的IRaster數據接口。
IRasterCursorPtr ipRaCur;
ipRa2->CreateCursorEx(ipPnt,&ipRaCur);// ipPnt相鄰圖幅重疊范圍IPnt接口。
IPixelBlockPtr ipPB;
ipRaCur->get_PixelBlock(&ipPB);
VARIANT vtValue;
ipPB->GetVal(0,m,n,&vtValue);//m,n為行列序號。
4.4.5 不同格網間距DEM的接邊檢驗
由于格網間距不同,鄰接圖幅的重疊區域就沒有明確的同名格網點,需要將其中一個文件的DEM數據進行雙線性內插,以求取同名點高程值進行比較,確定拼接精度[4]。雙線性內插數學模型為:設Q11、Q12、Q21、Q22分別為內插點位最鄰近的4個格網點,其點位坐標分別為(X1,Y1) (X2,Y2)(X3,Y3)(X4,Y4),dx為橫坐標差,內插點坐標為(xc,yc)。則插值如下。

其中,H為最終內插值。
4.5 數據格式轉換
實際應用中可能存在不同數據格式的DEM,這就需要將數據轉換為滿足程序運行要求的目標格式。在ArcMap下可以實現多數常用格式的轉換,也可以在VC++下進行編程處理。
4.6 進度條的實現
由于DEM數據文件較大、數量多,程序往往會運行較長時間,需要設計進度條來為用戶顯示程序處理的進度情況,這就要求進度條顯示更新與數據處理運行在不同空間。因此,該程序把DEM數據處理放在線程去運行,進度條在主進程中實時顯示刷新。代碼中用“static UINT DEMGCJD_Thread(LPVOID param)”語句申明線程,“UINT DEMJBJD_Thread(LPVOID param)”語句定義線程。所有DEM數據處理的代碼將在此線程中運行,使進度條在程序運行期間始終處于活躍更新狀態,避免了程序卡死的現象。
4.7 項目檢驗應用
該程序在國家基礎測繪1∶50 000和地方基礎測繪1∶10 000、1∶5 000等比例尺的DEM數據接邊檢驗工作中得到了充分使用。例如,重慶1∶5 000數字地形測量項目中,其生產期持續三年多,成果包含數千幅分幅的DEM數據,檢驗工作量巨大,人工檢驗不可能在期限內完成。通過對該項目全部DEM數據進行接邊檢驗、裁切范圍檢驗,在0.5 h內程序運行結束,以100%的準確率輸出全部數據檢驗報表。
在諸多DEM數據接邊檢驗的方法中,借助現行的工具軟件進行手工操作的方式只能適應極少量的數據。在面對大量數據時,通過借助ArcGIS提供的組件,結合VC++(或者VB等)程序開發工具,編寫一套程序來完成大批量數據接邊檢驗(也可以包括數據范圍檢驗),快速、高效、準確、可靠,是地理信息數據檢驗必要的方法。
[1] 李志林,朱慶. 數字高程模型[M].武漢:武漢大學出版社,2001
[2] CH/T 1015.2-2007.基礎地理信息數字產品1∶10 000 1∶50 000生產技術規程 第二部分:數字高程模型(DEM)[S].
[3] 彭珊鸰,李軍,廖明,等.基于arcengine10的dem數據自動接邊檢查實現[J].科技資訊,2012,(18):9
[4] 王佩,呂志勇.DEM產品質量檢查標準研究與實現[J].測繪與空間地理信息,2011,34(5):280-283
[5] 李書濤.C語言程序設計教程[M].北京:北京理工大學出版社,1993
[6] 楊秀伶.數字高程模型DEM的構建與應用[J].綠色科技,2014(5):315-316
[7] 劉錦軍,鈕利平,孫穎.1∶50 000數字高程模型(DEM)的質量控制[J].東北測繪,2002,25(2):44
P208
B
1672-4623(2016)12-0084-03
10.3969/j.issn.1672-4623.2016.12.028
裴光菊,工程師,主要從事測繪地理信息產品檢驗工作。
2015-12-15。