CN102663396A - 一种水稻乳熟期自动检测的方法 - Google Patents

Abstract

本发明提供一种水稻乳熟期的自动检测方法，以稻田中实时采集到的前下视图像为对象，自动地较准确地对图像中的稻穗进行分割，进而在特定的检索区域内，计算稻穗的弯曲角度。该方法以表征水稻稻穗的颜色特征以及稻穗弯曲角度的形态特征作为判断依据，可以实时地对水稻稻穗进行分割，进而对水稻乳熟期进行检测，检测结果准确率高，误差为前后三天，对水稻相应的农事活动具有重要的指导意义。

Description

一种水稻乳熟期自动检测的方法

技术领域

本发明属于数字图像处理和农业气象观测交叉的领域，具体涉及到一种水稻乳熟期的自动检测方法，即以田间实时拍摄的水稻图像序列为对象，从图像特征上实时自动检测水稻足否进入乳熟期的方法。

背景技术

水稻是我国主要的粮食作物之一，在我国南方广泛种植。一直以来，对于水稻各个发育期的观测主要是通过人工观测，受观测员主观因素的影响较大；同时由于水稻种植地域广、生长周期长，利用人工观测显然不够经济。因此，通过每天所拍摄的稻田图像，借助图像处理的手段，对其生长发育进行检测显得十分必要。水稻乳熟期检测是田间水稻生长期自动检测的一个重要环节，水稻乳熟期是水稻的第二个水分临界期，此时期水稻对水分非常敏感，如果乳熟期水稻受旱，会严重影响籽粒灌浆，造成秕粒增多，千粒重下降进而影响产量。准确识别水稻乳熟期，可以便与水稻后期相关的田间作业，它是水稻农业气象观测的一个重要内容。

2003年孟亚利在《中国农业科学》发表论文“基于生长过程的水稻阶段发育与物候期模拟模型”构建了预测水稻顶端发育阶段与物候期的模拟模型。每日热效应和光周期效应分别选用Beta函数和二次曲线函数来描述，并引入5个遗传参数调节水稻发育的遗传差异，使不同类型品种到达特定发育阶段所需的生理发育时间保持恒定。利用不同环境下的播期试验数据对不同类型品种的出苗期、穗分化期、抽穗期和成熟期进行了预测，不同类型品种4个生育期的预测误差RMSE(根均方差)平均分别为1.47、5.10、4.58和3.37d。2008年孙华生在博士学位论文“利用多时相MODIS数据提取中国水稻种植面积和长势信息”和2009年在《遥感学报》上发表的论文“利用MODIS数据识别水稻关键生长发育期”中研究利用EOS-MODIS数据空间中覆盖中国范围的数据，对全国范围水稻种植面积和生长信息进行提取，同时利用时间序列对水稻的主要生长期进行了识别，生长期识别误差的时间在±16天左右；2009年赵永林在硕士论文“基于MODIS数据的双季稻实时监测研究”中于2007年4月至2008年10月，在江西省南昌设定了4个双季水稻观测样点，每8天对双季稻进行一次观测，并记录水稻的干鲜重、株高、叶面积指数LAI、稻田灌水深度的变化以及水稻管理措施等。分析了MODIS植被指数与实地采样数据之间的关系，建立了水稻实时监测模型和双季晚稻的估产模型。结果表明：(1)MODIS-EVI(增强型植被指数)在水稻分蘖孕穗期和灌浆乳熟期与LAI有很强的相关关系，相关系数分别为0.90，0.89，拥有很强的稳定性和抗干扰能力，更适合研究区高温、高湿、高植被覆盖率的实际状况。(2)LSWI在抽穗扬花期与LAI的相关系数达到0.8以上，将MODIS-EVI和LSWI相结合可以实现研究区双季稻从分蘖期到乳熟期的实施监测。以上方法中孟亚利构建了预测水稻顶端发育阶段与物候期的模拟模型，对水稻的四个生长期进行了预测，但没有给出水稻乳熟期的具体预测方法，且从其他四个生长期的预测结果来看，其模型预测误差较大，此种方法不适合水稻乳熟期的实时自动检测；孙华生和赵永林是利用遥感数据对水稻的种植面积进行确定，利用时间序列中水稻的叶面积指数，但由于遥感图像分辨率低，且容易受到云层、云阴影和气溶胶等影响的原因，MODIS不太适合研究区早稻的实时监测。然而，在农业气象观测领域，主要通过利用稻田的水稻图像实时的对水稻种植密度进行较准确的计算，以便及时指导后期给水灌溉、估产等农事活动，通过分析上述利用遥感图像等方法都不可行。

