GMT在繪制遼寧地形圖中的應用

孔祥雪，王萬寧

（遼寧省地震局，遼寧 沈陽 110034）

0 引言

GMT，全稱Generic Mapping Tools，一般翻譯成“通用制圖工具”。其是擁有80多種命令行工具的開源繪圖軟件，支持30多種地圖投影，自帶模塊包括海岸線、河流、國界等數據信息，用于繪制不同類型的專題地圖[1-2]。GMT是純命令行的輸入方式，所有的繪圖操作都是通過命令行代碼執行的。命令行相對于圖形界面的優勢在于內存占用更低；精確控制圖形的顯示，例如線條寬度、圓的大小；方便寫成腳本，具有可批量處理、可重復、可自動化的特點。GMT將繪圖及數據處理功能劃分到不同的模塊，各個模塊之間相互獨立且代碼量少，易于更新和維護；繪圖步驟之間以及涉及的數據類型獨立性強，因而應用領域較廣；可在腳本語句中調用一系列程序，進而繪制更為復雜的圖片。GMT輸出的圖件格式為矢量圖片格式，具有可以任意放大縮小而不失真的特點。GMT充分利用這種特性，可以生成高質量的矢量圖件，并可以很容易地轉換為其他圖片格式。基于上述優點，GMT在地球科學領域得到了廣泛的應用[3-4]。

目前GMT還提供了與MATLAB兼容的程序接口，可以在MATLAB中直接調用GMT命令；而GMT命令生成的結果，如格網數據（GRD）、表格數據 （TABLE）、顏色表 （CPT）、文本文件（TXT）、圖片（PSFILE）等，都可以作為MATLAB變量調用運算；此外，MATLAB中的矩陣變量可以直接作為GMT的輸入數據[5]。上述特性將MATLAB數學運算與GMT成果繪制的功能優勢有效地融合在一起，進一步促進了GMT在地球科學研究中的應用。

目前，遼寧省地形圖通常使用MapInfo和ArcGIS等軟件進行繪制。其中，MapInfo的含義是“Mapping+Information（地圖+信息）”即：地圖對象+屬性數據。MapInfo是一種操作簡便的桌面地圖信息系統，它具有圖形的輸入與編輯、圖形的查詢與顯示、數據庫操作、空間分析和圖形的輸出等基本功能，通過菜單條命令或工具條按鈕進行操作，提供地圖窗口、瀏覽窗口、統計窗口以及幫助輸出設計的布局窗口等查看表窗口，并將結果方便地輸出到打印機或繪圖儀。總體來說，MapInfo軟件更加側重于地圖信息的屬性查詢，在地圖繪制功能方面相對較弱。ArcGIS for Desktop是對地理信息進行編輯、創建以及分析的一種GIS軟件，提供了一系列的工具用于數據采集和管理、可視化、空間建模分析以及高級制圖。ArcGIS不僅支持單用戶和多用戶的編輯，還可以進行復雜的自動化工作流程，同時便捷地實現柵格圖像數字化輸入、地圖的編輯與制圖輸出、數據的編輯與管理，以及空間分析等功能，是專題圖繪制的較好選擇[6]。雖然MapInfo和ArcGIS軟件操作簡單，使用方便，但是在地形繪制和等值線繪制方面存在一定的局限性。以遼寧省地形圖繪制為例，遼寧省基本地貌可概括為“六山一水三分田”，地勢大致呈現自北向南、由東西兩側向中部傾斜的趨勢，山地丘陵分列遼東、遼西兩側，遼河平原位于中部，呈馬蹄形指向遼西渤海沿岸狹長的海濱平原[7]，境內山脈分別位于東西兩側。考慮到遼寧省內平原與山地丘陵地區比例均衡，且境內山地丘陵海拔相對較低，與平原地區海拔差較小，利用MapInfo和ArcGIS軟件等現有手段繪制的二維地形圖存在山體陰影效果不明顯的現象，且難以實現三維地形圖的繪制。為了更好地展現遼寧省地形地貌特征，本文擬采用GMT軟件實現遼寧省二維及三維地形圖的繪制。

1 GMT數據類型介紹

GMT作為一種繪圖工具，同時作為一種數據處理工具，模塊涉及多種輸入和輸出數據，但概括來講，GMT軟件的輸入輸出數據類型可劃分為以下幾類。

第一類為表數據，英文稱為table data，也稱為列數據或者多列數據，常用于表示點和線。表數據中有N個記錄，每個記錄都有M個字段。表數據可以有三種形式：ASCII表、二進制表和netCDF表。最常用的表數據形式是ASCII表，可以使用編輯器直接編輯。ASCII表中有N行M列，每行代表一個記錄，每列代表一個字段。一個記錄內的字段之間可以用空格、制表符、逗號或分號分隔。

