从DEM到三维地形：用PLY、OBJ、glTF构建GIS可视化模型

本文节选自新书《GIS基础原理与技术实践》第7章。很多人以为三维建模只能靠 3ds Max 或 Blender，但在 GIS 中，我们完全可以从 DEM 出发，用代码手动生成带颜色、带纹理、甚至符合现代 glTF 标准的三维地形模型。本文带你一步步实现 PLY 白模、OBJ 纹理贴图、glTF 资产封装，揭开三维 GIS 的底层逻辑。

导言

在前三章中，笔者详细论述了矢量、栅格和地形相关的知识，这三种数据也是GIS中最为基本的三种地理空间数据。不过，这三种传统的地理空间数据通常被认为是二维的，其处理方法也大多数基于二维平面空间。随着越来越复杂的地理空间信息的需求，GIS也在逐渐向三维方向发展，三维GIS成为了地理信息系统科学中炙手可热的研究方向。将三维模型作为基本的地理空间数据的观点目前可能还未进入学校的经典教材，但已经有这个趋势。在本章中，笔者也总结了一些与GIS相关的三维模型的知识，希望给读者以参考。

7.1 初识三维模型

7.1.1 三维模型的数据载体

随着计算机图形技术的发展，我们或多或少都会见过或者听说过三维模型。笔者始终记得小时候第一次在电视上看到三维动画《变形金刚：超能勇士》的震撼感受；而现在我们已经可以在手机上玩三维游戏《王者荣耀》，实时操作造型精美的英雄模型了。这些东西的背后都离不开三维模型数据，他们往往是通过像Autodesk 3D Max这样的三维建模软件制作出来的。

如同栅格数据和矢量数据一样，三维模型数据也有形形色色的数据格式。这些不同的数据格式有时来源于不同的三维建模软件；一些机构或者组织出于标准化的目的，也会定义某种通用的三维模型数据格式。这些不同格式的数据文件构成了三维模型的数据载体。常用的三维数据格式如下表7.1所示：

名称	全称	特点
PLY	Polygon File Format	描述最简单
OBJ	Wavefront .obj file	经典通用性高
3DS	3D Studio	广泛应用的经典格式
MAX	3D Studio Max	3DMax专属格式
glTF	Graphics Language Transmission Format	现代自由开放，适合OpenGL管线
FBX	Filmbox	现代而全面的格式，游戏引擎中常用

对这些三维数据格式我们可以做一个大致的认识，因为后面可能会直接用到：

1. PLY是一种最简单的三维数据格式，一般用其存储不带纹理的模型数据，可通过文本和二进制两种形式来描述。
2. OBJ是非常通用的三维模型数据格式，与PLY相比增加了对材质的描述，包括纹理信息。因此一个典型的OBJ格式的文件除了.obj文件，同时还会附带一个.mtl文件用于描述材质。而在材质文件中就可以指定纹理图片的地址。很长一段时间内由于其对三维模型文件的描述比较全面，通常用于不同三维建模软件的中转。
3. 3DS和MAX属于著名三维建模软件Autodesk 3ds Max的专属的文件数据格式。理论上来说，3DS和MAX都属于商业三维数据格式，但是Autodesk 3ds Max在三维建模上的使用非常广泛，所以这两种数据格式也很常见。不同的是，3DS现在已经基本能够各种三维软件所识别，一些开源工具也能解析识别；Max则是Autodesk 3ds Max所专用，其他三维软件或者开源工具一般都不支持。Max还有一个笔者认为不太好的特点，就是版本迭代太快，低版本的Autodesk 3ds Max无法打开高版本的Max格式数据。
4. glTF是一种自由开放，无专利限制的，适合传输和加载3D模型和场景的文件格式。相比较前面介绍三维数据格式来说，glTF诞生的时间较晚，可以采用了更为现代的图形技术来封装和组织这个格式。Khronos Group制定和维护了glTF数据格式的标准，同时由于其也是OpenGL接口标准的指定者和维护者，因此glTF特别适合OpenGL系列（OpenGL，OpenGL ES，WebGL）的图形渲染流水线所需要进行的处理。这个特点意味着glTF足够轻量化。目前glTF有1.0和2.0两个版本，其中glTF2.0已经成为了ISO国际标准。
5. FBX同样也是Autodesk公司开发的一种通用三维数据格式。与glTF一样，FBX也是更为现代的三维数据格式，比如支持显示效果更为真实的PBR材质。FBX广泛应用于游戏开发领域，目前最火的三维游戏引擎Unity和Unreal都支持直接导入这种格式。Autodesk官方为FBX提供了开发包，支持解析和修改该格式数据文件。除此之外，也有一些开源第三方组件使用自己的方式兼容它。

