Matlab下地形图绘图包m

您所在的位置:网站首页 国家地图地形图怎么画 Matlab下地形图绘图包m

Matlab下地形图绘图包m

2024-07-11 04:03| 来源: 网络整理| 查看: 265

〇、引言

m_map是Matlab下用于绘制地图(地形图)的工具箱,和GMT有些相似。

1. 文档

m_map官网:https://www.eoas.ubc.ca/~rich/mapug.html

用法可以参考百度文库中的官网翻译版:M_Map1.4用户指南 https://wenku.baidu.com/view/32b9c4c8d4d8d15abf234e06.html

也可以参考CSDN这位老兄的翻译版:m_map中文文档 https://blog.csdn.net/good_learner_1/article/details/88767974?utm_medium=distribute.pc_relevant.none-task-blog-baidujs_baidulandingword-3&spm=1001.2101.3001.4242

2. 推荐教程

再推荐本人写的几个详细教程:

m_map绘制晕渲(shaded relief)地形图m_map绘图包绘制高分辨率海岸线、国界线与河流m_map绘图之添加线段比例尺函数m_rulerm_map如何绘制矩形等距投影图m_map绘制点、线m_map绘图添加遥感图片m_map导入本地地形数据m_map绘制多波束数据 一、官网绘图展示

m_map官网绘图展示

二、下载与安装 1. 下载

进入m_map官网:https://www.eoas.ubc.ca/~rich/mapug.html 下载tar压缩文件或者zip压缩文件。 在这里插入图片描述

2. 安装

将m_map工具包解压后放在matlab的toolbox中,添加路径即可。 对matlab不了解的话,可参考:M_Map绘图笔记——安装(一) https://blog.csdn.net/wokaowokaowokao12345/article/details/88354118

三、高分辨率地形与海岸线数据加载

m_map自带的底图数据分辨率只有1°,较低,在绘制区域地形图时分辨率无法满足要求,需要额外加载高分辨率的地形数据。如下:

1. 高分辨率海岸线数据文件下载 从NOAA官网下载1分(1’)分辨率的海岸线数据:

https://www.ngdc.noaa.gov/mgg/shorelines/gshhs.html 在这里插入图片描述

下载二进制(bin)格式的zip压缩文件: 在这里插入图片描述 2. 高分辨率地形数据文件下载 从NOAA官网下载1分(1’)分辨率的地形数据

https://www.ngdc.noaa.gov/mgg/global/relief/ETOPO1/data/ice_surface/grid_registered/binary/ 在这里插入图片描述

有32位浮点型(4字节),即f4,和16位整型(2字节),即i2两种格式。推荐选第二种:i2。 3. 数据解压 下载好这两个数据后分别解压,建议放在matlab的安装目录下,我在matlab的安装目录下新建了一个work文件夹,再把下载的m_map文件解压后放在work目录下。然后,在m_map目录下建立了etopo1和coastline文件夹(文件名随便起),本别放置地形数据和海岸线数据。 在这里插入图片描述 4. 路径设置

然后需要设置二者(1分的地形与海岸线数据)的路径,才可以自动调用。打开以下m文件m_etopo2.m:

把其中的

PATHNAME='/ocean/rich/more/mmapbase/etopo1/';

改成自己的地形数据放置路径,比如我的是/m_map/etopo1/:

PATHNAME='D:/install/MATLAB2019a/work/m_map/etopo1/';

这样,地形数据路径就设置好了。

同理,打开m_gshhs.m,将路径改为自己的海岸线数据放置路径,比如我的是/m_map/coastline/:

FILNAME='D:/install/MATLAB2019a/work/m_map/coastline/';

到此安装就大功告成!!!

四、用法 1. 设置投影方式 m_proj(‘mercator’); 2. 查看投影方式类别 m_proj(‘set’); 3. 如何选择投影方式

全球:

m_proj(‘mercator’); m_proj(‘miller’); m_proj(‘hammer’); m_proj(‘Mollweide’);

国家:

m_proj(‘robinson’);

纬度跨度大,经度跨度小:

m_proj(‘sinusoidal’); m_proj(‘Mollweide’);

小区域(如地中海,中国):

m_proj(‘conic’);

极地:

m_proj(‘stereographic’); 4. 网格与轴设置

m_grid是关于绘制网格线和边框的设置。 查看功能参数:

m_grid get;

或者

m_grid(‘get’);

几个常用参数: box是边框设置:on有,off无,fancy为黑白框 xtick:轴显示刻度线数 xticklabels:[数组]显示特定的标签 xlabeldir: 标签显示方向:middle正常横向,end竖向,ytickdir与之相反 ticklen:刻度线长度或者黑白框宽度 tickdir:刻度线方向。对于fancy来说in表示黑白框out表示线条框 tickstyle:坐标带不带NSEW方向color边框颜色 gridcolor:网格线颜色 backgroundcolor:背景色 linewidth:线宽 fontsize:轴字号 fontname :轴字体 zaxislocation:轴位置 linest:网格线型,有-- -. : - 四种

5. 添加海岸线、国界线与河流

打开函数文件可以查看说明:

