引用本文
李辉, 刘姁升, 蒋金豹, 陈绪慧, 张帅, 唐珂, 赵新伟, 杜兴强, 玉龙飞雪. 基于冬小麦高光谱图像的天然气微泄漏胁迫区域提取[J]. 光谱学与光谱分析, 2024,44(3): 770-776.
LI Hui, LIU Xu-sheng, JIANG Jin-bao, CHEN Xu-hui, ZHANG Shuai, TANG Ke, ZHAO Xin-wei, DU Xing-qiang, YU LONG Fei-xue. Extraction of Natural Gas Microleakage Stress Regions Based on Hyperspectral Images of Winter Wheat[J]. Spectroscopy and Spectral Analysis, 2024,44(3): 770-776.
Doi:10.3964/j.issn.1000-0593(2024)03-0770-07  
Permissions
Copyright©2024, 《光谱学与光谱分析》期刊社
《光谱学与光谱分析》期刊社 所有
基于冬小麦高光谱图像的天然气微泄漏胁迫区域提取
李辉1, 刘姁升2, 蒋金豹3,*, 陈绪慧4, 张帅5, 唐珂1, 赵新伟1, 杜兴强1, 玉龙飞雪1
1.中国四维测绘技术有限公司, 北京 100086
2.鄂尔多斯应用技术学院, 内蒙古 鄂尔多斯 017000
3.中国矿业大学(北京)地球科学与测绘工程学院, 北京 100083
4.生态环境部卫星环境应用中心, 北京 100094
5.中国资源卫星应用中心, 北京 100094
*通讯作者 e-mail: jjb@cumtb.edu.cn

作者简介: 李 辉, 1994年生, 中国四维测绘技术有限公司工程师 e-mail: 18303015502@163.com

收稿日期: 2022-10-19 修回日期: 2023-04-27 接受日期: 2023-04-27
基金: 国家自然科学基金项目(42271389, 41571412)资助
摘要

天然气在能源结构中逐渐占据重要地位。 由于天然气管道和储气库常年埋藏于地下, 无氧腐蚀、 自然灾害、 注入井和管道口松懈等因素会导致气体发生泄漏。 在地下储存天然气大规模泄漏之前确定泄漏点的位置并做出早期判断和预警十分必要。 以冬小麦为研究对象, 采集4期高光谱图像数据, 融合其图像、 空间、 时相特征, 探寻天然气胁迫下冬小麦胁迫范围半径和时长之间的关系, 并间接探测天然气微泄露点。 一方面对连续统去除后的冬小麦高光谱图像数据进行连续小波变换并利用 CWTmexh指数[ CWTmexh= CW7702 /(1- CW487)× CW550]对高光谱胁迫组图像进行微泄露点信息提取, 另一方面提取高光谱图像数据的PCA特征, 基于SVM分类器提取天然气胁迫区域。 最后将天然气微泄露识别结果进行数学形态学分析, 并利用最小二乘对胁迫区域进行圆曲线拟合, 探索天然气泄漏胁迫半径与胁迫天数的关系。 结果表明: (1) CWTmexh指数应用到成像高光谱数据表现出较好的识别性能; (2) SVM分类器可基于光谱差异性特征识别冬小麦胁迫区域, 分类精度较好, 最大分类精度可以达到99.25%, Kappa系数为0.97, 且识别精度随天然气胁迫持续而增加; (3) 冬小麦受胁迫区域半径和通气天数呈现强烈的线性相关。 因此, 结果表明, 在冠层尺度和低空尺度通过高光谱遥感监测地表植被间接识别天然气微泄漏点具有可行性, 可以预测地下天然气微泄漏随着时间变化引起的胁迫区域变化。 该工作为星载高光谱遥感监测地下储存天然气泄漏点提供科学依据, 为以后的工程应用提供技术支持。