综合来说，PLY、OBJ和3DS都属于比较早期的三维模型数据格式，受限于当年的图形技术的认知；而Max、glTF和FBX则设计得更为现代，文件组织结构更为合理，能提供更为强大的可视化效果。例如，现代三维模型数据格式已经不仅仅是像早期三维模型数据格式那样只包含模型数据本身，还会包括材质、动画、灯光甚至相机等，其描述的对象可以是整个三维场景。

7.1.2 从地形来认识三维模型（PLY格式）

如果没有三维图形的基础知识，上一小节的论述可能会让有的读者一头雾水。那么我们可以从GIS中的地形开始说起——在第6.3节中我们就已经使用过PLY格式的三维数据，将其表达成不规则三角网地形。但是，如图6.6所示的地形实在过于简陋，有没有办法给这个白模赋予着色信息，使其有更好的可视化效果呢？

一种最简单的可视化优化方案是，可以结合第6.5节中晕渲图的实现，创建一个带颜色信息的地形三维模型数据。如下例7.1所示：

//例7.1 DEM数据转换PLY三维模型
#include <gdal_priv.h>#include <algorithm>
#include <array>
#include <fstream>
#include <iostream>
#include <vector>using namespace std;struct VertexProperty {double x;double y;double z;uint8_t red;uint8_t green;uint8_t blue;
};size_t vertexCount;
vector<VertexProperty> vertexData;
size_t faceCount;
vector<int> indices;//颜色查找表
using F_RGB = std::array<double, 3>;
vector<F_RGB> tableRGB(256);//生成渐变色
void Gradient(F_RGB& start, F_RGB& end, vector<F_RGB>& RGBList) {F_RGB d;for (int i = 0; i < 3; i++) {d[i] = (end[i] - start[i]) / RGBList.size();}for (size_t i = 0; i < RGBList.size(); i++) {for (int j = 0; j < 3; j++) {RGBList[i][j] = start[j] + d[j] * i;}}
}//初始化颜色查找表
void InitColorTable() {F_RGB blue({17, 60, 235});   //蓝色F_RGB green({17, 235, 86});  //绿色vector<F_RGB> RGBList(60);Gradient(blue, green, RGBList);for (int i = 0; i < 60; i++) {tableRGB[i] = RGBList[i];}F_RGB yellow({235, 173, 17});  //黄色RGBList.clear();RGBList.resize(60);Gradient(green, yellow, RGBList);for (int i = 0; i < 60; i++) {tableRGB[i + 60] = RGBList[i];}F_RGB red({235, 60, 17});  //红色RGBList.clear();RGBList.resize(60);Gradient(yellow, red, RGBList);for (int i = 0; i < 60; i++) {tableRGB[i + 120] = RGBList[i];}F_RGB white({235, 17, 235});  //紫色RGBList.clear();RGBList.resize(76);Gradient(red, white, RGBList);for (int i = 0; i < 76; i++) {tableRGB[i + 180] = RGBList[i];}
}//根据高程选颜色
inline int GetColorIndex(double z, double min_z, double max_z) {int temp = (int)floor((z - min_z) * 255 / (max_z - min_z) + 0.6);return temp;
}void ReadDem() {string workDir = getenv("GISBasic");string demPath = workDir + "/../Data/Model/dem.tif";GDALDataset* dem = (GDALDataset*)GDALOpen(demPath.c_str(), GA_ReadOnly);if (!dem) {cout << "Can't Open Image!" << endl;return;}int srcDemWidth = dem->GetRasterXSize();int srcDemHeight = dem->GetRasterYSize();//坐标信息double geoTransform[6] = {0};dem->GetGeoTransform(geoTransform);double srcDx = geoTransform[1];double srcDy = geoTransform[5];double startX = geoTransform[0] + 0.5 * srcDx;double startY = geoTransform[3] + 0.5 * srcDy;double endX = startX + (srcDemWidth - 1) * srcDx;double endY = startY + (srcDemHeight - 1) * srcDy;size_t demBufNum = (size_t)srcDemWidth * srcDemHeight;vector<float> srcDemBuf(demBufNum, 0);int depth = sizeof(float);dem->GetRasterBand(1)->RasterIO(GF_Read, 0, 0, srcDemWidth, srcDemHeight,srcDemBuf.data(), srcDemWidth, srcDemHeight,GDT_Float32, depth, srcDemWidth * depth);GDALClose(dem);double minZ = *(std::min_element(srcDemBuf.begin(), srcDemBuf.end()));double maxZ = *(std::max_element(srcDemBuf.begin(), srcDemBuf.end()));vertexCount = (size_t)srcDemWidth * srcDemHeight;vertexData.resize(vertexCount);for (int yi = 0; yi < srcDemHeight; yi++) {for (int xi = 0; xi < srcDemWidth; xi++) {size_t m = (size_t)srcDemWidth * yi + xi;vertexData[m].x = startX + xi * srcDx;vertexData[m].y = startY + yi * srcDy;vertexData[m].z = srcDemBuf[m];int index = GetColorIndex(srcDemBuf[m], minZ, maxZ);vertexData[m].red = (uint8_t)(tableRGB[index][0] + 0.5);vertexData[m].green = (uint8_t)(tableRGB[index][1] + 0.5);vertexData[m].blue = (uint8_t)(tableRGB[index][2] + 0.5);}}faceCount = (size_t)(srcDemHeight - 1) * (srcDemWidth - 1) * 2;// indices.resize(faceCount);for (int yi = 0; yi < srcDemHeight - 1; yi++) {for (int xi = 0; xi < srcDemWidth - 1; xi++) {size_t m = (size_t)srcDemWidth * yi + xi;indices.push_back(m);indices.push_back(m + srcDemWidth);indices.push_back(m + srcDemWidth + 1);indices.push_back(m + srcDemWidth + 1);indices.push_back(m + 1);indices.push_back(m);}}
}void WriteDemModel() {string workDir = getenv("GISBasic");string demPath = workDir + "/../Data/Model/dst.ply";ofstream outfile(demPath);if (!outfile) {printf("write file error %s\n", demPath.c_str());return;}outfile << "ply\n";outfile << "format ascii 1.0\n";outfile << "comment CL generated\n";outfile << "element vertex " << to_string(vertexCount) << '\n';outfile << "property double x\n";outfile << "property double y\n";outfile << "property double z\n";outfile << "property uchar red\n";outfile << "property uchar green\n";outfile << "property uchar blue\n";outfile << "element face " << to_string(faceCount) << '\n';outfile << "property list uchar int vertex_indices\n";outfile << "end_header\n";outfile << fixed;for (int vi = 0; vi < vertexCount; vi++) {outfile << vertexData[vi].x << ' ';outfile << vertexData[vi].y << ' ';outfile << vertexData[vi].z << '\n';outfile << (int)vertexData[vi].red << ' ';outfile << (int)vertexData[vi].green << ' ';outfile << (int)vertexData[vi].blue << '\n';}for (size_t fi = 0; fi < faceCount; fi++) {outfile << 3;for (int ii = 0; ii < 3; ii++) {int id = indices[fi * 3 + ii];outfile << ' ' << id;}outfile << '\n';}
}int main() {GDALAllRegister();  //注册格式InitColorTable();ReadDem();WriteDemModel();return 0;
}

