基于长时序NDVI的土壤重金属胁迫甄别方法技术

本发明专利技术公开了一种基于长时序NDVI的土壤重金属胁迫甄别方法，包括步骤遥感影像预处理、NDVI长时间序列构建、EEMD分解、统计性描述指标计算、季节平均模型构建、土壤重金属胁迫序列甄别、土壤重金属胁迫稳定特征提取、地面数据实测、土壤重金属胁迫水平监测步骤。本发明专利技术利用NDVI长时间序列，将胁迫光谱特征与时间特征相结合，避免了单一时相缺乏时间特征的弊端；综合EEMD分解与季节平均模型构建，解决了噪声残留问题，更有利于捕捉土壤重金属胁迫的长期影响，实现对土壤重金属胁迫的甄别提取；综合作物物候期，在复杂的土壤环境中建立有效甄别土壤重金属胁迫的稳定特征指标，结合实测数据实现对大面积土壤重金属胁迫的监测。数据实现对大面积土壤重金属胁迫的监测。数据实现对大面积土壤重金属胁迫的监测。

【技术实现步骤摘要】
基于长时序NDVI的土壤重金属胁迫甄别方法


[0001]本专利技术涉及一种土壤重金属胁迫甄别方法，涉及一种基于长时序NDVI的土壤重金属胁迫甄别方法，属于遥感


技术介绍

[0002]传统的重金属监测方法是将实地采集的样本进行实验室化验，其结果精度高，但是速度慢、范围小、成本高，不能在空间分布上获取连续的信息。相较于传统方法，遥感技术具有大面积同步监测、传输速度快、信息丰富等特点，为其在重金属污染监测方面提供了借鉴。
[0003]利用遥感卫星影像数据对农田重金属污染大范围监测，逐渐成为重要的技术手段。近年来，国内外学者借助遥感技术，利用作物反射光谱响应特征以及光谱信息分析方面研究重金属污染状况。其中，光学遥感方法具备高空间分辨率的优势，光学植被指数数据与土壤重金属高度相关，是进行土壤重金属监测最常用的一种数据。归一化植被指数(Normalized Difference Vegetation Index，NDVI)定义为近红波段与红光波段的反射率差值除以近红外波段与红光波段的反射率和值，是反映农作物长势的重要参数之一。
[0004]自然农田生态系统中农作物生长环境的复杂性导致其在生长过程中受土壤水分胁迫、病虫害胁迫、重金属胁迫等多源复合胁迫因素的影响。在实现土壤重金属胁迫甄别时，需要对特定特征信号进行处理，从而提取出表征特定特征的分量。集合经验模态分解(Ensemble Empirical Mode Decomposition，EEMD)算法在经验模态分解(Empirical Mode Decomposition，EMD)算法的基础上进行了改进，能将原始长时间序列数据进行不断分解，获取不同本征模态函数(Intrinsic Mode Function，IMF)分量，最终分解为若干单一频率的序列以及一个残差组分的形式，为土壤重金属胁迫甄别提供技术支撑。
[0005]当前遥感监测作物重金属胁迫的主要局限在于采用的遥感指标的信息量匮乏、指标的稳定性不足，并且缺乏对农作物重金属胁迫光谱响应参数时间特征的系统分析。如果将胁迫光谱特征参数与作物生长阶段的时间特征相结合，将能够有效区分重金属胁迫和其他胁迫因子。因此，基于“时
‑
谱”分析是提取重金属胁迫特征的重要方法。

技术实现思路