关键词: 天然气微泄漏; 冬小麦; 高光谱图像; 空间特征
中图分类号:P237 文献标志码:A
Extraction of Natural Gas Microleakage Stress Regions Based on Hyperspectral Images of Winter Wheat
LI Hui1, LIU Xu-sheng2, JIANG Jin-bao3,*, CHEN Xu-hui4, ZHANG Shuai5, TANG Ke1, ZHAO Xin-wei1, DU Xing-qiang1, YU LONG Fei-xue1
1. China Siwei Surveying and Mapping Technology Co., Ltd., Beijing 100086, China
2. Ordos Institute of Technology, Ordos 017000, China
3. College of Geosciences and Surveying Engineering, China University of Mining and Technology, Beijing 100083, China
4. Satellite Application Center for Ecology and Environment, MEE, Beijing 100094, China
5. China Centre for Resources Satellite Data and Application, Beijing 100094, China
*Corresponding author
Abstract

Natural gas has gradually occupied an important position in the energy structure. As natural gas pipelines and gas storage are buried underground all year round, oxygen-free corrosion, natural disasters, looseness of injection wells and pipelines, and other factors will lead to gas leakage. So, it is necessary to determine the location of leakage points and make early judgments and warnings before large-scale leakage from underground natural gas storage. This paper collected four periods of hyperspectral image data of winter wheat. It integrated the spatial-temporal-spectral features of hyperspectral data to explore the relationship between the radius and duration of winter wheat stress under natural gas stress, thus indirectly detecting the microleakage point of natural gas. On the one hand, the index CWTmexh( CWTmexh= CW7702 /(1- CW487)× CW550), constructed by continuous wavelet transform of the canopy spectra after continuum removal, was used to classify pixels into non-stress and stress with threshold segmentation. On the other hand, PCA features of hyperspectral image data are extracted, and natural gas stress regions are identified with the SVM classifier. Finally, the results of both threshold segmentation and SVM classification are analyzed by mathematical morphology, and the stress area is fitted with a circular curve using the least square to explore the relationship between the stress radius of natural gas leakage and the stress days. The results show that: (1) The CWTmexh index can be applied to imaging hyperspectral data, showing good recognition performance; (2) SVM classifier can recognize winter wheat stress areas based on spectral difference characteristics with good classification accuracy (i.e., the maximum classification accuracy of 99.25% and kappa coefficient is 0.97) and the recognition accuracy increases with the continuation of natural gas stress; (3) There is a strong linear correlation between the radius of the stressed area and the ventilation days of winter wheat. Results of this study showed that it is feasible to indirectly identify natural gas micro-leakage points through hyperspectral remote sensing by monitoring surface vegetation at the canopy and low altitude scales and can predict time-dependent changes associated with underground natural gas micro leakage stress. The results can provide a theoretical basis for monitoring the leakage points of underground natural gas storage by spaceborne hyperspectral remote sensing and provide technical support for future engineering applications.

Keyword: Natural gas micro-leakage; Winter wheat; Hyperspectral image; Spatial characteristics
引言

天然气在能源结构中逐渐占据重要地位, 近年来我国天然气消费量逐年增加, 2020年天然气消费量已超越3 600亿m3, 在一次能源消费比例中占据8.3%~10%的高比[1]。 天然气储气库储存量大、 经济合理、 安全系数高, 在优化供气系统和保障供气安全上建设储气库颇为重要。 2015年我国天然气管道建设的长度约7.2万km[2], 已建成的地下储气库据不完全统计为25座[3]。 由于地下管道和地下储气库常年埋藏于地下, 无氧腐蚀、 自然灾害、 注入井和管道口松懈等因素会导致气体发生泄漏[4]。 天然气管线或储气库泄漏引起的爆炸事故, 带来的灾害、 伤亡严重威胁人类生命、 财产、 环境安全。

地下储存的油气等烃类化合物轻微泄漏会导致土地表面的岩石、 土壤发生裂变和腐蚀, 研究表明利用高光谱遥感研究岩石、 土壤的裂变腐蚀来探寻油气等烃类化合物的泄漏点可行[5]。 而埋藏在地下的天然气管道或是储气库轻微泄漏初期, 腐蚀能力较弱, 难以使地面的岩石、 土壤产生裂变和腐蚀。 但天然气微泄漏会占据土壤中O2含量限制地表植被光合和呼吸作用, 进而导致地面植被生长异常[6]。 高光谱遥感技术以植被为媒介可间接探测地下天然气管道泄漏点。

