CN107192671B - 一种基于光谱曲面匹配的土壤类型识别方法 - Google Patents

一种基于光谱曲面匹配的土壤类型识别方法 Download PDF

Info

Publication number
CN107192671B
CN107192671B CN201710192544.6A CN201710192544A CN107192671B CN 107192671 B CN107192671 B CN 107192671B CN 201710192544 A CN201710192544 A CN 201710192544A CN 107192671 B CN107192671 B CN 107192671B
Authority
CN
China
Prior art keywords
soil
spectrum
curved surface
soil profile
profile
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN201710192544.6A
Other languages
English (en)
Other versions
CN107192671A (zh
Inventor
解宪丽
李安波
李德成
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nanjing Normal University
Institute of Soil Science of CAS
Original Assignee
Nanjing Normal University
Institute of Soil Science of CAS
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nanjing Normal University, Institute of Soil Science of CAS filed Critical Nanjing Normal University
Priority to CN201710192544.6A priority Critical patent/CN107192671B/zh
Publication of CN107192671A publication Critical patent/CN107192671A/zh
Application granted granted Critical
Publication of CN107192671B publication Critical patent/CN107192671B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00—Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/17—Systems in which incident light is modified in accordance with the properties of the material investigated
    • G01N21/25—Colour; Spectral properties, i.e. comparison of effect of material on the light at two or more different wavelengths or wavelength bands
    • G—PHYSICS
    • G01—MEASURING; TESTING
    • G01N—INVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N2201/00—Features of devices classified in G01N21/00
    • G01N2201/12—Circuits of general importance; Signal processing

Landscapes

  • Physics & Mathematics (AREA)
  • Spectroscopy & Molecular Physics (AREA)
  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Chemical & Material Sciences (AREA)
  • Analytical Chemistry (AREA)
  • Biochemistry (AREA)
  • General Health & Medical Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Immunology (AREA)
  • Pathology (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Abstract

本发明公开了一种基于光谱曲面匹配的土壤类型识别方法，首先，针对待识别土壤剖面各发生层的光谱反射率进行剖面深度插值处理和一阶求导处理，生成待识别土壤剖面的光谱曲面；然后，分别针对土壤样本库中的各个土壤剖面，针对土壤剖面各发生层的光谱反射率，进行剖面深度插值处理和一阶求导处理，生成各个土壤剖面的光谱曲面；最后进行光谱曲面匹配，识别出待识别土壤剖面的土壤类型。如此，基于土壤剖面的光谱反射率进行深度插值处理、求导处理、光谱曲面生成，以及光谱曲面相似性判别，整合利用土壤剖面中不同深度土壤光谱数据进行匹配，实现最终土壤类型的识别。

Description

一种基于光谱曲面匹配的土壤类型识别方法
技术领域
本发明涉及一种土壤类型识别方法，具体说是一种基于光谱曲面匹配的土壤类型识别方法，属于土壤地理和模式识别技术领域。
背景技术
土壤分类是依据土壤实体性状的综合差异，划分出不同级别的土壤类型，以便因地制宜的推广农业技术和改良土壤。对于土壤类型的识别需要获取大量的土体形态、物理、化学乃至生物等诊断信息，其中部分信息(如形态、容重等)可以通过肉眼观察或简单测量得到，但大多数物理、化学信息在传统上需要经过实验室测试分析方可获得。所以，对土壤类型的识别和鉴定通常成本较高，而且需要土壤分类专家的参与。
目前，为了满足快速发展的大量高精度土壤信息需求，快速获取和持续更新土壤信息是土壤资源研究领域的核心研究内容之一。传统方法由于其成本高、效率低，通常存在采样尺度偏大、采样密度偏稀疏、调查频率偏低等问题，难于满足对土壤信息进行动态、快速、低成本获取和更新的需求。
发明内容
本发明所要解决的技术问题是提出了一种基于光谱曲面匹配的土壤类型识别方法，弥补了现有技术空白，利用土壤剖面中不同深度土壤光谱数据进行匹配，实现土壤类型的识别。
本发明为了解决上述技术问题采用以下技术方案：本发明设计了一种基于光谱曲面匹配的土壤类型识别方法，包括如下步骤：
步骤一.针对待识别土壤剖面各发生层的光谱反射率，分别进行剖面深度插值处理和一阶求导处理，获得待识别土壤剖面的光谱曲面；
步骤二.分别针对土壤样本库中的各个土壤剖面，针对土壤剖面各发生层的光谱反射率，进行剖面深度插值处理和一阶求导处理，分别获得各个土壤剖面的光谱曲面；
步骤三.根据各个土壤剖面的光谱曲面，针对待识别土壤剖面的光谱曲面进行光谱曲面匹配，识别出待识别土壤剖面的土壤类型。
作为本发明的一种优选技术方案：所述步骤一中和步骤二中，针对土壤剖面各发生层的光谱反射率，按如下步骤，进行剖面深度插值处理和一阶求导处理，获得土壤剖面的光谱曲面；
步骤1.针对土壤剖面，获得各发生层的光谱反射率的集合P＝{pij|i＝1,…,n；j＝1,…,m}，其中，n为光谱波段数，m为土壤剖面所对应发生层的数量,pij表示土壤剖面第j个发生层在波段i处的反射率值；
步骤2.采用等积二次样条插值方法，针对土壤剖面各发生层的反射率值pij，进行深度插值处理，获得预设距离间隔、不同深度处的土壤光谱反射率P'＝{pa'ij'|i＝1,…,n；j'＝1,…,m'}，其中，m'表示插值后的土壤光谱曲线数量，pai'j'表示插值所生成的第j'条光谱曲线在波段i处的反射率值；
步骤3.针对经深度插值处理所获得的土壤光谱反射率，首先进行11点移动平均平滑处理，然后采用Savitzky-Golay卷积求导法计算21点二次多项式平滑的光谱反射率一阶导数；
步骤4.针对光谱反射率一阶导数，以波段为横坐标，以剖面深度为纵坐标，生成光谱曲面。
作为本发明的一种优选技术方案：所述步骤三中，根据土壤样本库中各个土壤剖面的光谱曲面，针对待识别土壤剖面的光谱曲面，基于平均Hausdorff距离进行光谱曲面匹配，识别出待识别土壤剖面的土壤类型。
作为本发明的一种优选技术方案，所述步骤三，包括如下步骤：
步骤3-1.分别针对土壤样本库中各个土壤剖面的光谱曲面，计算获得各土壤剖面光谱曲面分别与待识别土壤剖面光谱曲面之间的平均Hausdorff距离；
步骤3-2.获得最小平均Hausdorff距离所对应的土壤剖面的土壤类型，即为待识别土壤剖面的土壤类型。
作为本发明的一种优选技术方案，所述步骤3-1中，分别针对土壤样本库中各个土壤剖面的光谱曲面，按如下步骤，计算获得各土壤剖面光谱曲面与待识别土壤剖面光谱曲面之间的平均Hausdorff距离:
步骤3-1-1.初始化点集A、点集B分别为长度为n的空集，以及初始化参数l＝1，n表示光谱波段数；
步骤3-1-2.从待识别土壤剖面Ct光谱曲面中抽取深度值l处的光谱线加入到点集A当中；同时，从土壤样本库中土壤剖面光谱曲面CSg中抽取深度值l处的光谱线加入到点集B当中；其中，g＝{1、…、G}，G表示土壤样本库中土壤剖面光谱曲面的个数；
步骤3-1-3.根据如下公式：
计算获得点集A和点集B之间的平均Hausdorff距离其中，dB(a)表示点集A中元素a到点集B距离的最小值，p为点集A中元素的个数；
步骤3-1-4.判断l是否大于光谱曲面Ct或CSg的最大深度值，是则进入步骤3-1-5；否则，l所对应的值进行加一更新，并返回步骤3-1-2；
步骤3-1-5.针对待识别剖面光谱曲面Ct中各条光谱线分别与土壤样本库土壤剖面光谱曲面CSg中各条光谱线之间的平均Hausdorff距离，采用求平均值方法，获得待识别剖面光谱曲面Ct与土壤样本库中土壤剖面光谱曲面CSg间的平均Hausdorff距离。
本发明所述一种基于光谱曲面匹配的土壤类型识别方法采用以上技术方案与现有技术相比，具有以下技术效果：本发明所设计基于光谱曲面匹配的土壤类型识别方法，针对每一土壤剖面，通过对不同深度采样点的可见-近红外反射光谱数据，进行剖面深度插值处理和求导处理后，生成三维曲面，进而通过基于平均Hausdorff距离的光谱曲面相似度进行判别，达到土壤类型自动识别的目的，综合利用土壤剖面中不同深度土壤光谱数据进行匹配，进而实现土壤类型的识别。
附图说明
图1是本发明所设计基于光谱曲面匹配的土壤类型识别方法的流程图；
图2是实施例待识别土壤剖面的光谱曲面示意图；
图3是实施例样本剖面的光谱曲面示意图。
具体实施方式
下面结合说明书附图对本发明的具体实施方式作进一步详细的说明。
可见-近红外漫反射光谱技术是快速测定土壤各类理化参数使用最频繁的光谱区域，可见-近红外漫反射光谱不仅可以定量分析具有明确光谱吸收特征的水分、粘粒、铁氧化物、有机碳和全氮等基本土壤参数，还可以较好地预测那些不具有明确吸收特征、与基本土壤参数存在相关关系的其它土壤参数(如阳离子交换量)。
土壤的基本物质组成变化既是土壤发生发育过程的反映，是确定土壤诊断指标和类型划分的依据，同时也是影响土壤反射光谱的主要因素。在土壤调查与制图中，利用光谱分析技术实现土壤类型的快速、准确、自动识别在理论上是可行的。实现上述应用的关键在于土壤光谱分析模型的建立和土壤光谱先验知识的积累。为此，需要广泛采集土壤样品，在相对标准的测定条件下获取土壤的光谱数据，建立完善的土壤光谱库，在此基础上，研发有效的模式识别方法以达到土壤光谱识别和分类的目标。
如图1所示，本发明设计了一种基于光谱曲面匹配的土壤类型识别方法，在实际应用过程当中，针对土壤剖面各发生层，通过可见-近红外漫反射光谱技术去获得土壤剖面各发生层的光谱反射率，土壤类型识别方法具体包括如下步骤：
步骤一.针对待识别土壤剖面各发生层的光谱反射率，分别进行剖面深度插值处理和一阶求导处理，获得待识别土壤剖面的光谱曲面。
步骤二.分别针对土壤样本库中的各个土壤剖面，针对土壤剖面各发生层的光谱反射率，进行剖面深度插值处理和一阶求导处理，分别获得各个土壤剖面的光谱曲面。
实际应用中，上述步骤一中和步骤二中，针对土壤剖面各发生层的光谱反射率，具体按如下步骤，进行剖面深度插值处理和一阶求导处理，获得土壤剖面的光谱曲面；
步骤1.针对土壤剖面，获得各发生层的光谱反射率的集合P＝{pij|i＝1,…,n；j＝1,…,m}，其中，n为光谱波段数，m为土壤剖面所对应发生层的数量,pij表示土壤剖面第j个发生层在波段i处的反射率值。
步骤2.采用等积二次样条插值方法，针对土壤剖面各发生层的反射率值pij，进行深度插值处理，获得预设距离间隔、不同深度处的土壤光谱反射率P'＝{pa'ij'|i＝1,…,n；j'＝1,…,m'}，其中，m'表示插值后的土壤光谱曲线数量，pa'ij'表示插值所生成的第j'条光谱曲线在波段i处的反射率值。
所述等积二次样条插值方法，公开在下列文件中：1.Bishop,T.F.A.,McBratney,A.B.,Laslett,G.M.,1999.Modelling soil attribute depth functions with equal-area quadratic smoothing splines.Geoderma 91,27–45.；2.Malone,B.P.,McBratney,A.B.,Minasny,B.,Laslett,G.M.,2009.Mapping continuous depthfunctions of soilcarbon storage and available water capacity.Geoderma 154,138–152.
步骤3.针对经深度插值处理所获得的土壤光谱反射率，首先进行11点移动平均平滑处理，然后采用Savitzky-Golay卷积求导法计算21点二次多项式平滑的光谱反射率一阶导数。
所述Savitzky-Golay卷积求导法，公开在下列文件中：Savitzky,A.,Golay,M.J.E.,1964.Smoothing and Differentiation of Data by Simplified Least SquaresProcedures.Analytical Chemistry 36(8),1627–1639.
步骤4.针对光谱反射率一阶导数，以波段为横坐标，以剖面深度为纵坐标，生成光谱曲面。
步骤三.根据土壤样本库中各个土壤剖面的光谱曲面，针对待识别土壤剖面的光谱曲面，基于平均Hausdorff距离进行光谱曲面匹配，识别出待识别土壤剖面的土壤类型。
实际应用中，上述步骤三，具体包括如下步骤：
步骤3-1.分别针对土壤样本库中各个土壤剖面的光谱曲面，按如下步骤，计算获得各土壤剖面光谱曲面与待识别土壤剖面光谱曲面之间的平均Hausdorff距离:
步骤3-1-1.初始化点集A、点集B分别为长度为n的空集，以及初始化参数l＝1，n表示光谱波段数。
步骤3-1-2.从待识别土壤剖面Ct光谱曲面中抽取深度值l处的光谱线加入到点集A当中；同时，从土壤样本库中土壤剖面光谱曲面CSg中抽取深度值l处的光谱线加入到点集B当中；其中，g＝{1、…、G}，G表示土壤样本库中土壤剖面光谱曲面的个数。
步骤3-1-3.根据如下公式：
计算获得点集A和点集B之间的平均Hausdorff距离其中，dB(a)表示点集A中元素a到点集B距离的最小值，p为点集A中元素的个数。
步骤3-1-4.判断l是否大于光谱曲面Ct或CSg的最大深度值，是则进入步骤3-1-5；否则，l所对应的值进行加一更新，并返回步骤3-1-2。
步骤3-1-5.针对待识别剖面光谱曲面Ct中各条光谱线分别与土壤样本库土壤剖面光谱曲面CSg中各条光谱线之间的平均Hausdorff距离，采用求平均值方法，获得待识别剖面光谱曲面Ct与土壤样本库中土壤剖面光谱曲面CSg间的平均Hausdorff距离。
步骤3-2.获得最小平均Hausdorff距离所对应的土壤剖面的土壤类型，即为待识别土壤剖面的土壤类型。
上述技术方案所设计基于光谱曲面匹配的土壤类型识别方法，针对每一土壤剖面，通过对不同深度采样点的可见-近红外反射光谱数据，进行剖面深度插值处理和求导处理后，生成三维曲面，进而通过基于平均Hausdorff距离的光谱曲面相似度进行判别，达到土壤类型自动识别的目的，综合利用土壤剖面中不同深度土壤光谱数据进行匹配，进而实现土壤类型的识别。
将上述所设计基于光谱曲面匹配的土壤类型识别方法，结合附图并通过描述一个土壤类型匹配识别实例，来进一步说明本发明的实际应用效果。
本实施例使用的土壤样本为中国某省土系调查的存档样本，采自104个剖面的474个发生层。这些剖面按照中国土壤系统分类划分为5个土纲、6个亚纲、16个土类及29个亚类。采用Cary 5000分光光度计(Agilent Technologies)测量样本的350-2500nm漫反射光谱，光谱采样间隔为1nm。由于350-400nm、2450-2500nm的测量值噪声较大，本实施例中仅选用400-2450nm的光谱进行土壤类型的匹配识别。
基于图1所示，具体实施过程如下：
(一)生成待识别土壤剖面的光谱曲面。待识别土壤剖面TestProfile的光谱曲面具体生成步骤，示例如下：
步骤一.针对一土壤剖面TestProfile，读取相应土壤样品的光谱反射率到集合P中，其中，光谱波段数n为2051(400nm-2450nm)，光谱曲线数m为6；其部分土壤剖面光谱反射率数据内容，如下表1所示。
Layer UpperDeptr LowerDenth 400nm 401nm 402nm 403nm 404nm 405nm 406nm 407nm 408nm 409nm 410nm
1 0 35 13.41364 13.49586 13.60696 13.67073 13.78241 13.86285 13.94667 14.01074 14.10604 14.15471 14.2464
2 35 50 14.50707 14.5954 14.67242 14.78604 14.85784 14.95909 15.03139 15.13364 15.21278 15.26049 15.34988
3 50 60 14.2031 14.28368 14.34982 14.43736 14.51813 14.5899 14.71283 14.78467 14.86383 14.91999 14.97972
4 60 90 14.39339 14.49908 14.58461 14.63667 14.71747 14.82547 14.90479 14.99014 15.06256 15.14267 15.24865
5 90 120 14.78277 14.87434 14.95373 15.00892 15.11342 15.1949 15.31476 15.36927 15.50681 15.55252 15.63557
6 120 140 13.81063 13.90401 13.99262 14.04254 14.13786 14.23267 14.27512 14.3706 14.49029 14.51679 14.58542
表1
步骤二.采用等积二次样条插值方法，对集合P中每一波段的土壤光谱反射率进行深度插值处理，得到0-140cm(1cm间隔)不同土壤深度处的土壤光谱反射率数据，插值后的光谱曲线数为140，其部分土壤剖面光谱反射率深度插值结果数据，如下表2所示。
UpperDepth LowerDepth 400nm 401nm 402nm 403nm 404nm 405nm 406nm 407nm 408nm 409nm 410nm
0 1 13.12648 13.20665 13.32772 13.37531 13.50007 13.57264 13.66358 13.71452 13.8152 13.86467 13.9554
1 2 13.12743 13.20761 13.32865 13.37628 13.50101 13.5736 13.66453 13.7155 13.81616 13.86563 13.95637
2 3 13.13027 13.21047 13.33144 13.37921 13.50383 13.57648 13.66736 13.71844 13.81906 13.86853 13.95927
3 4 13.13501 13.21525 13.33609 13.38407 13.50853 13.58129 13.67207 13.72334 13.8239 13.87335 13.96411
4 5 13.14165 13.22193 13.34259 13.39089 13.5151 13.58803 13.67868 13.7302 13.83066 13.8801 13.97089
5 6 13.15019 13.23053 13.35095 13.39965 13.52356 13.59669 13.68717 13.73902 13.83936 13.88879 13.97961
6 7 13.16062 13.24103 13.36117 13.41037 13.53389 13.60727 13.69754 13.7498 13.84999 13.8994 13.99026
7 8 13.17294 13.25344 13.37325 13.42302 13.54611 13.61977 13.70981 13.76254 13.86256 13.91194 14.00285
8 9 13.18717 13.26777 13.38718 13.43763 13.5602 13.6342 13.72396 13.77724 13.87705 13.92641 14.01737
9 10 13.20329 13.284 13.40297 13.45418 13.57617 13.65055 13.74 13.7939 13.89348 13.94281 14.03383
10 11 13.2213 13.30214 13.42062 13.47268 13.59402 13.66883 13.75792 13.81252 13.91185 13.96114 14.05223
11 12 13.24122 13.3222 13.44013 13.49313 13.61375 13.68903 13.77773 13.8331 13.93214 139814 14.07257
12 13 13.26303 13.34416 13.4615 13.51553 13.63535 1371116 13.79943 13.85564 13.95437 14.00358 14.09484
13 14 13.28673 13.36803 13.48472 13.53987 13.65884 13.73521 13.82302 13.88014 13.97853 14.0277 14.11905
14 15 13.31233 13.39381 13.50981 13.56616 13.6842 13.76118 13.84849 13.90661 14.00463 14.05375 14.1452
15 16 13.33983 13.42151 13.53675 13.5944 13.71145 13.78908 13.87585 13.93503 14.03266 14.08173 14.17328
16 17 13.36923 13.45111 13.56555 13.62459 13.74057 13.8189 13.9051 13.96541 14.06262 14.11163 14.2033
17 18 13.40052 13.48262 13.5962 13.65672 13.77157 13.85064 13.93623 13.99775 14.09451 14.14347 14.23525
18 19 13.43371 13.51604 13.62872 13.6908 13.80445 13.88431 13.96925 14.03205 14.12834 14.17723 14.26915
19 20 13.46879 13.55137 13.66309 13.72683 13.83921 13.9199 14.00415 14.06831 14.1641 14.21292 14.30498
20 21 1350578 13.58861 13.69932 13.76481 13.87585 13.95742 14.04095 14.10654 14.2018 14.25055 14.34274
21 22 13.54465 13.62776 13.73741 13.80473 13.91437 13.99686 14.07963 14.14672 14.24142 14.2901 14.38244
22 23 13.58543 13.66882 13.77736 13.8466 13.95477 14.03822 14.12019 14.18886 14.28298 14.33158 14.42408
23 24 13.6281 13.7118 13.81916 13.89042 13.99704 14.08151 14.16265 14.23296 14.32647 14.37499 14.46766
24 25 13.67266 13.75668 13.86283 13.93619 14.0412 14.12672 14.20699 14.27902 14.3719 14.42034 14.51317
25 26 13.71913 13.80347 13.90835 13.9839 14.08723 14.17385 14.25321 14.32705 14.41926 14.46761 14.56062
26 27 13.76749 13.85217 13.95573 14.03356 14.13514 14.22291 14.3.133 14.37703 14.46855 14.51681 14.61001
27 28 13.81775 13.90278 14.00496 14.08517 14.18493 14.2739 14.35133 14.42897 14.51977 14.55793 14.66133
28 29 13.8699 13.9553 14.05606 14.13873 14.2366 14.3268 14.40322 14.48287 14.57293 14.62099 14.71459
29 30 13.92395 14.00972 14.10901 14.19423 14.29015 14.38163 14.45699 14.53874 14.62802 14.67598 14.76979
30 31 13.97989 14.06606 14.16382 14.25168 14.34558 14.43839 14.51265 14.59656 14.68504 14.7329 14.82692
31 32 14.03774 14.12431 14.22049 14.31108 14.40288 14.49707 14.5702 14.65634 14.744 14.79174 14.88599
32 33 14.09748 14.18447 14.27902 14.37242 14.46207 14.55767 14.62963 14.71809 14.50489 14.85252 14.947
33 34 14.15911 14.24654 14.3394 14.43572 14.52313 14.6202 14.69096 14.78179 14.86771 14.91523 15.00994
34 35 14.22264 14.31052 14.40164 14.50096 14.58608 14.68465 14.75416 14.84745 14.93247 14.97986 15.07482
35 36 14.28807 14.37641 14.46575 14.56815 14.6509 14.75102 14.81926 14.91508 14.99916 15.04643 15.14164
36 37 14.3507 14.43945 14.527 14.63233 14.71281 14.81433 14.88151 14.97968 15.06285 15.11002 15.20536
37 38 1440583 14.49488 14.58071 14.68854 14.76701 14.86956 14.9362 15.03627 15.11863 15.16576 15.26097
38 39 14.45346 14.54271 14.62686 14.73679 14.8135 14.91674 14.98331 15.08486 15.1665 15.21363 15.30845
39 40 14.49359 14.58293 14.66546 14.77707 14.85229 14.95584 15.02285 15.12545 15.20645 15.25364 15.34781
表2
光谱反射率的一阶求导处理。首先把光谱反射率进行11点移动平均平滑处理，其次采用Savitzky-Golay卷积求导法计算21点二次多项式平滑的光谱一阶导数，求导处理后部分土壤剖面光谱反射率求导结果数据，如下表3所示。
UpperDepth LowerDepth 410nm 411nm 412nm 413nm 414nm 415nm 416nm 417nm 418nm 419nm 420nm
0 1 0.08 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09
1 2 0.08 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09
2 3 0.08 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09
3 4 0.08 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09
4 5 0.08 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09
5 6 0.08 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09
6 7 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
7 8 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
8 9 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
9 10 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
10 11 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
11 12 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
12 13 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
13 14 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
14 15 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
15 16 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
16 17 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09
17 18 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10
18 19 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10
19 20 0.08 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10
20 21 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.10
21 22 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.10
22 23 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.10
23 24 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.09 0.10
24 25 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10 0.10
25 26 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10 0.10
26 27 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10 0.10
27 28 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10 0.10
28 29 0.08 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10 0.10
29 30 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10
30 31 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10
31 32 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10
32 33 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10
33 34 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10 0.10
34 35 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10 0.10
35 36 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10 0.10
36 37 0.08 0.08 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10 0.10
37 38 0.08 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10 0.10 0.10
38 39 0.08 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10 0.10 0.10
39 40 0.08 0.09 0.09 0.09 0.09 0.09 0.10 0.10 0.10 0.10 0.10
表3
步骤三.根据深度插值处理和一阶求导处理后的光谱反射率数据，以波段为横坐标，以剖面深度为纵坐标，生成光谱曲面Ct。本实施例中，土壤剖面TestProfile的光谱曲面如图2所示。
(二)生成土壤样本库中各个土壤剖面的光谱曲面，具体步骤与步骤(一)一致，在此不再重复叙述；
(三)基于Hausdorff距离的光谱曲面匹配，具体步骤包括：
步骤3-1.分别针对土壤样本库中各个土壤剖面的光谱曲面，按如下步骤，计算获得各土壤剖面光谱曲面与待识别土壤剖面光谱曲面之间的平均Hausdorff距离：
步骤3-1-1.初始化点集A、点集B分别为长度为n的空集，以及初始化参数l＝1，n表示光谱波段数。
步骤3-1-2.从待识别土壤剖面Ct光谱曲面中抽取深度值l处的光谱线加入到点集A当中；同时，从土壤样本库中土壤剖面光谱曲面CSg中抽取深度值l处的光谱线加入到点集B当中；其中，g＝{1、...、G}，G表示土壤样本库中土壤剖面光谱曲面的个数。本实例中土壤剖面光谱曲面为DY-014，记为CS1，如图3所示，并且从光谱曲面Ct、CS1中各取一光谱曲线，记为la,lb,并分别将其读入集合A、B中。本实施例中，集合A为{0.08，0.08，0.08，0.08，0.08，0.09，0.09，0.09，0.09，0.09，…,-0.04,-0.04,-0.04,-0.04,-0.04,-0.04,-0.04,-0.04,-0.03,-0.03}，集合B为{0.08，0.09，0.09，0.09，0.09，0.10，0.10，0.10，0.10，0,10，…,-0.05,-0.04,-0.04,-0.04,-0.04,-0.04,-0.04,-0.04,-0.04,-0.04}。
步骤3-1-3.根据如下公式：
计算获得点集A和点集B之间的平均Hausdorff距离其中，dB(a)表示点集A中元素a到点集B距离的最小值，p为点集A中元素的个数，本实施例中，值为0。
步骤3-1-4.判断l是否大于光谱曲面Ct或CSg的最大深度值，是则进入步骤3-1-5；否则，l所对应的值进行加一更新，并返回步骤3-1-2，实际中，由于光谱曲面Ct、CS1的深度值分别为140cm、120cm，因此，取120为循环结束条件。
步骤3-1-5.针对待识别剖面光谱曲面Ct中各条光谱线分别与土壤样本库土壤剖面光谱曲面CSg中各条光谱线之间的平均Hausdorff距离，采用求平均值方法，获得待识别剖面光谱曲面Ct与土壤样本库中土壤剖面光谱曲面CSg间的平均Hausdorff距离。在本实施例中，所获光谱曲面Ct、CS1之间的平均Hausdorff距离为0.00897810465074511。
如此按上述步骤，计算出待识别剖面光谱曲面Ct与土壤样本库中每一土壤剖面光谱曲面CSg的平均Hausdorff距离，并将两光谱曲面相应的剖面编号、平均Hausdorff距离值记入集合S中。本实施例中，集合S的部分数据内容如表4所示。
待识别剖面 光谱库剖面 平均Hansdorff距离
TestProfile 34-002 0.005878722
TestProfile 34-006 0.00575133
TestProfile 34-013 0.012536485
TestProfile 34-014 0.002078553
TestProfile 34-015 0.007466477
TestProfile 34-016 0.003929145
TestProfile 34-017 0.013489106
TestProfile 34-023 0.004996959
TestProfile 34-024 0.013477257
TestProfile 34-025 0.003082999
TestProfile 34-026 0.011447511
TestProfile 34-027 0.004113817
TestProfile 34-028 0.005078846
TestProfile 34-030 0.016959213
TestProfile 34-032 0.014600401
TestProfile 34-033 0.011552242
TestProfile 34-035 0.00949656
TestProfile 34-036 0.010540808
TestProfile 34-038 0.012839614
TestProfile 34-039 0.01281967
TestProfile 34-040 0.016568236
TestProfile 34-042 0.010634925
TestProfile 34-043 0.008930264
TestProfile 34-044 0.007563428
TestProfile 34-045 0.010989086
TestProfile 34-054 0.008824906
TestProfile 34-055 0.009318296
TestProfile 34-056 0.010143285
TestProfile 34-057 0.009469448
TestProfile 34-059 0.015201773
TestProfile 34-061 0.010058566
表4
步骤3-2.遍历集合S，找出平均Hausdorff距离的最小值和对应的剖面编号，即获得最小平均Hausdorff距离所对应的土壤剖面的土壤类型，即为待识别土壤剖面的土壤类型，本实施例中，平均Hausdorff距离的最小值为0.00207855266793084，其对应土壤剖面编号为34-014。光谱中34-014的土类为淡色潮湿雏形土。则待测试土壤剖面TestProfile的土类为淡色潮湿雏形土，与实际土壤类型一致。
上述应用实施例中，仅进行了土壤土类的识别。该方法同样可适用于土纲、亚纲、亚类等不同级别土壤类型的识别；并且本实施例中，仅基于可见-近红外漫反射光谱数据进行土壤类型的识别，该方法同样可适用于中红外漫反射光谱等其它类型的土壤光谱数据。
上面结合附图对本发明的实施方式作了详细说明，但是本发明并不限于上述实施方式，在本领域普通技术人员所具备的知识范围内，还可以在不脱离本发明宗旨的前提下做出各种变化。

Claims (2)

1.一种基于光谱曲面匹配的土壤类型识别方法，其特征在于，包括如下步骤：
步骤一.针对待识别土壤剖面各发生层的光谱反射率，分别进行剖面深度插值处理和一阶求导处理，获得待识别土壤剖面的光谱曲面，然后进入步骤二；
步骤二.分别针对土壤样本库中的各个土壤剖面，针对土壤剖面各发生层的光谱反射率，进行剖面深度插值处理和一阶求导处理，分别获得各个土壤剖面的光谱曲面，然后进入步骤三；
其中，上述步骤一中和步骤二中，针对土壤剖面各发生层的光谱反射率，按如下步骤1至步骤4，进行剖面深度插值处理和一阶求导处理，获得土壤剖面的光谱曲面；
步骤1.针对土壤剖面，获得各发生层的光谱反射率的集合P＝{pij|i＝1,…,n；j＝1,…,m}，其中，n为光谱波段数，m为土壤剖面所对应发生层的数量,pij表示土壤剖面第j个发生层在波段i处的反射率值；
步骤2.采用等积二次样条插值方法，针对土壤剖面各发生层的反射率值pij，进行深度插值处理，获得预设距离间隔、不同深度处的土壤光谱反射率P'＝{pa′ij'|i＝1,…,n；j'＝1,…,m'}，其中，m'表示插值后的土壤光谱曲线数量，pa′ij'表示插值所生成的第j'条光谱曲线在波段i处的反射率值；
步骤3.针对经深度插值处理所获得的土壤光谱反射率，首先进行11点移动平均平滑处理，然后采用Savitzky-Golay卷积求导法计算21点二次多项式平滑的光谱反射率一阶导数；
步骤4.针对光谱反射率一阶导数，以波段为横坐标，以剖面深度为纵坐标，生成光谱曲面；
步骤三.按如下步骤3-1至步骤3-2，根据土壤样本库中各个土壤剖面的光谱曲面，针对待识别土壤剖面的光谱曲面，基于平均Hausdorff距离进行光谱曲面匹配，识别出待识别土壤剖面的土壤类型；
步骤3-1.分别针对土壤样本库中各个土壤剖面的光谱曲面，计算获得各土壤剖面光谱曲面分别与待识别土壤剖面光谱曲面之间的平均Hausdorff距离；
步骤3-2.获得最小平均Hausdorff距离所对应的土壤剖面的土壤类型，即为待识别土壤剖面的土壤类型。
2.根据权利要求1所述一种基于光谱曲面匹配的土壤类型识别方法，其特征在于，所述步骤3-1中，分别针对土壤样本库中各个土壤剖面的光谱曲面，按如下步骤，计算获得各土壤剖面光谱曲面与待识别土壤剖面光谱曲面之间的平均Hausdorff距离:
步骤3-1-1.初始化点集A、点集B分别为长度为n的空集，以及初始化参数l＝1，n表示光谱波段数；
步骤3-1-2.从待识别土壤剖面Ct光谱曲面中抽取深度值l处的光谱线加入到点集A当中；同时，从土壤样本库中土壤剖面光谱曲面CSg中抽取深度值l处的光谱线加入到点集B当中；其中，g＝{1、…、G}，G表示土壤样本库中土壤剖面光谱曲面的个数；
步骤3-1-3.根据如下公式：
计算获得点集A和点集B之间的平均Hausdorff距离其中，dB(a)表示点集A中元素a到点集B距离的最小值，p为点集A中元素的个数；
步骤3-1-4.判断l是否大于光谱曲面Ct或CSg的最大深度值，是则进入步骤3-1-5；否则，l所对应的值进行加一更新，并返回步骤3-1-2；
步骤3-1-5.针对待识别剖面光谱曲面Ct中各条光谱线分别与土壤样本库土壤剖面光谱曲面CSg中各条光谱线之间的平均Hausdorff距离，采用求平均值方法，获得待识别剖面光谱曲面Ct与土壤样本库中土壤剖面光谱曲面CSg间的平均Hausdorff距离。
CN201710192544.6A 2017-03-28 2017-03-28 一种基于光谱曲面匹配的土壤类型识别方法 Active CN107192671B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710192544.6A CN107192671B (zh) 2017-03-28 2017-03-28 一种基于光谱曲面匹配的土壤类型识别方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710192544.6A CN107192671B (zh) 2017-03-28 2017-03-28 一种基于光谱曲面匹配的土壤类型识别方法

Publications (2)

Publication Number Publication Date
CN107192671A CN107192671A (zh) 2017-09-22
CN107192671B true CN107192671B (zh) 2018-09-18

Family

ID=59871787

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710192544.6A Active CN107192671B (zh) 2017-03-28 2017-03-28 一种基于光谱曲面匹配的土壤类型识别方法

Country Status (1)

Country Link
CN (1) CN107192671B (zh)

Families Citing this family (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109211814B (zh) * 2018-10-29 2019-08-06 中国科学院南京土壤研究所 一种基于三维光谱曲面分区特征的土壤剖面类型识别方法
CN114216866B (zh) * 2021-11-05 2022-12-06 中国科学院南京土壤研究所 一种基于高光谱成像技术的土壤剖面颜色获取方法
CN114324216B (zh) * 2022-01-06 2023-08-01 中国科学院南京土壤研究所 一种基于土层组合特征的土壤数值分类方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101210875A (zh) * 2007-12-25 2008-07-02 浙江大学 基于近红外光谱技术的无损测量土壤养分含量的方法
CN103234922A (zh) * 2013-03-29 2013-08-07 浙江大学 一种基于大样本土壤可见-近红外光谱分类的土壤有机质快速检测方法
CN203132531U (zh) * 2013-03-06 2013-08-14 四川农业大学 一种土层详细划分与土壤剖面摄影的土壤标尺
CN105548031A (zh) * 2015-12-18 2016-05-04 北京农业智能装备技术研究中心 基于移动终端的土壤类型识别装置

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101210875A (zh) * 2007-12-25 2008-07-02 浙江大学 基于近红外光谱技术的无损测量土壤养分含量的方法
CN203132531U (zh) * 2013-03-06 2013-08-14 四川农业大学 一种土层详细划分与土壤剖面摄影的土壤标尺
CN103234922A (zh) * 2013-03-29 2013-08-07 浙江大学 一种基于大样本土壤可见-近红外光谱分类的土壤有机质快速检测方法
CN105548031A (zh) * 2015-12-18 2016-05-04 北京农业智能装备技术研究中心 基于移动终端的土壤类型识别装置

Non-Patent Citations (7)

* Cited by examiner, † Cited by third party
Title
A global spectral library to characterize the world"s soil;R.A. Viscarra Rossel 等;《Earth-Science Reviews》;20160206;全文 *
Empirical estimates of uncertainty for mapping continuous depth functions of soil attributes;B.P.Malone 等;《Geoderma》;20101218;全文 *
G.M. Vasques 等.Soil classification using visible/near-infrared diffuse reflectance spectra from multiple depths.《Geoderma》.2014,全文. *
三维自由表面匹配及其应用;冯义从 等;《测绘工程》;20050930;第14卷(第3期);全文 *
典型地物实测光谱的相似性测度与实验分析;施蓓琦 等;《同济大学学报（自然科学版）》;20110215;第39卷(第2期);全文 *
基于Savizky-Golay多项式的三维荧光光谱的曲面平滑方法;杜树新 等;《光谱学与光谱分析》;20110215;第31卷(第2期);全文 *
基于高光谱非线性向量空间的光谱曲线特征差异性分析;戴晓爱 等;《遥感技术与应用》;20151215;第30卷(第6期);全文 *

Also Published As

Publication number Publication date
CN107192671A (zh) 2017-09-22

Similar Documents

Publication Publication Date Title
CN107192671B (zh) 一种基于光谱曲面匹配的土壤类型识别方法
CN103854305B (zh) 一种基于多尺度建模的模型传递方法
US11231519B2 (en) Method and device for simulating discharge, and computer device
CN109344550A (zh) 一种基于国产高分卫星遥感数据的森林蓄积量反演方法及系统
CN107729651B (zh) 基于多维度的坡面细沟发育形态特征综合量化方法
CN103020935A (zh) 一种自适应在线字典学习的图像超分辨率方法
Ye et al. Analysing the potential of UAV point cloud as input in quantitative structure modelling for assessment of woody biomass of single trees
CN113420412B (zh) 基于成像光谱的土壤有机碳含量连续深度分布提取方法
CN113436153B (zh) 一种基于高光谱成像和支持向量机技术的原状土壤剖面碳组分预测方法
CN102798607A (zh) 一种利用中红外光谱技术估测土壤有机碳含量的方法
CN206741554U (zh) 基于深度摄像头的室内3d扫描设备的户型3d建模系统
CN117576583A (zh) 一种基于无人机的冬小麦生物量快速高精准估测方法
Madonna et al. Verification of North Atlantic warm conveyor belt outflows in ECMWF forecasts
CN106771071B (zh) 一种基于油水相渗的密闭取心饱和度校正方法
CN113076645B (zh) 一种结合遥感数据的宇宙射线中子仪空间特征模拟方法
CN107144222B (zh) 标准放大图测量方法、测量装置及测量系统
CN113887064A (zh) 大尺度地下水储量遥感动态监测及驱动因子定量拆分方法
CN104182953B (zh) 像元解混逆过程：规格化多端元分解的高光谱重构方法
CN113570538B (zh) 一种叶片rgb图像偏态分布参数信息采集及分析方法
CN116403048B (zh) 一种基于多模态数据融合的农作物生长估计模型构建方法
CN115238327B (zh) 基于卫星监测的三维空间甲烷分布格局模拟方法
CN112613152B (zh) 一种套管式地埋管岩土热物性参数估算方法
CN114189313B (zh) 一种电表数据重构方法及装置
CN111768101B (zh) 一种顾及物候特征的遥感耕地变化检测方法及系统
CN110096743B (zh) 基于遥感数据和高程信息估算地表水汽压方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant