(GSHHG) Shapefile格式转换海岸线数据
这个例子展示了如何:
提取海岸线的一个子集数据从全球有层次高分辨率地理(GSHHG)数据集
控制多边形特性添加湖泊和其他内陆水体内部多边形戒指(“洞”)
保存已修改的数据集在MATLAB®shapefile供以后使用,或用于出口的地理信息系统
全球有层次高分辨率地理(GSHHG;以前全球有层次高分辨率海岸线或GSHHS)的数据集,由保罗·韦塞尔和沃尔特·h·f·史密斯提供了一组一致的层级结构安排封闭的多边形。可以用它们来构建基础地图,或者在应用程序或分析涉及到操作,比如地理搜索或海岸线的统计特性。
步骤1:定义一个工作目录
这个例子中创建一些临时文件并使用变量workingFolder来表示他们的位置。这里使用的价值取决于的输出tempdir命令,但是您可以轻松地定制。
workingFolder = tempdir;
步骤2:GNU®解压和索引的粗糙的分辨率GSHHG层
GSHHG可在大范围的空间分辨率。金宝搏官方网站这个示例使用情况数据,从二进制文件gshhs_c.b。GNU压缩副本中包含该文件的映射工具箱™数据文件夹,在MATLAB的道路。
使用MATLABgunzip函数来减压gshhs_c.b.gz和创建文件gshhs_c.b在位置所示workingFolder。然后创建一个索引文件,gshhs_c.i,在同一文件夹中。一般来说,拥有一个索引文件有助于加速后调用gshhs函数。请注意,当您使用“方法createindex”选项,gshhs不提取数据。
文件= gunzip (“gshhs_c.b.gz”,workingFolder);文件名={1}文件;indexfile = gshhs(文件名,“方法createindex”);
步骤3:为南美GSHHG导入数据
选择数据为一个特定的经度四边形和导入映射工具箱“geostruct”数组:
latlim = 15 [-60];lonlim = (-90 - -30);S = gshhs(文件名、latlim lonlim);
如果你已经完成了提取数据,您可以删除解压GSHHS文件和索引文件。
删除(文件名)删除(indexfile)
步骤4:检查数据集
检查geostruct数组的第一个元素年代。除了纬度和朗坐标数组,注意存在的各种属性字段。
(1)
ans =结构体字段::边界框(几何:“多边形”大小)[2 x2双]纬度:[-53.0004 -53 -52.5017 -52.7963 -52.0434 -52.0838 -52.3638 -51.9692 -52.4763 -51.7925 -51.5217 -51.1569 -51.4650 -51.7722 -51.5550 -52.0364 -51.6342 -52.0427 -51.6283 -51.1117 -50.9233 -50.6942 -50.3871 -50.8513 -50.3750……朗:[-73.3617 -73.3626 -72.8909 -73.7034 -73.7392 -72.9983 -72.9884 -72.7708 -72.9184 -72.4792 -73.1442 -72.9335 -73.4259 -72.6031 -73.3626 -73.1710 -73.4026 -73.3898 -73.9309 -73.6609 -74.2375 -73.3392 -73.5700 -74.0367 -73.9624……西方]南:-53.9004北:71.9942:191.8947:325.2054:3.7652 e + 07级:1 LevelString:“土地”NumPoints: 971 FormatVersion: 3资料来源:”“全球价值调查”主要根据“CrossesGreenwich: 0 GSHHS_ID: 1
的海岸线GSHHS包括四个层次:
一级-“土地”
二级-“湖”
三级-“湖岛”
四级-“池塘在湖岛”
查看数据导入包括水平。的水平数字字段包含数值水平。
水平= [S.Level];独特的(水平)
ans =1×31 2 3
的LevelString场提供了解释。例如,
年代(104).LevelString
ans =“湖”
显示功能104是一个湖(2级特性)。
在这个例子中,由于低分辨率或者空间构造子集,没有四级功能。
第五步:前两个层次上提取到单独的Geostruct数组
这个例子中操纵GSHHS顶部两个级别的层次结构,每一个“湖”插入周围的土地质量。
提取GSHHS一级(大陆和海洋岛屿的海岸线外):
L1 = S(水平= = 1);
海岸线提取2级(一级多边形内的湖泊和海洋):
L2 = S(水平= = 2);
看到他们的空间关系,可以映射一级边缘如红线蓝线和2级边缘:
图axesm (“墨卡托”,“MapLatLimit”latlim,“MapLonLimit”,lonlim) gridm;mlabel;plabel geoshow ([L1。Lat], [L1.Lon],“颜色”,“蓝”(L2) geoshow (。纬度], [L2.Lon],“颜色”,“红色”)tightmap
第六步:将2级多边形合并到1级
定义一个匿名谓词函数来检测限定框十字路口(返回true,如果一对边界盒相交和假否则)。输入一个和B2×2限定框矩阵的形式吗
[min(朗)最小(lat)最大(朗)马克斯(lat)]。
boxesIntersect =…@ (A, B)(~(任何(一个(2:)< B(1:)) | |任何(B (2:) < (1,:))));
为了方便在循环,将一级边界框复制到一个三维数组:
L1boxes =重塑([L1。BoundingBox],[2 2 numel(L1)]);
检查每个级别1 - 2条特性可能的交点。看看polybool返回任何输出,但避免调用polybool除非先检测到边界框的十字路口:
为k = 1:元素个数(L2)为j = 1:元素个数(L1)%看看边界盒相交如果boxesIntersect (L2 (k)。BoundingBox, L1boxes(:,:,j))%看看实际特征相交如果~ isempty (polybool (“十字路口”,…L2 (k)。朗,L2 (k)。纬度,L1 (j)。朗,L1 (j)。纬度))%逆二级顶点顺序合并之前%正确东方内部环L1 (j)。朗= [L1(j).Lon fliplr(L2(k).Lon) NaN]; L1(j).Lat = [L1(j).Lat fliplr(L2(k).Lat) NaN];结束结束结束结束
第七步:在一个Shapefile保存结果
与一个调用shapewrite,您可以创建三个文件,
gshhs_c_SouthAmerica。轴马力gshhs_c_SouthAmerica。shx gshhs_c_SouthAmerica.dbf
在你的工作目录。
shapepath = fullfile (workingFolder,“gshhs_c_SouthAmerica.shp”);shapepath shapewrite (L1)
第八步:验证Shapefile
验证的结果shapewrite阅读新shapefile到地理空间表中南美洲:
南美洲= readgeotable (shapepath,“CoordinateSystemType”,“地理”)
南美洲=79×13表形状南部西北东部地区水平LevelString NumPoints FormatVersi源CrossesGree GSHHS_ID _______ _________交_____ _____ ___________ _____………………________ geopolyshape -53.9 71.994 191.89 325.21 3.7652 e + 07年1“土地“971 3”“全球价值调查”主要根据“0 1 geopolyshape -55.055 -52.469 287.98 294.89 46269 1“土地”19日3”“全球价值调查”主要根据“0 42 geopolyshape -52.36 -51.235 300.26 302.29 4402.4 1“土地”11 85 geopolyshape”“全球价值调查”主要根据“0 -55.675 -54.933 289.99 292.02 2669.1 1“土地”10 97 geopolyshape”“全球价值调查”主要根据“0 -50.065 -48.67 284.53 285.62 4021.9 1“土地”9”“全球价值调查”主要根据“0 115 geopolyshape -52.26 -51.336 298.92 300.81 4326.2 1“土地”9”“全球价值调查”主要根据“0 121 geopolyshape -54.135 -53.38 286.33 287.85 3809.7 1“土地”8 143 geopolyshape”“全球价值调查”主要根据“0 -43.448 -41.768 285.55 286.65 10188 1“土地”8 144 geopolyshape”“全球价值调查”主要根据“0 -53.478 -52.73 287 288.46 3006.8 1“土地”7 151 geopolyshape”“全球价值调查”主要根据“0 10.041 10.841 298.07 299.09 5174 1“土地”7 174 geopolyshape”“全球价值调查”主要根据“0 -54.276 -53.566 289.1 289.8 1420.3 1“土地”5 292 geopolyshape”“全球价值调查”主要根据“0 -44.937 -44.364 286.54 287.31 955.5 1“土地”5 329 geopolyshape”“全球价值调查”主要根据“0 -54.401 -53.936 288.29 289.04 704.5 1“土地”5 330 geopolyshape”“全球价值调查”主要根据“0 -51.125 -50.71 285.04 285.64 462 1“土地”5 339 geopolyshape”“全球价值调查”主要根据“0 -50.458 -50.001 284.54 285.24 612.9 1“土地”5 341 geopolyshape”“全球价值调查”主要根据“0 -48.711 -48.137 284.98 285.53 1039.8 1“土地”5 3”“全球价值调查”主要根据“0 343⋮
注意,这两个字段名最长,“FormatVersion”和“CrossesGreenwich”,已被截断为11个字符。这发生在调用shapewrite是不可避免的,因为严格的xBASE表(11-character限制。在shapefile DBF格式)用于存储属性。编写shapefile时(一般情况下,您可能希望重新定义字段名超过11个字符,以避免或控制自动截断的影响。)
可选地,把新shapefile从你的工作文件夹。(本例中需要清理;在真实的应用程序中你可能想省略此步骤)。
删除([shapepath”。*’])
从新的shapefile显示数据。注意各种“洞”在南美多边形指示湖泊和海岸线的其他扩展水体内部的大陆。
图ax = axesm (“墨卡托”,“MapLatLimit”latlim,“MapLonLimit”,lonlim);斧子。颜色=“青色”;gridm;mlabel;plabel geoshow(南美洲,“FaceColor”tightmap [0.15 0.8 0.15])
参考
韦塞尔,P。,和W. H. F. Smith, 1996, A global self-consistent, hierarchical, high-resolution shoreline database, Journal of Geophysical Research, Vol. 101, pp. 8741-8743.
额外的数据
下载完整的GSHHG数据集可能来自美国国家海洋和大气管理局(NOAA)的网站。发现数据集的更多信息,请参阅发现地理空间矢量数据。
学分
GSHHG数据映射工具箱中提供的文件是由夏威夷大学的保罗·韦塞尔博士沃尔特·h·f·史密斯博士和美国国家海洋和大气管理局。
有关更多信息,运行:
> >类型gshhs_c.txt