目前国内外已有许多学者致力于利用高光谱遥感技术分析植被胁迫症状并间接识别油气微泄漏点[7, 8, 9, 10, 11, 12, 13]。 崔鑫利用航空高光谱数据, 结合野外实测光谱, 提取准噶尔盆地东北缘油气微渗漏引起的烃类及相关蚀变矿物信息, 结果显示, 研究区的油气渗漏异常区具有环带状分布特征[11]。 Noomen等设计盆栽玉米试验, 分别通入天然气、 CH4、 C2H6进行胁迫, 并利用波段深度(band depth, BD)指数来识别不同气体泄漏胁迫下的玉米; 在实验中发现, 三种气体的泄漏会导致土壤中氧气含量的减少[12]。 赵欣梅利用EO-1 Hyprion成像光谱数据以植被为媒介, 观察植被光谱变化进而识别油气渗漏点。 但这些研究都只局限于室内, 未在野外开展实验。 van der Werff等利用机载高光谱数据, 通过多次霍夫变换(Hough transform)结合天然气微泄漏的行为规律, 检测出近似圆形的天然气微泄漏点[8]; van der Werff等依据植被光谱反射率一阶微分最大值在受到外界胁迫后会“ 红边蓝移” , 利用Hymap成像高光谱数据, 对野外管道微泄漏点进行了成功探寻[13]; 李辉[7]等模拟天然气地下管道和储气库微泄漏对地表冬小麦的胁迫实验, 采集冬小麦冠层光谱信息, 利用连续小波变换构建CWTmexh指数模型, 可以较好识别受胁迫和健康冬小麦。 但现有研究大多利用单一的植被冠层光谱或成像光谱数据, 将冠层尺度构建的指数模型应用到成像尺度的案例不多。

本工作拟将冬小麦冠层光谱构建的CWTmexh指数模型应用到成像高光谱数据, 分析CWTmexh指数模型的可识别性、 稳定性, 探寻天然气微泄漏胁迫下的空间特征。

1 实验部分
1.1 实验设计

实验区位于北京市大兴区长子营镇, 实验场地长40 m, 宽20 m, 设计8个冬小麦实验小区, 分别为4个天然气胁迫实验区和4个对照区, 实验区与对照区间隔排列, 每个区大小为2.5 m× 2.5 m, 区之间的间隔均为0.5 m。 天然气泄漏点位于实验区的中心下方60 cm处(图1中红色圆点即微泄漏位置)。 实验时间为2016.10-2017.06。 天然气泄漏速率为1 L· min-1, 于2017年4月11日在小麦返青期前开始持续通气至6月1日实验结束。 种植的冬小麦品种为京冬14号, 该品种成穗率较高, 抗倒伏、 抗病性较好。

图1 实验区空间分布图Fig.1 Spatial distribution map of the experimental area

1.2 高光谱图像采集与预处理

实验采用SOC710-VP高光谱成像仪, 其波长范围为400~1 000 nm, 镜头焦距为8 mm, 视场角为21° 。 采集高光谱成像数据时, 为完全覆盖实验田2.5 m× 2.5 m的范围, 将成像光谱仪搭载在5 m高的升降平台上, 并保持成像光谱仪的镜头垂直向下。 选择天气晴朗, 风力小于三级, 无云或少云的时段进行采集, 观测时间为10:00-14:00之间, 此时太阳高度角足够大。 共选取4期图像数据, 分别为5月1日, 5月11日, 5月21日, 5月31日, 其中5月1日的数据是在天然气胁迫第21天之后采集。 数据采集后用SOC710-VP配套的软件SRAnal710进行反射率转换, 得到反射率图像。 反射率图像数据选用MNF正逆变换和5点平滑相结合的方法进行降噪平滑处理。

图2 (a) 实验田概况; (b)高光谱成像数据采集Fig.2 (a) Experimental field overview; (b) Hyperspectral imaging data acquisition

1.3 图像处理及分析方法

处理及分析方法有连续统去除、 连续小波变换、 Kapur最大熵阈值分割、 PCA及SVM分类等。 具体处理流程如图3所示。

图3 基于冬小麦高光谱图像天然气微泄漏胁迫区域提取流程Fig.3 Extraction process of natural gas micro-leakage stress region based on winter wheat hyperspectral image

(1) CWTmexh指数模型

利用CWTmexh指数模型提取天然气微泄漏胁迫区域。 CWTmexh指数是基于墨西哥帽(Mexihat)母小波(尺度参数32)对连续统去除后的冠层光谱进行连续小波变换的特征, 具体请参照文献[7]。