第二類為網格文件，GMT可以繪制2D網格數據。通常，2D網格文件的X方向代表經度、Y方向代表緯度，Z值可以表示高程、重力值、溫度、速度等。也可以將XYZ格式的1D表數據網格化得到2D網格數據。GMT可以使用的2D網格數據有三種格式：netCDF格式、Sun光柵文件和自定義的二進制數據格式。其中最常見也較為推薦的網格數據格式是netCDF格式，一般以.nc或.grd作為文件后綴。

第三類為CPT文件，CPT全稱是color palette table，也稱為調色板或色標文件。CPT文件可以在grdimage、psxy、psxyz等命令中使用。一般情況下，不同的數據范圍顏色要求也不盡相同，所以可以使用makecpt或grd2cpt命令對已有的CPT（例如GMT內置的CPT）文件進行自定義編輯。當然，也可以根據不同程序中的顏色配比需求手寫CPT文件，或者使用awk、perl之類的文本處理工具自動生成CPT文件。

第四類為PS文件，GMT生成的圖片為PS文件格式，全稱是PostScript。PostScript是一種用于描述矢量圖形的頁面描述語言。PS是矢量圖形格式，即用點、線段或多邊形等數學方程的幾何元素來表示圖像。因而可以任意旋轉與縮放而不會出現圖像失真的情況。GMT提供了psconvert命令，可以很方便地將PS文件以任意精度轉換為jpeg、eps、png、pdf等圖片格式，以滿足不同情形下的需求。

2 地形圖的繪制

GMT繪圖命令較為復雜，從該軟件網站可以獲取較為詳盡的使用手冊和模塊手冊介紹，本文以遼寧省及周邊作為繪圖區域，利用手冊總結出二維和三維地形圖的繪制方法。

2.1 2D地形圖的繪制

利用GMT繪制2D地形圖的步驟概括如下：

（1）繪制地形底圖，加載1′×1′格網分辨率的earth_relief_01m.grd全球地形模型文件。

gmt grdcut earth_relief_01m.grd-R＄R-G%grd%

gmt grdgradient%grd%-Ne0.7-A50-G%grd%

grdcut命令指從大區域地形中裁剪出自定義范圍的數據（本例中經度為118°～126°E，緯度38°～44°N），該操作可降低后續梯度計算的復雜度，提高繪圖效率。grdgradient命令為制作梯度信息的網格文件，得到某一個方向的方向導數。

（2）生成并編輯調試CPT文件。

gmt grd2cpt%grd%-Cglobe-S-10000/10000/200-Z-D＞%cpt%

grd2cpt命令為制作顏色表文件。-C選項指定生產的cpt文件類型，-S選項設置顏色表最小值和最大值和間隔，-Z選項生產連續顏色的CPT文件。

（3）繪制2D地形圖。

gmt grdimage%grd%-R-J-I%grd%-C%cpt%-Q-O-K ＞＞＄PS

grdimage命令為在地圖上繪制網格數據文件，-I選項接grdgradient計算得到的梯度文件，-C選項接顏色表文件。

（4）加載行政區邊界等信息。

gmt psxy%CN-border-La% -J-R＄R-W0.5p-O-K ＞＞＄PS

gmt psxy%prov%-J-R＄R-W1.5p/red-O-K ＞＞＄PS

gmt psxy%city%-J＄J-R＄R-Sc0.08c-Gblack-K-O＞＞＄PS

gmt psxy%city%-J＄J-R＄R-Sc0.15c-W0.2p,black-K-O ＞＞＄PS

gmt pstext%city%-J＄J-R＄R-Dj0.15c/0.15c-P-K-O ＞＞＄PS

psxy為繪制線段，符號命令%CN-border-La%，%prov%，%city%分別為國界、省界、市政府所在地，-W選項設置線段或符號輪廓的畫筆屬性，-G選項設置符號或多邊形的填充色，-S選項表示要繪制符號類型。pstext在圖上寫文本，-D選項表示文本在指定坐標的基礎上相對x軸和y軸的偏移。pstext命令中加載的市政府所在地標注信息，其數據格式是經度、緯度、字符大小、旋轉角度、字體、文本。

（5）編輯圖像的標題、圖例等信息。

gmt psscale-DjCB+w18c/0.3c+o0/-1.7c+h-R-J-C%cpt%-BWSEN-Bxa100f50-G0/400-O-K＞＞＄PS

psscale命令為繪制色標，-D選項表示色標在圖中的位置，+w18c/0.3c表示色標長18厘米，寬0.3厘米，-C選項接顏色表文件，-B選項表示色標間隔，-G選項表示截斷范圍為0～400 m。

（6）將繪制的地形圖輸出為可移植網絡圖形格式（PNG）文件。

gmt psconvert＄PS-A-Tg-P-Z

-A選項表示對PS文件進行裁剪，僅保留其中有繪圖的部分。-T選項指定要轉換的圖片格式。-Z選項轉換完成后刪除輸入的PS文件。

基于上述步驟，利用GMT軟件繪制的遼寧省2D地形圖如圖1（左）所示。

