CN109856064A - 基于光子传输模拟的苹果高光谱品质检测方法 - Google Patents

Abstract

本发明公开了一种基于光子传输模拟的苹果高光谱品质检测方法。本发明分析了光子入射最佳位置和源探距离，并用点光源高光谱仪实际拍摄红富士苹果进行验证。分析表明，光子在苹果赤道位置入射，具有73.12％概率到达更深的深度。源探位置与苹果的光学参数有关，形状为圆环，其源探内外径位置半径为1.5mm‑10.15mm。点光源高光谱仪采集的红富士苹果入射位置为赤道，源探位置为半径为3.6mm‑10.8mm的圆环，与模拟数据分析结果基本一致。蒙特卡洛光子传输模拟方法为研究高光谱苹果品质无损检测开辟了新思路，分析结果可以为研究高光谱品质检测实验设计和苹果便携式品质检测光学仪器设计提供理论基础。

Description

基于光子传输模拟的苹果高光谱品质检测方法

技术领域

本发明涉及检测领域，具体是一种基于光子传输模拟的苹果高光谱品质检测方法。

背景技术

中国苹果资源丰富，年产量占世界总产量的1/3，居世界第一。高光谱图像是将二维图像技术和光谱技术融合在一起的先进技术，图像反映水果的外部特征，光谱检测水果的物理结构、化学成分等，而组织结构、物理特征以及化学成分等决定了水果的光学参数。很多学者研究了利用高光谱或近红外光谱检测水果品质、成熟度等指标的可能性，通过改进实验过程和数据处理模型，获得了较好的检测效果。随着光谱仪设备的改进更新，便携式光谱仪使用灵活轻便，更适合现场作业，对便携式光谱仪的开发和将便携光谱仪应用到水果品质检测的研究受到越来越多学者的重视。目前，利用高光谱进行果蔬品质检测的相关研究思路为：利用光谱仪采集样本光谱数据，通过改进波段选择算法和拟合模型来提高拟合精度。