CWTmexh=CW7702(1-CW487)×CW550(1)

式(1)中, CW为连续小波能量系数值, 其下标为波长。

(2) Kapur最大熵阈值分割算法

利用Kapur最大熵阈值分割区分图像的目标和背景。 Kapur最大熵阈值熵越大, 分布越均匀。 具体原理为: 使选择的阈值分割图像目标、 背景两部分灰度统计的信息量为最大。

(3) 主成分分析(PCA)

利用主成分分析(PCA)对高光谱数据有效降维并减小数据冗余。 PCA前几个主成分基本能够概括所有波段的95%以上信息, 且每个主成分能反应高光谱图像的不同信息和特征[10]。

(4) 支持向量机(SVM)

利用SVM对图像数据的PCA特征进行监督分类, 采用的SVM核函数为径向基(RBF)核函数。 其核心思想是在线性可分的数据当中, 寻找最优分隔面来对数据进行分隔, 能够最大程度地将待分样本分离, 且保证分隔距离最大[10]。

(5) 图像形态学分析

利用数学形态学开运算达到图像增强。 数学形态学分析的基本方法有腐蚀(erosion)、 膨胀(dilation)、 开运算(open)、 闭运算(close)。 开运算在操作上具有优势, 可以做到不明显改变目标面积的同时平滑目标的边缘。

(6) 最小二乘拟合圆

最小二乘在回归问题解算中用于估计和预测输入的量和输出变量两者之间存在的关联。 如图4所示为最小二乘拟合圆的示意图, 点(Xi, Yi)到圆心的距离平方和半径的平方差为σ i, 则

σi=di2-R2=(Xi-A)2+(Yi-B)2-R2=Xi2+Yi2-2AXi-2BYi+A2+B2-R2(2)

图4 最小二乘拟合圆示意图Fig.4 The diagram of the fitting circle used by the least squares

令Q为σ i的平方和, 则Q为关于A, B, R的函数, A, B, R使得Q值最小时的最优解即为所求的圆心坐标和半径。

Q(A, B, R)=∑1nσi2(3)

在最小二乘拟合前, 先对形态学分析开运算处理后的冬小麦二值化图像进行连通区域标记, 利用8邻接连通区域分析, 获取冬小麦在天然气胁迫区域下的轮廓和边界。

2 结果与讨论
2.1 基于CWTmexh指数的天然气微泄漏胁迫区域提取

CWTmexh指数对高光谱图像进行单一波段计算, 运算得到的结果如图5所示, 经过CWTmexh指数运算后, 土壤等出现异常值, 从灰度图上明显看到胁迫区域小麦较亮, 随着天然气胁迫的持续, 胁迫范围的晕圈在不断扩大。

图5 冬小麦高光谱指数运算结果
(a): 2017.5.1; (b): 2017.5.11; (c): 2017.5.21; (d): 2017.5.31Fig.5 Winter wheat hyperspectral index calculation results
(a): May 1, 2017; (b): May 11, 2017; (c): May 21, 2017; (d): May 31, 2017

Kapur最大熵阈值分割对指数运算结果进行处理, 得到二值化图像, 在Kapur最大熵阈值分割算法前先对灰度图像进行中值滤波和归一化处理。 定义白色区域(1值)为天然气泄漏胁迫区域, 黑色区域(0值)为其他, 如图6所示, 可有效区分出胁迫和健康小麦区域。 整体上胁迫区域都呈现出晕圈状空间特征, 随天然气泄漏胁迫持续, 晕圈状逐渐变大, 且逐渐趋于圆形。

图6 Kapur最大熵阈值分割后的二值化图像
白色区域代表胁迫区域, 黑色代表其他
(a): 2017.5.1; (b): 2017.5.11; (c): 2017.5.21; (d): 2017.5.31Fig.6 Binary images after Kapur maximum entropy threshold segmentation
The white area represents the stress area and the black represents the other
(a): May 1, 2017; (b): May 11, 2017; (c): May 21, 2017; (d): May 31, 2017

2.2 基于SVM分类器的天然气微泄漏胁迫区域提取

