[轉載]GMT地形資料總結

原文地址：GMT地形資料總結作者：Jason

轉載：http://seisman.info/gmt-topo-grid-datas.html

 

目前接觸到的地形資料主要有三種：etopo、GTOPO30以及SRTM。按照精度從低到高排序如下：

etopo5：GMT網格資料之etopo5

etopo2：GMT網格資料之ETOPO2v2

etopo1：GMT網格資料之ETOPO1

GTOPO30：GMT網格資料之GTOPO30

SRTM30_PLUS: 

全球地形起伏資料：SRTM30_PLUS

SRTM3：高精度地形網格資料SRTM

ASTER GDEM：全球數字高程資料：ASTER GDEM

其中etopo資料包含了全球地形和海水深度資料；GTOPO30資料包含了全球地形資料；SRTM30_PLUS精度：30弧秒，網格精度約 1 km,包括全球地形和海洋深度資料，大部分割槽域使用SRTM30中的全球地形資料，高緯度區域使用GTOPO30資料（無SRTM資料），南緯81度到北緯81度的海底地形資料來自於Smith和Sandwell的全球1分資料等;SRTM3資料包含了全球大部分割槽域的地形資料。ASTER GDEM資料覆蓋範圍為北緯83°到南緯83°之間的所有陸地區域，達到了地球陸地表面的99%。

 

 

對於全球尺度或者大尺度地形圖來說，可以使用etopo1資料；

對於3度*3度-10度*10度的區域，可以使用GTOPO30資料；

對於3度*3度以下的區域，可以使用SRTM3資料。

 

    etopo5是全球地形及水深資料，其網格取樣間隔為5弧分，目前已經被更高精度的etopo1所取代，但是在繪製大區域地形（比如全球）時，5弧分的精度也夠用；

    etopo2v2與ETOPO5類似，也是高程+海底地形資料，其取樣間隔為2弧分；

    etopoO1包含了全球地形和海洋深度，取樣間隔為1弧分，是目前精度最高的global relief資料；

    GTOPO30是USGS公佈的全球數字高程模型，其取樣間隔為30弧秒；

    SRTM，全稱為Shuttle Radar Topography Mission，該專案獲取了北緯60度至南緯60度之間的雷達影像資料，進一步處理製成了數字地形高程模型（DEM），即現在的SRTM地形產品資料。

 

ETOPO1：

其分為兩個版本，Ice Surface和Bedrock，兩個版本基本一致。不同之處在於在處理南極洲和Greenland地形時，前者給出的是加上冰蓋層之後的高程，後者給出的是岩床的高程。

對於每個版本又細分為 grid-registered和cell-registered，其中grid-registered是權威版本，cell-registered是衍生版本，因而推薦下載使用grid-registered版本。

在每個子版本下又有多種資料格式， netCDF，binary， xyz，tiff。

我選擇的是grid-registered版本的binary格式。

    下載：

http://www.ngdc.noaa.gov/mgg/global/relief/ETOPO1/data/ice_surface/grid_registered/binary/etopo1_ice_g_i2.zip

解壓得.bin為二進位制網路檔案，.hdr為頭段檔案

將.bin拷貝至/usr/local/GMT-4.5.9/share/dbase：

sudo cp etopo1_ice_g_i2.bin /usr/local/GMT-4.5.9/share/dbase

修改grdraster.info：

9 "ETOPO1 Ice Surface"  "m"  -R-180/180/-90/90  -I1m  GG i 1  0 -32768 etopo1_ice_g_i2.bin L

同理對於bedrock版本的網格資料，其grdraster.info為：

10 "ETOPO1 Bedrock" "m"  -R-180/180/-90/90 -I1m GG i 1 0 -32768  etopo1_bed_g_i2.bin L

如果下載的是netCDF格式的網格檔案，需要利用如下命令將資料轉換為binary格式：

grdreformat ETOPO1_Ice_g_gmt4.grd etopo1_ice_g_i2.bin=bs -N -V

畫圖測試：

#!/bin/bash

verbose=-V

#verbose=

grdraster 9 -Rg -I2m -Gout.grd $verbose

makecpt -Cglobe -T-10500/8000/1000 -Z $verbose > colors.cpt