我们将生成PLY格式的三维数据导入到开源三维软件MeshLab中，其显示的效果如下图7.1所示。可以看到虽然我们更换了一个地形数据，但是其展示的效果与第6.5节中晕渲图的效果比较类似。其实准确来说，二维晕渲图的可视化效果正是来自于三维渲染一定光照条件下的实现。

例6.4使用PLY数据格式来表达不规则三角网地形，而本例表达的则是规则格网地形。但是只要是保存为三维数据格式，其存储的数据信息都是相同的，都包含顶点信息和索引信息。其中顶点信息不再只包含位置信息了，还包含了每个顶点的颜色信息RGB，因此本例封装了一个顶点属性的结构体来表达一个顶点：

struct VertexProperty {double x;double y;double z;uint8_t red;uint8_t green;uint8_t blue;
};

我们生成PLY格式是文本格式，通过记事本打开也可以直接看到位置信息和颜色信息，如下图 7.2所示：

通过这个简单的例子，就可以知道为什么需要三维模型数据分成顶点信息和索引信息来保存。在本例中，是将DEM每个方形格网转换成两个三角形，如果格网DEM为m行n列，这意味着存在(m-1)⋅(n-1)个格网，即2⋅(m-1)⋅(n-1)个三角形。如果我们以一个三角形顶点接着一个三角形顶点来描述三维模型文件，那么就需要6⋅(m-1)⋅(n-1)个顶点。但其实DEM中的顶点个数很明确，就是 m ⋅n个——这意味着至少存在这4到5倍的数据冗余。