发明内容

本发明目的在于提供一种水稻乳熟期的自动检测方法，能够利用田间实时获取的水稻数字图像准确地检测出水稻是否进入乳熟期。

一种水稻乳熟期的自动检测方法，首先根据水稻历史图像序列离线训练得到稻穗像素点的水稻稻穗颜色信息表，水稻稻穗颜色信息表包含信息有稻穗像素点在各亮度i＝1，…，255下对应的(a，b)颜色值的分布情况，然后按照如下方式在线判断水稻是否进入水稻乳熟期：

(1)将当前水稻图像t转化到L-a-b颜色空间，将图像t中各像素点的L-a-b颜色空间色彩分量与所述水稻稻穗颜色信息表进行比对，确定满足颜色值分布情况的色彩分量对应的像素点为稻穗像素点；

(2)从由稻穗像素点构成的稻穗区域截取最能表征稻穗完整形态结构的子区域作为稻穗角度检测区域；

(3)计算稻穗角度检测区域的稻穗弯曲角度；

(4)依据稻穗弯曲角度判断水稻是否进入乳熟期。

进一步地，所述根据水稻穗像历史图像序列离线训练得到稻穗像素点的水稻稻穗颜色信息表的具体实现方式为：

(01)根据水稻穗像历史图像序列生成稻穗样本图像集；

(02)获取稻穗样本图像集中所有图像的稻穗像素点在L-a-b颜色空间的色彩分量(L，a，b)，其中L表示亮度，a表示红或绿的对应值，b表示黄或蓝的对应值，其中L，a，b为8位无符号整形；

(03)统计在每一亮度i＝1，…，255下(a，b)颜色值相同的稻穗像素点个数；

(04)在每一亮度i＝1，…，255下，按照(a，b)颜色值相同的稻穗像素点个数越多，则稻穗像素点在亮度i下表现为(a，b)颜色值的可能性越大的原则选取稻穗(a，b)颜色的标定值；

(05)定义三维ColorMap全零矩阵，第一维表示a颜色值，第二维表示b颜色值，第三维表示L亮度，令步骤(04)在亮度i下选取的稻穗(a，b)颜色标定值表示为(s，t，i)，将ColorMap全零矩阵中的(s，t，i)元素值ColorMap(s，t，i)赋值为1，从而得到在每一亮度i下表征稻穗颜色的二值图像。

(06)对ColorMap中每一亮度i下的二值图像依次进行膨胀和腐蚀，得到水稻稻穗颜色信息表。

进一步地，所述膨胀和腐蚀的具体实现方式为：采用结构元素大小为5的圆形结构元素对ColorMap中每一亮度i下的二值图像进行膨胀，再采用结构元素大小为3的圆形结构元素对ColorMap中每一亮度i的二值图像进行腐蚀。

进一步地，所述步骤(3)计算稻穗角度检测区域的稻穗弯曲角度的具体实现方式为：

(31)对所述稻穗角度检测区域作连通域标记；

(32)计算第k个稻穗的弯曲角度 $\mathrm{\θ}\left(k\right)=\arctan \left(\frac{\mathrm{RegionLength}\left(k\right)}{\mathrm{Regionwidth}\left(k\right)}\right),$ k＝1，2，…，regionum，RegionLength(k)和Regionwidth(k)分别为第k个连通域的长和宽，regionum为连通域个数；

(33)计算稻穗角度检测区域的稻穗弯曲角度RiceAngle₁＝mean(θ(k))。

进一步地，所述步骤(1)还对得到的稻穗区域作去噪处理。

本发明的技术效果体现在：本发明自动对所采集的实时前下视水稻田间图像进行特征提取，并利用颜色信息对相机场景中的稻穗进行图像分割，对稻穗分割图像的特定区域，计算稻穗的弯曲角度信息，进而判断该块稻田中的水稻是否进入乳熟期。该方法以表征水稻生长状况的重要形态参数作为判断依据，实时地对水稻生长图像进行计算，检测结果准确率高，对乳熟期以及后续的相关农事活动具有重要的指导意义。