为了弥补以上几个缺点，提高光谱测量的准确性，拟利用模拟方法获得大量无噪声数据。光子在苹果组织中的运动属于随机过程，因此，可以利用蒙特卡洛方法对苹果中的光子运动状态进行模拟。蒙特卡洛方法由Ulam和Von Neumann首先提出。Wang(Wang L,Jacques S L,Zheng L.MCML--Monte Carlo modeling of light transport in multi-layered tissues.[M]//Optical-Thermal Response of Laser-IrradiatedTissue.Springer US,1995:131-46.)综合前辈的理论知识将其编写成一个完整的蒙特卡洛模拟MCML程序。在此之后，有学者将蒙特卡洛方法应用到生物组织模拟中，并针对研究对象的不同提出了不同的改进方法。在苹果光子传输模拟上，Qin(Qin J,Lu R.Monte Carlosimulation for quantification of light transport features in apples[J].Computers&Electronics in Agriculture,2009,68(1):44-51.)测量了苹果的高光谱图像和光学参数，并对光的传输规律进行了蒙特卡洛模拟，主要研究了影响光子入射深度的因素。Askoura M L(Askoura M L,Vaudelle F,L'Huillier J P.Numerical Study ofLight Transport in Apple Models Based on Monte Carlo Simulations[J].Photonics,2015,3(11):26.)0研究证明了苹果果皮组织对于模拟模型的重要性。石舒宁(水果组织的光学描述与光传输规律研究[D].华中农业大学,2016.)用蒙特卡洛方法模拟双层苹果模型在808nm波段下的光子传输过程，其苹果模型为双层标准圆模型。

利用光谱仪实际测量存在如下问题：1、高光谱仪采集光谱图像质量受采集人操作和采集环境影响，光源和探头摆放位置固定，影响采集数据的质量。2、普通光谱相机体积大，不便于携带，便携式光谱仪检测精度受限。3、实验使用的苹果样品数量不足，品种单一会导致分析数据覆盖面不全，建立的模型不具有普适性。4、实际测量环节耗费大量人力物力。

发明内容

本发明针对背景技术中存在的问题，提出了一种基于光子传输模拟的苹果高光谱品质检测方法。

技术方案：

本发明公开了一种基于光子传输模拟的苹果高光谱品质检测方法，它包括以下步骤：

(1)基于椭球曲面方程构建苹果模型；

(2)将点光源对准苹果赤道，源探距离选择1.5mm～10.15mm，基于苹果模型进行蒙特卡洛光子传输模拟；

(2.1)光子初始化；

(2.2)计算光子运动方向和步长，同时进行越界和出界判断；

(2.21a)越界判断，越界则进行步骤(2.21b)，否则进行步骤(2.22a)；

(2.21b)出界判断，出界则进行出界处理并转步骤(2.5)，否则进行步骤(2.22b)；

(2.22a)光子吸收与散射计算，转入步骤(2.3)；

(2.22b)光子折射、反射计算，转入步骤(2.3)；

(2.3)判断光子权值是否过小，是则进行步骤(2.4)，否则返回步骤(2.2)；

(2.4)判断光子是否存亡，是则进行步骤(2.5)，否则返回步骤(2.2)；

(2.5)判断是否为最后一个光子，是则结束，否则返回步骤(2.1)；

(3)取大量光子和不同光学参数组合重复步骤(2)，获得苹果模型的无噪声光亮度分布图。

具体的，步骤(1)中，基于椭球曲面方程构建苹果模型，苹果果皮和果肉的椭球曲面方程为：

其中：a、b、c，p1，p2为常数，共同控制苹果的长短轴半径和上下凹陷程度；a∈[3.1,5]，b∈[3.1,5]，c∈[4,6]，p₁∈[3,35]，p₂∈[5,670]；-π/2≤u≤π/2,0≤v≤2π；

苹果果核的椭球曲面方程为：

其中：a、b、c，p3为常数，共同控制苹果的长短轴半径和上下凹陷程度；k＝p₃sin5v，p3＝0.5；a∈[3.1,5]，b∈[3.1,5]，c∈[4,6]，-π/2≤u≤π/2,0≤v≤2π。

具体的，步骤(2.2)中，计算光子运动方向和步长：光子的运动方向由方位角ψ和散射偏转角θ确定，其中：

方位角ψ采用下式计算：

ψ＝2πξ

式中，ξ是区间(0，1)上的随机数；

散射偏转角θ采用下式计算：

式中，θ为散射偏转角，g为组织的各向异性，ξ是区间(0，1)上的随机数；

光子步长的计算公式为：

其中：

μ_t＝μ_a+μ_s

ξ是随机数，μ_a为物质的吸收系数，μ_s为物质的散射系数。

具体的，步骤(2.2)中，光子入射点的计算步骤为：

A、取越层运动位移的中点c，判断光子是否处于越层过程：若处于越层过程，则计算光子前一步位置坐标a和当前位置坐标b并进行步骤B；若未处于越层过程，则进入步骤C；

B、判断c点是否在模型交界点上，若是则返回该点的值，进入步骤E，否则进入步骤C；

C、判断a与c是否在一个位层：若是则将c赋给a；否则，将c赋给b；

D、判断b是否在交界点，否则返回A，是则进入步骤E；

E输出c。

具体的，光子位层的计算方法为：光子位于点P(x,y,z)，通过无误差转换计算得到参数(u,v)，将得到的(u,v)带入苹果果皮椭球曲面方程得到在该参数下果皮点p1(x1,y1,z1)，将p1和p点分别与原点的距离相比较：

-若p1点较远，则p点在果皮内部；

-若p1点较近，则p点在果皮之外；

-若p1于p点重合，则说明点在果皮上；

重复带入苹果果肉椭球曲面方程、苹果果核椭球曲面方程，得到果肉和果核位层的判断。

具体的，光子位层的计算方法为：光子位于点P(x,y,z)，通过无误差转换计算(极坐标系和直角坐标系)得到参数u，将得到的u带入苹果果皮椭球曲面方程得到在该参数下果皮点Z轴上的分量z1，通过对比z与z1的长短得到点P的位置：

-若z1较长，则p点在果皮内部；

-若z1较长，则p点在果皮之外；

-若z1与z相等，则说明点在果皮上；

重复带入苹果果肉椭球曲面方程、苹果果核椭球曲面方程，得到果肉和果核位层的判断。

具体的，所述无误差转换的计算公式为：

u＝tan^-1(az/(cx cos v))

v＝tan^-1(ay/bx)

式中，x、y、z分别为点P的直角坐标，a、b、c分别为系数，a∈[3.1,5]，b∈[3.1,5]，c∈[4,6]。

具体的，步骤(2.2)中，引入模型厚度，苹果模型包括：果皮外交界层、果皮层、果皮果肉交界层、果肉层、果肉果核交界层、果核层；光子的运动采取越层约束算法，算法步骤为：

A、计算光子移动前后的层位，若光子移动前层位a和移动后层位b不属于同一个层位则进行步骤B判断；

B、判断层位a与层位b是否为相邻层位：相邻则进入步骤F，不相邻则进入步骤C；

C、计算出光子运动方向；

D、计算出与光子运动方向同方向的、并且与光子位移前位层a相邻的位层，称为位层c；

E、暂认为光子由位层a位移至位层c并带入继续计算；

F、计算交界点，进入反射折射计算。

具体的，步骤(3)中，所述光学参数包括果皮吸收系数μ_a1、果肉吸收系数μ_a2、果皮散射系数μ_s1、果肉散射系数μ_s2。

具体的，果皮吸收系数μ_a1分为5类，果肉吸收系数μ_a2分为4类，果皮散射系数μ_s1分为5类，果肉散射系数μ_s2分为5类，各取中间值进行组合获得500组光学参数组合。

本发明的有益效果

本发明在标准算法的基础上对算法进行了并行和优化，使得算法运行速度大大提高，在分析模拟数据的基础上用实际光谱测量的方法验证结论，证明了模拟结果的正确性。模拟方法解决了实际测量中的误差问题，且模拟范围不受样品种类数量限制，降低了人力物力成本。源探位置的分析结果可以提高便携光谱仪的检测精度。为高光谱水果品质检测提供新的研究思路，也为便携式水果品质检测光谱仪的设计提供理论基础。具体的：

本发明将点光源设置在苹果赤道位置：光子有更大的可能性进入到更深的深度。

本发明将点光源对准苹果赤道，源探距离选择1.5mm～10.15mm，即检测时探头在距离入射光源1.5mm-10.15mm的位置检测最好。

通过蒙特卡洛模拟可以快速得到不同光学参数组合下的无噪声光亮度分布图，弥补了实际测量噪声大、采集数据覆盖面小、耗费人力物力的缺点。

由于苹果的果皮和果肉均存在厚度，标准蒙特卡洛模拟忽略了介质厚度这一点，本发明引入模型厚度，因此，光子入射点的确定不再是两条线的交点问题而是一条线与一个平面的交点求解问题。确定交点所在位置，可以用解析几何的方法联立方程得到精确解，但是解方程效率低，采用二分法的方法找到光子的入射点。一般循环两次即可找到入射位置。

蒙特卡洛具有随机性，由于模型具有一定厚度，因此，光子可以在3层组织中任意穿梭，也可能从果皮层溢出。当光子越层时，光子的运动状态受到两层的影响，而苹果的果皮层薄，故对越层做约束处理避免跨越几层的现象。

附图说明

图1为基于椭球曲面的苹果组织模型

图2为基于椭球曲面的苹果多层组织模型的蒙特卡洛模拟

图3为苹果组织多层模型蒙特卡洛模拟图

图4a为内径与几何、光学参数的相关度图

图4b为外径与几何、光学参数的相关度图

图5a为神经网络拟和内径大小测试图

图5b为神经网络拟和外径大小测试图

图6为源探距离与糖度水分相关度分布图

图7a为神经网络预测实测数据内径分布图

图7b为神经网络预测实测数据外径分布图

具体实施方式

下面结合实施例对本发明作进一步说明，但本发明的保护范围不限于此：

1基于多层苹果组织模型的蒙特卡洛光子传输模拟

1.1基于椭球曲面方程的苹果模型构建

三维苹果几何模型如图1所示。苹果为三层组织模型，三层分别为果皮层，果肉层和果核层。果皮层和果肉层在形状上近乎相似，都是上下略有凹陷的椭球形。苹果果皮和果肉的椭球曲面方程为公式1。

其中，苹果的果核层是一个类似杨桃的形状。构建的方程如公式1，其中k发生变化，k＝p₃sin5v。苹果的几何形状可以根据椭球曲面方程中的参数改变而调整。

1.2基于椭球曲面苹果模型的蒙特卡洛光子传输模拟

图2是基于椭球曲面苹果模型的蒙特卡洛模拟流程图。光子的运动状态主要分为初始化、光子运动方向和步长计算、越界判断、出界判断和光子消亡五个主要步骤。基于椭球曲面的苹果模型蒙特卡洛模拟与标准蒙特卡洛模拟的不同主要在于介质形状为不规则曲面，而且曲面的厚度是有限的，光子在越界时需要约束，发生折射和反射时需要求解该点在曲面上的法向量，再计算折射角和反射角。

1.2.1交界点计算

由于苹果的果皮和果肉均存在厚度，标准蒙特卡洛模拟忽略了介质厚度这一点，在这里引入模型厚度，因此，光子入射点的确定不再是两条线的交点问题而是一条线与一个平面的交点求解问题。确定交点所在位置，可以用解析几何的方法联立方程得到精确解，但是解方程效率低，采用二分法的方法找到光子的入射点。一般循环两次即可找到入射位置。算法可用伪代码表示如下：

1.2.2光子越层约束算法

蒙特卡洛具有随机性，在1.2.1中引入了模型厚度概念，光子可以在3层组织中任意穿梭，也可能从果皮层溢出。当光子越层时，光子的运动状态受到两层的影响，而苹果的果皮层薄，如不作约束容易发生一次步长过大，跨越几层的现象，这里需要对这种越层做约束处理。即光子越多层时，找到光子运动方向，强制光子只越一层，约束算法如下：

1.2.3曲面上光子的反射和折射

基于曲面的光线折射反射计算首先需要进行曲面上法向量求解。

(1)光线在曲面上的反射。假设入射向量为p，光子入射轨迹与曲面相交于M点，法向量则出射向量q计算可通过向量计算。假设所有向量都是单位向量，入射角余弦值应该是根据向量计算三角形法则：反射角与入射角相等，可得到该三角形为等腰三角形，所以系数k为两倍的的模与余弦值的乘积，在这里的模为1，系数为2cosα，即得到公式2。

(2)光线在曲面上的折射。射角与折射角分居于表面两侧。由于两个向量的线性运算得到的向量仍然属于这两个向量所在的平面认为存在系数k和t，使得折射向量满足等式3。通过确定系数k和t，即可得到折射向量因为只有线性计算，大大降低了计算成本。

同样的，三个向量同样存在向量的三角法则，为方便计算，将其表示在二维平面上，在这里认为为入射向量坐标为(cosα，sinα)，为出射向量坐标为(cosβ，sinβ)，法向量坐标(0，-1)，可以得到cosβ＝kcosα，通过向量关系可得到sinβ＝k·sinα-t。分解等式并带入计算可以得到结果如公式4,5。

k＝cosβ/cosα#(4)

t＝sinαcosβ/cosα-sinβ#(5)

再由二维平面推广到三维空间中。由折射定理可知，入射光纤，折射光线与法线三线共面，可以理解为只要确定一个平面，折射现象就是在二维平面上完成的。而使用入射向量和法向量线性加减求得的向量，一定与这两个向量共面。所以可以得到，二维平面中的公式适用于三维空间中折射向量的求解。

1.2.4算法并行化

若想在短时间内获得大量模拟结果，使用单核运行模拟程序是不现实的，因此，在实际模拟中使用了并行化技术。用于并行化计算的机器为超微Super Server Main ServerChassis。采用matlab软件的并行化工具箱进行并行计算。利用24万入射光子，一个模型一组光学参数组合进行测试，使用单核运行，运行时间为1927.588s，8核运行时间为249.3s，24核运行时间73.0s。由对比可以说明并行化可以大大提高运行效率。对比实际光谱测量，可以节省大量人力物力成本。

1.3蒙特卡洛模拟结果

随机挑选一个苹果模型和一组光学参数，设定光子数为十万，用蒙特卡洛方法进行模拟，可以获得了光子运动轨迹图如图3所示。

从图中可以看出，绝大多数的光子的运动轨迹都集中在入射点周围，并呈发散状，以入射点为圆心均匀的向四周散射。观察其垂直入射方向，可以看出大多数光子轨迹集中在苹果表层之下，越靠近苹果内部的果核，光子轨迹越少，只有极少部分的光子能到达果核，这与现实情况基本相符。通过蒙特卡洛模拟可以快速得到不同光学参数组合下的无噪声光亮度分布图，弥补了实际测量噪声大、采集数据覆盖面小、耗费人力物力的缺点。

2基于模拟图像的光子入射最佳位置及源探距离分析

基于高光谱图像的苹果品质检测研究中，存在两个问题，1、针对光子入射位置，研究者一般选择苹果赤道位置入射，原因在于苹果两端为果柄与部分果核，而人们更关注果肉，但光子入射位置选择仅凭经验，缺乏相关论证。2、模拟实验中，光亮度值随着距离光斑中心的加大而减小，在同样光照强度下，随着各参数的变化，光斑大小和光亮度也会发生变化，因此猜测存在一个区域，在该区域内光亮度变化明显，探头在该区域会探测到更多的信息，从而获得更好的检测结果。针对以上两个问题和相关猜测，利用模拟方法证明光子入射最佳位置、确定源探距离及影响源探距离的因素。

2.1基于模拟图像的光子入射最佳位置分析

光子进入苹果的深度越深，光子携带的果肉信息越多，越有利于苹果品质分析，光子的入射深度与光源放置位置有关，因此，找到合适的光源放置位置对于苹果的品质分析非常重要。

由实际测量得到了200个红富士苹果的光学参数范围，将果皮吸收系数μ_a1分为5类，果肉吸收系数μ_a2分为4类，果皮散射系数μ_s1分为5类，果肉散射系数μ_s2分为5类，取中间值进行组合，共500个参数组合。选取造型差异较大的38个几何模型，100万个入射光子，分别在纬度π/4，π/2，3π/4位置入射(为便于后续说明，分别标号为1，2，3)，用第一章中的蒙特卡洛方法进行计算。得到19000组光子亮度分布图。

入射深度的计算方法为：对每一列高亮度像素点进行统计，记录连续分布点的个数，取数值最高的值作为该列的深度。对于每一组数据，分别统计1，2，3号入射位置的深度，当某个位置的深度同时大于等于其余两个位置，则认为该入射点深度最深。几个入射深度对比时可能出现深度相近的情况，当较大值与较小值的差值与较大值之比小于2.5％则认为近似相等。对19000组数据的统计结果如表1。由表1可以得出，在纬度为π/2位置入射更可能到达更深的深度，光子中含有更多的果肉信息。

表1模拟最佳入射位置对比表

2.2基于模拟图像的源探距离分析

2.2.1有效光子比

对于苹果来说，需要探究的苹果品质信息更多集中在苹果果肉中，在这里可以认为，果皮是干扰层，果肉是有效层。因此，苹果的有效光子是指在果肉层中运动的路程多于在果皮中运动路程的光子。理论上说，光子在有效层内运动的距离越长，所携带的有效层光学信息就越多，也越能反映有效层的光学特性。一个光子在运动结束时同时达到以下条件就可被认为是有效光子：

1)从苹果模型表面溢出。

2)光子运动的所有路径中，走过果肉的路径大于走过果皮的路径(由于苹果的高散射性，光子很难到达果核)。

光斑上的光子溢出存在有效光子溢出Ph_e和无效光子溢出Ph_i两种，每个光子在溢出时所携带的光亮度I由光子在组织中行走的路径长短决定。自定义苹果组织表面一点溢出的有效光子的光亮度总和I_e与全部光子光亮度总和I_a的比率为有效光子比。如公式6。

理论上来说距离入射点越近，反射光越强，有效光子比越低；随着距离的增加，有效光子比会逐渐变大趋近于1，反映出的光学信息也随之增加。但由于光强度逐渐减弱，噪声越来越大，可以被探测到的有用信息也会逐渐变少。因此，在光源附近存在一个范围，探测器在该范围能探测到更多关于苹果果肉的信息，这个范围被称为源探位置。为了找到有效区域，分别使用蒙特卡洛模拟方法，通过分析模型相关度和光学参数相关度，寻找源探位置的分布规律。

由于苹果果皮较薄，认为在苹果表面的有效光子与总光子之比大于0.5，即有效光子数量大于无效光子则认为该点具备探测价值。同时，在模拟环境中当总光子权值和小于0.01则认为，光亮度过小，在实际应用中，探测效果会由于亮度太低，相对噪声较大导致探测效果不佳。在这里通过两个限制条件将有效探测区域限制在以光斑为圆心的圆环中，圆环内径由有效光子比限制，圆环外径由相对光子权值和限制。由以上规则可计算得，模拟数据的内径半径大小为1.5mm-3mm，外径半径大小为7.65mm-10.15mm，因此模拟数据的源探距离为1.5mm-10.15mm。

