投影和显示栅格数据

投影栅格数据

要将与地理栅格参考对象相关联的数据进行投影，首先需为栅格中的每个点创建一个经纬度坐标网格。然后，将地理坐标投影到 x-y 地图坐标系中。

例如，将高程栅格数据作为数组和地理像元引用对象导入。使用 geographicGrid 函数获取栅格中每个点的经纬度坐标。

[Z,R] = readgeoraster('n39_w106_3arc_v2.dt1');
[lat,lon] = geographicGrid(R);

既然已经有了网格，请选择一种地图投影法，用于对坐标进行投影。在此示例中，请为北半球的 UTM 第 13 区创建一个 projcrs 对象。然后，将经纬度坐标投影到 x-y 坐标系中。

p = projcrs(32613);
[x,y] = projfwd(p,lat,lon);

通过调用 mapshow 并指定 x-y 坐标和高程数组，将投影栅格显示为曲面。添加坐标轴标签，并应用适合高程数据的颜色图。

figure
mapshow(x,y,Z,'DisplayType','surface')
xlabel('x (meters)')
ylabel('y (meters)')
demcmap(Z)

如果经纬度坐标的地理坐标参考系统 (CRS) 与投影坐标参考系统 (CRS) 不一致，则投影坐标可能会不准确。您可以通过查询 projcrs 对象或地理栅格参考对象的 GeographicCRS 属性，来获取其地理坐标参考系统 (CRS)。

p.GeographicCRS.Name

ans = 
"WGS 84"

R.GeographicCRS.Name

ans = 
"WGS 84"

本示例中使用的 DTED 文件由美国地质调查局提供。

对栅格数据进行反投影

要对与地图栅格参考对象关联的数据进行反投影，首先需为栅格中的每个点创建一个由 x-y 坐标组成的网格。然后，将 x-y 地图坐标反投影为地理坐标。

例如，将一张波士顿的图片作为数组和一个单元格引用映射对象导入。通过查询引用对象的 ProjectedCRS 属性，获取作为 projcrs 对象的地图投影信息。

[Z,R] = readgeoraster('boston.tif');
p = R.ProjectedCRS;

使用 worldGrid 函数获取栅格中每个点的 x-y 坐标。

[x,y] = worldGrid(R);

使用 projinv 函数，并指定 projcrs 对象和坐标网格，将 x-y 坐标反投影为经纬度坐标。

[lat,lon] = projinv(p,x,y);

通过调用 geoshow 并指定经纬度坐标和图像数组，显示未反投影的图像。默认情况下，geoshow 使用正方形投影显示坐标。然后，添加坐标轴标签。

figure
geoshow(lat,lon,Z)
xlabel('Longitude (degrees)')
ylabel('Latitude (degrees)')