CN109635432A - 地块三维形貌评价方法及实现所述方法的系统 - Google Patents

Abstract

本发明提供一种地块三维形貌评价方法，所述方法包括对地块内各点进行采样生成三维形貌数据；去除三维形貌数据中的趋势项并进行平滑处理；利用小波分析进行地表三维形貌分析，计算粗糙度评定基准面与表面粗糙度、计算地表三维参数，通过三维表面轮廓算术平均偏差、三维表面轮廓均方根偏差、三维平面轮廓的偏斜度和三维表面高度分布的陡度判断地表三维形貌的粗糙度。本发明还提供实现所述方法的系统。本发明的方法与系统能够可以直观的看到土壤地表三维形貌图以及评价参数值，通过分析可以直接得到测量地块三维表面粗糙度的状态，结果简单易懂。

Description

地块三维形貌评价方法及实现所述方法的系统

【技术领域】

本发明涉及农业信息技术及分析领域，具体涉及一种对耕作土地进 行地块三维形貌评价方法，还涉及实现所述方法的系统。

【背景技术】

物体表面存在的连续凹凸不平的微几何形状统称为物体表面形貌， 因此就耕作土壤地表而言，其三维形貌主要是在不同的农业机械作业后 产生的土壤与空气之间的分界面。

耕作土壤地表三维形貌不仅在农业方面有着重要影响作用，而且对 于土壤含水量的预测，土壤风蚀预测、预防，天气和气候的预测等众多 领域的研究也都起着非常重要的影响作用：

(1)耕作土壤地表三维形貌直接影响土壤的灌溉效率和灌溉均匀 度，状态良好的耕作土壤地表三维形貌可以提高水资源、人力资源利用 率。

(2)耕作土壤地表三维形貌会导致农业机械作业时产生振动，从 而影响机具生产性能和效率，而且长时间振动会导致机身零部件过早地 出现磨损、疲劳和损坏现象，不仅会缩短农业机械使用寿命，还会导致 农业机械维修次数增加，耗时费力，不利于农业生产。以移栽机、播种 机为例，在作业过程中由于耕作土壤表面的凹凸不平而引起的振动，容 易将秧苗、种子抖落掉，从而产生漏苗现象，特别是会影响栽植、播种 深度，造成秧苗或种子深浅不一，从而达不到移栽机、播种机在作业过 程中对秧苗、种子深度均匀性的要求，导致移栽、播种效果降低。

(3)对于播种、地面喷药、施肥、收获机械而言耕作土壤地表三 维形貌对于种子的利用率、发芽率、成活率，农药、化肥的喷洒均匀性， 农作物的收获率都有直接的影响。

由此可见，耕作土壤地表三维形貌是实现精细农业技术体系中所要 研究的一个重要指标，其对于灌溉均匀、进行机械化耕作、发挥机具效 率等都有着至关重要的作用，将直接影响到农作物高产稳产。因此，对 耕作土壤地表三维形貌测量与分析是极为必要的。

耕作土壤地表三维形貌测量的方法常分为两类：接触式和非接触式 测量。接触式测量会对耕作土壤表面产生一定程度的破坏，导致不能精 确地测量到地表粗糙度；非接触测量方法测量精度高、效率高，不会使 土壤表面变形，便于重复测量。在非接触测量方法中，测量方法有立体 摄像法、超声波测量法、图像阴影法、激光扫描法等，其中立体摄像法在后期处理数据时所需时间较长；超声波测量法测量原理是声波发射， 具有声波的扇形发射特性，不适合精准的定点测量，并且其测量精度只 能精确到厘米级；图像阴影法测量所得的地表粗糙度受测量时的光照强 度影响较大，具有不确定性；激光扫描法是由一束极小的激光进行测量， 适合精确测量，且其测量精度可以达到毫米级，因此被广泛应用。

目前耕作土壤地表三维形貌的评价分析主要是对自相关长度、功率 谱密度、分形维数等参数的分析，大多停留在二维层面上，局限性较强， 随着农业技术的不断发展，对耕作土壤地表三维形貌的要求不断提高， 原有的二维信息显然已经不能满足要求，所以研究耕作土壤地表三维形 貌空间点阵测量和建立了耕作土壤地表三维形貌评价系统是极具有必 要性的。

【发明内容】

本发明的目的是提供一种分析精度高、评价准确度好的地块三维形 貌评价方法。

为了实现上述目的，本发明提供一种地块三维形貌评价方法，所述 方法包括以下步骤：

(1)通过距离传感器对地块内各点至标准水平面的距离进行采样， 所得采样数据生成地块的三维形貌数据；

(2)通过最小二乘法去除步骤(1)的三维形貌数据中的趋势项， 再对去除了趋势项的三维形貌数据进行平滑处理，得到预处理的三维形 貌数据；

(3)利用小波分析进行地表三维形貌分析，包括：

