2023年亚太杯APMCM数学建模大赛
A题 水果采摘机器人的图像识别
原题再现
中国是世界上最大的苹果生产国,年产量约3500万吨。同时,中国也是世界上最大的苹果出口国,世界上每两个苹果中就有一个是中国出口的,世界上超过六分之一的苹果是中国出口的。中国提出了“一带一路”倡议,这是建设未来共享的国际社会的重要支柱。由于这一举措,越南、孟加拉国、菲律宾、印度尼西亚等沿线国家已成为中国苹果的主要出口目的地。
苹果采摘主要依靠手工采摘。苹果成熟后,几天内就需要大量的采摘工人。但大多数当地农民在自己的果园里种植苹果。此外,农业工人的老龄化和年轻人离开村庄工作的现象导致了苹果采摘季节的劳动力短缺。为了解决这一问题,中国从2011年左右开始研究能够摘苹果的机器人,并取得了重大进展。
然而,由于果园环境与受控实验环境的差异,各种苹果采摘机器人在世界范围内的推广应用还不够理想。在复杂、非结构化的果园环境中,现有的机器人大多无法准确识别障碍物,如“树叶遮挡”、“树枝遮挡”、“水果遮挡”、“混合遮挡”等,如果直接采摘苹果而不根据实际场景做出准确判断,则损坏水果的风险很高,甚至对采摘手和机械臂造成伤害。这对收获效率和果实质量产生不利影响,导致更大的损失。此外,对不同采收果实的识别和分类也非常重要,如分类、加工、包装、运输等程序。然而,许多水果的颜色、形状和大小与苹果非常相似,这给苹果的采后识别带来了很大的困难。
本次竞赛的目的是通过分析和提取标记水果图像的特征,建立一个识别率高、速度快、准确度高的苹果图像识别模型,并对图像进行数据分析,如自动计算图像中苹果的数量、位置、成熟度,估计图像中苹果的质量等。具体任务如下:
问题1:数苹果
基于附件1中提供的已采收苹果图像数据集,提取图像特征,建立数学模型,统计每个图像中的苹果数,绘制附件1中所有苹果分布的直方图。
问题2:估计苹果的位置
根据附件1中提供的已收获苹果的图像数据集,以图像左下角为坐标原点,识别每个图像中苹果的位置,并绘制附件1中所有苹果几何坐标的二维散点图。
问题3:估计苹果的成熟状态
根据附件1提供的已采收苹果图像数据集,建立数学模型,计算每幅图像中苹果的成熟度,并绘制附件1中所有苹果成熟度分布的直方图。
问题4:估计苹果的质量
根据附件1中提供的已收获苹果的图像数据集,以图像左下角为坐标原点,计算每幅图像中苹果的二维面积,估计苹果的质量,并绘制附件1中所有苹果的质量分布直方图。
问题5:苹果的识别
基于附件2提供的采集水果图像数据集,提取图像特征,训练苹果识别模型,识别附件3中的苹果,并绘制附件3中所有苹果图像ID号的分布直方图。
整体求解过程概述(摘要)
为了应对苹果采摘机器人识别率低、速度慢、图像识别精度低等问题。利用改进的Hough圆变换方法建立了苹果识别的第一类模型。为了在复杂环境下更准确地识别苹果,通过CNN卷积网络提取苹果的特征值信息,对提供的苹果图像进行深度学习,建立了基于YOLOv4和VGG19的第二类苹果识别模型。另外,在计算图像中苹果的质量时,我们自己发明了一种亮度估计方法,使计算出的苹果质量更符合实际情况。
问题1要求我们识别图像中的苹果数量。采用Lab颜色空间变换和高斯滤波方法对图像进行预处理。利用成熟苹果为红色的特点,在Lab颜色空间中绘制a*通道图并进行二值化处理,然后利用苹果可以近似视为二维平面上的圆的特点,采用改进的Hough圆变换方法,得到附录1中确定的苹果数为3094个。
问题2和3要求我们确定图像中苹果的位置和成熟度。利用CNN卷积网络中的YOLOv4网络进行深度学习,输入人工帧选择和识别的图像数据进行模型训练,并用训练后的模型进行识别,发现苹果位置的结果与第一个模型的结果相似,表明该模型训练良好。在计算问题3中的成熟度时,我们还定义了成熟度函数H来量化苹果的成熟度。
问题4要求我们计算苹果的质量。进一步考虑了“大近、小远”的成像原理,分析了焦距和苹果与镜头的距离对成像结果的影响,认为如果仅仅根据苹果的二维图像面积来计算苹果的质量,可能会与实际情况有很大的不同。在此基础上,提出了一种亮度估计的校正方法,并利用前三个问题得到的数据结果计算校正后的苹果品质。
问题5要求我们鉴别不同水果中的苹果。我们注意到yolov4是一种目标检测模型,更偏向于图像分类,附录2中的样本数量过多,人工划分训练集的效率较低。因此,我们使用具有更多卷积层的VGG19网络来训练模型。与YOLOv4相比,它具有更高的效率和学习深度。
本文的创新之处在于运用几何方法、深度学习方法,从不同角度建立模型,并对所得结果进行分析比较,使模型具有较高的可信度。同时,我们也考虑了许多实际问题,这些问题中可能没有考虑到,并给出了我们的解决办法。
模型假设:
假设1
在图像的背景中往往有许多小而多的苹果,难以计数。为了保证图像识别的速度和准确性,我们将忽略图像中的小苹果。
说明:
这一假设是合理和充分的。对于假设的情况,在现实生活中,可以通过走近缩小距离、拍摄局部图像、改变拍摄角度来更好地识别苹果,而不是将其强加在背景中数量众多、难以计数的苹果上。检测,不能保证识别的速度和准确性。
假设2
在果实成熟期,果实体积和密度基本保持不变。如果忽略微小误差,则可将其视为固定值。
说明:
一般来说,果实在结果期发育,整体形状不断扩大。当果实形状基本稳定时,开始进入果实成熟期。果实成熟期的外在表现主要是果皮颜色的变化。其形状和质量基本不变,其密度也可视为不变。
假设3
所有的苹果都是红苹果。青苹果被认为是未成熟的红苹果,我们认为它们的成熟度与苹果从绿到红的梯度成线性关系。
说明:
青苹果和未熟红苹果的外观非常相似,由于我们所能获得的图像数据只是像素较低的二维图像,因此无法准确区分两者。考虑到大多数图像数据类型都是红苹果,我们可以认为所有的苹果类型都是红苹果,而苹果的成熟度只能通过颜色来区分,因此假设苹果的颜色变化和成熟度变化是线性的,相关的变化使我们更容易建立模型,与实际情况的差距很小。
假设4
在给定的图像数据中,所有照片都是在均匀柔和的自然光下拍摄的。拍摄对象均为静物,拍摄过程中未产生鬼影。所用拍摄设备的焦距相同。
说明:
采摘苹果的时间一般在白天,因此不需要另外的强光源,强光源会对成像产生很大影响,影响识别结果。我们认为光是均匀柔和的自然光,容易消除干扰因素。在正常情况下,受外界环境的影响,苹果的位置变化不大。至少可以认为拍摄时苹果是静止的,消除了成像过程中鬼影图像对识别结果的干扰。由于采摘机器人的成像设备是相同的,在实际应用中可以认为所有照片的拍摄焦距是相同的。
问题分析:
问题一分析
在问题1中,基于附件1中提供的苹果图像数据,我们旨在识别能够帮助我们识别苹果相关信息的特征。据观察,成熟的苹果表现出红色和具有近似球形外观(在二维图像中可以被视为圆形)的特征。基于这两个特征,可以对图像中的苹果进行识别。在建立模型之前,需要对图像数据进行预处理。
在这里,我们将原始RGB图像转换为Lab图像。然后,在Lab颜色空间中对a*分量(代表由绿色到红色的光谱变化)进行滤波,得到a通道图像。我们对这幅A通道图像进行二值化,并应用高斯滤波进行降噪。然后,利用圆形特征,采用改进的Hough圆变换对图像中的圆形进行检测和计数。此外,我们合并高重叠的圆。最后,我们对附件1中的数据进行求解并绘制直方图。
问题二分析
在问题2中,我们需要识别每个图像中苹果的位置。通过对问题1中的所有图像进行预处理,得到每幅图像中的苹果数,利用问题1中建立的模型对苹果进行编号后得到的质心,可以粗略地确定苹果的位置。但经过分析比较,发现苹果的几何中心与由Hough圆得到的质心仍有一定的偏差。
因此,我们采用了一种新的基于深度学习的识别方法,利用YOLOv4网络确定每个苹果的坐标。我们首先使用特征清晰的图像训练模型,然后将训练后的模型应用于附件1中的所有图像。最后,将识别结果与训练结果和问题1的结果进行比较,对模型进行灵敏度分析,绘制出所有苹果的几何坐标散点图。
问题三分析
在问题3中,我们可以遵循问题2中使用的方法,利用YOLOv4网络进行模型训练。但是,在训练过程中,需要用苹果的成熟状态(成熟、半成熟和未成熟)对手动添加的绑定框进行注释。然后对模型进行训练,利用训练后的模型识别图像中每个苹果的成熟状态。我们获得有关成熟状态和每个识别的苹果置信水平的信息。基于这些信息,我们可以设计一个函数来测量苹果的成熟度,从而计算出每个苹果的成熟度。最后,生成描述所有苹果成熟度分布的直方图以供进一步分析。
问题四分析
在问题4中,我们要解决的问题是计算二维图像的面积和估计苹果的质量。对于面积问题,我们不妨使用我们在第一个问题中获得的半径rad来计算它作为一个简单的估计。正如我们在假设中所说,我们假设苹果的密度在发育周期中没有显著变化,这是一个固定值,可以记录为ρ。此外,我们的团队注意到附录1中的图像存在明显的差异,即“大-近-远-小”问题导致的误差结果是,对于图像中的某些苹果,实际苹果较大,而图像中的苹果较小。为此,我们的研究小组查阅了通过二维平面图像计算物体实际距离的相关文献。遗憾的是,文献中回顾的方法不适用于我们团队面临的困境。经过讨论,最后提出了一种新的近似估计苹果实际大小的方法——亮度估计法。在这种方法中,我们需要知道拍摄时相机的焦距,最后绘制出所有苹果的质量分布直方图。
问题五分析
这个问题要求我们在图像中区分苹果和其他水果,这就要求我们提取更多关于苹果的信息特征,而不是像前面的问题那样计算苹果的位置、成熟度、质量等。如果继续使用YOLOv4网络,当训练次数较多时,时间复杂度会过高。因此,我们应该选择一种卷积层数更多、计算速度更快的深度学习算法,以提高算法的精度。有效识别苹果。同时,考虑到苹果ID号过多,我们选择每1000个连续的ID进行聚类,统计1000个ID内的苹果图像数,然后绘制直方图。
模型的建立与求解整体论文缩略图
全部论文及程序请见下方“ 只会建模 QQ名片” 点击QQ名片即可
部分程序代码:
clear
clc
% Generate file path, 200 is the number of images in the folder Fpath="xxx"; % Fill in the file path at xxx
Fspath="xxx"; % Fill in the file path at xxx
for i=1:200
path(i,1)=Fpath+num2str(i)+".jpg";
spath(i,1)=Fspath+num2str(i)+".jpg";
end
% Counting matrix, circle center matrix and radius matrix count=[];
center=[];
rad=[];
for num =1:200
% Get image data
data = imread(path(num,1));
% Convert to lab space
ldata = rgb2lab(data);
% only a* component is required
gdata = ldata(:,:,2);
% Set the number of filters
n=4;
% Start filtering
for M=1:n
mean_gdata=mean(gdata(:));
var_gdata=var(gdata(:));
for i=1:size(gdata,1)
for j=1:size(gdata,2)
if gdata(i,j)<mean_gdata+sqrt(var_gdata)
gdata(i,j)=0;
end
end
end
end
% Calculate the threshold k of Gaussian filter
k = graythresh(gdata);
% Gaussian filter
W = fspecial('gaussian',[5,5],2/3);
G = imfilter(gdata, W, 'replicate');
% OTSU split
bw=im2bw(G,k);
% Hough circle transformation
[centers,radii] = imfindcircles(bw,[10 70],"ObjectPolarity","bright",
"Sensitivity",0.88,"EdgeThreshold",0.10,"Method","twostage");
% draw circle
%imshow(data)
% viscircles(centers,radii);
% Generate distance matrix
for i=1:size(centers,1)
for j=i:size(centers,1)
d(i,j)=sqrt((centers(i,1)-centers(j,1))^2+(centers(i,2)-centers(j,2))^2);
end
end
% Generate position vector
p=[];
row=1;
% Marker for overlapping circle position
for i=1:size(centers,1)
for j=i+1:size(centers,1)
if max(radii(i),radii(j))>d(i,j)+0.4*min(radii(i),radii(j)) && d(i,j)~=Inf if radii(i)<radii(j)
d(i,:)=Inf;
d(:,i)=Inf;
p(row,1)=i;
row=row+1;
elseif radii(j)<radii(i) && d(i,j)~=Inf
d(j,:)=Inf;
d(:,j)=Inf;
p(row,1)=j;
row=row+1;
end
end
end
end
% Sort the position vectors
p=sort(p);
% Eliminate overlapping circles
for i=row-1:-1:1
centers(p(i,1),:)=[];
radii(p(i,1),:)=[];
end
% generate image
imshow(data)
viscircles(centers,radii)
f=gcf;
exportgraphics(gcf,spath(num,1))
close
% count
count=size(radii,1);
center{num,1}=centers;
rad{num,1}=radii;
end
% The code of Problem 4
clear
clc
k=1;
rho=0.75;
f=0.6;
load('bright.mat');
load('center.mat')
load('rad.mat')
for i =1:200
centers=round(center{i,1});
for j = 1:size(rad{i,1},1)
alpha=ones(185,1);beta=ones(1,270);
alpha(max(centers(j,1)-rad{i,1}(j,1),1):min(centers(j,1)+rad{i,1}(j,1),185),1)=1; beta(1,max(centers(j,2)-rad{i,1}(j,1),1):min(centers(j,2)+rad{i,1}(j,1),270))=1; R=alpha*beta;
lin=bright{i};
H{i,j}=sum(sum(lin.*R))/(4*rad{i,1}(j,1)^2);
m(k,1)=pi*rho/(6*f^3)*(rad{i,1}(j,1)^3)*sqrt(pi^3*H{i,j}^3/827^3);
k=k+1;
end
end
k=1;
for i=1:200
for j=1:size(rad{i,1},1)
stand(k,1)=log(H{i,j}^3/rad{i,1}(j,1)^6);
stand(k,2)=i;
stand(k,3)=j;
k=k+1;
end
end