宋佳音，池志祥，张晓鹏，朱庆林，宋文龙，杨柳松

(1. 东北林业大学机电工程学院，哈尔滨 150040； 2. 国网黑龙江省送变电工程有限公司，哈尔滨 150040)

森林是陆地生态系统的重要组成部分，对保持水土和维持气候具有重要作用。及时掌握林分蓄积量及森林生物量对生态系统研究具有重要价值。树高作为当前衡量森林系统稳定最广泛的指标，其高度的精确测量具有重要意义[1-2]。目前常用的测量树高方法主要有布鲁莱斯测高器[3]测量和超声波测高器[4]测量，布鲁莱斯测高器结构简单，但测量时需要先用皮尺测出仪器与待测树木的距离，存在读数误差及测量精度低、应用不便等缺点；超声波测高器是以超声波完成距离测量，再辅助角度数据计算树木高度，该方法测量时较为烦琐，且超声波测距时存在发散角，因此在林区树木较密的情况下，信号容易被其他树木遮挡反射，造成测距的不准确而影响高度的测量。

为提高测量精度，减小测量复杂度，近景摄影测量技术[5-6]发挥了重要作用。陈相武等[7]提出了一种基于单目视觉的立木高度测量方法，该方法具有较高的精度，但测量视角小，不能在近距离获取树木的完整图像，不适用于密集的森林测量。对于解决视角问题的有效方法是采用双目视觉测量，佘银海等[8]提出了一种双目视觉测量方法，该方法有效解决了测量视角的问题，但因为需要对采集的左右两幅图像进行特征点提取，匹配和三维像素点重建，存在特征点匹配困难的缺点。鱼眼镜头具有短焦距、大视角的特点，且采集的鱼眼图像特征点提取简单，宋佳音等[9]提出了一种基于半球图像的立木高度测量方法，解决了测量视角和特征点匹配的问题，但是忽略了鱼眼镜头虚拟距离，导致测量结果精度不足。

在鱼眼镜头测量树高的基础上，现引入Scaramuzza鱼眼镜头标定模型和鱼眼镜头虚拟距离参数，提高测量的精确度。首先根据鱼眼镜头成像模型确立世界坐标系与图像像素坐标系中对应点的关系；然后以畸变系数的等距投影模型为基础，建立测量系统模型，测量系统所需参数用Scaramuzza鱼眼镜头标定方法求取；最后利用多算子融合算法对鱼眼图像进行轮廓提取，用改进四叉树算法进行图像分割，Hough变换直线检测提取极值点，实现树高的测量。利用手机搭配鱼眼镜头，对黑白棋盘格进行鱼眼图像获取，利用棋盘格内角点距离去验证测量系统模型的准确性；对待测树木进行鱼眼图像采集，用该测量系统模型对树高进行反演，用提出的计算树木高度方法与文献[9]测量树高方法进行对比，检验该测量方法的优越性。

1 鱼眼相机成像模型

鱼眼镜头成像模型能确定世界坐标系与图像像素坐标系之间的坐标变换关系，将世界坐标系中的目标点投影到图像像素坐标系中，鱼眼镜头成像模型如图1所示。

图1中坐标系OwXwYwZw是世界坐标系，它是一个在适当空间中选取的参考坐标系，P(xw，yw，zw)表示世界坐标系中目标点P的坐标。坐标系OXYZ是鱼眼镜头的参考坐标系。坐标系OcXcYcZc是相机坐标系，原点Oc位于摄像机光学中心，坐标轴Xc和Yc分别平行于X轴和Y轴，Zc轴是相机的主光轴实际上点Oc和O重合，P′(xc，yc，zc)是相机坐标系中点P对应的成像点，r′是点P′到光轴的距离。坐标系ouv是图像像素坐标系，原点o是光学中心点，坐标轴u和v分别平行于Xc轴和Yc轴，P″(u，v)是图像像素坐标系中点P对应的成像点，r″是点P″到光轴的距离。

根据等距投影定理[10]得

r′=fω

(1)

ω=tan-1(r/L)=tan-1[(xω2+yω2)1/2/L]