2.2.2基于模拟图像的源探距离与几何模型、光学参数相关性分析

随机选取了15个形状各异的苹果模型，分别对每个模型使用相同的三组光学参数模拟，选取实测范围内，差异较大的三组光学参数(μ_a1、μ_a2、μ_s1、μ_s2)，光学参数的数值如下所示。模拟光子量为50万个，统计圆环内径大小如图4。

光学参数数值分别为(按照μ_a1、μ_a2、μ_s1、μ_s2的顺序)：第一组：0.7mm^-1、0.5mm^-1、30mm^-1、12mm^-1；第二组：2.6mm^-1、5mm^-1、30mm^-1、56mm^-1；第三组：4.5mm^-1、5mm^-1、190mm^-1、56mm^-1。在图4中分别标注为1，2，3。

由统计图可以分析，所有模型在相同的光学参数组合下内径保持在一个较小的范围内浮动，而对于不同的光学参数，内径和外径差异十分明显。由此可以初步推论，源探位置的内径和外径与模型参数相关性十分微小；内径和外径大小对于光学参数比较敏感。

2.2.3基于模拟图像的源探距离与光学参数相关度分析

用模拟方法获得了645组数据，随机选取517条数据作为训练，128条数据作为测试。利用神经网络对四个光学参数和内径大小做拟合，选取的神经网络结构为输入层为4，隐藏层为20，输出层大小为1，训练迭代次数为1000次，学习率为0.0001。