(3.1)建立地表三维形貌粗糙度评定的数学模型

用f(x,y)表示步骤(2)得到的预处理的三维形貌数据，其为表面形 状误差、表面波纹度、表面粗糙度等多种信号之和；用f₁(x,y)表示小波 提取的低频信号粗糙度基准面，为表面形状误差、表面波纹度等低频信 号之和；用f₂(x,y)表示表面粗糙度，为高频信号，综上得到地表三维形 貌粗糙度评定的数学模型为：

f₂(x,y)＝f(x,y)-f₁(x,y)

(3.2)计算粗糙度评定基准面与表面粗糙度

设ω₀为地表三维形貌粗糙度与基准面的分界频率，则f₁(x,y)与 f₂(x,y)的频谱分别满足：

式中，分别是f(x,y)、f₁(x,y)、f₂(x,y)的 二维傅里叶变换；

对所述地表三维形貌的原始数据信号f(x,y)进行小波分解得：

式中：Ca_j,k和Cd_j,k分别代表低频和高频成分 在尺度2^-j上的系数；

信号的低频成分A_j(x,y)表示为：

信号的高频成分D_j(x,y)表示为：

所以：

式中：

其中：η(x,y)、η'(x,y)、η″(,y)和分别表示表面粗糙度、表面波 度、表面形状误差和多尺度成分；和分别表示分解和重构的 尺度函数；和ψ(x,y)分别表示分解和重构的小波函数；为软 阈值估计，为硬阈值估计；

即：

得到：

f₂(x,y)＝η(x,y)

其中：f₁(x,y)表示粗糙度评定基准面；f₂(x,y)表示表面粗糙度；

(3.3)计算地表三维参数，包括：

(3.3.1)计算三维表面轮廓算术平均偏差和三维表面轮廓均方根偏 差

以小波分析提取的基准面计算耕作土壤地表三维参数，其中耕作土 壤地表上的任易一点到对应小波提取的基准面的偏距为l(x,y,z)；采样 区域在基准面上的正投影为S；三维表面轮廓算术平均偏差Sa，表示被 测表面上各点到基准表面偏距的平均值；三维表面轮廓均方根偏差S_q， 表示轮廓偏离基准面的程度；关系式如下：

(3.3.2)计算三维平面轮廓的偏斜度Ssk：Ssk用于度量轮廓凹凸幅 值分布曲线的形状；

S_sk的关系式如下：

(3.3.3)计算三维表面高度分布的陡度S_ku：S_ku用于描述形貌高 度分布的形状，是形貌高度分布的峰度和陡度的度量；

S_ku关系式如下：

(4)判断地表三维形貌的粗糙度

通过S_a、S_q评价地表三维形貌的粗糙度：S_a、S_q越小，耕作土壤地 表越平整，S_a、S_q越大，耕作土壤地表相对越粗糙；

通过S_sk、S_ku评价地表的偏斜度、陡度：，当S_sk为负值时，表示在 基准面下的三维形貌中有较大的“深谷”；当S_sk为正值时，表示在基准 面下的三维形貌“尖峰”形状；当S_sk为零时，表示的是深谷上出现尖 锋的复合形状；当S_ku＝3时，服从于标准正态分布的轮廓凸峰；当S_ku>3 时，形貌高度集中在中心，具有很尖峭凸峰的轮廓，称为尖峰态；当 S_ku<3时，形貌高度分散，表面轮廓凸峰呈现坦峰态。

另一方面，本发明还提供一种地块三维形貌测量系统，所述系统包 括支架1、设置于支架1上部的一对y轴导轨2、分别设置在y轴导轨 2上的一对y轴滑块4，在y轴滑块4上设置一根x轴导轨6，x轴导轨 6上设一x轴滑块5，将距离传感器9与x轴滑块5机械连接；所述系统还包括与x轴滑块5相连、用于控制x轴滑块5位移的x轴电机3 和与y轴滑块4相连、用于控制y轴滑块4位移的y轴电机8，所述系 统还包括与x轴电机3和y轴电机8电连接的开关7；

所述系统还包括对x轴电机3和y轴电机8进行控制并对距离传感 器9获取的数据进行转换的控制台；

所述系统还包括与所述控制台连接的上位机，所述上位机具有对控 制台输出的数据进行计算和分析的装置。

本发明中，所述距离传感器9用于获取地表至传感器所在的距离， 例如采用激光位移传感器。

可选地，所述开关7是电感式接近开关。

可选地，所述支架1是方通支架。

根据一种优选的实施方式，所述系统还包括y轴导轨2和x轴导轨 6的水平校准装置。

在本发明中，所述对控制台输出的数据进行计算和分析的装置包 括：

(1)根据控制台输出的数据生成三维形貌数据的装置；