help m_gshhs

摘取如下:

% RES: A one-char string (optionally 2 or 3) % % First char: resolution - one of % 'c' crude % 'l' low % 'i' intermediate % 'h' high % 'f' full % % Second char: type - one of % 'c' GSHHS coastline (default) % 'b' WDB Border % 'r' WDB River % % Third char - if 2nd char is 'b': % '1' Country borders % '2' State/Province and Country borders % - if 2nd char is 'r': '1','2','3','4' % add successively more tributaries

c、b、r分别表示海岸线,国界线与河流 c、l、i、h、f 分别表示粗糙/低/中/高/满分辨率 具体使用时,c/b/r与 c/l/i/h/f 二者组合使用 例子:

m_gshhs('lc','patch','k'); % 低分辨率海岸线,陆地填充:黑色 m_gshhs('ib’,’color','k'); % 中分辨率国界线,黑色 m_gshhs('fr','color','b'); % 满分辨率河流,蓝色

当添加国界线b时,可以选择第三个参数进一步显示省界/州界:1表示国界,2表示省界/州界。 同理,当添加河流r时,可以选择第三个参数进一步显示河流分支:1, 2, 3, 4分别为更加细化的河流分支。

6. 绘制各种类型的图形 m_plot(LONG,LAT,...line properties) m_line(LONG,LAT,...line properties) % 线条 m_text(LONG,LAT,‘string’) % 文本 m_quiver(LONG,LAT,U,V,S) % 矢量箭头 m_patch(LONG,LAT,..patch properties) % 色块 m_annotation(‘line’,LON,LAT) % 注释 m_contour(LONG,LAT,VALUES) % 等值线 m_contourf(LONG,LAT,VALUES) % 等值线填充 m_image(LON,LAT,DATA) % 影像 m_pcolor(LON,LAT,DATA) % 色块 [IM,X,Y]=m_image(LON,LAT,DATA); % 将经纬度坐标转为XY坐标 m_shadedrelief(X,Y,IM,'coords','map’) % 地形渲染 m_etopo2(OPTION) % 画地形图 m_elev(OPTION) % 生成高程数据 m_ruler([.5 .8],.9,'tickdir','out','ticklen',[.007 .007]); % 添加距离比例尺 m_northarrow(-123-4.5/60,49+19.5/60,1/60,'type',4,'aspect',1.5); % 添加指北针 五、例子 5.1 本人写的例子 m_map绘制晕渲(shaded relief)地形图m_map绘图包绘制高分辨率海岸线、国界线与河流m_map绘图之添加线段比例尺函数m_rulerm_map如何绘制矩形等距投影图m_map绘制点、线m_map绘图添加遥感图片m_map导入本地地形数据m_map绘制多波束数据 5.2 官网例子

下面看几个官网例子。

例1:One Ocean Projection m_proj('azimuthal equal-area','radius',156,'lat',-46,'long',-95,'rot',30); ax1=subplot(2,2,1,'align'); m_coast('patch','r'); m_grid('xticklabel',[],'yticklabel',[],'linestyle','-','ytick',[-60:30:60]); ax2=subplot(2,2,2,'align'); m_elev('contourf',[-7000:1000:0 500:500:3000],'edgecolor','none'); colormap(ax2,[m_colmap('blues',70);m_colmap('gland',30)]); caxis(ax2,[-7000 3000]); m_grid('xticklabel',[],'yticklabel',[],'linestyle','-','ytick',[-60:30:60]); ax3=subplot(2,2,3,'align'); colormap(ax3,[m_colmap('blues',70);m_colmap('gland',30)]); caxis(ax3,[-7000 3000]); m_elev('image'); m_grid('xticklabel',[],'yticklabel',[],'linestyle','-','ytick',[-60:30:60]); ax4=subplot(2,2,4,'align'); colormap(ax4,[m_colmap('blues')]); caxis(ax4,[-8000 000]); m_elev('shadedrelief','gradient',.5); m_coast('patch',[.7 .7 .7],'edgecolor','none'); m_grid('xticklabel',[],'yticklabel',[],'linestyle','-','ytick',[-60:30:60]); ha = axes('Position',[0 0 1 1],'Xlim',[0 1],'Ylim',[0 1],'Box','off',... 'Visible','off','Units','normalized', 'clipping' , 'off'); text(0.5, 0.98,'This projection shows all oceans connected to each other',... 'horiz','center','fontsize',20);

在这里插入图片描述

例2:Bathymetry load /ocean/rich/home/dens14/VENTS lp=load('/ocean/rich/home/dens14/Linep_201402'); m_proj('lambert','long',[-160 -115],'lat',[32 60]); [CS,CH]=m_etopo2('contourf',[-7000:1000:-1000 -500 -200 0 ],'edgecolor','none'); m_gshhs_f('patch',[.7 .7 .7],'edgecolor','none'); h1=m_line(vents.lon,vents.lat,'marker','s','color',[0 .5 0],... 'linest','none','markerfacecolor','w','clip','point'); h2=m_line(lp.POS(:,2),lp.POS(:,1),'marker','o','color','r','linewi',2,... 'linest','none','markersize',8,'markerfacecolor','w'); m_grid('linest','none','tickdir','out','box','fancy','fontsize',16); legend([h1(1),h2(1)],'Known Hydrothermal vents','Line-P Stations','location','southwest'); colormap(m_colmap('blues')); caxis([-7000 000]); [ax,h]=m_contfbar([.55 .75],.8,CS,CH,'endpiece','no','axfrac',.05); title(ax,'meters') set(gcf,'color','w'); % otherwise 'print' turns lakes black

Bathymetry

例3:Shaded Relief (Example 1) m_proj('lambert','lat',[5 24],'long',[105 125]); set(gcf,'color','w') % Set background colour before m_image call caxis([-6000 0]); colormap(flipud([flipud(m_colmap('blues',10));m_colmap('jet',118)])); m_etopo2('shadedrelief','gradient',3); m_gshhs_i('patch',[.8 .8 .8]); m_grid('box','fancy'); ax=m_contfbar(.97,[.5 .9],[-6000 0],[-6000:100:000],'edgecolor','none','endpiece','no'); xlabel(ax,'meters','color','k');

South China Sea

例4:Shaded Relief (Example 2) 注意:这个需要开发者的多波束数据 m_proj('utm','ellipse','grs80','zone',10,'lat',[49+15.7/60 49+21/60],... 'long',[-123-15/60 -123-3/60]); % Uses multibeam bathymetry with 10m horizontal resolution % Already regularly gridded in UTM coords with vector x2/y2, and % matrix Z2. caxis([-150 0]); colormap([m_colmap('water',128)]); m_shadedrelief(x2,y2,-Z2,'lightangle',-45,'gradient',8,'coord','z'); % Add some contours hold on; [cs,h]=contour(x2,y2,Z2,[0:20:150],'color','k'); clabel(cs,h,'fontsize',6); hold off; % Land parts from a previously saved high-resolution coastline col=[255 214 140]/255; % CHS chart land colour m_usercoast('/ocean/rich/more/mmapbase/bcgeo/PNW.mat','patch',col); m_usercoast('/ocean/rich/more/mmapbase/bcgeo/PNWrivers.mat','patch',col); % Lat/long AND a UTM grid m_grid('tickdir','out','fontsize',12,'linest','none','xaxisloc','top','yaxisloc','right'); m_utmgrid('xcolor','b','ycolor','b','linest','-'); m_ruler([.5 .8],.9,'tickdir','out','ticklen',[.007 .007]); m_northarrow(-123-4.5/60,49+19.5/60,1/60,'type',4,'aspect',1.5); xlabel('Vancouver Harbour','color','k');

在这里插入图片描述

例5:SAR image of internal waves (HDF-5 format) 注意:这个需要HDF格式的SAR数据 fname='SAR_IMP_1PNESA19920724_190439_00000018C086_00199_05354_0000.h5'; %h5disp(fname) % See the structure titlestr=h5readatt(fname,'/metadata/MPH','STATE_VECTOR_TIME'); datsize=double([ h5readatt(fname,'/bands/Amplitude','raster_width') ... h5readatt(fname,'/bands/Amplitude','raster_height')]); tielat=h5read(fname,'/tie_point_grids/latitude'); tielon=h5read(fname,'/tie_point_grids/longitude'); stp=[h5readatt(fname,'/tie_point_grids/latitude','sub_sampling_x') ... h5readatt(fname,'/tie_point_grids/latitude','sub_sampling_y') ]; % Pull out a subsection in by [2000 2500] from one corner and in % [600 2000] from the opposite corner istart=[2000 2500]; strd=[3 3]; cnt=fix([(datsize(1)-istart(1)-600)/strd(1) (datsize(2)-istart(2)-2000)/strd(2)]); % ...and read. subimg=h5read(fname,'/bands/Amplitude',istart,cnt,strd); % ....smooth it a bit... subf=filter2(ones(3,3)/9,subimg); % Now generate lat/lon for all pixels by interpolating from % the tie points. Ty=[0:size(tielat,2)-1]*stp(2)+1; Tx=[0:size(tielat,1)-1]*stp(1)+1; Iy=istart(2)+[0:size(subimg,2)-1]*strd(2); Ix=istart(1)+[0:size(subimg,1)-1]*strd(1); sublat=interp2(Ty',Tx,tielat,Iy',Ix); sublon=interp2(Ty',Tx,tielon,Iy',Ix); % Now make the map m_proj('lambert','lon',[-123-25/60 -122-40/60],'lat',[48+42/60 49+9/60]); m_pcolor(sublon,sublat,subf);shading flat; m_grid('box','fancy','tickdir','out'); m_ruler(1.03,[.15 .5],'ticklen',[.01]); caxis([0 450]); colormap(gray); title(titlestr)

SAR image of internal waves (HDF-5 format)

注:

内容部分参考了本人的新浪博客文章:m_map简介 http://blog.sina.com.cn/s/blog_15c32b1470102z6f2.html



【本文地址】


今日新闻


推荐新闻


CopyRight 2018-2019 办公设备维修网 版权所有 豫ICP备15022753号-3