对内径的拟合结果效果如图5a均方跟误差0.3287，相关系数为0.9968。用同样结构的神经网络对光学参数和外径大小作拟合，对外径的拟合结果效果如图5b，均方跟误差3.6203，相关系数为0.9888。

由两组相关系数可得，通过有效光子比计算的有效测量区域的内、外径大小与光学参数有密不可分的联系。在获得样本的一组光学参数后，利用神经网络进行计算，可以得到在模拟情况下的源探位置，这为实际测量源探位置的寻找提供一个良好的指引。

3基于高光谱图像的光子入射最佳位置及源探距离分析

第二章用模拟方法确定了光入射最佳位置和源探距离，本章用高光谱实际测量的方法验证第二章的结论。

3.1基于高光谱图像的光子入射最佳位置分析

3.1.1点光源高光谱图像采集

试验选用了上海五铃光电科技有限公司的HSI-VNIR-0001高光谱成像系统。软件由Spectral-image取像软件和HISAnalyzer分析软件组成。采集波段为373.54-1033.87nm。图像采集过程描述如下:将苹果对半切开，光纤光源放在经度π/4，纬度为高中低三个位置分别平行照射苹果。

3.1.2基于苹果点光源高光谱图像的光子入射最佳位置分析

在3.3.1中，分别采集了30个样本，光源分别位于高中低三个纬度入射的光谱数据，分别标号为1，2，3。选择450-900nm波段的图像进行分析，共计13530组数据。