(2)通过最小二乘法去除三维形貌数据中的趋势项、再对去除了 趋势项的三维形貌数据进行平滑处理、获得平滑处理后的三维形貌数据 的装置；

(3)利用小波分析进行地表三维形貌分析，获取三维表面轮廓算 术平均偏差、三维表面轮廓均方根偏差、偏斜度和陡度的装置；

(4)根据三维表面轮廓算术平均偏差、三维表面轮廓均方根偏差、 偏斜度和陡度判断地表三维形貌的装置。

其中，所述利用小波分析进行地表三维形貌分析，获取三维表面轮 廓算术平均偏差、三维表面轮廓均方根偏差、偏斜度和陡度的装置的计 算方法包括：

(1)建立地表三维形貌粗糙度评定的数学模型

用f(x,y)表示经过平滑处理的地表三维形貌数据，其为表面形状误 差、表面波纹度、表面粗糙度等多种信号之和；用f₁(x,y)表示小波提取 的低频信号粗糙度基准面，为表面形状误差、表面波纹度等低频信号之 和；用f₂(x,y)表示表面粗糙度，为高频信号，综上得到地表三维形貌粗 糙度评定的数学模型为：

f₂(x,y)＝f(x,y)-f₁(x,y)

(2)计算粗糙度评定基准面与和表面粗糙度

设ω₀为地表三维形貌粗糙度与基准面的分界频率，则f₁(x,y)与 f₂(x,y)的频谱分别满足：

式中，分别是f(x,y)、f₁(x,y)、f₂(x,y)的 二维傅里叶变换；

对所述地表三维形貌的原始数据信号f(x,y)进行小波分解得：

式中：Ca_j,k和Cd_j,k分别代表低频和高频成分 在尺度2^-j上的系数；

信号的低频成分A_j(x,y)表示为：

信号的高频成分D_j(x,y)表示为：

所以：

式中：

其中：η(x,y)、η'(x,y)、η″(,y)和分别表示表面粗糙度、表面波 度、表面形状误差和多尺度成分；和分别表示分解和重构的 尺度函数；和ψ(x,y)分别表示分解和重构的小波函数；为软 阈值估计，为硬阈值估计；

即：

得到：

f₂(x,y)＝η(x,y)

其中：f₁(x,y)表示粗糙度评定基准面；f₂(x,y)表示表面粗糙度；

(3)计算地表三维参数，包括：

(3.1)以小波分析提取的基准面计算耕作土壤地表三维参数，其 中耕作土壤地表上的任易一点到对应小波提取的基准面的偏距为 l(x,y,z)；采样区域在基准面上的正投影为S；三维表面轮廓算术平均偏 差S_a，表示被测表面上各点到基准表面偏距的平均值；三维表面轮廓 均方根偏差S_q，表示轮廓偏离基准面的程度；关系式如下：

(3.2)计算三维平面轮廓的偏斜度S_sk，S_sk用于度量轮廓凹凸幅值 分布曲线的形状；S_sk的关系式如下：

(3.3)计算三维表面高度分布的陡度S_ku：S_ku用于描述形貌高度 分布的形状，是形貌高度分布的峰度和陡度的度量；

S_ku关系式如下：

本发明的方法与系统能够可以直观的看到测量耕作土壤地表三维 形貌图以及评价参数值，通过分析可以直接得到测量地块三维表面粗糙 度的状态，结果简单易懂。除此之外，所述评价系统还可以通过改变输 入参数，处理不同面积试验田，应用范围灵活。

【附图说明】

图1为激光位移传感器的测量原理示意图；

图2为地块三维形貌测量系统的立体图；

图3为通过双轴控制器控制的激光传感器移动路径示意图；

图4为地块三维形貌测量系统的上位机的交互示意图；

图5为旋耕刀三维地表形貌；

图6为双铧犁三维地表形貌；

图7为旋耕刀和平地耙三维地表形貌；

图8为双铧犁和平地耙三维地表形貌；

图9为旋耕刀旋耕后的原始二维曲线图；

图10为图9去除趋势项后的二维曲线图；

图11为图10平滑处理后的二维曲线图；

图12为MATLAB绘制的旋耕刀旋耕后地表三维形貌图；

图13为图12经过小波方法提取地表三维基准面；

图14为经过预处理的旋耕刀三维地表形貌图；

图15为经过预处理的双铧犁三维地表形貌图；

图16为经过预处理的旋耕刀和平地耙三维地表形貌；

图17为经过预处理的双铧犁和平地耙三维地表形貌；

图18为旋耕刀旋耕后的地表状况；

图19为双铧犁处理后的地表状况；

图20为旋耕刀和平地耙处理后的地表状况；

图21为双铧犁和平地耙旋耕刀旋耕后的地表状况。

【具体实施方式】

以下实施例用于非限制性地解释本发明的技术方案。

设置一套非接触式的地表三维形貌测量系统，测量范围为1×1m²， 距离分辨率为1mm，空间分辨率为0.001mm。