对四期图像数据进行PCA降维处理, 并保留其前四个主成分PC1, PC2, PC3, PC4, 图7所示为5月1日图像数据的前四个主成分的结果。 前四个主成分基本上能占据图像99.82%的信息量, 能很好的反映原始高光谱图像。

图7 5月1日数据的第四个主成分: (a)PC1, (b)PC2, (c)PC3, (d)PC4Fig.7 First four principal components of the data on May 1st: (a) PC1, (b)PC2, (c)PC3, (d)PC4

对图像数据进行SVM分类。 先对每期图像进行训练样本选取, 每期图像分布均匀的选取受胁迫、 健康小麦训练样本各40个, 验证样本各20个, 训练样本和验证样本不重叠, 且验证样本尽量只选择纯净像元。 在选取训练样本前, 为防止土壤对分类结果产生影响, 利用掩膜滤除部分土壤。 基于样本选取的结果, 计算两类样本的J-M距离, 区分两类样本的可分离性, J-M踞离可分离性原理请参照参考文献[7]。 最终选取的训练样本特征分离度如表1所示。

表1 选取的训练样本特征分离度 Table 1 Selected training sample feature separation

对PCA降维后的图像数据进行SVM分类, 生成分类图像, 红色代表受胁迫区域的小麦, 如图8所示。 从分类结果来看, 随着天然气的持续通入, 冬小麦受天然气胁迫区域逐渐增大, 且在空间特征上逐渐接近圆形, 与指数模型的分类结果保持一致。 使用验证样本建立混淆矩阵, 得到每期图像的分类精度和Kappa系数, 如表2所示。 4期图像数据分类精度都优于93%, Kappa系数大于0.83, 5.31日数据分类精度达到99.25%, Kappa系数为0.97。 由表2也可以看到, 冬小麦随着天然气胁迫持续, 受胁迫小麦的分类精度在增大, Kappa系数在增大, 表明可分性逐渐增加。

图8 SVM分类结果
(a): 2017.5.1; (b): 2017.5.11; (c): 2017.5.21; (d): 2017.5.31Fig.8 SVM classification results
(a): May 1, 2017; (b): May 11, 2017; (c): May 21, 2017; (d): May 31, 2017

表2 4期数据的分类精度 Table 2 Classification accuracy of 4 groups of data
2.3 数学形态学图像处理和最小二乘拟合圆分析

天然气小孔泄漏胁迫区域呈晕圈状空间特征, 随着胁迫时间的持续, 黑色晕圈在扩大, 逐渐呈圆形分布且具有连续性, 可以通过数学形态学和连通区域分析来提取天然气泄漏胁迫区域。 对二值化图像进行腐蚀、 膨胀开运算操作, 之后对其进行填充、 连续边界提取, 处理结果如图9所示。 利用CWT和SVM提取的冬小麦受天然气胁迫的天数与半径关系如表3所示。 结合天然气泄漏行为及其空间信息, 天然气胁迫天数和胁迫区域半径呈一元线性回归关系, 结果如图10所示。 CWTmexh指数模型下冬小麦受天然气泄漏胁迫天数与胁迫半径的关系为: y=0.013x+0.492(x≥ 21), R2=0.97, 基于SVM分类下冬小麦受天然气泄漏胁迫天数与半径的关系为: y=0.015x+0.439(x≥ 21), R2=0.94, 其中x均表示胁迫天数(单位: d), y均表示天然气胁迫半径(单位: m)。 该结果表明冬小麦在天然气泄漏胁迫下的半径与泄漏时间成近似线性关系。

图9 冬小麦受天然气胁迫范围拟合结果
白色区域代表胁迫区域, 黑色区域代表其他, 红色线代表拟合边界
(a): 2017.5.1; (b): 2017.5.11; (c): 2017.5.21; (d): 2017.5.31Fig.9 Fitting areas of winter wheat under natural gas stress
The white area represents the stress area, the black area represents the other, and the red line represents the fit boundary
(a): May 1, 2017; (b): May 11, 2017; (c): May 21, 2017; (d): May 31, 2017

表3 不同模型下冬小麦受天然气胁迫的天数与半径关系 Table 3 Relationship between days and radius of winter wheat under natural gas stress

