如何用GMT绘制GTOPO30的DEM地形图.docx

转载如何用GMT绘制GTOPO30的DEM地形图 (2012-05-17 19:06:57)转载标签：转载多谢分享！原文地址：如何用GMT绘制GTOPO30的DEM地形图作者：forever上次提到了GTOPO30的命名法则以及下载地址，弄完，当然就要开始绘图了。一、解压缩把下载好的*.tar.gz文件放到一个文件夹里，然后在终端进入，输入ls *.tar.gz | gawk print tar zxvf,$1; | sh，就会在当前目录下把各个压缩包里的东东解压出来。一般来说，每个包都含有8个文件，后缀分别是DEM, DMW, GIF, HDR, PRJ, SCH, SRC, STX，在使用GMT绘图时只需用到DEM文件，但建议其它格式文件也都保留，另外的一些绘图工具要用到（比如MATLAB好像要用到HDR文件）。二、定义文件环境之所以这样做是因为要用到grdraster，它从DEM文件中抽取需求范围内的数据并形成grd文件，而使用该工具的前提是要对DEM数据的文件代码、别称、及数据类型进行描述，这些描述必须存储在文件中。假设GMT的安装目录为/usr/local/gmt，那么会在/usr/local/gmt/share/dbase里。进入到这个目录，并把解压缩好了的DEM数据拷贝进该目录，然后编辑文件，按要求加入DEM文件的描述。对于中国来说，E060N90, E100N90, E060N40, E100N40这4个DEM文件就足够了。我的是这么写的：29 GTOPO30 E060N90m-R60/100/40/90-I0.5mP i 10-9999E060N90.DEM B30 GTOPO30 E060N40m-R60/100/-10/40-I0.5mP i 10-9999E060N40.DEM B32 GTOPO30 E100N90m-R100/140/40/90-I0.5mP i 10-9999E100N90.DEM B33 GTOPO30 E100N40m-R100/140/-10/40-I0.5mP i 10-9999E100N40.DEM B分别是文件代码、文件别称、单位、范围等等，倒数第二列是文件名，最后一个字符“B”在我的系统上必须要加，不然出的图就像是电视没信号，都是雪花（似乎见过一些朋友没加但也用得挺好的）。写完后保存退出。三、绘图假设绘一个四川的地图，范围为97/110/25/35，这个范围涉及到两个DEM数据文件E060N40, E100N40，可先用grdraster分块取出然后用grdpaste合并之。步骤是grdraster 30 -R97/100/25/35 -I0.5m -Gsc1.grdgrdraster 33 -R100/110/25/35 -I0.5m -Gsc2.grd30和33是文件中对应的文件代码，-R中所指定的范围不要超过DEM文件的范围，且二者需要有共同的边界，-I指定采样间隔，-G输出文件grdpaste sc1.grd sc2.grd -Gsc.grd嗯，将两个文件合并在一起，输出为sc.grd。然后制作调色板（CPT）文件，如果已经有了grd文件，则可使用grd2cpt工具来制作grd2ctp sc.grd -Ctopo -S-200/6000/200 -Z sc.cpt-C指定生产的cpt文件类型，GMT里内置了很多的，可以查看帮助文档，-S指定Z（高程）的开始和结束的范围以及间隔，-Z生产连续颜色的CPT文件。然后绘图：grdimage -JM6i -R97/110/25/35 -B2/2 -Csc.cpt sc.grd -K -P sc.ps-J指定投影方式，-R指定绘图范围，-B地图旁的刻度，-C调色板文件然后加入色标，用psscale来加：pssca