所述系统包括下位机硬件部分和上位机部分。下位机硬件部分由 运动测试台和控制箱两部分组成，而上位机软件部分采用美国国家仪器 公司开发的LabVIEW软件平台进行编写。

控制箱中具有双轴控制器，用于控制运动测试台上激光位移传感 器的运动；而控制箱中的采集模块则负责将采集到的表示地表高程的电 信号上传到上位机中，上位机(包括软件)将采集到的电信号进行数据 的存储，转化为空间点阵数据，实时显示地表三维形貌。

本实施例的地表三维形貌测量系统需要具备以下功能：

(1)通过双轴控制器能够控制激光位移传感器按照设计好的运动 轨迹进行空间点阵测量；

(2)采集模块能够将电信号通过串口上传到上位机软件；

(3)上位机能够实时显示地表三维形貌、空间点阵数据并储存。

此外，系统技术要求如下：

(1)测量精度不低于±1mm；

(2)最大测量范围需要达到1×1m²；

(3)测量系统在采集数据过程中采样间隔需保持一致，采集到的 高程数据能够形成空间点阵数据，并可以通过不同方式的提取数据，以 计算不同方向上耕作土壤地表粗糙度(例如：垂直于耕作方向、平行于 耕作方向等)；

(4)测量装置装配、拆卸方便、便于运输，适于在田间作业环境 下工作。

其中，测量系统中采用激光测距技术实现距离测量，具体选用北 京通测意达科技有限公司的TC–01型激光位移传感器。利用激光位移 传感器测量传感器到地面的距离，根据田间试验要求，采用测量精度及 环境要求较为合适的相位法对地表三维形貌进行测量。

激光测距的工作原理是：对激光束进行幅度调制并测定调制光往 返一次所产生的相位延迟，再根据调制光的波长，换算此相位延迟所代 表的距离。测量原理示意图如图1所示。

若调制光角频率为ω，在待测量距离L_AB上往返一次产生的相位延 迟为则对应时间t可表示为：

其中

m表示激光往返L_AB所经历的整数个波长；

Δm表示不足一个波长的分量；

则待测距离L_AB可表示为：

激光位移传感器的标定方法是通过等间隔改变传感器到被测物的 距离，获得相应输出电压值，得到测量距离与输出电压信号的对应关系。

为此，设置如图2所示的运动测试台。包括铝合金支架1、平行设 置于支架1上部的一对y轴同步导轨2、分别设置在y轴导轨2上的一 对y轴滑块4，在y轴滑块4上设置一根x轴导轨6，x轴导轨6上设 一x轴滑块5，将距离传感器9与x轴滑块5机械连接；所述系统还包 括与x轴滑块5相连、用于控制x轴滑块5位移的x轴电机3和与y 轴滑块4相连、用于控制y轴滑块4位移的y轴电机8，所述系统还包 括与x轴电机3和y轴电机8电连接的开关7。通过控制箱内的双轴控 制器控制步进驱动器间歇驱动X、Y轴电机带动滑块上的激光位移传感 器按编辑好的运动轨迹完成对待测区域的测量。

为了实现测量平面的水准基准，将一对水平尺安装在支架相互垂直 的两边，对测量系统的水平位置、垂直位置和45度角倾斜位置偏离程 度进行测量、调试、纠偏。

控制箱包括电源、全数字高性能两相步进驱动器ZD-2HD542、 KBM-30D单通道数据采集模块、MGR-1D220D10单相固态继电器、串口 转换器、CM4OL智能双轴步进控制器、开关电源等，用于接收上位机指 令，控制测试台上的滑块按测试要求进行移动，并利用串口转换器将采 集模块采集到的数据传送到上位机软件。

测量运动轨迹如图3所示，图中滑块所在的位置为程序零点，将起 始点所在位置作为机械零点，由电感式接近开关的位置设置机械零点。 测量时，首先，控制器驱动滑块做回机械零运动到达A点；其次控制器 间接驱动X、Y轴电机带动安装在滑块上的激光位移传感器沿着 A→B→C→D→E轨迹运动，与此同时，当滑块到达A点时，采集系统开 始采集；随后按实际测量范围要求循环上一步操作，直至滑块到达M 点，再次驱动X轴电机使滑块带动传感器完成M→N的运动，当滑块到 达N点时采集系统停止采集，此时完成对整个测量区域的空间点阵测 量；最后，控制器驱动滑块做回程序零运动，以待下次测量。

地表三维形貌测量系统的上位机软件部分与上述下位机硬件部分 采用串口通信方式。上位机数据采集系统一般由传感器、信号调理器、 数据采集卡、计算机和测量采集软件组成，如图4所示，数据采集系统 由LabVIEW 2013编写，程序的编写属于本领域技术人员已经掌握的知 识，在此不做赘述。