(2)

式中：f为光学系统的物方焦距；ω为点P相对于光轴的入射角；L为世界坐标系中的点到鱼眼镜头中心的水平距离。

由于鱼眼镜头存在畸变，为保证图像的均匀性，引入失真系数λ得

r′=λfω

(3)

相机平面中心点为Oc(x0，y0)，P′点坐标为(xc，yc)，P点坐标为(xw，yw，zw)，令Xc和Yc轴方向的畸变系数分量为λx和λy，则

(4)

(5)

式中：θ为点P的方位角，也是相机坐标系中点P′的方位角。

图像像素坐标系中心点o的坐标为(u0，v0)，P′通过等距投影得到P″，相机坐标系与图像像素坐标系中对应点关系为

(6)

式中：mx、my为放大系数。

令kx=λxmxf，ky=λymyf，联立式(2)～式(6)得

(7)

式中：kx、ky为畸变系数。式(7)是引入畸变系数的等距投影模型，建立世界坐标系与图像像素坐标系中对应点的转换关系。该模型需要求取的参数为光学中心点o(u0，v0)，畸变系数kx和ky，水平距离L。根据图1可得

L=h+l

(8)

式(8)中：l为鱼眼镜头的虚拟成像距离；h为世界坐标系中的点到鱼眼镜头顶部的水平距离。

2 测量模型的建立

以畸变系数的鱼眼镜头等距投影模型为基础，构建测量系统模型。测量系统模型由鱼眼镜头、测距仪和智能手机组成，如图2所示。

图2 测量系统模型Fig.2 Measurement system model

2.1 测量系统数学模型

在图2中，A′和B′为成像平面中的参考点，处于图像像素坐标系中，设A′坐标为(uA′，vA′)，B′坐标为(uB′，vB′)。A、B为被测物中两点，处于世界坐标系中。A′和A存在对应关系，B′和B存在对应关系，由式(7)可得，世界坐标系与图像像素坐标系中对应点的关系为

联立式(7)～式(9)可得

(10)

令zw=0，根据式(10)，可得坐标A(xA，yA)，B(xB，yB)的值。计算AB之间的距离H，H为测量系统模型所求的结果。

H=[(xA-xB)2+(yA-yB)2]1/2

(11)

2.2 光学中心点获取

为获取光学中心点，引入Scaramuzza鱼眼镜头标定模型[11]。在Scaramuzza模型中，世界坐标系与相机坐标系中对应点的关系为

(12)

式(12)中：R为旋转矩阵；T为平移矩阵。

相机坐标系中的点与理想图像像素点的关系为

(13)

式(13)中： (u″，v″)为没有发生畸变的理想图像像素点；μ为标量因子；a0、a2、a3、a4为Scaramuzza模型的多项式系数；p为关于u″、v″的函数。

(14)

实际图像像素坐标系中的点与理想图像像素点的关系为

(15)

式(15)中： (u0，v0)为光学中心点； (u，v)为发生畸变的实际图像像素坐标系中的点。

2.3 虚拟成像距离的获取

世界坐标系是一个参考坐标系，可以人为选择。根据式(12)，当zw=0时，zc=t3，由图1得，世界坐标系到相机坐标系的距离为

t3=L

(16)

联立式(8)、式(16)得

l=t3-h

(17)

2.4 畸变系数的获取

引入畸变系数kx和ky，为求取畸变系数，确定标定板平面与成像平面的对应关系，建立标定畸变系数模型，标定板平面位于世界坐标系中，成像平面位于图像像素坐标系中，标定畸变系数模型如图3 所示。

图3 标定畸变系数模型Fig.3 Calibration distortion coefficient model

o为光学中心点，坐标为(u0，v0)，r″为点P″到光学中心点o的距离，r1″为r″在u轴方向的分量，r2″为r″在v轴方向的分量。ow为光轴与标定板平面的交点，坐标为(xow，yow)，r为点P到ow的距离，r1为r在Xw轴方向的分量，r2为r在Yw轴上的分量。