附图说明

图1是训练阶段的流程图；

图2是水稻稻穗图像分割的流程图；

图3是自动检测水稻抽穗期的流程图；

图4是待检测序列中随机挑选的一张图像；

图5是图4中图像的最终分割结果；

图6是图5中稻穗连通域标记后的结果，图中的白色方框部分为进行水稻稻穗角度计算的检索区域；

图7是一个特定稻穗角度检测的示意图，角度θ为此稻穗的弯曲角度；

图8是2010年晚稻相机未乳熟至乳熟期区间内稻穗角度原始数据曲线，横轴为时间轴，纵轴为稻穗弯曲角度，图中标记出位置为水稻进入乳熟期对应时间；

图9是2011年早稻未乳熟至乳熟期区间内稻穗角度原始数据曲线，横轴为时间轴，纵轴为稻穗弯曲角度，图中标记出位置为水稻进入乳熟期对应时间；

图10是2011年晚稻未乳熟至乳熟期区间内稻穗角度原始数据曲线，横轴为时间轴，纵轴为稻穗弯曲角度。

具体实施方式

下面结合附图来详细说明本发明一实施例。

实施例中使用图像序列采集与离地面高5米的相机，镜头焦距为14毫米，水平视场角为46度，垂直视场角为35度，相机分辨率不低于400万像素。实施例以每一天为一检测时段，每一检测时段内拍摄w张水稻图像(w＝5)。每天为一检测阶段，有利于识别水稻的主要关键生长期。此发明旨在自动检测水稻抽穗期。整个方法分为训练阶段、分割阶段和检测阶段。

1.训练阶段：通过前一年拍摄的稻田水稻的历史图像数据，统计图像中稻穗部分在L-a-b颜色空间的颜色信息，其流程如图1所示，具体步骤如下：

(1)生成稻穗样本图像。

从往年历史图像数据中，随机选取各种光照条件下拍摄的包含稻穗的水稻图像80幅左右，手工扣取上述图像中包含稻穗的图像部分图，部分图大小为200×200像素(图像块大小不需要固定)，将上述稻穗图像块转换到L-a-b颜色空间，在L-a-b颜色空间下利用k-means聚类将部分图分为5类，每一类对应像素位置保留部分图的颜色信息，其它像素各通道颜色值置为0。然后手工选出仅包含有稻穗的图像块作为样本图像集。

上述图像部分图可以扣取仅包含稻穗的区域，也可扣取包含稻穗和水稻叶的区域，因为稻穗所占像素小，因此扣取仅包含稻穗的区域难度比较大，优选后者。本步骤也不局限于聚类分割方法，还可使用图割、分水岭和区域生长等分割方法。

(2)获取稻穗样本图像集中所有图像的稻穗像素点在L-a-b颜色空间的色彩分量(L，a，b)，其中L表示亮度，a表示红或绿的对应值，b表示黄或蓝的对应值，其中L，a，b为8位无符号整形。

设稻穗样本图像集为S_m，m＝1…150，图像中的第i行与第j列像素点相对应的在L-a-b颜色空间的色彩分量L(i，j)、a(i，j)、b(i，j)和在R-G-B颜色空间的色彩分量R(i，j)、G(i，j)、B(i，j)，依次从每张样本图像的每个像素点查找，记num＝1，若R(i，j)≠0(背景区域为纯黑色)，则获取此水稻稻穗像素点特征为[a_num，b_num，L_num]，且令num＝num+1。依次查找直到终止，将最后的num-1记为TotalNum。于是上述过程得到一个TotalNum×3的水稻稻穗颜色信息矩阵F^TotalNum×3。

(3)统计在不同亮度i＝1，…，255下(a，b)颜色值相同的稻穗像素点个数。