其中，为了能够实时显示测量耕作土壤地表的三维形貌图，在程序 中引用了MathScript节点，并调用mesh函数得出耕作土壤地表三维形 貌图。

以一具体地块为例，利用测量装置进行田间试验之前，需要对测量 装置进行精度检验试验。试验选择精度等级为IT2级的标准量块为检测 对象，其中1–3mm的量块其误差范围为0–1.2×10^-3mm,4–10mm的量 块其误差范围为0–1.5×10^-3mm。试验时用所述水平仪找正平台，将厚 度为1–10mm的10块量块(T1-T10)分别摆放在平台上，调平所述测 试台。由激光位移传感器测量相对平台的距离，通过编辑控制系统驱动 所述激光传感器对所述平台上的量块进行测量，测量结果如表1所示：

表1测量装置测量量块数据

对比分析测量均值数据与量块数据，拟合得到拟合曲线为 y＝1.0003x+0.0205,相关系数R²为0.9997，测量数据的均方根误差为 0.017mm，证明测量装置的测量精度较高，可以较准确的测量出耕作土 壤的地表粗糙度。

为了能够高效完成对整个平面的测量，在所述精度检验试验中采集 频率均选择4Hz，同步带运动精确到0.001mm，采样间隔可以根据不同 的试验要求设定，当设置采样间隔为1、2、5、10mm等间隔时(这里X、Y轴选择相同的采样间隔)，在采样频率确定的情况下对应的同步带速 度分别为240、480、1200、2400mm/min。完成1×1m²的测量需要时间 为69.444、17.361、2.778、0.694h。通过所述检验试验，发现采样间 隔为10mm的时候也能较好的反映出物体的表面轮廓，结合测量效率， 设置田间试验采样间隔为10mm。

在本发明的实施例中，通过田间试验测量不同整地机具下耕作土壤 地表三维形貌并进行分析评价，试验时，利用微耕机旱地旋耕刀、四轮 拖拉机用双铧犁、平地耙三种农用机具组合成4种不同形式的平整土地 地段，每种地段选5个地块样本进行测量，共计20次试验。试验所用 到的主要设备有：

1、耕作土壤地表三维形貌空间点阵测量装置

2、22kW拖拉机

3、微耕机旱地旋耕刀；

4、四轮拖拉机用双铧犁；

5、平地耙。

4种形式平整土地地段包括：拖拉机带动旋耕机耕后的地段、旋耕 机耕后再经平地耙耙后的地段、双铧犁犁后的地段、双铧犁犁后再经平 地耙耙后的地段，每种地段长15m，宽2m。

试验过程中，上位机采集到的数据表示的是激光位移传感器到地表 的距离，当某一点数据较其他数据大时，表示此点是凹处，反之，则是 凸起处，如图5-8所示为所述测量装置对所述4种地段测量所得三维形 貌图，从图中可以看出图中都有一些突起，四幅图在Z方向都有一种向 上的趋势，因此在对这些数据进行利用之前，需要对其进行分析和处理，如剔除误差较大或明显不正确的点，以提高数据的准确性。

在本发明的实施例中，对数据进行分析时，如果不处理趋势项将会 对分析结果产生巨大的影响，因此消除趋势项对数据处理极为重要。

如图9所示是以旋耕刀旋耕后所测得的原始数据，图中的曲线有一 种明显向上的趋势，在本发明的实施例中，用最小二乘法去除趋势项。 图10所示是以旋耕刀旋耕后所测得的数据利用最小二乘法消除趋势项 的结果。

在本发明的实施例中，去除趋势项后的数据还需经过平滑处理，得 到经过预处理的有效数据，以旋耕刀旋耕后所测得的数据为例，平滑处 理后如图11所示。

为更好地反应地表形貌部分特性，在进行田间数据分析时，利用小 波分析对耕作土壤地表三维形貌进行三维粗糙度以及基准面的提取。提 取时，用f(x,y)表示地表三维形貌的原始数据信号，其为表面形状误差、 表面波纹度、表面粗糙度等多种信号之和；用f₁(x,y)表示小波提取的 低频信号粗糙度基准面，为表面形状误差、表面波纹度等低频信号之和； 用f₂(x,y)表示表面粗糙度，为高频信号，综上得到地表三维形貌粗糙度 评定的数学模型为：

f₂(x,y)＝f(x,y)-f₁(x,y)

设ω₀为地表三维形貌粗糙度与基准面的分界频率，则f₁(x,y)与 f₂(x,y)的频谱分别满足：

式中分别是f(x,y)、f₁(x,y)、f₂(x,y)的 二维傅里叶变换。

对所述地表三维形貌的原始数据信号f(x,y)进行小波分解得：

式中：Ca_j,k和Cd_j,k分别代表低频和高频 成分在尺度2^-j上的系数。

信号的低频成分A_j(x,y)可表示为：

