在GMT中，我們可以加載nc或者grd等格網數(shù)據(jù)，而遙感影像一般是geotiff格式，因此無法在GMT中直接加載遙感影像。本人也一直被這個問題困擾，直到最近看到的一篇微信推送才解決了這個問題。這里我主要介紹兩類方法：多波段組合和單波段加載。本節(jié)主要介紹多波段的加載。

【1】多波段組合（已經實現(xiàn)）

所需數(shù)據(jù)：經過前期處理好，經過波段合成的真彩色影像。

下一步需要將其進行坐標轉換，即將其轉換成wgs84地理坐標。有很多的方法可以實現(xiàn)這一轉換，但是我在這里像介紹利用gdal_translate工具直接轉換。我在這里使用時gmt5版本，其中沒有包含這一函數(shù)。而gmt6新版本中包含了這一函數(shù)，可以直接調用

坐標轉換的命令：

gdalwarp -t_srs EPSG:4326 test.tif test_wgs.tif

test.tif為轉換前圖像名稱，test_wgs.tif是轉換后的影像名稱，直接修改即可。

另外如果需要進行小區(qū)域的裁剪，可以使用一下命令：

gdalwarp -t_srs EPSG:4326 -te 116 40 117 41.2 test_wgs.tif teat_wgs_out.tif

需要修改的是范圍大小，輸入范圍框的左下角和右上角經緯度坐標。

然后利用一下的代碼即可完成遙感影像的加載：

set R=116/117/40/41.2

set J=M12c

set PS=map2.ps

gmt gmtset FORMAT_GEO_MAP=ddd:mm:ss

gmt psbasemap -R%R% -J%J% -Bpxa50mf -Bpya50mf -BnWSe -Xc -Yc? -K >%PS%?

gmt pscoast -R%R% -J%J% -O -A500 -W -I1 -K? -Lx1.0i/0.15i+c120/36+w5k >>%ps%

gmt grdimage? teat_wgs_out.tif+b0? teat_wgs_out.tif+b1 teat_wgs_out.tif+b2 -J%J% -R%R% -K -O >> %ps%

gmt psxy ucas.txt -J%J% -R%R% -Sc0.05c -W1,red -O -K>> %ps%

REM # Add jpg

REM gmt psimage app.png -Dx0.3/6.5+w11.4c -O -K >> %ps%

REM # Add small map in Up-Right corner

REM gmt pscoast -Rg -JG120/41/2c? -Dc -A5000 -Gpink -Sblue? -O -X2.5c -Y6.5c -K >>%ps%

REM echo? 120:26:30 36:02:57 | gmt psxy -Sa0.1c -Gred -R -J -O? >> %ps%

REM gmt psconvert %ps% -A -P -Tg

創(chuàng)作不易，期待大家的評論收藏分享！！