实验选用美国Headwall公司产的Hyperspec VNIR-A型号高光谱成像光谱仪。 成像系统由1 600× 2 500高分辨率CCD镜头、 全反射式光栅、 狭缝及镜头组成。 本实验通过以2.5 nm(平均间隔)波段间隔, 对实验样品在380~1 000 nm共234个波段进行光谱成像。 光谱仪具体参数如表1所示。

表1

Table 1

表1(Table 1)

表1 高光谱成像相机参数 Table 1 Hyperspectral imaging camera parameters Hyperspec VNIR-A高光谱成像光谱相机 光谱范围 380~1 000 nm 光谱分辨率 2~3 nm 光谱通道数 837 空间通道数 1 004 狭缝宽度 25 μ m 数值孔径 F/2.0 光谱采样值 0.74 nm· pixel-1 Smile & Keystone 2× 2 3× 3 484 探测器类型 CCD 像元间距 7.4 μ m 最大帧频 90 fps 相机接口 Cameralink 重量 ~2.8 kg

表1 高光谱成像相机参数 Table 1 Hyperspectral imaging camera parameters

高光谱成像系统如图1(a, b, c)所示, 实验样本图像的采集在暗室中进行, 系统光源为卤素灯模拟的太阳光, 通过光源控制调节系统光源亮度并稳压, 避免图像过饱和; 设置相机积分时间和凝视时间防止图像发生畸变; 调节镜头焦距, 使待测目标成像清晰; 测量待测样本和相机狭缝的距离, 并计算设置移动载物平台的推扫速度, 使采集到的样本图像清晰且不失帧, 保证样本光谱信息和空间信息的完整性。

图1

Fig.1

	Figure Option ViewDownloadNew Window
	图1 高光谱成像系统 (a): 原理图; (b): 实物图; (c): 暗室环境图Fig.1 Hyperspectral imaging system (a): Schematic diagram; (b): Scene diagram; (c): Darkroom environment diagram

(1)模拟破片散布

以铁质破片为研究对象, 选取从爆炸场收集的铁质破片两枚、 石子两个、 树叶一片放在载有沙土的沙盘上。 采集其高光谱图像, 得到模拟的空间图像及其光谱数据, 空间图形如图2所示。

图2

Fig.2

	Figure Option ViewDownloadNew Window
	图2 模拟破片散布空间图Fig.2 Spatial diagram of simulated fragment distribution

(2)建立破片识别模型数据库

根据模拟破片散布场的特点, 由于沙土与破片、 石子、 树叶, 相对无规则且复杂, 将其视为第一背景, 通过空间分割将其去除。 然后对沙土上的目标用谱域信息进行分类识别, 因此需要建立包含石子和树叶的破片识别数据库。 另取铁质破片、 石头、 树叶若干样本, 采集它们的高光谱图像, 感兴趣区域提取像素点作为样本库, 为了使提取的样本点具有代表性, 采集数据时要进行样本多角度采集, 选取样本点也要广泛全面。 图3(a, b, c)为采集的破片及其背景样本图例。

图3

Fig.3

	Figure Option ViewDownloadNew Window
	图3 破片及其背景样本图列 (a): 铁质破片; (b): 石子; (c): 树叶Fig.3 Fragments and their background samples (a): Iron fragments; (b): Stones; (c): Leaves

1.3.1 图像数据预处理

为了消除相机暗电流和光线带来的图像噪声, 需要对采集的高光谱图像进行黑白校正。 分别采集反射率为99.9%的白板得到白板标定图像R_white和盖上相机盖后的全暗参考图像R_dark, 然后对采集的原始图像R_raw, 采用校正公式(1)处理得到校正图像R_norm, 用于后续分析。

${{R}_{\text{norm}}}=\frac{{{R}_{\text{raw}}}\ \ \ -\ \ \ {{R}_{\text{dark}}}}{{{R}_{\text{white}}}\ \ \ -\ \ \ {{R}_{\text{dark}}}}\ \ \ $(1)