[0006]本专利技术的目的是提供一种基于长时序NDVI的土壤重金属胁迫甄别方法。
[0007]为解决上述技术问题，本专利技术采用的技术方案是：
[0008]一种基于长时序NDVI的土壤重金属胁迫甄别方法，包括以下步骤：
[0009]步骤1：遥感影像预处理：对N幅遥感影像数据进行辐射定标、大气校正以及几何校正预处理，N＞1；
[0010]步骤2：NDVI长时间序列构建：基于预处理后的遥感影像数据计算NDVI长时间序列x(n)，1≤n≤N：
[0011][0012]其中ρ
NIR
(n)为第n幅遥感影像近红外波段的反射率，ρ
RED
(n)为第n幅遥感影像红光波段的反射率；
[0013]步骤3：EEMD分解：基于EEMD算法对NDVI长时间序列x(n)进行分解，其分解流程如下：
[0014]将具有标准正态分布的白噪声信号z
i
(n)加入到NDVI长时间序列x(n)上，产生加噪NDVI长时间序列x
i
(n)：
[0015]x
i
(n)＝x(n)+z
i
(n)，1≤i≤M
[0016]式中，M为加入标准正态分布的白噪声信号的次数；
[0017]对加噪NDVI长时间序列x
i
(n)分别进行EMD分解，得到各自IMF和的形式：
[0018][0019]式中，c
i，j
(n)为加噪NDVI长时间序列x
i
(n)的第j个IMF分量，r
i，j
(n)是加噪NDVI长时间序列的EMD残余函数，J是IMF分量的数量；
[0020]计算NDVI长时间序列x(n)的EEMD分解的第j个IMF分量：
[0021][0022]步骤4：计算统计性描述指标：计算NDVI长时间序列x(n)的EEMD分解的第一个至第J个IMF分量的极值点个数K
j
、波动周期P
j
：
[0023][0024]步骤5：季节平均模型构建：对NDVI长时间序列各年同月数据求取平均值，得到NDVI长时间序列的季节平均序列；对NDVI长时间序列的季节平均序列使用所述步骤3进行EEMD分解，得到NDVI长时间序列的季节平均序列的IMF分量imf
r
，计算NDVI长时间序列的季节平均序列的IMF分量imf
r
的周期，r＝1，2，...R；
[0025]步骤6：土壤重金属胁迫序列甄别：按照预设甄别条件选择符合长周期土壤重金属胁迫特征的NDVI长时间序列x(n)的EEMD分解的各IMF分量累加合成为土壤重金属胁迫序列Y
d
(n)；
[0026]步骤7：土壤重金属胁迫稳定特征提取：计算土壤重金属胁迫序列Y
d
的重金属胁迫稳定特征：
[0027][0028]式中：Y
dt
(n)为土壤重金属胁迫序列Y
d
(n)的前一时相的数值，Y
dt+1
(n)为土壤重金属胁迫序列Y
d
(n)的后一时相的数值，ΔDAY为前一时相对于后一时相对的日期间隔；
[0029]步骤8：地面数据实测：使用便携式XRF分析仪对土壤中的预设种类的重金属元素含量测定；
[0030]步骤9：构建重金属预测模型：拟合重金属胁迫稳定特征Y
df
(n)与土壤的预设种类的重金属元素含量间的关系，作为重金属预测模型；
[0031]步骤10：监测土壤重金属胁迫程度：采集遥感影像，进行所述步骤1
‑
步骤7的处理，得到对应的重金属胁迫稳定特征，使用所述重金属预测模型预测土壤重金属胁迫程度。
[0032]进一步，所述步骤6的预设甄别条件为NDVI长时间序列的IMF分量的波动周期P
r
为大于T，T为NDVI长时间序列的IMF分量与其对应的NDVI长时间序列的季节平均序列的IMF分量的共同波动周期。
[0033]进一步，所述步骤6的预设甄别条件为NDVI长时间序列的IMF分量的波动周期P
r
为大于6。
[0034]进一步，所述步骤6的预设甄别条件为NDVI长时间序列的IMF分量的波动周期P
r
大于12。
[0035]进一步，所述步骤9中使用二次方曲线拟合重金属胁迫稳定特征Y
df
(n)与土壤的预设种类的重金属元素含量间的关系。
[0036]采用上述技术方案，本专利技术的有益效果是：
[0037]1.本专利技术利用长时本文档来自技高网...

【技术保护点】

【技术特征摘要】
1.一种基于长时序NDVI的土壤重金属胁迫甄别方法，包括以下步骤：步骤1：遥感影像预处理：对N幅遥感影像数据进行辐射定标、大气校正以及几何校正预处理，N>1；步骤2：NDVI长时间序列构建：基于预处理后的遥感影像数据计算NDVI长时间序列x(n),1≤n≤N：其中ρ
NIR
(n)为第n幅遥感影像近红外波段的反射率，ρ
RED
(n)为第n幅遥感影像红光波段的反射率；步骤3：EEMD分解：基于EEMD算法对NDVI长时间序列x(n)进行分解，其分解流程如下：将具有标准正态分布的白噪声信号z
i
(n)加入到NDVI长时间序列x(n)上,产生加噪NDVI长时间序列x
i
(n)：x
i
(n)＝x(n)+z
i
(n)，1≤i≤M式中，M为加入标准正态分布的白噪声信号的次数；对加噪NDVI长时间序列x
i
(n)分别进行EMD分解，得到各自IMF和的形式:式中，c
i,j
(n)为加噪NDVI长时间序列x
i
(n)的第j个IMF分量，r
i,j
(n)是加噪NDVI长时间序列的EMD残余函数，J是IMF分量的数量；计算NDVI长时间序列x(n)的EEMD分解的第j个IMF分量：步骤4：计算统计性描述指标：计算NDVI长时间序列x(n)的EEMD分解的第一个至第J个IMF分量的极值点个数K
j
、波动周期P
j
：步骤5：季节平均模型构建：对NDVI长时间序列各年同月数据求取平均值，得到NDVI长时间序列的季节平均序列；对NDVI长时间序列的季节平均序列使用所述步骤3进行EEMD分解，得到NDVI长时间序列的季节平均序列的IMF分量imf
r
，计算NDVI长时间序列的季节平均序列的IMF分量imf
r
的周期，r＝1,2,
…
R；步骤6：土壤重金属胁迫序列甄别：按照预设甄别条件选...

【专利技术属性】
技术研发人员：李旭青，卢成乾，丁雪瑶，顾会涛，金永涛，何跃君，张文豪，吴艳萍，王延仓，刘天娇，吴琼，孙肖，曹启运，李凌飞，王玉博，王小丹，张文龙，赵辰雨，陈文博，吴冬雪，
申请(专利权)人：北华航天工业学院，
类型：发明
国别省市：