由于光谱仪采集图像噪声较大，对光子入射深度的统计影响较大，因此，采用了如下处理过程：

1)使用SG平滑算法将高光谱原始数据信息平滑处理。将相对较低的光强信息直接作为噪声屏蔽掉。

2)循环对每一列高亮度像素点研究，计算连续分布点的个数，取连续数最高的值作为该列的深度。

3)为避免偶然性，对所有列分析后，得到的100个深度值按大小排序，取最大的前十五个值的平均值为该波长下入射深度。

深度比较标准与2.1节相同，统计结果如表2。

表2实测最佳入射位置对比表

3.2基于点光源高光谱图像的源探距离分析

3.2.1苹果点光源高光谱图像采集及品质参数采集

苹果的点光源采集过程采集了200个苹果样本的点光源高光谱图像，由于3.1.2已经做过证明，因此，光源入射位置选择在赤道。

糖度测量：本实验中用可溶性固形物代替糖度，可溶性固形物是指液体或流体食品中所有溶解于水的化合物的总称，苹果的可溶性固形物主要是糖。测量仪器是ATAGOPAL-1数显糖度计。每个样本测量三次，取三次平均值为该样本的糖度值。

水分测量：苹果水分的方法测定方法为直接干燥法。使用的设备有干燥箱、电子天平、蒸发皿等。其中干燥箱为上海精宏DHG-9070A，天平为奥豪斯CP423C电子天平。水分含量的测定方法：首先称量蒸发皿的质量为m1，放入重量范围在40-50g之间的苹果块，用保鲜膜封闭，测出重量m2，取下保鲜膜放入到干燥箱中进行干燥，每干燥3小时秤重一次，直到连续两次测得质量差值小于0.3毫克为止，此时重量记为m3，样本水分含量计算公式7。

3.2.2基于苹果点光源高光谱图像的源探距离分析

点光源照射苹果会在苹果表面形成一个光斑，以光斑中最亮点为中点，以步长为1个像素点，对两边区域求平均值作为感兴趣区域，每个样本选择15个感兴趣区域。对200个样本的15个感兴趣区域与其对应的糖度和水分做相关性分析，分析方法过程如下：(1)对光谱数据进行一阶导处理，消除噪声。(2)使用随机划分数据算法对数据集进行划分。(3)使用偏最小二乘法对品质数据进行拟合分析。分析结果如图6所示。

