高光谱成像（二）光谱角映射 SAM

光谱角映射（Spectral Angle Mapping，简称 SAM）是高光谱图像分析中的一种经典分类与匹配方法。它基于像素光谱与参考光谱之间的角度相似性，实现对地物类型的识别。相比其他方法，SAM 对光照强度和地形起伏的干扰具有鲁棒性，因此在遥感领域得到广泛应用。本文将从光谱角映射的基本概念出发，逐步深入讲解其数学原理、实现步骤，并通过代码示例演示如何从零实现 SAM 分类。内容设计由浅入深，适合初学者到进阶开发者。—## 1. 光谱角映射的基本概念### 1.1 什么是光谱角？在高光谱图像中，每个像素对应一条光谱曲线，该曲线记录了该像素在不同波长下的反射率或辐射值。如果我们将 n 个波段视为 n 维空间，那么每个像素的光谱值就是一个 n 维向量。两个光谱向量之间的“角度”可以衡量它们的形状相似性：角度越小，说明光谱形状越接近；角度越大，则说明光谱差异越大。光谱角映射正是利用这一特性，计算每个像素光谱与参考光谱之间的夹角，从而判断像素属于哪一类地物。### 1.2 SAM 的核心优势-对光照强度不敏感：因为角度只取决于光谱形状，而不受整体亮度缩放的影响。-简单高效：无需复杂的训练过程，只需提供参考光谱即可。-物理意义明确：角度越小，光谱越相似，直观易懂。—## 2. SAM 的数学原理### 2.1 向量点积与角度设两个光谱向量为a= [a₁, a₂, …, aₙ] 和b= [b₁, b₂, …, bₙ]，则它们之间的夹角 θ 由以下公式计算：[\cos(\theta) = \frac{\mathbf{a} \cdot \mathbf{b}}{|\mathbf{a}| \cdot |\mathbf{b}|}]其中：- (\mathbf{a} \cdot \mathbf{b} = \sum_{i=1}^n a_i b_i) 是点积- (|\mathbf{a}| = \sqrt{\sum_{i=1}^n a_i^2}) 是向量模长最终的角度为：[\theta = \arccos\left( \frac{\mathbf{a} \cdot \mathbf{b}}{|\mathbf{a}| \cdot |\mathbf{b}|} \right)]### 2.2 分类决策对于每个像素，计算其与所有参考光谱的角度，选择角度最小的参考类别作为该像素的类别。如果最小角度大于某个阈值（例如 0.1 弧度），则该像素可能不属于任何已知类别，视为“未知”。—## 3. SAM 的实现步骤1.准备数据：高光谱图像（三维数组：行 × 列 × 波段）和参考光谱（二维数组：类别数 × 波段）。2.归一化处理：对每个光谱向量进行归一化，只保留方向信息。3.计算角度：对每个像素，计算其与所有参考光谱的夹角。4.分类：为每个像素选择夹角最小的参考类别。5.阈值过滤：排除角度过大的像素。—## 4. 代码示例一：从零实现 SAM 分类以下代码演示了如何在模拟的高光谱数据上实现 SAM 分类。我们使用 numpy 进行向量化计算，以提高效率。pythonimport numpy as npdef sam_classify(image, ref_spectra, angle_threshold=0.1): """ 使用光谱角映射对高光谱图像进行分类 参数: image: 高光谱图像，形状为 (rows, cols, bands) ref_spectra: 参考光谱矩阵，形状为 (num_classes, bands) angle_threshold: 角度阈值（弧度），超过此值的像素标记为 -1 返回: class_map: 分类结果，形状为 (rows, cols) """ # 获取图像尺寸 rows, cols, bands = image.shape num_classes = ref_spectra.shape[0] # 归一化：对每个光谱向量除以模长，只保留方向 # 处理图像：避免除以零 image_norm = image / (np.linalg.norm(image, axis=2, keepdims=True) + 1e-10) # 处理参考光谱 ref_norm = ref_spectra / (np.linalg.norm(ref_spectra, axis=1, keepdims=True) + 1e-10) # 计算所有像素与所有参考光谱的点积 -> 得到余弦相似度矩阵 # 将图像展平为 (rows*cols, bands) 形状 flat_image = image_norm.reshape(-1, bands) # 形状: (N, bands) # 点积: (N, bands) @ (bands, num_classes) = (N, num_classes) cos_sim = flat_image @ ref_norm.T # 形状: (N, num_classes) # 将余弦相似度限制在 [-1, 1] 内，避免数值误差导致 arccos 出错 cos_sim = np.clip(cos_sim, -1.0, 1.0) # 计算角度（弧度） angles = np.arccos(cos_sim) # 形状: (N, num_classes) # 对每个像素，找到最小角度及其对应类别 min_angles = np.min(angles, axis=1) # 形状: (N,) class_indices = np.argmin(angles, axis=1) # 形状: (N,) # 应用阈值：如果最小角度大于阈值，标记为 -1（未知） class_indices[min_angles > angle_threshold] = -1 # 恢复为原始图像形状 class_map = class_indices.reshape(rows, cols) return class_map# 示例：生成模拟数据np.random.seed(42)rows, cols, bands = 10, 10, 20# 模拟图像：随机生成一些像素image = np.random.rand(rows, cols, bands)# 模拟参考光谱：假设有 3 个类别ref_spectra = np.array([ [0.5, 0.6, 0.7, 0.8, 0.9, 0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8], [0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, 0.1, 0.2, 0.3, 0.4, 0.5], [0.8, 0.7, 0.6, 0.5, 0.4, 0.3, 0.2, 0.1, 0.2, 0.3, 0.4, 0.5, 0.6, 0.7, 0.8, 0.9, 0.8, 0.7, 0.6, 0.5]])# 执行分类class_map = sam_classify(image, ref_spectra, angle_threshold=0.3)print("分类结果（10x10 网格）：")print(class_map)# 输出示例：每个数字代表类别索引（0,1,2）或 -1（未知）代码说明：- 使用np.linalg.norm进行向量归一化。- 通过矩阵乘法一次性计算所有像素与所有参考光谱的余弦相似度。- 使用np.arccos将余弦转换为角度。- 通过阈值过滤，将不相似像素标记为 -1。—## 5. 进阶应用与优化### 5.1 并行加速对于大规模高光谱图像（例如 1000×1000 像素，200 波段），上述方法仍然依赖 CPU 串行计算。我们可以利用 GPU 加速（如 CuPy）或分块处理来提升性能。### 5.2 自适应阈值固定阈值可能不适用于复杂场景。一种改进方法是根据数据分布动态计算阈值，例如使用所有像素角度的统计值（如 90% 分位数）作为阈值。### 5.3 多尺度 SAM在某些应用中，需要同时考虑局部和全局光谱特征。多尺度 SAM 通过在不同空间尺度上计算光谱角，提高分类的鲁棒性。—## 6. 代码示例二：使用光谱角进行图像匹配与可视化以下示例展示了如何在真实高光谱数据（使用开源的 Indian Pines 数据集模拟）上应用 SAM，并可视化分类结果。pythonimport numpy as npimport matplotlib.pyplot as pltdef plot_sam_classification(image, ref_spectra, angle_threshold=0.15): """ 运行 SAM 分类并可视化结果 """ class_map = sam_classify(image, ref_spectra, angle_threshold) # 创建一个彩色分类图 # 假设有 3 个类别，分别用红、绿、蓝表示 color_map = np.zeros((*class_map.shape, 3), dtype=np.uint8) # 类别 0 -> 红色 color_map[class_map == 0] = [255, 0, 0] # 类别 1 -> 绿色 color_map[class_map == 1] = [0, 255, 0] # 类别 2 -> 蓝色 color_map[class_map == 2] = [0, 0, 255] # 未知类别 -> 黑色 color_map[class_map == -1] = [0, 0, 0] # 显示原始图像（选取三个波段模拟 RGB） rgb_image = image[:, :, [5, 10, 15]] # 假彩色 rgb_image = (rgb_image - rgb_image.min()) / (rgb_image.max() - rgb_image.min()) # 归一化到 0-1 plt.figure(figsize=(12, 5)) plt.subplot(1, 2, 1) plt.imshow(rgb_image) plt.title("原始图像（假彩色）") plt.axis('off') plt.subplot(1, 2, 2) plt.imshow(color_map) plt.title("SAM 分类结果") plt.axis('off') plt.tight_layout() plt.show()# 生成模拟的高光谱图像（类似真实数据）np.random.seed(123)rows, cols, bands = 50, 50, 30image = np.random.rand(rows, cols, bands)# 定义三个参考光谱（模拟不同地物）ref_spectra = np.array([ np.sin(np.linspace(0, np.pi, bands)) + 0.5, # 类别 0 np.cos(np.linspace(0, np.pi, bands)) + 0.5, # 类别 1 np.exp(-np.linspace(0, 3, bands)) * 0.8 + 0.2 # 类别 2])# 运行并显示结果plot_sam_classification(image, ref_spectra, angle_threshold=0.2)代码说明：- 使用假彩色显示原始图像（选取三个波段模拟 RGB）。- 将分类结果映射为彩色图像，便于视觉分析。- 阈值设为 0.2 弧度，大于该值的像素视为未知。—## 7. 总结本文从光谱角映射的基本概念出发，逐步讲解了其数学原理、实现方法和优化方向。通过两个完整的代码示例，我们展示了如何从零构建 SAM 分类器，并将其应用于模拟高光谱数据。### 关键要点回顾：1.核心思想：通过计算光谱向量之间的夹角来衡量相似性，角度越小越相似。2.数学基础：基于向量点积和 arccos 函数，计算简单且物理直观。3.实现技巧：利用矩阵运算和归一化避免循环，提高效率；阈值过滤可处理未知类别。4.应用价值：SAM 特别适合光照变化剧烈的遥感场景，如卫星影像、无人机高光谱数据分类。### 进一步学习建议：- 尝试将 SAM 与机器学习方法（如 SVM）结合，提升分类精度。- 研究光谱角与混合像元分解的关系，探索亚像素级分类。- 使用真实高光谱数据集（如 AVIRIS、HyMap）进行实验，熟悉数据预处理流程。希望本文能帮助你掌握 SAM 的核心思想与实现方法，为高光谱图像分析打下坚实基础。