信号的高频成分D_j(x,y)可表示为：

所以：

式中：

其中：η(x,y)、η'(x,y)、η″(,y)和分别表示表面粗糙度、表面波 度、表面形状误差和多尺度成分；和分别表示分解和重构的 尺度函数；和ψ(x,y)分别表示分解和重构的小波函数。为软 阈值估计，为硬阈值估计。

即：

得：

f₂(x,y)＝η(x,y)

其中：f₁(x,y)表示粗糙度评定基准面；f2(x,y)表示表面粗糙度。

在本发明的实施例中，将所述地表三维形貌数据信号f(x,y)进行预 处理后在MATLAB画出的旋耕刀旋耕后地表三维形貌如图12。

如图13是利用所述小波方法提取的旋耕刀旋耕后的地表三维基准 面。

在本发明的实施例中，以小波分析提取的基准面对耕作土壤地表三 维参数进行计算，其中耕作土壤地表上的任易一点到对应小波提取的基 准面的偏距为l(x,y,z)；采样区域在基准面上的正投影为S；三维表面 轮廓算术平均偏差Sa，为被测表面上各点到基准表面偏距的平均值； 三维表面轮廓均方根偏差Sq，反应轮廓偏离基准面的程度。

三维平面轮廓的偏斜度Ssk，用来度量轮廓凹凸幅值分布曲线的形 状，功能与二维中偏态系数Rsk相同。当三维平面轮廓的偏斜度为负值 时，表示在基准面下的三维形貌中有较大的“深谷”；当三维平面轮廓 的偏斜度为正值时，表示在基准面下的三维形貌“尖峰”形状；当三维 平面轮廓的偏斜度为零时，表示的是深谷上出现尖锋的复合形状。

S_ku─三维表面高度分布的陡度：描述形貌高度分布的形状，是形 貌高度分布的峰度和陡度的度量，功能与二维中偏态系数Rku相同。

当S_ku＝3时，服从于标准正态分布的轮廓凸峰；

当S_ku>3时，形貌高度集中在中心，具有很尖峭凸峰的轮廓，称 为尖峰态；

当S_ku<3时，形貌高度分散，表面轮廓凸峰呈现坦峰态。

如表2所示对四种不同耕具下20块试验田求得三维参数。

表2不同耕具下耕作土壤地表三维参数

在本发明的实施例中，Sa、Sq被普遍运用于描述物体三维表面的 粗糙度特征，用于描述轮廓凸峰的相对高度差，从表2可以看出四种耕 作土壤表面轮廓凸峰的相对高度差从小到大依次是旋耕刀和平地耙<双 铧犁和平地耙<旋耕刀<双铧犁,即旋耕刀和平地耙的耕作表面相对较平 坦；旋耕刀&平地耙、双铧犁&平地耙、旋耕刀耕作三维平面轮廓的偏 斜度S_sk均大于0，在基准面上的耕作土壤地表呈现尖峰形状，双铧犁 耕作三维平面轮廓的偏斜度S_sk小于0，在基准面下的耕作土壤地表轮 廓凸峰呈现深谷形状；三维表面高度分布的陡度S_ku的趋势与Sa、Sq参 数的成正相关关系，旋耕刀和平地耙耕作后的表面高度分布的陡度S_ku大于3轮廓的凸峰呈现尖峰态，旋耕刀、双铧犁耕作后的维表面高度分 布的陡度S_ku小于3轮廓的凸峰呈现坦峰态。

针对耕作土壤地表粗糙度而言，主要由Sa、Sq的大小来评价，Sa、 Sq越小，耕作土壤地表相对越平整，Sa、Sq越大，耕作土壤地表相对 越粗糙，S_sk、S_ku主要用来描述耕作土壤地表的偏斜度、陡度。由上可 得，平地耙工作后，耕作土壤表面的三维参数Sa、Sq会变小，耕作土 壤表面变得更加平坦，产生尖峰形状的凸峰。

根据已采集到的所述耕作土壤地表三维参数，在本发明的实施例 中，在MATLAB软件中实现数据预处理、三维形貌图的显示、三维评价 参数的计算，并利用MATLABGUI编辑耕作土壤地表三维形貌评价系统， 该评价系统主要由输入参数、图形显示、参数计算三部分组成，系统评 价界面，在运行之前根据需要处理试验田的面积大小，输入其测量长度、 测量宽度参数，输入数据采集系统采集时设置好的采样间隔；参数设置 完成后，点击“开始处理”按钮，会跳出选择文件位置的选项，根据本 次需要处理的数据选择对应的文件，系统进行处理；处理完成后，在图 形显示区域会显示耕作土壤地表三维形貌原始图、小波三维基准面、小 波三维粗糙面、二维曲线图，在参数统计区域会显示此次处理试验田的 三维评价参数Sa、Sq、Ssk、Sku；点击“保存处理结果”会将图形显 示区域、参数统计区域内的图像、参数根据个人需要进行保存；处理完 成后，按“退出”按钮退出评价系统。