在15个感兴趣区域里面，糖度R值最高达到0.70,最低为0.54。水分R值最高达到0.75，最低为0.26。

分别选择与糖度和水分相关度最高的6个区域，分别为：与糖度相关区域为12,10,4,11,8,5。与水分相关区域为8,6,13,4,11,12。重合区域有：4,8,11,12。选择糖分与水分相关度最低的6个区域分别为：糖度：14,15,1,7,2水分：15,14,9,2,10。重合区域有2,14,15。由以上数据的对比和图像可发现，探头位置选择对品质检测是有影响的，要选择距离光斑中心点一定距离的位置进行检测，由于光子是均匀向四周分布的，所以可以认为最佳检测区域位置在距离光斑中心点为圆心，4-12个像素点为半径的一个环形区域内。

4综合分析

4.1光子入射最佳位置分析

由于光子进入苹果的深度越深，光子携带的果肉信息越多，因此，可以用光子进入苹果的深度来确定最佳光子入射位置，即最佳光源位置。根据2.1和3.1.2的结果分析表明，无论模拟或实测，光子从苹果赤道位置入射都有更大的概率到达更深的位置，携带更多的果肉信息，因此，光子入射的最佳位置为赤道。

4.2源探距离分析

模拟数据的源探距离为1.5mm-10.15mm，实测感兴趣区域的源探距离为4-12像素点，经换算即3.6mm-10.8mm，模拟范围基本覆盖实测范围，有部分偏差的原因是实测过程中会存在仪器测量误差和其他光线干扰等因素。将实际测量得到的光学参数输入到2.2.3节构建好的两个神经网络中，得到的内径外径分布预测图如图7。

对108240组数据内外径分别统计可以看出，内径分布比较密集，集中在3mm-3.46mm处，与实际测量结果基本相同；而外径部分相对较分散，大多分布在15.3mm-20.3mm附近，导致的原因是光子到达外圈的数量较少，有较大的偶然性。根据上述分析证明模拟数据有效，模拟分析结论可以用于实际测量源探位置的确定。

结论

为了弥补光谱测量过程的缺点，提出了基于苹果多层组织模型的蒙特卡洛模拟算法，并在matlab平台上实现了算法的并行计算,开辟了光谱检测水果品质的新思路。在模拟数据的基础上进行了光子入射最佳位置和探测器放置最佳位置研究，并用高光谱相机实际测量做出了验证，证明了模拟结果的正确性，同时，得到以下结论：1)光子在苹果赤道位置入射效果更好；2)探头放置最佳位置与光学参数有关，红富士苹果的最佳探测距离为距离光源1.5mm-10.15mm的环形区域为，可以根据苹果光学参数不同，调整探头位置。以上两点结论可以为光谱实际测量实验设计和便携式水果品质检测光谱仪的开发提供理论基础。

本文中所描述的具体实施例仅仅是对本发明精神做举例说明。本发明所属技术领域的技术人员可以对所描述的具体实施例做各种各样的修改或补充或采用类似的方式替代，但并不会偏离本发明的精神或者超越所附权利要求书所定义的范围。

Claims (10)

1.一种基于光子传输模拟的苹果高光谱品质检测方法，其特征在于它包括以下步骤：

(1)基于椭球曲面方程构建苹果模型；

(2)将点光源对准苹果赤道，源探距离选择1.5mm～10.15mm，基于苹果模型进行蒙特卡洛光子传输模拟；

(2.1)光子初始化；

(2.2)计算光子运动方向和步长，同时进行越界和出界判断；

(2.21a)越界判断，越界则进行步骤(2.21b)，否则进行步骤(2.22a)；

(2.21b)出界判断，出界则进行出界处理并转步骤(2.5)，否则进行步骤(2.22b)；

(2.22a)光子吸收与散射计算，转入步骤(2.3)；

(2.22b)光子折射、反射计算，转入步骤(2.3)；

(2.3)判断光子权值是否过小，是则进行步骤(2.4)，否则返回步骤(2.2)；

(2.4)判断光子是否存亡，是则进行步骤(2.5)，否则返回步骤(2.2)；

(2.5)判断是否为最后一个光子，是则结束，否则返回步骤(2.1)；

(3)取大量光子和不同光学参数组合重复步骤(2)，获得苹果模型的无噪声光亮度分布图。

2.根据权利要求1所述的方法，其特征在于步骤(1)中，基于椭球曲面方程构建苹果模型，苹果果皮和果肉的椭球曲面方程为：