为充分利用水稻稻穗在各种光照条件下的颜色特征，将F^TotalNum×3矩阵按行随机排列，再截取其中的前RowNum行特征，记为f^RowNum×3。定义z^{255×255×255}为全零矩阵，依次从i＝1，2，…，RowNum，按照f^RowNum×3(也即[a_RowNum，b_RowNum，L_RowNum])行的顺序，令z(f(i，1)，f(i，2)，f(i，3))＝z(f(i，1)，f(i，2)，f(i，3))+1。于是获得稻穗颜色信息的分布矩阵z，下面对矩阵z进一步处理。RowNum的取值区间为[10000，50000]，此处我们选择RowNum＝50000，获得了较好的结果。

(4)在第i＝1，…，255亮度下，按照(a，b)颜色值相同的稻穗像素点个数越多，则稻穗像素点在亮度i下表现为(a，b)颜色值的可能性越大的原则选取稻穗(a，b)颜色的标定值。

在像素点亮度为L＝i，i＝1，…，255的情况下，将z(a，b，i)矩阵转化为向量形式，并将其按照元素由大到小的顺序进行排列得到向量mapLine，并计算mapLine的向量和mapLineSum，向量长度mapLineLength，依次从j＝1，…，mapLineLength，令：

mapLineSumTemp＝mapLineSumTemp+mapLine(j)；

直至mapLineSumTemp＞0.95×mapLineSum，记此时的j为TagLength，再依次从j＝1，…，TagLength，令z(s，t，i)＝mapLine(j)，由上式求得的(s，t)即为水稻稻穗在亮度i下(a，b)颜色的标定值。

(5)定义尺寸为255×255×255的三维ColorMap全零矩阵，第一维表示a颜色值，第二维表示b颜色值，第三维表示L亮度，令步骤(04)在亮度i下选取的稻穗(a，b)颜色的标定值表示为(s，t，i)，将ColorMap全零矩阵中的(s，t，i)元素值ColorMap(s，t，i)赋值为1，从而得到的ColorMap中，第i层表示在亮度i下表征稻穗颜色信息的二值图像。

(6)膨胀和腐蚀处理。

在样本数量有限情况下，ColorMap中本来应该标记为1的元素位置未能进行标记仍然为0，这就需要对ColorMap中的每一层进行后续处理，以便填补稻穗样本有限进而带来的“空洞”。

ColorMap中的每层为255×255的二值图像，定义结构元素大小为5的圆形结构元素，对ColorMap中的每层二值图像进行依次膨胀(van den Boomgard，R，andR.van Balen，″Methods for Fast Morphological Image Transforms Using Bitmapped Images，″Computer Vision，Graphics，and Image Processing：Graphical Models and Image Processing，Vol.54，Number 3，pp.254-258，May 1992.)，然后定义结构元素大小为3的圆形结构元素，对ColorMap中的每层二值图像进行依次腐蚀(van den Boomgard，R，and R.van Balen，″Methods for Fast Morphological Image Transforms Using Bitmapped Images，″Computer Vision，Graphics，and Image Processing：Graphical Models and Image Processing，Vol.54，Number 3，pp.254-258，May 1992.)，于是最终我们获取稻穗的颜色信息统计表ColorMap。

2.分割阶段，利用训练阶段生成稻穗的颜色信息统计表ColorMap，对待检测实时前下视稻田水稻图像进行水稻稻穗分割，流程图如图2所示，具体操作步骤如下：

(1)利用生成的稻穗颜色信息统计表，对实时前下视稻田图像进行水稻稻穗区域进行粗选。将待分割图像序列记为S_t，t＝1，2，…n，文中所有下标t均代表待测图像序列中的每张图像S。将图像S转化到L-a-b颜色空间，此时图像中的第p行与第q列像素点相对应的在L-a-b颜色空间的色彩分量L_t(p，q)、a_t(p，q)、b_t(p，q)，定义行列数与S相同的全零矩阵逻辑矩阵SL_t。

若

ColorMap(a_t(p，q)，b_t(p，q)，L_t(p，q))＝1

则此像素点为水稻稻穗对应像素点，将SL_t此位置标记为1。由此得到粗选后水稻稻穗区域的分割结果逻辑矩阵SL_t。

作为优化，对稻穗区域进行去噪处理。利用标记算法，对二值图像SL_t进行8邻域连通域标记，第k个连通域记为SL_t(k)，连通域的面积即像素个数记为num_t(k)，根据如下的公式：