同时为了更加清楚直观地了解空间图形特征, 对模拟破片散布的空域图像进行去噪和增强等预处理, 方便后期的区域分割提取。

1.3.2 破片识别模型

(1) ReliefF特征选择

由于实验样本比较少, 光谱维度较高, 为避免信息冗余和出现过拟合, 采用有监督的特征选择方法ReliefF对样本数据进行降维, 经过训练选择对分类识别任务权重高的的波段来代替原始高维数据, 去除冗余。

经过改进的ReliefF算法主要处理多分类问题。 通过每次从训练样本集中随机取出一个样本R, 然后从与R同类的样本集中找出R的k个近邻样本(near Hits), 从每个R的不同类的样本集中寻找出k个近邻样本(near Misses), 然后更新每个特征的权重^[7]。 计算公式如式(2)所示

$\begin{array}{*{35}{l}} \text{W}\left( \text{A} \right)=\text{W}\left( \text{A} \right)-\overset{k}{\mathop{\underset{j-1}{\mathop \sum }\,}}\,\text{diff}\left( A,R,{{H}_{j}} \right)/\left( mk \right)+ \\ \underset{C\ne \text{class}\left( R \right)}{\mathop \sum }\,\left\{ \frac{p\left( C \right)}{1-p\left[ \text{Class}\left( R \right) \right]}\overset{k}{\mathop{\underset{j-1}{\mathop \sum }\,}}\,\text{diff}\left[ A,R,{{M}_{j}}\left( C \right) \right] \right\}/\left( mk \right) \\ \end{array}$(2)

式(2)中, m为设定的迭代次数, class(R)为样本R所述类别, p(C)为类别C的先验概率, p[Class(R)]为样本R所属类别的先验概率, M_j(C)为C中的第j个最近邻样本, diff(A, R₁, R₂)为样本R₁和样本R₂在特征上A的差, 其计算公式如式(3)所示

diff(A, R1, R2)=0R1[A]=R2[A]1R1[A]≠R2[A](3)

经过训练集样本的训练, 输出每个波段在分类中所占权重, 权重越高表示该波段分类能力越强, 作为特征波段输出。 在分类任务中根据每增加一个特征, 训练样本需要增加10的原则, 设置样本的抽样次数为100, 设定输出特征波段的权重阈值, 最近邻样本个数设为4个; 从234个波段中输出63个特征子集, 在63个波段下铁质破片、 石子和树叶的参考光谱图分别为如图4(a, b, c)所示。

图4

Fig.4

	Figure Option ViewDownloadNew Window
	图4 (a)、 (b) 和(c) 分别为破片、 石子和树叶的参考光谱Fig.4 (a), (b) and (c) are the reference spectra of fragments, stones and leaves

从图4可知, 铁质破片、 石子、 树叶在这63个波段中均存在不同位置的波峰(波谷), 且波形趋势与反射率明显不同, 表明特征选择提取的63个特征波段可以代表原始样本用于分类。

(2)决策树分类模型

决策树作为一种有监督的机器学习方法, 能够有效处理光谱反射率的二分类或多分类问题。 利用信息论中的信息增益寻找样本集中具有最大信息量的波段, 将其作为决策树的一个结点, 再根据该波段的不同取值建立树的分支, 再在每个分支子集中重复建立下层结点和分支, 从而完成分类任务^[8]。 其中决策树的每一个节点代表的是对其中某个特征的测试, 分支则代表这个特征的测试结果, 决策树的每一个叶子节点代表了一类样本。

构造决策树的具体过程^[9]: 首先寻找初始分裂, 整个训练集作为产生决策树的集合, 训练集每个样本必须是带有标签的, 以此来决定哪个属性能够做出最准确的判断, 做法是对所有属性域的分裂指标做出量化, 计算每个分裂的多样性指标来求出最好的一个分裂, 然后, 重复上述步骤, 直至每个叶节点的判断结果都属于同一类且增长到一棵完整的树。 然后从树的底层开始向根节点计算, 判断是否需要剪枝, 通过剪枝操作后得到一个泛化能力更强的决策树, 用于分类任务。