a、b、c、d分别对应图 5中的图5(a)、图 5(b)、图 5(c)和图 5(d)； Fmax表示对两幅图像取合集； Fmin表示对两幅图像取交集图4 多算子融合算法示意图 Fig.4 Schematic diagram of multi operator fusion algorithm

(18)

ω=tan-1(r/h+l)

(19)

联立式(12)、式(13)、式(15)得

(20)

通过式(20)求得ow，根据图3可得

(21)

(22)

联立式(18)～式(22)，畸变系数的计算公式为

(23)

3 图像处理与特征点提取

使用多算子融合算法对树木进行轮廓提取[12]，算法步骤如图4所示。为增强图像的处理能力，将采集到的彩色图像转换成灰度图像，如图5(a)所示；对灰度图像进行灰度变换增强，如图5(b)所示；由于相机自身缺陷和外界因素的限制，使获取到的图像含有噪声，对图像进行自适应滤波去除噪声[13]，如图5(c)所示；多算子融合算法的轮廓提取效果如图5(d)所示；用改进四叉树算法对图像进行分割[14]，如图5(e)所示；最后运用Hough变换直线提取树木极值点[15]，通过Hough变换检测提取全幅图像直线信息，标识的直线中包含着树木轮廓的边缘，直线的交点可以选择性地获取树木轮廓极值点信息，进而得到其像素坐标值，如图5(f)为Hough变换检测提取的树木图像，选出直线相交的红色点A即为树木冠层的最高点，红色交点B为树木的最低点。

图5 图像处理过程Fig.5 Processing of the images

4 实验

4.1 实验准备

实验设备包括手机、鱼眼镜头、黑白棋盘标定板和激光测距仪，鱼眼图像的有效像素为3156×3156，棋盘格的间隔为50mm。选取10幅平行于标定板但距离不同的鱼眼图像。实验在Windows10系统平台上进行，处理软件为MATLAB2018b。

图 6 鱼眼图像Fig.6 Fisheye image

4.2 实验结果分析

求光学中心点用Scaramuzza模型，直接调用MATLAB2018b自带鱼眼镜头标定工具箱对鱼眼图进行处理。图6(a)是棋盘格图，图6(b)是棋盘格角点提取图。求得光学中心点o(1578.50，1579.60)。

表 1 l和ow的值

表 2 AB、AC、AD测量结果

表 3 树木高度测量结果

图 7 kx与ky标定结果Fig.7 kx and ky calibration results

图 8 A、B、C、D位置图Fig.8 Location map of A，B，C and D

为验证模型的精度，拍摄5组不同距离的棋盘格图片，图8为拍摄距离612mm和1502mm下的棋盘格图。在棋盘格上取3组距离来检验测量模型精度。表2为AB、AC和AD之间距离的测量结果。根据表2中棋盘格的计算值与实测值的误差分析可知，测量距离时，最大测量误差为8.283mm，最大相对误差为1.274%，最小相对误差为0.323%，平均相对误差为0.823%，证明了模型的有效性和准确性。利用模型对树高进行测量，实际树高用全站仪获取，测量结果如表3所示。

表 3 可得利用本文方法测量树高，测量误差在-0.196～0.195m，测量最高相对误差为2.75%，平均相对误差为1.505%；用文献[9]测树高方法测量，最高相对误差为6.13%，平均相对误差为3.509%。通过对比可以得出所提出的测量树高方法更加准确。

5 结论

提出一种基于鱼眼镜头的树高测量方法。通过建立测量模型、求解参数来反演树高，用黑白棋盘格间的角点验证测量模型，并通过该模型对树木进行实验测量，得到如下结论。

(1)棋盘格内角点间距测量实验结果表明，棋盘格的实际距离与测量距离最大误差为8.283mm，最大相对误差为1.274%，最小相对误差为0.323%，平均相对误差为0.823%，说明了该测量模型的有效性。

(2)用该测量模型对树木高度进行测量，测量误差介于-0.196～0.195m，树木测量最高相对误差为2.75%，平均相对误差为1.505%，与文献[9]提出的测量方法进行对比，测量误差更小，测量结果更加精确。