制作三维地形是一个结合数据、算法和可视化的综合任务。根据地形数据来源的不同，实现方法主要分为两类：**基于真实地形数据**和**基于过程化生成算法**。两种方法各有侧重，前者适用于地理信息系统、规划仿真等需要真实地理信息的场景[ref_5]；后者则广泛应用于游戏开发、虚拟现实等需要高效生成可控、风格化地形的领域[ref_3][ref_6]。 以下是核心方法对比及具体实现： | 方法类别 | 核心思路 | 数据/算法 | 优点 | 适用场景 | | :--- | :--- | :--- | :--- | :--- | | **基于真实数据** | 从外部加载数字高程模型等地理数据，并转换为可渲染的高程矩阵[ref_5]。 | 使用 `gdal`, `rasterio` 等库读取 `GeoTIFF`, `DEM` 文件。 | 结果真实、准确，忠实反映现实地形。 | 地理分析、工程规划、科学可视化。 | | **基于过程化生成** | 利用数学算法（如噪声、分形）生成随机但可控的地形高度场[ref_2][ref_3]。 | Perlin噪声、Simplex噪声、Diamond-Square算法等。 | 灵活可控、无需外部数据、可无限生成。 | 游戏地形、虚拟世界、艺术创作。 | ### 一、基于真实地形数据的实现 此方法的关键步骤是**读取数据**和**三维可视化**。以读取 `GeoTIFF` 格式的 DEM 数据并使用 Matplotlib 进行可视化为例[ref_5]： ```python import numpy as np import matplotlib.pyplot as plt from mpl_toolkits.mplot3d import Axes3D from osgeo import gdal # 或使用 rasterio def load_dem_with_gdal(filepath): """ 使用GDAL库读取DEM文件，返回高程数据和地理信息。 """ dataset = gdal.Open(filepath, gdal.GA_ReadOnly) if dataset is None: raise FileNotFoundError(f"无法打开文件: {filepath}") band = dataset.GetRasterBand(1) elevation = band.ReadAsArray() transform = dataset.GetGeoTransform() dataset = None # 关闭数据集 return elevation, transform def create_3d_surface(elevation, stride=5): """ 根据高程矩阵创建三维曲面。 stride参数用于对数据进行下采样，避免网格点过多导致渲染缓慢。 """ # 获取地形网格的X, Y坐标（此处简化为索引） x = np.arange(0, elevation.shape[1]) y = np.arange(0, elevation.shape[0]) X, Y = np.meshgrid(x, y) # 下采样以加速渲染 X_small = X[::stride, ::stride] Y_small = Y[::stride, ::stride] Z_small = elevation[::stride, ::stride] # 创建三维图形 fig = plt.figure(figsize=(12, 8)) ax = fig.add_subplot(111, projection='3d') # 绘制曲面图 surf = ax.plot_surface(X_small, Y_small, Z_small, cmap='terrain', # 使用地形色图 linewidth=0, antialiased=True, alpha=0.8) # 添加颜色条 fig.colorbar(surf, ax=ax, shrink=0.5, aspect=5, label='Elevation (m)') ax.set_xlabel('X (Column)') ax.set_ylabel('Y (Row)') ax.set_zlabel('Elevation (m)') ax.set_title('3D Terrain from DEM Data') # 调整视角以获得更好的立体感 ax.view_init(elev=30, azim=120) plt.tight_layout() plt.show() # 使用示例 if __name__ == "__main__": # 注意：需要替换为实际的DEM文件路径 dem_file = "path/to/your/dem.tif" try: elevation_data, _ = load_dem_with_gdal(dem_file) create_3d_surface(elevation_data, stride=10) except Exception as e: print(f"处理过程中发生错误: {e}") ``` 此代码的核心在于 `gdal` 读取栅格数据并转换为 `numpy` 数组，然后利用 Matplotlib 的 `plot_surface` 函数进行三维渲染。为提升大文件处理性能，代码中引入了 `stride` 参数进行下采样[ref_5]。 ### 二、基于过程化生成算法的实现 过程化生成的核心是创建一个具有自然地貌特征的高度图。这里介绍两种主流算法。 **1. 使用 Perlin 噪声生成地形** Perlin 噪声能生成连续、平滑的随机值，非常适合模拟自然地形[ref_6]。 ```python import numpy as np import matplotlib.pyplot as plt from mpl_toolkits.mplot3d import Axes3D def generate_perlin_noise_2d(shape, scale=100.0, octaves=6, persistence=0.5, lacunarity=2.0): """ 生成二维Perlin噪声高度场。 参数: shape: 输出噪声图的大小 (height, width)。 scale: 噪声采样频率，越小噪声越平滑。 octaves: 噪声层数，增加细节。 persistence: 振幅衰减因子 (0~1)，决定每层噪声对最终结果的贡献度。 lacunarity: 频率倍增因子 (>1)，决定每层噪声的频率增加倍数。 """ height, width = shape noise = np.zeros(shape) # 生成多个倍频程（Octave）的噪声并叠加 for i in range(octaves): freq = lacunarity ** i amplitude = persistence ** i # 为每个倍频程创建网格点 x = np.linspace(0, scale * freq, width, endpoint=False) y = np.linspace(0, scale * freq, height, endpoint=False) X, Y = np.meshgrid(x, y) # 简化版：使用插值模拟噪声（实际应用中推荐使用`noise`库） # 此处生成一个基础噪声网格 base_noise = np.random.randn(height, width) # 使用简单的双线性插值平滑（此为示意，真实的Perlin噪声更复杂） from scipy.ndimage import zoom if i == 0: octave_noise = base_noise else: # 通过缩放和下采样模拟不同频率 zoom_factor = 2 ** i resampled = zoom(base_noise, zoom_factor, order=1) # 裁剪到原始尺寸 crop_h = min(resampled.shape[0], height) crop_w = min(resampled.shape[1], width) octave_noise = resampled[:crop_h, :crop_w] # 如果裁剪后尺寸不匹配，填充零 if octave_noise.shape != shape: pad_h = shape[0] - octave_noise.shape[0] pad_w = shape[1] - octave_noise.shape[1] octave_noise = np.pad(octave_noise, ((0, pad_h), (0, pad_w)), mode='constant') noise += octave_noise * amplitude # 将噪声值归一化到[0,1]区间 noise_min, noise_max = noise.min(), noise.max() if noise_max > noise_min: noise = (noise - noise_min) / (noise_max - noise_min) return noise # 生成并绘制Perlin噪声地形 terrain_size = (100, 150) height_map = generate_perlin_noise_2d(terrain_size, scale=50.0, octaves=5) height_map = height_map * 100 # 缩放高程值 # 创建网格 x = np.arange(0, terrain_size[1]) y = np.arange(0, terrain_size[0]) X, Y = np.meshgrid(x, y) fig = plt.figure(figsize=(10, 7)) ax = fig.add_subplot(111, projection='3d') surf = ax.plot_surface(X, Y, height_map, cmap='terrain', linewidth=0, antialiased=True) fig.colorbar(surf, ax=ax, shrink=0.5, aspect=5) ax.set_title('Procedural Terrain using Perlin Noise') ax.set_xlabel('X') ax.set_ylabel('Y') ax.set_zlabel('Height') ax.view_init(elev=45, azim=60) plt.show() ``` **2. 使用 Diamond-Square 分形算法生成地形** Diamond-Square 算法是一种经典的分形地形生成算法，能快速产生具有自相似特征的崎岖地形[ref_2]。 ```python def diamond_square(size, roughness=0.7): """ 实现Diamond-Square算法生成分形地形高度图。 参数: size: 必须是2^n + 1 (如 129, 257)。 roughness: 粗糙度系数，值越大地形越崎岖。 """ n = int(np.log2(size - 1)) if 2 ** n + 1 != size: raise ValueError("Size must be 2^n + 1 (e.g., 129, 257, 513)") height_map = np.zeros((size, size), dtype=float) # 初始化四个角点 height_map[0, 0] = np.random.rand() * 2 - 1 height_map[0, -1] = np.random.rand() * 2 - 1 height_map[-1, 0] = np.random.rand() * 2 - 1 height_map[-1, -1] = np.random.rand() * 2 - 1 step = size - 1 scale = roughness while step > 1: half = step // 2 # Diamond步骤：为每个正方形的中心点赋值 for y in range(half, size, step): for x in range(half, size, step): avg = (height_map[y - half, x - half] + height_map[y - half, x + half] + height_map[y + half, x - half] + height_map[y + half, x + half]) / 4.0 height_map[y, x] = avg + np.random.rand() * scale * 2 - scale # Square步骤：为每个菱形的中心点赋值 for y in range(0, size, half): # 处理偏移的行 row_offset = half if y % step == 0 else 0 for x in range(row_offset, size, step): total, count = 0.0, 0 # 收集上、下、左、右四个邻居的值（如果存在） if y - half >= 0: total += height_map[y - half, x] count += 1 if y + half < size: total += height_map[y + half, x] count += 1 if x - half >= 0: total += height_map[y, x - half] count += 1 if x + half < size: total += height_map[y, x + half] count += 1 if count > 0: height_map[y, x] = total / count + np.random.rand() * scale * 2 - scale step //= 2 scale *= roughness # 随着迭代，扰动幅度减小 return height_map # 生成并绘制Diamond-Square地形 ds_size = 129 # 2^7 + 1 = 129 ds_height_map = diamond_square(ds_size, roughness=0.6) ds_height_map = (ds_height_map - ds_height_map.min()) / (ds_height_map.max() - ds_height_map.min()) * 50 x = np.arange(0, ds_size) y = np.arange(0, ds_size) X, Y = np.meshgrid(x, y) fig, (ax1, ax2) = plt.subplots(1, 2, figsize=(14, 5), subplot_kw={'projection': '3d'}) # 3D曲面图 surf1 = ax1.plot_surface(X, Y, ds_height_map, cmap='terrain', linewidth=0, antialiased=True) fig.colorbar(surf1, ax=ax1, shrink=0.6) ax1.set_title('Diamond-Square Terrain (3D Surface)') ax1.view_init(elev=30, azim=45) # 3D线框图 ax2.plot_wireframe(X, Y, ds_height_map, color='darkgreen', linewidth=0.5, rstride=5, cstride=5) ax2.set_title('Diamond-Square Terrain (3D Wireframe)') ax2.view_init(elev=20, azim=75) plt.tight_layout() plt.show() ``` ### 三、进阶技术与优化 生成基础地形后，通常需要进一步处理以增强真实感或提升性能。 1. **纹理映射与着色**：使用 `cmap` 参数根据高度进行着色（如 `'terrain'`, `'gist_earth'`）是最简单的方法。更高级的做法是将高度图与卫星影像、岩层贴图等进行映射[ref_2][ref_6]。 2. **地形平滑与侵蚀模拟**：对生成的高度图应用高斯模糊或特定的水侵蚀、热侵蚀算法，可以使地形更自然[ref_6]。 3. **LOD（细节层次）技术**：对于大规模地形渲染，距离相机远的地块使用低精度网格，近处使用高精度网格，这是平衡画面质量与渲染性能的关键技术[ref_3]。 4. **光照计算**：添加方向光源，并根据网格法线计算光照强度，可以极大地增强地形的立体感。法线向量可以通过计算高度图的梯度来近似[ref_6]。 综上所述，选择哪种方法取决于具体需求：若追求地理真实性，应采用基于真实 DEM 数据的方法[ref_5]；若追求开发效率和可控性，过程化生成算法（如 Perlin 噪声和 Diamond-Square 算法）是更佳选择[ref_2][ref_6]。实际项目中，常将两者结合，例如用真实数据作为基础，再用噪声添加细节，或通过过程化算法生成基础地形，再根据真实生物群系规则进行修饰。