作為對比，利用ArcGIS軟件繪制2D地形圖的步驟如下：

（1）使用ArcToolbox→Spatial Analyst工具→表面分析→山體陰影工具，輸入dem圖層，生成新的柵格數據圖層；

（2）將柵格數據圖層透明度設置為35%；

（3）將柵格數據圖層移至dem圖層上方，選擇dem圖層中符號系統選項中適合的色帶。

根據上述步驟，利用ArcGIS軟件繪制的遼寧省2D地形圖如圖1（右）所示。

圖1 基于GMT（左）和ArcGIS（右）繪制的遼寧省2D地形圖Fig.1 The 2D topographic map of Liaoning Province by using GMT （left） and ArcGIS （right） respectively

通過根據上述命令制作的地形圖，可以大致上了解遼寧地形地貌，但從圖1（左）中可以看出在地勢比較平緩的中部地區顏色的區分差別不大，東西兩側山體陰影效果不是十分顯著。通過左右兩幅圖對比，利用ArcGIS繪制的地形圖山體陰影效果相對較弱，僅僅在顏色上面區分高程較為明顯。

業務部門在使用ArcGIS軟件繪制專題圖時通常只繪制二維圖形，在三維展示模塊具有一定的局限性。因此GMT的靈活、開源等優勢更加得到凸顯，下面將進一步討論GMT的三維成圖效果，并與圖1進行對比。

2.2 3D地形圖的繪制

利用GMT繪制3D地形圖的步驟概括如下：

（1）網格化數據，生成grd文件。

surface%dat%-Gtmp.grd-R118/126/38/44-I2m

surface命令是數據網格化的一種方法。%dat%選項為網格數據，-G選項為網格化后的數據文件名，-R選項為網格化經緯度范圍，-I選項為指定x方向和y方向的增量。

（2）生成并編輯調試CPT文件。

grd2cpt tmp.grd-Cglobe-S-10000/10000/200-Z-V＞%cpt%

grd2cpt命令制作顏色表文件。同二維地形圖。

（3）制作梯度信息的網格文件。

grdgradient tmp.grd-A270-Gtmp.grad grdhisteq tmp.grad-Gtmp.his-N grdmath tmp.his 0.4 x=tmp.int

grdgradient命令是制作梯度信息的網格文件，得到某一個方向的方向導數，輸入文件為tmp.grd，-A選項設置方向導數的方向為270°，輸出的文件為tmp.grad。grdhisteq命令表示網格文件柱狀圖均化。grdmath命令是對grd文件進行數字計算。

（4）繪制3D地形圖。

grdview-K-O-S-P tmp.grd-JM10-JZ2.5-Ba1/a1/800:"Altitude(m)":SEwnZ-Qs-N-1/235-Itmp.int-X2-C%cpt%-R-E150/45-Y13＞＞＄psf

grdcontour tmp.grd-J-P-O-K-E230/30-RC100-A1000+p+g255+s5-Wa3/0-X-1-Y-11.5＞＞＄psf

grdview-O-K-S-P tmp.grd-J-JZ2.5i-Ba1/a1/800SeWnz-Qs-N-100-Itmp.int-C%cpt%-RE230/30＞＞＄psf

grdview命令為基于二維網格數據文件繪制三維灰度或彩色圖像。-B選項設置地圖的邊框注記以及刻度間隔，-Qs選項表示采用表面繪制，-I選項表示給出一個含有照明強度的文件名，-E選項設置觀察點的方位角和高度角。

grdcontour命令為繪制等值線。-C選項表示指定等值線繪制方式，-A選項表示標注標簽信息，-W選項表示指定畫筆屬性。

（5）編輯圖例等信息。

psscale-E-D15/5/7/0.3-C%cpt%-Ba1000g 500/:"(m)":-O-P-V＞＞＄psf

（6）將繪制的地形圖輸出為可移植網絡圖形格式（JPEG）文件。

ps2raster-A-Tj＄psf

通過上述步驟，繪制遼寧省3D地形圖，如圖2所示。

圖2 遼寧省3D地形圖Fig.2 3D topographic map of Liaoning Province

圖2上圖為俯視的遼寧省三維地形圖，下圖為含有等高線的三維地形圖。圖2相比于圖1來說，在展現山體陰影方面有了很大的提高，同時在中部地區的平地也可以看出地形的凹凸起伏。相比二維地形圖，三維地形圖在地形地貌的表達上具有較大的優勢。

3 結論

通過上述例子可以看出，GMT可以繪制出精美直觀的二維和三維的地形圖，并且在此基礎上可添加行政邊界和政府所在地等信息，結果表明GMT在繪制地形圖上具有操作方便、圖形美觀等優勢。此外，GMT軟件同樣實現了地形圖由平面到立體，由抽象到直觀，打破了二維平面展示和二維可視的限制，更好的展現了遼寧省地形，為后續基于地形圖繪制地震應急專題圖等方面提供積極性的參考。