先描述顶点信息，再描述索引信息，这样可以兼容一些共顶点的情况，因而可以最大化减少数据量，毕竟一个索引比一个顶点的数据量更少。如下图7.3所示，是本例生成的PLY格式文件中的索引信息。其中每一行代表一个面，3表示绘制的是三角形，后面三个数则表示顶点数据中每个顶点的索引编号。当然，PLY格式也可以将每个面描述成四边形或者多边形。但是目前大多数图形API或渲染引擎都将三角形作为绘制的基本图元，以三角形作为最小的绘制单位是效率最高的。

7.1.3 地形和影像组成三维模型（OBJ格式）

在上一节中展示了基于地形的三维模型的可视化（图7.1），但这种效果其实是一种风格化的效果。所谓风格化效果，就是不一定写实，但是由于抓住了事物对象主要特征，我们可以很容易确信展示就是渲染的就是该事物对象，例如卡通风格就是一种典型的风格化效果。

与风格化效果相对应的就是写实效果，写实效果能够让用户有更为真实的感受。在这里，如果要让这个地形三维模型数据得到这种写实的效果，那么就可以利用我们在第5章介绍过的栅格影像（DOM），将其铺在地形三维模型的表面而不是使用顶点着色。问题在于，如何将这个影像铺在三维模型上呢？

这个时候我们就要用到除了位置和颜色之外的，另一种顶点属性信息：纹理坐标。在计算机图形中，影像/图片数据在被传输到GPU后，就被封装成一种名为“纹理”的数据对象。顶点的纹理坐标就是该顶点对应于这张纹理图片的位置，一个三角形面片有三个顶点，也就对应了纹理图片上的三个位置，从而可以让我们取得纹理上的颜色值。而三角形内部的顶点的颜色值，就直接从纹理上的三角形面片的区域取值内插得到。

在这里说的内插过程，有点像我们之前第5章中介绍的影像进行图像内插过程，但并不完全准确。这涉及到计算机图像渲染流水线中光栅化的过程，是一个很复杂的过程。本章我们只用关心三维模型数据本身，可视化的问题我们后面再介绍。在这里我们只需要知道，顶点信息可以附带纹理坐标信息，从而将纹理图片的颜色值映射到模型上。

如果我们要在一个三维模型中附带纹理和纹理坐标信息，那么使用PLY格式的三维模型数据就不是很方便了。如下例7.2所示，我们使用OBJ格式的三维模型数据，来表达带纹理和纹理坐标信息的地形：