将经过去噪声、 黑白校正、 反演光谱反射率信息等预处理后的样本高光谱数据, 感兴趣区域随机提取铁质破片样本点250个, 石子样本点250个, 树叶样本点250个, 三类样本像素点共750个样本点, 随机选取600个点作为训练集其余作为测试集, 通过训练后得到破片、 石子、 树叶的预测准确度分别为88%、 88%、 94%, 平均预测准确度为90%的决策树分类模型。

1.3.3 空间分割

为了简化识别模型、 提高识别精度、 美化识别结果, 采用边缘检测结合区域生长的方法进行空间分割。 具体为通过Sobel边缘检测得到的边缘来限制区域生长, 从而提取到破片散布场中相对有规则、 有形状的目标。

(1) Sobel边缘检测

Sobel算子是一种用于边缘检测蹬离散微分算子, 其结合了高斯平滑和微分求导, 通过计算破片散布图像中明暗程度近似值, 根据图像边缘旁边明暗程度把该区域内超过某个数的特定点记为边缘。 对于噪声较夺, 灰度渐变的高光谱空域图像, Sobel算子的边缘定位更准确。 其算法模板^[10]如式(4)所示, 其中d_x表示水平方向, d_y表示垂直方向。

dx=-101-202-101dy=-1-2-1000121(4)

(2) 区域生长

将在破片散布场中, 破片、 树叶、 石子各自区域内像素集像素相似。 对每个需要分割的目标区域设定一个或者多个种子像素作为生长的起点, 就可以将种子像素周围邻域中与种子像素具有相同或相似性质的像素合并到种子像素所在的区域中, 将这些新像素当做新的种子像素继续进行上面的过程, 直到再没有满足条件的像素可被包括进来。 这样就可以将破片散布场中的待分割目标区域提取出来。

1.3.4 空谱联合

针对破片高光谱图像特点, 在前期研究的基础上, 采用前后级联空谱融合方法。 在空域分割信息基础上, 对各分割区域进行谱域的类型识别, 以获得准确的目标和背景信息。 具体流程如图5所示。

图5

Fig.5

	Figure Option ViewDownloadNew Window
	图5 空谱联合识别流程图Fig.5 Flow chart of joint identification of space spectrum

将高光谱相机采集的模拟破片散落在背景中的空间图像经过增强、 去噪等预处理后, 采用边缘检测结合区域生长的前后级联空谱联合法进行空间分割, 然后对分割的结果进行形态学处理后得到空域分割结果图。 其中灰度值255代表了提取到的待分类目标。 采用matlab里自带的image labeler对空间原始图进行逐像素标记生成标定图像, 用于评价该空间分割算法的分割效果, 主要结果如图6所示。

图6

Fig.6

	Figure Option ViewDownloadNew Window
	图6 空域分割结果图 (a): 空间增强; (b): 边缘检测结合区域生长结果图; (c): 空间分割二值图; (d): Image labeler标定二值图Fig.6 The result of spatial segmentation (a): Spatial enhancement; (b): The result of edge detection combined with region growth; (c): Binary map of spatial segmentation; (d): Binary map of image labeler calibration

从图6(a)到图6(c)呈现了分割的过程, 经过图像增强的图像(a)中的待分割样本轮廓更加清晰, 经过边缘检测后在边缘检测的基础上进行区域生长得到图6(b), 在图6(b)中存在过分割和欠分割问题, 通过形态学处理补充区域内的欠分割和减少过分割得到最终分割图6(c)。 图6(d)展现了手动标定的结果, 用于和算法分割图像(c)进行对比评价。

采用交并比(IOU)和真阳率(TPR)两种指标对空间分割结果进行评价, 计算公式如式(5)

$\text{IOU}=\frac{TP}{FN+TP+FP}~\text{TPR}=\frac{TP}{TP+FN}$(5)