图10 CWTmexh指数、 SVM分类提取出的冬小麦受胁迫范围和胁迫时间的一元回归分析
(a): CWTmexh指数; (b): SVM分类Fig.10 One-dimensional regression analysis of stress range and stress time of winter wheat extracted by CWTmexh index or SVM classification
(a): CWTmexh index; (b): SVM classification

3 结论

通过建立野外试验场, 模拟天然气地下管道和储气库微泄露对地表冬小麦的胁迫实验, 基于时间序列的冬小麦高光谱图像数据, 融合其图像、 空间、 时相特征, 探寻被天然气胁迫的冬小麦胁迫范围半径和胁迫时长之间的关系, 以提高并验证通过植被光谱信息间接检测天然气泄漏点的可识别性、 稳定性。 主要结论如下:

基于冬小麦实验数据, CWTmexh指数模型在图像尺度识别天然气泄漏特征上同样适用;

基于SVM分类器可以把受胁迫的冬小麦区域提取出来, 分类精度较好, 最大分类精度可以达到99.25%, Kappa系数为0.97。 冬小麦随着天然气胁迫持续, 受胁迫小麦的分类精度在增大;

小麦受胁迫区域半径和通气时间呈现强烈的线性相关, 可以预测地下天然气微泄漏随着时间变化引起的胁迫区域变化。

将CWTmexh指数应用到同期成像高光谱数据, 表现出较好的识别性能。 虽然该实验研究在探索利用高光谱识别天然气微泄漏上取得了一些进展, 但本实验的高光谱图像数据仅采用的是5 m高的地面升降平台系统, 而随着卫星数据源的不断丰富, 星载高光谱遥感大范围对天然气管道和储气库的监测会逐渐普及, 如何利用星载高光谱图像数据来识别天然气微泄漏点还有待探究。

参考文献
[1] ZHANG Ting-shou(张廷首). Chemical Engineering Design Communications(化工设计通讯), 2020, 46(3): 2. [本文引用:1]
[2] Tubb Rita. Pipeline and Gas Journal, 2016, 243(1): 238. [本文引用:1]
[3] CAI Xiao-long, LIU Yong-jun, DONG Xiao-qi(蔡晓龙, 刘永军, 董晓琪). Journal of Oil and Gas Technology(石油天然气学报), 2021, 43(2): 6. [本文引用:1]
[4] Khanna S, Santos M J, Koltunov A, et al. Remote Sensing, 2017, 9(2): 169. [本文引用:1]
[5] YOU Jin-feng, XING Li-xin, PAN Jun, et al(尤金凤, 邢立新, 潘军, 等). Journal of Jilin University: Earth Science Edition(吉林大学学报: 地球科学版), 2016, 46(5): 1589. [本文引用:1]
[6] Jiang J B, Steven M D, He R Y, et al. International Journal of Greenhouse Gas Control, 2015, 37: 1. [本文引用:1]
[7] LI Hui, JIANG Jin-bao, CHEN Xu-hui, et al(李辉, 蒋金豹, 陈绪慧). Spectroscopy and Spectral Analysis(光谱学与光谱分析), 2019, 39(12): 3743. [本文引用:4]
[8] van der Werff H M A, Bakker W H, van der Meer F D, et al. Computers & Geosciences, 2006, 32(5): 1334. [本文引用:2]
[9] Sanches I D, Filho C R S, Kokaly R F. ISPRS Journal of Photogrammetry and Remote Sensing, 2014, 97: 111. [本文引用:1]
[10] TONG Qing-xi, ZHANG Bing, ZHENG Lan-fen(童庆喜, 张兵, 郑兰芬). Hyperspectral Remote Sensing: Principle, Technology and Applications(高光谱遥感: 原理技术与应用). Beijing: Higher Education Press(北京: 高等教育出版社), 2006. [本文引用:3]
[11] CUI Xin, ZHAO Ying-jun, TIAN Feng, et al(崔鑫, 赵英俊, 田丰, 等). Acta Geologica Sinica(地质学报), 2019, 93(4): 928. [本文引用:2]
[12] Noomen Marleen F, Skidmore Andrew K, van der Meer Freek D, et al. Remote Sensing of Environment, 2006, 105(3): 262. [本文引用:2]
[13] van der Werff Harald, van der Meijde Mark, Jansma Fokke, et al. Sensors, 2008, 8(6): 3733. [本文引用:2]