//例7.2 DEM数据转换OBJ三维模型
#include <gdal_priv.h>#include <algorithm>
#include <array>
#include <fstream>
#include <iostream>
#include <vector>using namespace std;struct VertexProperty {double x;double y;double z;double texCoordX;double texCoordY;
};size_t vertexCount;
vector<VertexProperty> vertexData;
size_t faceCount;
vector<int> indices;void ReadDem() {string workDir = getenv("GISBasic");string demPath = workDir + "/../Data/Model/dem.tif";GDALDataset* dem = (GDALDataset*)GDALOpen(demPath.c_str(), GA_ReadOnly);if (!dem) {cout << "Can't Open Image!" << endl;return;}int srcDemWidth = dem->GetRasterXSize();int srcDemHeight = dem->GetRasterYSize();//坐标信息double geoTransform[6] = {0};dem->GetGeoTransform(geoTransform);double srcDx = geoTransform[1];double srcDy = geoTransform[5];double startX = geoTransform[0] + 0.5 * srcDx;double startY = geoTransform[3] + 0.5 * srcDy;double endX = startX + (srcDemWidth - 1) * srcDx;double endY = startY + (srcDemHeight - 1) * srcDy;size_t demBufNum = (size_t)srcDemWidth * srcDemHeight;vector<float> srcDemBuf(demBufNum, 0);int depth = sizeof(float);dem->GetRasterBand(1)->RasterIO(GF_Read, 0, 0, srcDemWidth, srcDemHeight,srcDemBuf.data(), srcDemWidth, srcDemHeight,GDT_Float32, depth, srcDemWidth * depth);GDALClose(dem);double minZ = *(std::min_element(srcDemBuf.begin(), srcDemBuf.end()));double maxZ = *(std::max_element(srcDemBuf.begin(), srcDemBuf.end()));vertexCount = (size_t)srcDemWidth * srcDemHeight;vertexData.resize(vertexCount);for (int yi = 0; yi < srcDemHeight; yi++) {for (int xi = 0; xi < srcDemWidth; xi++) {size_t m = (size_t)srcDemWidth * yi + xi;vertexData[m].x = startX + xi * srcDx;vertexData[m].y = startY + yi * srcDy;vertexData[m].z = srcDemBuf[m];vertexData[m].texCoordX = (double)xi / (srcDemWidth - 1);vertexData[m].texCoordY = (double)yi / (srcDemHeight - 1);}}faceCount = (size_t)(srcDemHeight - 1) * (srcDemWidth - 1) * 2;// indices.resize(faceCount);for (int yi = 0; yi < srcDemHeight - 1; yi++) {for (int xi = 0; xi < srcDemWidth - 1; xi++) {size_t m = (size_t)srcDemWidth * yi + xi;indices.push_back(m);indices.push_back(m + srcDemWidth);indices.push_back(m + srcDemWidth + 1);indices.push_back(m + srcDemWidth + 1);indices.push_back(m + 1);indices.push_back(m);}}
}void WriteDemModel() {string workDir = getenv("GISBasic");string demPath = workDir + "/../Data/Model/dst.obj";ofstream outfile(demPath);if (!outfile) {printf("write file error %s\n", demPath.c_str());return;}outfile << "mtllib dst.mtl\n";outfile << fixed;for (int vi = 0; vi < vertexCount; vi++) {outfile << "v" << ' ';outfile << vertexData[vi].x << ' ';outfile << vertexData[vi].y << ' ';outfile << vertexData[vi].z << '\n';}for (int vi = 0; vi < vertexCount; vi++) {outfile << "vt" << ' ';outfile << vertexData[vi].texCoordX << ' ';outfile << vertexData[vi].texCoordY << '\n';}outfile << "usemtl dst\n";for (size_t fi = 0; fi < faceCount; fi++) {outfile << "f";for (int ii = 0; ii < 3; ii++) {int id = indices[fi * 3 + ii] + 1;outfile << ' ' << id << '/' << id;}outfile << '\n';}string mtlPath = workDir + "/../Data/Model/dst.mtl";ofstream mtlfile(mtlPath);if (!mtlfile) {printf("write file error %s\n", mtlPath.c_str());return;}mtlfile << "newmtl dst\n";mtlfile << "illum 2\n";mtlfile << "map_Ka tex.jpg\n";mtlfile << "map_Kd tex.jpg\n";mtlfile << "map_Ks tex.jpg\n";mtlfile << "Ns 10.000\n";
}int main() {GDALAllRegister();  //注册格式//设置Proj数据std::string projDataPath = getenv("GISBasic");projDataPath += "/share/proj";CPLSetConfigOption("PROJ_LIB", projDataPath.c_str());ReadDem();WriteDemModel();return 0;
}

程序运行完成之后，生成.obj格式后缀的三维模型文件，其数据内容与.ply格式后缀的三维模型文件差不多，都是由顶点数据和索引数据组成的，只不过两者的数据组织形式不同。OBJ格式的数据描述形式是先描述顶点位置信息，如下图7.4所示：

接着描述顶点的纹理坐标信息，如下图7.5所示：

最后是索引数据信息，如下图7.6所示。OBJ格式的索引数据的设计稍微复杂了一点，将索引划分成顶点位置数据的索引，以及顶点纹理坐标的索引。如果有顶点法向量数据，还可以加上顶点法向量数据的索引。但这里不用进行区分，直接都使用同一个索引值：