式(5)中, TP为分割结果和标定结果的交集, FN为欠分割部分, FP为过分割部分。 对空域分割评价如表2。

表2

Table 2

表2(Table 2)

表2 空域分割评价表 Table 2 Airspace segmentation evaluation table 评价指标 IOU TPR 结果 0.935 0.974

表2 空域分割评价表 Table 2 Airspace segmentation evaluation table

交并比(IOU)是算法分割和标定结果的交叠率, 即它们的交集与并集的比值, 比值越高说明分割效果越好, 当IOU为1的时为完全分割。 真阳率(TPR)代表了分割结果与标定结果的交集与标定结果的比值, 即在标定结果上实际分割了多少, 当TPR为1时也代表完全分割。 本次实验结果0.935的IOU和0.974的TPR均大于0.9, 表明该算法取得了很好的分割效果。

在空间分割的二值图上, 设计相关算法读取目标区域每一个像素的光谱信息, 然后输入到训练好的决策树模型判断得到该点的类别, 对不同类别的像素点赋予不同的颜色, 将分类结果可视化得到图7(a)。 然后设计了区域数目最大点代表该区域类别算法, 通过计数比较后赋值得到空谱联合最终识别结果见图7(b)。

图7

Fig.7

	Figure Option ViewDownloadNew Window
	图7 空谱联合结果图 (a): 分类结果图; (b): 最终识别结果图Fig.7 Combined result of empty spectrum (a): Classification result; (b): Final recognition result

在图7(a)中, 黄色代表石子的像素点、 银白色代表铁质破片、 绿色代表树叶。 对谱域破片分类识别图和前后级联空谱联合图进行像素点计数分析。 对于决策树分类模型的预测准确度, 采用算法分割图像和标定图像交集评价谱域的识别准确度, 避免过分割引起的误差, 分类结果评价如表3。 对空谱联合的识别结果与标定图像的交集作为最终识别率进行评价即真阳率(TPR), 空谱联合最终识别结果评价如表4所示。

表3

Table 3

表3(Table 3)

表3 光谱分类结果统计表 Table 3 Statistics of spectral classification results 样本 验证点数 样本点归类 准确度 /% 石子 铁质破片 树叶 石子 146 172 127 170 16 326 2 676 87 铁质破片 50 484 43 211 6 158 1 115 86 树叶 213 438 204 900 4 629 3 909 96

表3 光谱分类结果统计表 Table 3 Statistics of spectral classification results

表4

Table 4

表4(Table 4)

表4 空谱联合识别结果统计表 Table 4 Statistical table of joint identification results of space spectrum 样本 样本点数 识别点数 识别率/% 石子 155 502 146 172 94 铁质破片 52 045 50 484 97 树叶 217 794 213 438 98

表4 空谱联合识别结果统计表 Table 4 Statistical table of joint identification results of space spectrum

从空间分割的评价结果来看, 边缘检测结合区域生长的空间分割效果很好, 达到了预期效果, 可以将模拟破片散布场中的有形态目标分割出来, IOU和TPR均在90%以上, 但存在不足, 内部的边缘影响到区域生长, 从而影响到分割精度, 后续研究需弥补此不足。

上述对谱域识别的统计结果表明, 标定图像中三类样本原本点数分别为155 502、 52 045、 217 794, 识别到的点数分别为146 172、 50 484、 213 438, 基于前后级联方法的空谱联合识别方法以先聚类后识别的方法将识别率提高到94%、 97%、 98%。 基于决策树识别分类模型训练得到的预测准确度分别为88%、 88%、 94%, 满足设计的空谱联合分类算法识别的要求, 并且成功地将模拟破片散布场中的目标分类表示, 测试样本点的分类准确度分别为87%、 86%、 96%。 影响分类准确度的因素主要是提取的样本点太少, 以至于分类识别模型准确度比较低, 另一方面是石子与铁质破片的光谱反射率特征曲线有一定相似性, 这一点在表3和表4中铁质破片与石头的相互误判比较高体现出来。