用matlab处理数字高程模型（DEM）之点云数据（point cloud data）

        数字高程模型（Digital Elevation Model，简称DEM)是通过有限的地形高程数据实现对地面地形的数字化模拟，它是用一组有序数值阵列形式表示地面高程的一种实体地面模型，是数字地形模型(Digital Terrain Model，简称DTM)的一个分支。DTM是描述包括高程在内的各种地貌因子，如坡度、坡向、坡度变化率等因子在内的线性和非线性组合的空间分布，其中DEM是零阶单纯的单项数字地貌模型，其他如坡度、坡向及坡度变化率等地貌特性可在DEM的基础上派生。

        航空摄影测量一直是地形图测绘和更新最有效也是最重要的手段，其获取的影像是高精度大范围DEM生产最优价值的数据源。由于测绘所采取的方式和手段的不同，DEM有不同的数据样式，其中SRTM最为有名。SRTM DEM有多个版本(V1,V2,V4)，多种格式(hgt/Geotiff/Bil/Arc Grid)，多种精度（SRTM1/ SRTM3/ SRTM30）其中V1为原始版本，V2为利用现有水体数据库在V1基础上进行修正的版本，V4版是在V2版缺失数据区域进行插值和修补。SRTM1是以地球等角坐标系的1角秒作为采样间隔(约30m)，SRTM3和SRTM30分别是以3角秒和30角秒为采样间隔(约90m和900m)。可以从https://dds.cr.usgs.gov/srtm/version2_1/SRTM3/下载V2版本的SRTM3数据，其数据是以一个经度一个纬度存储在格式为hgt的文件中，覆盖了全球的大部分地区。

       hgt格式的文件可以直接用Global Mapper打开，并生成相应区域的等高线，也可以通过Global Mapper进行格式转化成ArcGIS可以处理的文件。从该网站下载覆盖28°~32°E和15°~17°N区域的数据，拖入到Global Mapper中处理。截取我们需要研究的区域，并利用Global Mapper强大的格式转换能力将其生成点云数据（point cloud data），文件后缀为xyz，可以直接用load函数进行读取数据。最终生成一个三列的矩阵，第一列存储经度，第二列存储纬度，第三列存储海拔，长度均为9554000（2000*4777）（具体数据视分析的区域而定，下同）。将经度、纬度、海拔分别用列向量x,y,z进行存储，经过matlab的分析，我们发现在列向量x中，从第一个元素开始，每2000个元素的数值都相等；在列向量y中，数据的周期为2000；所以得出在点云文件中，是以纬度为主序，将海拔信息存储在一列里，相当于一块区域被分隔为9554000个小网格，每个小网格的海拔近似是相等的。所以我们将存储着海拔的向量z转化为2000*4777的矩阵进行存储，行对应着纬度，列对应着经度。在进行着一步之前，我们需要将所有不同的经度和纬度提取出来，分别存储在1*4777以及1*2000的行向量中。由此，便可用matlab的contour或contourf画出该地区的等高线图。

A=load('data.xyz');%load point cloud data and store it in matrix A
x=A(:,1);%x is column vector presenting longitude
y=A(:,2);%y is column vector presenting latitude
z=A(:,3);%z is column vector presenting altitude
x1=x(1:2000:9554000);
X=x1';%X is row vector presenting the range of the longitude in the area
y1=y(1:2000);
Y=y1';%Y is row vector presenting the range of the latitude in the area
Z=reshape(z,2000,4777);%shape z as a 2000*4777 matrix Z

%plot the contour map of the area
contour(X,Y,Z,50);
set(gca,'xtick',[28:0.23:32]);
set(gca,'ytick',[-17:0.1:-15]);
title('Contour map of the area');
xlabel('East(Degree)');
ylabel('South(Degree)');
colorbar;
set(get(colorbar,'Title'),'string','Altitude(m)');
grid on;

最终结果如图所示：