grdimage out.grd -Ba60g30 -Rg -Yc -Xc -JN0/25c -Ccolors.cpt -K $verbose > etopo5.ps

psscale -Ba2500f500::/:"m": -Ccolors.cpt -D12.5c/-2c/15c/.35ch -O $verbose >> etopo5.ps

rm out.grd colors.cpt

 

GTOPO30：

    下載：

ftp://edcftp.cr.usgs.gov/data/gtopo30/global/

解壓複製：將所有DEM檔案拷貝至/share/dbase

修改grdraster.info：

# GTOPO30 N90-N40

11 "GTOPO30 W180N90" "m" -R-180/-140/40/90 -I0.5m PG i 1 0 -9999 W180N90.DEM B

12 "GTOPO30 W140N90" "m" -R-140/-100/40/90 -I0.5m PG i 1 0 -9999 W140N90.DEM B

13 "GTOPO30 W100N90" "m" -R-100/-60/40/90 -I0.5m PG i 1 0 -9999 W100N90.DEM B

14 "GTOPO30 W060N90" "m" -R-60/-20/40/90  -I0.5m PG i 1 0 -9999 W060N90.DEM B

15 "GTOPO30 W020N90" "m" -R-20/20/40/90 -I0.5m PG i 1 0 -9999 W020N90.DEM B

16 "GTOPO30 E020N90" "m" -R20/60/40/90 -I0.5m PG i 1 0 -9999 E020N90.DEM B

17 "GTOPO30 E060N90" "m" -R60/100/40/90 -I0.5m PG i 1 0 -9999 E060N90.DEM B

18 "GTOPO30 E100N90" "m" -R100/140/40/90  -I0.5m PG i 1 0 -9999 E100N90.DEM B

19 "GTOPO30 E140N90" "m" -R140/180/40/90  -I0.5m PG i 1 0 -9999 E140N90.DEM B

 

# GTOPO30 N40-S10

20 "GTOPO30 W180N40" "m" -R-180/-140/-10/40 -I0.5m PG i 1 0 -9999 W180N40.DEM B

21 "GTOPO30 W140N40" "m" -R-140/-100/-10/40 -I0.5m PG i 1 0 -9999 W140N40.DEM B

22 "GTOPO30 W100N40" "m" -R-100/-60/-10/40 -I0.5m PG i 1 0 -9999 W100N40.DEM B

23 "GTOPO30 W060N40" "m" -R-60/-20/-10/40 -I0.5m PG i 1 0 -9999 W060N40.DEM B

24 "GTOPO30 W020N40" "m" -R-20/20/-10/40  -I0.5m PG i 1 0 -9999 W020N40.DEM B

25 "GTOPO30 E020N40" "m" -R20/60/-10/40 -I0.5m PG i 1 0 -9999 E020N40.DEM B

26 "GTOPO30 E060N40" "m" -R60/100/-10/40  -I0.5m PG i 1 0 -9999 E060N40.DEM B

27 "GTOPO30 E100N40" "m" -R100/140/-10/40 -I0.5m PG i 1 0 -9999 E100N40.DEM B

28 "GTOPO30 E140N40" "m" -R140/180/-10/40 -I0.5m PG i 1 0 -9999 E140N40.DEM B

 

# GTOPO30 S10-S60

29 "GTOPO30 W180S10" "m" -R-180/-140/-60/-10 -I0.5m PG i 1 0 -9999 W180S10.DEM B

30 "GTOPO30 W140S10" "m" -R-140/-100/-60/-10 -I0.5m PG i 1 0 -9999 W140S10.DEM B

31 "GTOPO30 W100S10" "m" -R-100/-60/-60/-10 -I0.5m PG i 1 0 -9999 W100S10.DEM B

32 "GTOPO30 W060S10" "m" -R-60/-20/-60/-10 -I0.5m PG i 1 0 -9999 W060S10.DEM B

33 "GTOPO30 W020S10" "m" -R-20/20/-60/-10 -I0.5m PG i 1 0 -9999 W020S10.DEM B

34 "GTOPO30 E020S10" "m" -R20/60/-60/-10  -I0.5m PG i 1 0 -9999 E020S10.DEM B

35 "GTOPO30 E060S10" "m" -R60/100/-60/-10 -I0.5m PG i 1 0 -9999 E060S10.DEM B