程序在生成.obj文件的同时，生成了一个后缀名为.mtl的文件。mtl是英文单词material（材质）的缩写，这个文件定义了OBJ格式文件的材质信息。在图7.4中我们可以看到，.obj文件在第一行描述信息就引用了这个.mtl：

mtllib dst.mtl

这一行描述信息表示dst.mtl文件是该.obj文件的材质文件。在这个材质文件中，描述了一些材质参数，其中就包括我们前面讲到的纹理图片文件，如下图7.7所示，tex.jpg文件就是我们使用的纹理：

材质是三维可视化最为重要的概念之一，决定了物体对象以什么样的可视化效果渲染展示出来。例如，同一个物体的质感体现是金属、木头还是塑料，需要通过材质来进行定义。不过这个问题很复杂，目前我们只需要知道通常在材质中使用纹理，是材质的关键参数。

三维物体不会仅仅只包含一个材质，材质文件中可能会包含多个材质。每个材质通常与一段顶点索引数据相关联，表示这一段图元是通过该材质进行渲染的。如图7.6所示在描述顶点索引信息之前，使用了如下描述语句：

usemtl dst

这表示以下三角面图元是通过材质文件dst.mtl中的dst材质来进行渲染的。

最后，将生成OBJ格式的三维数据导入到开源三维软件MeshLab中，显示的三维渲染效果如下图7.8所示。可以看到相比例7.1的结果来说，具有更好的写实效果，可以看到突起的山峰以及峡谷的河流。这是因为使用了栅格影像（DOM）来作为纹理图片，而DOM影像通常由光学影像拍摄真实的地形拍摄而来，真实感效果更好。

7.1.4 认识现代三维模型数据（glTF格式）

三维模型数据是进行三维可视化的初始载体。在进行图形渲染的起始阶段，会将硬盘中的三维模型数据读取到内存中，然后再传入显存做进一步处理。因此，三维数据格式总是会随着计算机图形技术的发展而发展，要么会出现更新的三维模型数据的格式，要么会在已有的三维模型数据上作扩展。

相比较前面介绍的PLY格式或者OBJ格式，glTF是一种更为现代的三维模型数据格式。这种现代性不仅仅是体现在时间上，更是体现在多方面的：

• 现代三维数据格式包含的数据内容越来越多，大多数与三维场景渲染相关的信息都可以进行定义和保存。
• 现代三维数据格式数据定义的概念与三维渲染流程越来越适配，很多早期的三维数据格式并没有考虑到三维渲染。
• 现代三维数据格式往往会借用一些其他已经定义好的数据格式，是一个数据文件的复合体。
• 现代三维模型数据已经不单纯是一个数据资源，更是一个数据资产（asset）。

glTF的英文全称是GL Transmission Format，从这个名称就可以看出这个三维数据格式的设计目的就是为了最大化数据传输的效率。这个数据传输不仅仅是指网络端到本地端数据传输，也包括本地数据与内存传输，以及内存到显存数据传输。对于数据处理这一类程序来说，CPU或者GPU的运算速度已经足够快了，数据传输反而是程序的性能瓶颈。因此，glTF的设计思路尽可能轻量化，最大程度减少3D资产的大小，节约解析和使用这些资产所需的运行时处理的时间。

glTF有glTF1.0和glTF2.0两个版本，本书描述的内容以glTF2.0规范为准。通常来说，glTF格式数据包含以下几个部分：

1. 三维场景数据描述：使用JSON来描述。JSON（JavaScript Object Notation）是一种轻量级的数据交换格式，易于人类阅读和编写，也易于机器解析和生成，提升网络传输效率。另外，JSON尤其适配三维场景这种树形结构数据的表达。
2. 缓冲区数据：保存为二进制文件。缓冲区数据指的就是前面提到的顶点数据和顶点索引数据，其数据量通常比较大，以二进制的形式进行一次或者少数几次传输，可以最大化减少数据预处理以及数据传输的性能损耗。在OpenGL的图形渲染中，缓冲区数据读取后可以直接被其API接口调用。
3. 资源文件：例如纹理数据可以使用jpg格式图片来表达，jpg格式图片压缩比较高，利于进行数据传输。另外还有一些三维图形专用的纹理格式，例如DXT和KTX2，也可以作为单独的文件被glTF使用。

作为对照，这里还是使用将DEM数据转换成三维模型的例子，如下例7.3所示：