其中：a、b、c，p1，p2为常数，共同控制苹果的长短轴半径和上下凹陷程度；a∈[3.1,5]，b∈[3.1,5]，c∈[4,6]，p₁∈[3,35]，p₂∈[5,670]；-π/2≤u≤π/2,0≤v≤2π；

苹果果核的椭球曲面方程为：

其中：a、b、c，p3为常数，共同控制苹果的长短轴半径和上下凹陷程度；k＝p₃sin5v，p3＝0.5；a∈[3.1,5]，b∈[3.1,5]，c∈[4,6]，-π/2≤u≤π/2,0≤v≤2π。

3.根据权利要求1所述的方法，其特征在于步骤(2.2)中，计算光子运动方向和步长：光子的运动方向由方位角ψ和散射偏转角θ确定，其中：

方位角ψ采用下式计算：

ψ＝2πξ

式中，ξ是区间(0，1)上的随机数；

散射偏转角θ采用下式计算：

式中，θ为散射偏转角，g为组织的各向异性，ξ是区间(0，1)上的随机数；

光子步长的计算公式为：

其中：

μ_t＝μ_a+μ_s

ξ是随机数，μ_a为物质的吸收系数，μ_s为物质的散射系数。

4.根据权利要求1所述的方法，其特征在于步骤(2.2)中，光子入射点的计算步骤为：

A、取越层运动位移的中点c，判断光子是否处于越层过程：若处于越层过程，则计算光子前一步位置坐标a和当前位置坐标b并进行步骤B；若未处于越层过程，则进入步骤C；

B、判断c点是否在模型交界点上，若是则返回该点的值，进入步骤E，否则进入步骤C；

C、判断a与c是否在一个位层：若是则将c赋给a；否则，将c赋给b；

D、判断b是否在交界点，否则返回A，是则进入步骤E；

E输出c。

5.根据权利要求4所述的方法，其特征在于光子位层的计算方法为：光子位于点P(x,y,z)，通过无误差转换计算得到参数(u,v)，将得到的(u,v)带入苹果果皮椭球曲面方程得到在该参数下果皮点p1(x1,y1,z1)，将p1和p点分别与原点的距离相比较：

-若p1点较远，则p点在果皮内部；

-若p1点较近，则p点在果皮之外；

-若p1于p点重合，则说明点在果皮上；

重复带入苹果果肉椭球曲面方程、苹果果核椭球曲面方程，得到果肉和果核位层的判断。

6.根据权利要求4所述的方法，其特征在于光子位层的计算方法为：光子位于点P(x,y,z)，通过无误差转换计算得到参数u，将得到的u带入苹果果皮椭球曲面方程得到在该参数下果皮点Z轴上的分量z1，通过对比z与z1的长短得到点P的位置：

-若z1较长，则p点在果皮内部；

-若z1较长，则p点在果皮之外；

-若z1与z相等，则说明点在果皮上；

重复带入苹果果肉椭球曲面方程、苹果果核椭球曲面方程，得到果肉和果核位层的判断。

7.根据权利要求5或6所述的方法，其特征在于所述无误差转换的计算公式为：

u＝tan^-1(az/(cxcosv))

v＝tan^-1(ay/bx)

式中，x、y、z分别为点P的直角坐标，a、b、c分别为系数，a∈[3.1,5]，b∈[3.1,5]，c∈[4,6]。

8.根据权利要求1所述的方法，其特征在于步骤(2.2)中，引入模型厚度，苹果模型包括：果皮外交界层、果皮层、果皮果肉交界层、果肉层、果肉果核交界层、果核层；光子的运动采取越层约束算法，算法步骤为：

A、计算光子移动前后的层位，若光子移动前层位a和移动后层位b不属于同一个层位则进行步骤B判断；

B、判断层位a与层位b是否为相邻层位：相邻则进入步骤F，不相邻则进入步骤C；

C、计算出光子运动方向；

D、计算出与光子运动方向同方向的、并且与光子位移前位层a相邻的位层，称为位层c；

E、暂认为光子由位层a位移至位层c并带入继续计算；

F、计算交界点，进入反射折射计算。

9.根据根据权利要求1所述的方法，其特征在于步骤(3)中，所述光学参数包括果皮吸收系数μ_a1、果肉吸收系数μ_a2、果皮散射系数μ_s1、果肉散射系数μ_s2。

10.根据根据权利要求9所述的方法，其特征在于果皮吸收系数μ_a1分为5类，果肉吸收系数μ_a2分为4类，果皮散射系数μ_s1分为5类，果肉散射系数μ_s2分为5类，各取中间值进行组合获得500组光学参数组合。