36 "GTOPO30 E100S10" "m" -R100/140/-60/-10 -I0.5m PG i 1 0 -9999 E100S10.DEM B

37 "GTOPO30 E140S10" "m" -R140/180/-60/-10 -I0.5m PG i 1 0 -9999 E140S10.DEM B

 

# GTOPO30 S60-S90

38 "GTOPO30 W180S60" "m" -R-180/-120/-90/-60 -I0.5m PG i 1 0 -9999 W180S60.DEM B

39 "GTOPO30 W120S60" "m" -R-120/-60/-90/-60 -I0.5m PG i 1 0 -9999 W120S60.DEM B

40 "GTOPO30 W060S60" "m" -R-60/0/-90/-60  -I0.5m PG i 1 0 -9999 W060S60.DEM B

41 "GTOPO30 W000S60" "m" -R0/60/-90/-60 -I0.5m PG i 1 0 -9999 W000S60.DEM B

42 "GTOPO30 E060S60" "m" -R60/120/-90/-60 -I0.5m PG i 1 0 -9999 E060S60.DEM B

43 "GTOPO30 E120S60" "m" -R120/180/-90/-60 -I0.5m PG i 1 0 -9999 E120S60.DEM B

說明

將全球分為33個區塊，如下圖所示，北緯90度到南緯60度分為3份，每份緯向50度，經向40度，共計3*9=27個區塊；南極附近為南緯60度到90度，經向60度，計6個區塊，共計33個區塊。

取樣間隔30弧秒，即0.008333333333333度，全球分為21600行、43200列。

座標參考為WGS84；

其值為垂直方向高出海平面的距離，即海拔，單位為m；

海拔範圍-407到8752，這裡不包含海洋深度資訊，負值應該是大陸架的海拔；

海洋處海拔標記為-9999，海岸線以為大陸海拔至少為1；

小於1平方千米的島不考慮；

檔案以該區塊的左上角座標為檔名；

檔案字尾：DEM為資料高程模型；其他的檔案都是一些輔助資訊；

維度50度*經度40度檔案大小：57600000 bytes

維度30度*經度60度檔案大小：28800000 bytes

其他說明:

這些DEM檔案為Pixel Registration格式

資料值儲存為有符號16位整型

這些資料都是Big-endian的

繪圖測試:

#!/bin/bash

ID=11

Range=-180/-140/40/90

B=a10g10

Inc=0.5m

PS=${ID}.ps

verbose=-V

#verbose=

grdraster $ID -R$Range -I$Inc -Gout.grd $verbose

makecpt -Ctopo -T1/8000/1000 -Z $verbose > colors.cpt

grdimage out.grd -B$B -R$Range -Yc -Xc -JQ15c -Ccolors.cpt $verbose > $PS 

rm out.grd colors.cpt

這裡只繪製一個第一個區塊的地形，若需要同時繪製多個區塊，需要多次呼叫grdraster從不同的網格檔案中提取網格，然後用grdpaste貼上成一個網格再繪圖。

南極區域

在南極區域還有一個特殊的高程資料ANTARCPS.DEM，其包含了S60-S90區域的六個區塊的資料，但是是經過極座標方位立體投影的。因而橫座標範圍為（-2700，2700），縱座標範圍為（-2700，2700）。這個處理起來要麻煩一點。

 

SRTM:

SRTM地形資料主要包含兩類資料：SRTM3和SRTM1。

    SRTM3精度為3弧秒，即90m一個點，包括非洲、北美、南美、歐亞、澳大利亞以及部分島嶼；

    SRTM1精度為1弧秒，即30m一個點，僅限美國地區；

所有的資料被分為1度*1度的單元，共計上萬個數據檔案。檔名代表了該單元左下角的位置；不同網站下載的資料檔名不同，但檔名類似N37W105的 資料檔案，其單元左下角的座標為(37N,105W)；每個SRTM3資料檔案包含了1201*1201個取樣點；高程資料單位為m，參考水準面為 WGS84大地水準面模型；網格劃分採用gridline registration方式，因而單元的南北邊及左右邊與相鄰單元的邊重合。

資料下載：

http://srtm.datamirror.csdb.cn/admin/datademMain.jsp