去除较小的连通域以降低噪声对抽穗分割带来的干扰。面积阈值Threshold根据需要进行选择，一般为[5，20]。一般可取Threshold为6。操作完成后得到最后的SL_t。查找SL_t中(i，j)位置为0的点，取S中此位置像素点为黑色，于是我们得到图像D_t中水稻稻穗的最后分割结果DResult_t，图像序列中随机选取的一张图像如图4，则它的分割结果如图5所示。

3.检测阶段，由分割阶段我们获得二值图像SL_t，使水稻稻穗对应像素点值为1，背景像素点值为0。检测阶段分为水稻稻穗角度检测区域截取，稻穗弯曲角度计算与水稻乳熟期综合判断两个部分，具体操作步骤如下：

(1)从稻穗区域截取最能表征稻穗完整形态结构的子区域作为稻穗角度检测区域。

对于稻穗分割结果SL_t，t＝1，2，…n我们依次截取SL_t中，起点坐标为(180，20)，长为440像素，宽为120像素的稻穗弯曲角度检索区域，检索区域如图6中的白色方框所示区域，选取次检测区域可以尽可能获取稻穗的完整形态结构，以便于稻穗弯曲角度的计算，将上述获取的图像块序列定义为SSL_t，t＝1，2，…n；

(2)稻穗弯曲角度计算。

稻穗弯曲角度计算可采用最小二乘法曲线拟合、Ransac曲线拟合等方法，本实施例给出另外一种较为简单的计算方式：

对于上述获取的图像序列SSL_t，t＝1，2，…，n。SSL_t为二值图像序列，稻穗对应像素值为1，非稻穗区域对应像素值为0，依次检索SSL_t中的每一个稻穗连通域，设在SSL₁中检测到有regionum个连通域，每个连通域的长宽分别为RegionLength(k)，Regionwidth(k)，k＝1，2，…，regionum，于是定义每个稻穗的弯曲角度为θ(k)，单个稻穗的弯曲角度如图7所示，有：

$\mathrm{\θ}\left(k\right)=\arctan \left(\frac{\mathrm{RegionLength}\left(k\right)}{\mathrm{Regionwidth}\left(k\right)}\right),k=\mathrm{1,2},...,\mathrm{regionum};$

定义RiceAngle₁＝mean(θ(k)) k＝1，2，…，regionum，即以θ(k)的均值为SSL₁图像块检测到的稻穗弯曲角度。对序列图中的其它图像依次进行上述稻穗角度计算，于是可以得到表征稻穗弯曲角度的数据序列RiceAngle_t，t＝1，2，…，n；

(3)水稻乳熟期判断，根据2010年晚稻与2011年早稻的两个相机检测到的稻穗弯曲角度历史数据，以及人工对于水稻乳熟期的观测数据有：2010年晚稻相机检测到稻穗弯曲角度为52.26°时进入乳熟期，如图8所示；2011年早稻相机检测到稻穗弯曲角度为52.01°时进入乳熟期，如图9所示。于是给出当连续检测到稻穗弯曲角度低于53°累计达到五次时，判断水稻进入乳熟期。利用累计次数可以防止个别由于风等带来的干扰。

水稻乳熟期判断的具体步骤如下：

对于待检测图像序列SSL_t，t＝1，2，…，n处理后获取的稻穗弯曲角度的数据RiceAngle_t，t＝1，2，…，n，按照时间顺序依次判断RiceAngle_t，定义logical＝0，依次从t＝1，2，…，n，若RiceAngle_t＜53，则logical＝logical+1；直到logical＝＝5，则判读此时水稻进入乳熟期。

上述检测方法可行有一个的前提条件是在水稻乳熟检测过程中稻田没有遇到持续的大风天气。如果水稻乳熟检测过程中稻田遇到持续的大风天气，则上述算法会在可能发生乳熟的日期之前便检测到稻穗弯曲角度持续低于53°(如图10所示)，上述条件将造成算法失效。根据湖北、浙江、江西、湖南和广西五省于1992年至2011年二十年早稻与晚稻生长期的统计分析有：早稻平均于抽穗期12天后进入乳熟期，晚稻于抽穗期14天后进入乳熟期(如表1所示)。于是如果算法会在可能发生乳熟的日期之前(抽穗期9天后)便检测到稻穗弯曲角度持续低于53°，则判断遇到持续大风等天气现象的干扰，根据人工历史统计数据以及计算机自动检测出的水稻抽穗期日期，对水稻乳熟期时间进行估测。综上所述，水稻乳熟期判断的具体流程如图3所示。