//例7.3 DEM数据转换glTF三维模型
#include <gdal_priv.h>#include <fstream>
#include <iomanip>
#include <iostream>
#include <nlohmann/json.hpp>using namespace std;
using namespace nlohmann;size_t pointNum = 0;
size_t binBufNum = 0;
size_t indicesNum = 0;void CreateBinFile() {string workDir = getenv("GISBasic");string demPath = workDir + "/../Data/Model/dem.tif";GDALDataset *img = (GDALDataset *)GDALOpen(demPath.c_str(), GA_ReadOnly);if (!img) {printf("Can't Open Image!");return;}int bufWidth = img->GetRasterXSize();   //图像宽度int bufHeight = img->GetRasterYSize();  //图像高度int bandNum = img->GetRasterCount();    //波段数if (bandNum != 1) {printf("DEM波段数不为1");return;}int depth = GDALGetDataTypeSize(img->GetRasterBand(1)->GetRasterDataType()) /8;  //图像深度//获取地理坐标信息double padfTransform[6];if (img->GetGeoTransform(padfTransform) == CE_Failure) {printf("获取仿射变换参数失败");return;}double startX = padfTransform[0];double dX = padfTransform[1];double startY = padfTransform[3];double dY = padfTransform[5];//申请bufsize_t imgBufNum = (size_t)bufWidth * bufHeight * bandNum;float *imgBuf = new float[imgBufNum];//读取img->RasterIO(GF_Read, 0, 0, bufWidth, bufHeight, imgBuf, bufWidth, bufHeight,GDT_Float32, bandNum, nullptr, bandNum * depth,
bufWidth * bandNum * depth, depth);pointNum = (size_t)bufWidth * bufHeight;size_t position_texture_num = pointNum * 5;float *position_texture = new float[position_texture_num];for (int yi = 0; yi < bufHeight; yi++) {for (int xi = 0; xi < bufWidth; xi++) {size_t n = (size_t)(bufWidth * 5) * yi + 5 * xi;position_texture[n] = dX * xi;position_texture[n + 1] = dY * yi;size_t m = (size_t)(bufWidth * bandNum) * yi + bandNum * xi;position_texture[n + 2] = imgBuf[m];position_texture[n + 3] = float(xi) / (bufWidth - 1);position_texture[n + 4] = float(yi) / (bufHeight - 1);}}//释放delete[] imgBuf;imgBuf = nullptr;string binPath = workDir + "/../Data/Model/new.bin";ofstream binFile(binPath, std::ios::binary);binFile.write((char *)position_texture, position_texture_num * sizeof(float));size_t vertexBufNum = position_texture_num * sizeof(float);binBufNum = binBufNum + vertexBufNum;int mod = vertexBufNum % sizeof(uint16_t);if (mod != 0) {int spaceNum = sizeof(float) - mod;char *space = new char[spaceNum];binBufNum = binBufNum + sizeof(char) * spaceNum;memset(space, 0, sizeof(char) * spaceNum);binFile.write(space, sizeof(char) * spaceNum);delete[] space;space = nullptr;}indicesNum = (size_t)(bufWidth - 1) * (bufHeight - 1) * 2 * 3;uint16_t *indices = new uint16_t[indicesNum];for (int yi = 0; yi < bufHeight - 1; yi++) {for (int xi = 0; xi < bufWidth - 1; xi++) {uint16_t m00 = (uint16_t)(bufWidth * yi + xi);uint16_t m01 = (uint16_t)(bufWidth * (yi + 1) + xi);uint16_t m11 = (uint16_t)(bufWidth * (yi + 1) + xi + 1);uint16_t m10 = (uint16_t)(bufWidth * yi + xi + 1);size_t n = (size_t)(bufWidth - 1) * yi + xi;indices[n * 6] = m00;indices[n * 6 + 1] = m01;indices[n * 6 + 2] = m11;indices[n * 6 + 3] = m11;indices[n * 6 + 4] = m10;indices[n * 6 + 5] = m00;}}binFile.write((char *)indices, sizeof(uint16_t) * indicesNum);binBufNum = binBufNum + sizeof(uint16_t) * indicesNum;delete[] position_texture;position_texture = nullptr;delete[] indices;indices = nullptr;
}int main() {GDALAllRegister();CPLSetConfigOption("GDAL_FILENAME_IS_UTF8", "NO");  //支持中文路径//设置Proj数据std::string projDataPath = getenv("GISBasic");projDataPath += "/share/proj";CPLSetConfigOption("PROJ_LIB", projDataPath.c_str());ordered_json gltf;gltf["asset"] = {{"generator", "CL"}, {"version", "2.0"}};gltf["scene"] = 0;gltf["scenes"] = {{{"nodes", {0}}}};gltf["nodes"] = {{{"mesh", 0}}};ordered_json positionJson;positionJson["POSITION"] = 1;positionJson["TEXCOORD_0"] = 2;ordered_json primitivesJson;primitivesJson = {{{"attributes", positionJson}, {"indices", 0}, {"material", 0}}};gltf["meshes"] = {{{"primitives", primitivesJson}}};ordered_json pbrJson;pbrJson["baseColorTexture"]["index"] = 0;gltf["materials"] = {{{"pbrMetallicRoughness", pbrJson}}};CreateBinFile();gltf["textures"] = {{{"sampler", 0}, {"source", 0}}};gltf["images"] = {{{"uri", "tex.jpg"}}};gltf["samplers"] = {{{"magFilter", 9729},{"minFilter", 9987},{"wrapS", 33648},{"wrapT", 33648}}};gltf["buffers"] = {{{"uri", "new.bin"}, {"byteLength", binBufNum}}};ordered_json indicesBufferJson;indicesBufferJson["buffer"] = 0;indicesBufferJson["byteOffset"] = pointNum * 5 * 4;indicesBufferJson["byteLength"] = indicesNum * 2;indicesBufferJson["target"] = 34963;ordered_json positionBufferJson;positionBufferJson["buffer"] = 0;positionBufferJson["byteStride"] = sizeof(float) * 5;positionBufferJson["byteOffset"] = 0;positionBufferJson["byteLength"] = pointNum * 5 * 4;positionBufferJson["target"] = 34962;gltf["bufferViews"] = {indicesBufferJson, positionBufferJson};ordered_json indicesAccessors;indicesAccessors["bufferView"] = 0;indicesAccessors["byteOffset"] = 0;indicesAccessors["componentType"] = 5123;indicesAccessors["count"] = indicesNum;indicesAccessors["type"] = "SCALAR";indicesAccessors["max"] = {18719};indicesAccessors["min"] = {0};ordered_json positionAccessors;positionAccessors["bufferView"] = 1;positionAccessors["byteOffset"] = 0;positionAccessors["componentType"] = 5126;positionAccessors["count"] = pointNum;positionAccessors["type"] = "VEC3";positionAccessors["max"] = {770, 0.0, 1261.151611328125};positionAccessors["min"] = {0.0, -2390, 733.5555419921875};ordered_json textureAccessors;textureAccessors["bufferView"] = 1;textureAccessors["byteOffset"] = sizeof(float) * 3;textureAccessors["componentType"] = 5126;textureAccessors["count"] = pointNum;textureAccessors["type"] = "VEC2";textureAccessors["max"] = {1, 1};textureAccessors["min"] = {0, 0};gltf["accessors"] = {indicesAccessors, positionAccessors, textureAccessors};string workDir = getenv("GISBasic");string jsonFile = workDir + "/../Data/Model/new.gltf";ofstream binFile(jsonFile, std::ios::binary);std::ofstream outFile(jsonFile);outFile << std::setw(4) << gltf << std::endl;
}

例7.3最终得到的glTF格式使的三维模型文件如下图7.9所示（也可以直接从本书的在线主页中获取）。.gltf后缀的文件就是用于三维场景数据描述的JSON文件，.bin后缀的文件就是储存缓存区数据的二进制文件，.jpg文件就是三维模型用到的纹理图片。一些在线网站（如https://gltf-viewer.donmccurdy.com）提供了对glTF模型数据的浏览，可以得到如图7.8所示的渲染效果。

在这里，我们就不继续介绍例7.3的转换过程了，因为其关键就在于glTF格式数据解析。在下一节中，我们会根据这个成果数据详细介绍glTF数据格式规范，加深对三维模型数据的认识。

---

本文节选自作者新书《GIS基础原理与技术实践》第7章。书中系统讲解 GIS 核心理论与多语言实战，适合开发者与高校师生。

📚 配套资源开源：GitHub | GitCode
🛒 支持正版：京东｜当当