进一步地，利用所述耕作土壤地表三维形貌评价系统处理所述4 种不同整地机具下的地块，输入其测量长度为1000mm、测量宽度为 1000mm、采样间隔为10mm，点击“开始处理”，将处理结果进行保存， 如图14-17所示是4种不同整地机具下经过预处理后的耕作土壤地表三 维形貌图，其利用评价系统处理得到的三维评价参数如表3。

表3三维评价参数表

可见：旋耕刀和平地耙工作后的地表最为平整，双铧犁和平地耙、 旋耕刀较为平整，耕作土壤基准面上的地表形貌呈现的凸峰均为尖峰 状；双铧犁工作后的地表最为粗糙，在其耕作土壤基准面小的地表形貌 呈现“深谷”，其计算结果与实际情况(如图18-21所示)相符合。

因此，所述评价系统可以直观的看到测量耕作土壤地表三维形貌图 以及评价参数值，通过分析可以直接得到测量地块三维表面粗糙度的状 态，结果简单易懂。除此之外，所述评价系统还可以通过改变输入参数， 处理1×1m²内的不同面积试验田，应用范围灵活。

Claims (7)

1.地块三维形貌评价方法，所述方法包括以下步骤：

(1)通过距离传感器对地块内各点至标准水平面的距离进行采样，所得采样数据生成地块的三维形貌数据；

(2)通过最小二乘法去除步骤(1)的三维形貌数据中的趋势项，再对去除了趋势项的三维形貌数据进行平滑处理，得到预处理的三维形貌数据；

(3)利用小波分析进行地表三维形貌分析，包括：

(3.1)建立地表三维形貌粗糙度评定的数学模型

用f(x,y)表示步骤(2)得到的预处理的三维形貌数据，其为表面形状误差、表面波纹度、表面粗糙度等多种信号之和；用f₁(x,y)表示小波提取的低频信号粗糙度基准面，为表面形状误差、表面波纹度等低频信号之和；用f₂(x,y)表示表面粗糙度，为高频信号，综上得到地表三维形貌粗糙度评定的数学模型为：

f₂(x,y)＝f(x,y)-f₁(x,y)

(3.2)计算粗糙度评定基准面与表面粗糙度

设ω₀为地表三维形貌粗糙度与基准面的分界频率，则f₁(x,y)与f₂(x,y)的频谱分别满足：

式中，分别是f(x,y)、f₁(x,y)、f₂(x,y)的二维傅里叶变换；

对所述地表三维形貌的原始数据信号f(x,y)进行小波分解得：

式中：Ca_j,k和Cd_j,k分别代表低频和高频成分在尺度2^-j上的系数；

信号的低频成分A_j(x,y)表示为：

信号的高频成分D_j(x,y)表示为：

所以：

式中：

其中：η(x,y)、η'(x,y)、η”(,y)和分别表示表面粗糙度、表面波度、表面形状误差和多尺度成分；和分别表示分解和重构的尺度函数；和ψ(x,y)分别表示分解和重构的小波函数；为软阈值估计，为硬阈值估计；

即：

得到：

f₂(x,y)＝η(x,y)

其中：f₁(x,y)表示粗糙度评定基准面；f₂(x,y)表示表面粗糙度；

(3.3)计算地表三维参数，包括：

(3.3.1)计算三维表面轮廓算术平均偏差和三维表面轮廓均方根偏差

以小波分析提取的基准面计算耕作土壤地表三维参数，其中耕作土壤地表上的任易一点到对应小波提取的基准面的偏距为l(x,y,z)；采样区域在基准面上的正投影为S；三维表面轮廓算术平均偏差S_a，表示被测表面上各点到基准表面偏距的平均值；三维表面轮廓均方根偏差S_q，表示轮廓偏离基准面的程度；关系式如下：

(3.3.2)计算三维平面轮廓的偏斜度S_sk

S_sk用于度量轮廓凹凸幅值分布曲线的形状；S_sk的关系式如下：

(3.3.3)计算三维表面高度分布的陡度S_ku：S_ku用于描述形貌高度分布的形状，是形貌高度分布的峰度和陡度的度量；

S_ku关系式如下：

(4)判断地表三维形貌的粗糙度

通过S_a、S_q评价地表三维形貌的粗糙度：S_a、S_q越小，耕作土壤地表越平整，S_a、S_q越大，耕作土壤地表相对越粗糙；

通过S_sk、S_ku评价地表的偏斜度、陡度：当S_sk为负值时，表示在基准面下的三维形貌中有较大的“深谷”；当S_sk为正值时，表示在基准面下的三维形貌“尖峰”形状；当S_sk为零时，表示的是深谷上出现尖锋的复合形状；当S_ku＝3时，服从于标准正态分布的轮廓凸峰；当S_ku>3时，形貌高度集中在中心，具有很尖峭凸峰的轮廓，称为尖峰态；当S_ku<3时，形貌高度分散，表面轮廓凸峰呈现坦峰态。