表1 湖北、浙江、江西、湖南和广西五省于1992年至2011年早稻与晚稻二十年抽穗期至乳熟期人工观测历史数据表

Claims (5)

1.一种水稻乳熟期的自动检测方法，首先根据水稻历史图像序列离线训练得到稻穗像素点的水稻稻穗颜色信息表，水稻稻穗颜色信息表包含信息有稻穗像素点在每一亮度i＝1，…，255下对应的(a，b)颜色值的分布情况，然后按照如下方式在线判断水稻是否进入水稻乳熟期：

(1)将当前水稻图像t转化到L-a-b颜色空间，将图像t中各像素点的L-a-b颜色空间色彩分量与所述水稻稻穗颜色信息表进行比对，确定满足颜色值分布情况的色彩分量对应的像素点为稻穗像素点；

(2)从由稻穗像素点构成的稻穗区域截取最能表征稻穗完整形态结构的子区域作为稻穗角度检测区域；

(3)计算稻穗角度检测区域的稻穗弯曲角度；

(4)依据稻穗弯曲角度判断水稻是否进入乳熟期。

2.根据权利要求1所述的水稻乳熟期的自动检测方法，其特征在于，所述根据水稻穗像历史图像离线训练得到稻穗像素点的水稻稻穗颜色信息表的具体实现方式为：

(01)根据水稻穗像历史图像序列生成稻穗样本图像集；

(02)获取稻穗样本图像集中所有图像的稻穗像素点在L-a-b颜色空间的色彩分量(L，a，b)，其中L表示亮度，a表示红或绿的对应值，b表示黄或蓝的对应值，其中L，a，b为8位无符号整形；

(03)统计在每一亮度i＝1，…，255的(a，b)颜色值相同的稻穗像素点个数；

(04)在每一亮度i＝1，…，255下，按照(a，b)颜色值相同的稻穗像素点个数越多，则稻穗像素点在亮度i下表现为(a，b)颜色值的可能性越大的原则选取稻穗(a，b)颜色的标定值；

(05)定义三维ColorMap全零矩阵，第一维表示a颜色值，第二维表示b颜色值，第三维表示L亮度，令步骤(04)在亮度i下选取的稻穗(a，b)颜色标定值表示为(s，t，i)，将ColorMap全零矩阵中的(s，t，i)元素值ColorMap(s，t，i)赋值为1，从而得到在每一亮度i下表征稻穗颜色的二值图像。

(06)对ColorMap中每一亮度i下的二值图像依次进行膨胀和腐蚀，得到水稻稻穗颜色信息表。

3.根据权利要求2所述的水稻乳熟期的自动检测方法，其特征在于，所述膨胀和腐蚀的具体实现方式为：采用结构元素大小为5的圆形结构元素对ColorMap中每一亮度i下的二值图像进行膨胀，再采用结构元素大小为3的圆形结构元素对ColorMap中每一亮度i的二值图像进行腐蚀。

4.根据权利要求1所述的水稻乳熟期的自动检测方法，其特征在于，所述步骤(3)计算稻穗角度检测区域的稻穗弯曲角度的具体实现方式为：

(31)对所述稻穗角度检测区域作连通域标记；

(32)计算第k个稻穗的弯曲角度 $\mathrm{\θ}\left(k\right)=\arctan \left(\frac{\mathrm{RegionLength}\left(k\right)}{\mathrm{Regionwidth}\left(k\right)}\right),$ k＝1，2，…，regionum，RegionLength(k)和Regionwidth(k)分别为第k个连通域的长和宽，regionum为连通域个数；

(33)计算稻穗角度检测区域的稻穗弯曲角度RiceAngle₁＝mean(θ(k))。

5.根据权利要求1所述的水稻乳熟期的自动检测方法，其特征在于，所述步骤(1)还对得到的稻穗区域作去噪处理。

Images