中國科學院計算機網路資訊中心提供：Geotiff 和 ARC ASCII格式資料，資料是5度*5度一個單元，而不是標準的1度*1度，免去了使用者自己做網格合併的麻煩；GMT可以處理ARC ASCII格式，但無法處理Geotiff格式，所以下面只用ARC ASCII格式；在中國區域，還提供一些衍生資料，如坡度、坡向、坡位，資料格式為img格式，細節未知；

資料轉換：

將ARC ASCII格式轉換為GMT grid檔案：

xyz2grd srtm_63_02.asc -Gstrm_63_02.grd -E -V

繪圖例項：

R=130/135/50/55

B=a1f0.5WSen

J=Q4i

xyz2grd srtm_63_02.asc -Gstrm_63_02.grd -E -V

#xyz2grd N52E130.hgt -GN52E130.grd -I3c/3c -R$R -N-32768 -ZTLhw -V

grd2cpt N52E130.grd -Cnrwc.cpt > colors.cpt

grdimage N52E130.grd -R$R -Ccolors.cpt -B$B:.SRTM3: -J$J -K > test.ps

如果是USGS下載的標準資料，例子：較大的崑崙山區域，下載的資料為N37E085.hgt和N37E086.hgt 2個檔案。

#首先使用xyz2grd將二進位制的hgt檔案轉成GMT使用的grd檔案，或者nc字尾也可。

xyz2grd N37E085.hgt -Gg001.grd -I3c -R85/86/37/38 -N-32768 -ZTLhw

xyz2grd N37E086.hgt -Gg002.grd -I3c -R86/87/37/38 -N-32768 -ZTLhw

#合併grd檔案

grdpaste g001.grd g002.grd -Ga1.grd

makecpt -Cglobe -T-500/7000/500 -Z > mydata.cpt

grdgradient a1.grd -A30 -Gchina.grd.int -Ne0.8 -V -M

#使用grdimage繪製彩圖

grdimage a1.grd -R85/87/37/38 -Yc -Xc -JM10i -Cmydata.cpt  -Ichina.grd.int -K > world1.ps

pscoast -R -J -Di -N1 -A1000 -B0.5g0.5/0.5g0.5 -W0.25p  -O -K  >>world1.ps

psscale -Ba2500f500::/:"m": -Cmydata.cpt -D12.5c/-2c/15c/.35ch -O -U"author:YangLei">> world1.ps

del *.cpt *.int

 

ASTER GDEM

ASTER GDEM，即先進星載熱發射和反射輻射儀全球數字高程模型，與SRTM一樣為數字高程DEM，其全球空間解析度為30米。該資料是根據 NASA的新一代對地觀測衛星Terra的詳盡觀測結果製作完成的。其資料覆蓋範圍為北緯83°到南緯83°之間的所有陸地區域，達到了地球陸地表面的99%。（摘自百度百科）。號稱是“迄今最完整的全球地形資料”。

ASTER GDEM：空間解析度：1弧度秒 （約30 米），精度：垂直精度20米，水平精度30米，SRTM資料的緯度覆蓋範圍是[-60,60]，ASTER GDEM資料的緯度覆蓋範圍為[-83,83]；SRTM的空間解析度一般為90m，只有美國境記憶體在空間解析度為30m的資料；ASTER GDEM的空間解析度為30m；

網站：http://gdem.ersdac.jspacesystems.or.jp/

網站需要註冊，即可免費下載資料；

點選“Search”即可進入資料選擇介面；

資料被劃分為1度*1度的區塊，每個區塊檔名以區塊左下角的經緯度值來命名；

資料選擇介面有四種方式用來選擇區塊，介面的設計不太習慣，注意看英文；

資料可以單獨下載，也可以批量下載；批量下載貌似需要處理一段時間；

下載的檔案為zip壓縮檔案，包含了兩個檔案：dem和num，二者均為GeoTiff格式；

dem資料即為數字高程資料；

num資料也稱為QA檔案，即質量評估檔案，其包含了每一點的高程資料的來源等資訊，一般用不到；

使用GDAL提供的工具將其轉換為GMT可識別的netCDF格式：

gdal_translate -of GMT ASTGTM2_N36E111_dem.tif ASTGTM2_N36E111_dem.grd