2.地块三维形貌测量系统，其特征在于所述系统包括支架(1)、设置于支架(1)上部的一对y轴导轨(2)、分别设置在y轴导轨(2)上的一对y轴滑块(4)，在y轴滑块(4)上设置一根x轴导轨(6)，x轴导轨(6)上设一x轴滑块(5)，将距离传感器(9)与x轴滑块(5)机械连接；所述系统还包括与x轴滑块(5)相连、用于控制x轴滑块(5)位移的x轴电机(3)和与y轴滑块(4)相连、用于控制y轴滑块(4)位移的y轴电机(8)，所述系统还包括与x轴电机(3)和y轴电机(8)电连接的开关(7)；

所述系统还包括对x轴电机(3)和y轴电机(8)进行控制并对距离传感器(9)获取的数据进行转换的控制台；

所述系统还包括与所述控制台连接的上位机，所述上位机具有对控制台输出的数据进行计算和分析的装置。

3.根据权利要求1所述系统，其特征在于所述距离传感器(9) 是激光位移传感器。

4.根据权利要求2所述系统，其特征在于所述开关(7)是电感式接近开关。

5.根据权利要求2所述系统，其特征在于所述系统还包括y轴导轨(2)和x轴导轨(6)的水平校准装置。

6.根据权利要求2所述系统，其特征在于所述对控制台输出的数据进行计算和分析的装置包括：

(1)根据控制台输出的数据生成三维形貌数据的装置；

(2)通过最小二乘法去除三维形貌数据中的趋势项、再对去除了趋势项的三维形貌数据进行平滑处理、获得预处理的三维形貌数据的装置；

(3)利用小波分析进行地表三维形貌分析，获取三维表面轮廓算术平均偏差、三维表面轮廓均方根偏差、偏斜度和陡度的装置；

(4)根据三维表面轮廓算术平均偏差、三维表面轮廓均方根偏差、偏斜度和陡度判断地表三维形貌的装置。

7.根据权利要求6所述的系统，其特征在于所述利用小波分析进行地表三维形貌分析，获取三维表面轮廓算术平均偏差、三维表面轮廓均方根偏差、偏斜度和陡度的装置的计算方法包括：

(1)建立地表三维形貌粗糙度评定的数学模型

用f(x,y)表示预处理的三维形貌数据，其为表面形状误差、表面波纹度、表面粗糙度等多种信号之和；用f₁(x,y)表示小波提取的低频信号粗糙度基准面，为表面形状误差、表面波纹度等低频信号之和；用f₂(x,y)表示表面粗糙度，为高频信号，综上得到地表三维形貌粗糙度评定的数学模型为：

f₂(x,y)＝f(x,y)-f₁(x,y)

(2)计算粗糙度评定基准面和表面粗糙度

设ω₀为地表三维形貌粗糙度与基准面的分界频率，则f₁(x,y)与f₂(x,y)的频谱分别满足：

式中，分别是f(x,y)、f₁(x,y)、f₂(x,y)的二维傅里叶变换；

对所述地表三维形貌的原始数据信号f(x,y)进行小波分解得：

式中：Ca_j,k和Cd_j,k分别代表低频和高频成分在尺度2^-j上的系数；

信号的低频成分A_j(x,y)表示为：

信号的高频成分D_j(x,y)表示为：

所以：

式中：

其中：η(x,y)、η'(x,y)、η”(,y)和分别表示表面粗糙度、表面波度、表面形状误差和多尺度成分；和分别表示分解和重构的尺度函数；和ψ(x,y)分别表示分解和重构的小波函数；为软阈值估计，为硬阈值估计；

即：

得到：

f₂(x,y)＝η(x,y)

其中：f₁(x,y)表示粗糙度评定基准面；f₂(x,y)表示表面粗糙度；

(3)计算地表三维参数，包括：

(3.1)以小波分析提取的基准面计算耕作土壤地表三维参数，其中耕作土壤地表上的任易一点到对应小波提取的基准面的偏距为l(x,y,z)；采样区域在基准面上的正投影为S；三维表面轮廓算术平均偏差S_a，表示被测表面上各点到基准表面偏距的平均值；三维表面轮廓均方根偏差S_q，表示轮廓偏离基准面的程度；关系式如下：

(3.2)计算三维平面轮廓的偏斜度S_sk，S_sk用于度量轮廓凹凸幅值分布曲线的形状；S_sk的关系式如下：

(3.3)计算三维表面高度分布的陡度S_ku：S_ku用于描述形貌高度分布的形状，是形貌高度分布的峰度和陡度的度量；

S_ku关系式如下：

。