基于B样条曲线和粒子群算法的堆石料图像检测级配曲线矫正方法与流程

基于b样条曲线和粒子群算法的堆石料图像检测级配曲线矫正方法
技术领域
1.本发明涉及土石坝堆石料级配检测领域，具体涉及一种基于b样条曲线和粒子群算法的堆石料图像检测级配曲线矫正方法。


背景技术：

2.土石坝由于其对基础条件的良好适应性、能就地取材、能充分利用建筑物的开挖料、造价较低、结构稳定性好，在水电工程建设中得到了广泛应用。心墙堆石坝具有填筑方量大、填筑施工强度高、质量控制要求严等特点。
3.级配是影响土石料压实后力学性能的主要因素。为了保证土石坝工程的填方压实质量，提高填方工程的抗渗性能和抗变形能力，采用合理级配的填筑料尤为重要。而土石坝填筑过程中，由于填筑方量较大，势必要对大量的土石料样本做级配检测，以保证工程质量。如何快速、准确检测并判断填筑土石料的级配的合理性一直是工程界关注的焦点问题。
4.传统的土石坝堆石料级配检测方法是通过人工随机抽取检测样本，再采用人工筛分或摇筛机对土石料进行筛分，然后统计各筛核中剩余料的质量百分数。该方法虽然技术成熟，但是耗时、费力、代表性不强，与现行的高效机械化填筑施工不相配套，无法快速进行大量填筑土石料的级配检测。为此，不少学者进行了大量的探索，其中最主流的方向就是基于图像识别技术的级配检测研究。图像识别技术的相关理论在当下已经得到长足的发展。然而，堆石料具有检测方量大、粒径跨度广、色泽较暗等特点，故在使用图像识别技术进行级配检测时，不可避免的出现偏差。经大量实验表明，这种检测偏差主要出现在小粒径石料颗粒上(粒径小于5mm的石料颗粒)，表现为漏检、错检，即在级配曲线的起始段，传统筛分实验获得的级配曲线位于图像识别技术检测的级配曲线之上。
5.根据对大量的两类数据的分析比较，发现两类级配曲线在起始段存在一般性差异，且这种差异表现为：(1)筛分曲线大于图像检测曲线；(2)不同样本曲线存在差异的曲线段长度基本相同；(3)在存在差异曲线段内，筛分曲线表现为强劲的外凸性，而图像检测级配曲线则表现为弱外凸性或内凹性。
6.分析以上差异，发明人尝试求存在差异的曲线段两端点处的切线方程，并求两切线方程的交点，再利用该交点以及曲线段两端点作为控制点，拟合b样条曲线，发现拟合后的b样条曲线与筛分曲线可以很好的重合在一起。于是用拟合后的b样条曲线替换图像检测曲线中与筛分曲线存在差异的曲线段，则可提高图像检测曲线的准确率。现在问题变为：当仅知道一组堆石料样本的图像检测级配曲线，而不知道筛分曲线时，如何确定“存在差异的曲线段”，即如何确定存在差异的曲线段的上限坐标和起点位置处的切线方向角。


技术实现要素：

7.基于上述技术问题，本发明提供一种基于b样条曲线和粒子群算法的堆石料图像检测级配曲线矫正方法，在大量的历史筛分数据的基础上，用粒子群算法寻找存在差异的
曲线段的统计性上限坐标和统计性起点位置处的切线方向角，并保存该统计性上限坐标和方向角，将该统计性上限坐标和方向角直接应用于后面的堆石料图像检测级配曲线的b样条曲线矫正中，使图像检测技术在堆石料级配曲线的获取中得到实际应用，从而提高施工效率，降低生产成本。
8.本发明采取的技术方案为：
9.基于b样条曲线和粒子群算法的堆石料图像检测级配曲线矫正方法，包括以下步骤：
10.步骤1：利用三角函数对堆石料的传统筛分级配曲线l
s
和图像检测级配曲线l
p
进行拟合，拟合的目的是使由不同数量的散点坐标组成的曲线l
s
和l
p
可以参与微分运算和任意点的插值计算，拟合公式为：
11.f(x)＝a0 a1cos(wx) b1sin(wx) a2cos(2wx) b2sin(2wx) a3cos(3wx) b3sin(3wx)
ꢀꢀ
(1)
12.式(1)中，x为堆石料粒径，f(x)为堆石料小于粒径x的质量百分数，其余参数均为待拟合常数：a0、a1、a2、a3分别表示各余弦项的常系数，b1、b2、b3分别表示各正弦项的常系数；w为三角函数振幅。
13.曲线拟合时，使用最小二乘法进行拟合，最小二乘法计算公式如下：
[0014][0015]
式(2)中，e
r
为拟合曲线f(x)与原散点曲线l
s,p
的残差平方和，n为数据点数量，拟合后，传统筛分级配曲线的拟合方程为f
s
(x)，图像检测级配曲线的拟合方程为f
p
(x)；l
s,p,i
为原散点曲线l
s,p
上的第i个数据点。
[0016]
步骤2：确定拟合之后的传统筛分级配曲线与图像检测级配曲线的公共起点处的方向角，求起点处的切线方程：
[0017]
设公共起点处曲线的方向角为a
n
，起点坐标为a(x
a
,y
a
)，则起点处的切线方程ac为：
[0018]
y
ac
＝tan(a
n
)(x
ac
‑
x
a
) y
a
ꢀꢀ
(3)
[0019]
式(3)中，x
ac
和y
ac
分别为切线方程ac的自变量和因变量。
[0020]
步骤3：确定拟合之后的图像检测级配曲线待矫正段的上限位置，并求该位置处的切线方程：
[0021]
设图像检测级配曲线待矫正段的上限位置为b(x
b
,y
b
)，则曲线f
p
(x)在点b(x
b
,y
b
)处的切线方程bc为：
[0022][0023]
式(4)中，x
bc
和y
bc
分别为切线方程bc的自变量和因变量。
[0024]
步骤4：计算两切线方程的交点：
[0025]
联立求解切线方程ac和切线方程bc，解得两条切线的交点坐标为c(x
c
,y
c
)，求解方程组如下：
[0026]
[0027]
式(5)中，y
ac
(x)和y
bc
(x)分别为切线ac和切线bc的切线方程。
[0028]
步骤5：利用起点、待矫正段上限位置和交点作为控制点，拟合b样条曲线：
[0029]
其中，起点坐标为a(x
a
,y
a
)，待矫正段上限位置为b(x
b
,y
b
)，交点为c(x
c
,y
c
)，将三个控制点的集合{a,b,c}定义为p
j
，则b样条曲线定义如下：
[0030][0031]
式(6)中，m为b样条曲线上节点数量，β为b样条曲线阶数，其取值范围是大于等于1的正整数；x
j
为曲线上的节点，其取值范围满足：x
j
∈[x0,x
m
‑1],x0<x1<
…
<x
m
‑1；x0、x1…
x
m
‑1为曲线上的节点。
[0032]
b
j,β
(x)为β阶b样条函数，定义为：
[0033][0034]
式(7)中，x
j 1
为曲线上第j 1个控制点，x
j β
为曲线上第j β个控制点，x
j β 1
为曲线上第j β 1个控制点，b
j,β
‑1(x)为第j个控制点处的β
‑
1阶b样条函数值，b
j 1,β
‑1(x)为第j 1个控制点处的β
‑
1阶b样条函数值。
[0035]
步骤6：利用拟合后的b样条曲线s(x)替代图像检测级配曲线l
p
的待矫正段ab，以此获得矫正后的级配曲线l
bp
(x)；
[0036]
步骤7：计算拟合后的传统筛分级配曲线f
s
(x)与矫正后的图像检测级配曲线l
bp
(x)之间的相对误差e
i
，计算公式如下：
[0037][0038]
式(8)中，n为曲线上数据点数量。
[0039]
步骤8：以相对误差计算公式为适应度函数，利用粒子群算法更新公共起点处的方向角和待矫正段的上限位置，重新计算两条切线，并重复步骤4至步骤7；
[0040]
步骤9：判断步骤7所得相对误差是否满足要求或是否达到粒子群算法最大迭代次数，若不满足要求，则转至步骤8，直至满足最大允许相对误差或达到最大迭代次数为止，即：e
i
≤e
min
或i>n。
[0041]
所述步骤8中，设最大允许相对误差为e
min
，算法最大迭代次数为n，粒子数量为a，各粒子的初始位置的集合为x
i
，x
i
为a行2列的二维数值矩阵，第一列表示方向角a
n,i
，第二列表示待矫正段上限位置b点横坐标x
b,i
；每次迭代结束时，各粒子的最优位置为x
i,best
，整个种群的最优位置为x
best
；各粒子的初始运动速度的集合为v
i
，v
i
亦为a行2列的二维数值矩阵，第一列表示各粒子在方向角a
n,i
方向的运动速度，第二列表示各粒子在x
b,i
方向的运动速度；各粒子的最佳适应度函数值为e
i,best
，则各粒子的位置更新公式为：
[0042][0043]
式(9)中，为各粒子下一时刻的位置，满足x
min
,x
max
表示粒
子在各方向上的有效活动范围。各粒子的速度更新公式为：
[0044]
v
inext
＝v
i
c1r1(x
i,best
‑
x
i
) c2r2(x
best
‑
x
i
)
ꢀꢀ
(10)
[0045]
式(10)中，c1,c2表示学习因子，r1,r2为0和1的随机数。v
inext
是各粒子下一时刻的运动速度，满足v
inext
∈[v
min
,v
max
]，v
min
,v
max
表示粒子在各方向上的有效速度范围。
[0046]
本发明一种基于b样条曲线和粒子群算法的堆石料图像检测级配曲线矫正方法，技术效果如下：
[0047]
1)本发明所涉及方法，使矫正后的图像检测级配曲线与筛分级配曲线基本吻合，二者相对误差低于2％，满足设计误差要求；
[0048]
2)本发明所涉及方法，使图像识别技术在堆石料级配曲线快速检测中得以应用，为土石坝填筑施工过程中大量堆石料的级配曲线快速检测提供了技术支持，及大地提高了堆石料级配检测效率。
附图说明
[0049]
图1为本发明方法的流程图。
[0050]
图2为堆石料级配曲线的传统筛分法和图像检测法结果对比图。
[0051]
图3为两切线方程的交点c求解结果示意图。
[0052]
图4为粒子群算法参数寻优结果示意图。
[0053]
图5为b样条曲线对图像检测级配曲线的矫正效果图。
具体实施方式
[0054]
基于b样条曲线和粒子群算法的堆石料图像检测级配曲线矫正方法，包括以下步骤：
[0055]
步骤1：利用三角函数对堆石料的传统筛分级配曲线l
s
和图像检测级配曲线l
p
进行拟合，拟合的目的是使由不同数量的散点坐标组成的曲线l
s
和l
p
可以参与微分运算和任意点的插值计算，拟合公式为：
[0056]
f(x)＝a0 a1cos(wx) b1sin(wx) a2cos(2wx) b2sin(2wx) a3cos(3wx) b3sin(3wx)
ꢀꢀ
(1)
[0057]
式(1)中，x为堆石料粒径，f(x)为堆石料小于粒径x的质量百分数，其余参数均为待拟合常数：a0、a1、a2、a3分别表示各余弦项的常系数，b1、b2、b3分别表示各正弦项的常系数；w为三角函数振幅。
[0058]
曲线拟合时，使用最小二乘法进行拟合，最小二乘法计算公式如下：
[0059][0060]
式(2)中，e
r
为拟合曲线f(x)与原散点曲线l
s,p
的残差平方和，n为数据点数量，拟合后，传统筛分级配曲线的拟合方程为f
s
(x)，图像检测级配曲线的拟合方程为f
p
(x)；l
s,p,i
为原散点曲线l
s,p
上的第i个数据点。
[0061]
步骤2：确定拟合之后的传统筛分级配曲线与图像检测级配曲线的公共起点处的方向角，求起点处的切线方程：
[0062]
设公共起点处曲线的方向角为a
n
，起点坐标为a(x
a
,y
a
)，则起点处的切线方程ac
为：
[0063]
y
ac
＝tan(a
n
)(x
ac
‑
x
a
) y
a
ꢀꢀ
(3)
[0064]
式(3)中，x
ac
和y
ac
分别为切线方程ac的自变量和因变量。
[0065]
步骤3：确定拟合之后的图像检测级配曲线待矫正段的上限位置，并求该位置处的切线方程：
[0066]
设图像检测级配曲线待矫正段的上限位置为b(x
b
,y
b
)，则曲线f
p
(x)在点b(x
b
,y
b
)处的切线方程bc为：
[0067][0068]
式(4)中，x
bc
和y
bc
分别为切线方程bc的自变量和因变量。
[0069]
步骤4：计算两切线方程的交点：
[0070]
联立求解切线方程ac和切线方程bc，解得两条切线的交点坐标为c(x
c
,y
c
)，求解方程组如下：
[0071][0072]
式(5)中，y
ac
(x)和y
bc
(x)分别为切线ac和切线bc的切线方程。
[0073]
步骤5：利用起点、待矫正段上限位置和交点作为控制点，拟合b样条曲线：
[0074]
其中，起点坐标为a(x
a
,y
a
)，待矫正段上限位置为b(x
b
,y
b
)，交点为c(x
c
,y
c
)，将三个控制点的集合{a,b,c}定义为p
j
，则b样条曲线定义如下：
[0075][0076]
式(6)中，m为b样条曲线上节点数量，β为b样条曲线阶数，其取值范围是大于等于1的正整数；x
j
为曲线上的节点，其取值范围满足：x
j
∈[x0,x
m
‑1],x0<x1<
…
<x
m
‑1；x0、x1…
x
m
‑1为曲线上的节点。
[0077]
b
j,β
(x)为β阶b样条函数，定义为：
[0078][0079]
式(7)中，x
j 1
为曲线上第j 1个控制点，x
j β
为曲线上第j β个控制点，x
j β 1
为曲线上第j β 1个控制点，b
j,β
‑1(x)为第j个控制点处的β
‑
1阶b样条函数值，b
j 1,β
‑1(x)为第j 1个控制点处的β
‑
1阶b样条函数值。
[0080]
步骤6：利用拟合后的b样条曲线s(x)替代图像检测级配曲线l
p
的待矫正段ab，以此获得矫正后的级配曲线l
bp
(x)；
[0081]
步骤7：计算拟合后的传统筛分级配曲线f
s
(x)与矫正后的图像检测级配曲线l
bp
(x)之间的相对误差e
i
，计算公式如下：
[0082][0083]
式(8)中，n为曲线上数据点数量。
[0084]
步骤8：以相对误差计算公式为适应度函数，利用粒子群算法更新公共起点处的方向角和待矫正段的上限位置，重新计算两条切线，并重复步骤4至步骤7；
[0085]
步骤9：判断步骤7所得相对误差是否满足要求或是否达到粒子群算法最大迭代次数，若不满足要求，则转至步骤8，直至满足最大允许相对误差或达到最大迭代次数为止，即：e
i
≤e
min
或i>n。
[0086]
所述步骤8中，设最大允许相对误差为e
min
，算法最大迭代次数为n，粒子数量为a，各粒子的初始位置的集合为x
i
，x
i
为a行2列的二维数值矩阵，第一列表示方向角a
n,i
，第二列表示待矫正段上限位置b点横坐标x
b,i
；每次迭代结束时，各粒子的最优位置为x
i,best
，整个种群的最优位置为x
best
；各粒子的初始运动速度的集合为v
i
，v
i
亦为a行2列的二维数值矩阵，第一列表示各粒子在方向角a
n,i
方向的运动速度，第二列表示各粒子在x
b,i
方向的运动速度；各粒子的最佳适应度函数值为e
i,best
，则各粒子的位置更新公式为：
[0087][0088]
式(9)中，为各粒子下一时刻的位置，满足x
min
,x
max
表示粒子在各方向上的有效活动范围。各粒子的速度更新公式为：
[0089]
v
inext
＝v
i
c1r1(x
i,best
‑
x
i
) c2r2(x
best
‑
x
i
)
ꢀꢀ
(10)
[0090]
式(10)中，c1,c2表示学习因子，r1,r2为0和1的随机数。v
inext
是各粒子下一时刻的运动速度，满足v
inext
∈[v
min
,v
max
]，v
min
,v
max
表示粒子在各方向上的有效速度范围。
[0091]
实施例：
[0092]
某砾石心墙堆石坝的部分堆石料样品的传统筛分级配曲线和图像检测级配曲线如图2所示。由图2可知，当粒径小于400mm时，图像检测级配曲线与传统筛分曲线存在一般性差异，即在0～400mm粒径范围内，传统筛分实验获得的级配曲线与图像识别技术检测的级配曲线之间的差异主要表现在：(1)筛分曲线大于图像检测曲线；(2)不同样本的两类级配曲线存在差异的粒径范围基本相同；(3)在存在差异的级配曲线段内，筛分曲线表现为严格的外凸性，而图像检测级配曲线则表现为弱外凸性或弱内凹性。
[0093]
取一组堆石料样品的图像检测级配曲线进行矫正，首先根据步骤1进行筛分曲线和图像检测曲线拟合，根据步骤2和步骤3，初始化起点方向角an，令an为60
°
，初始化待矫正段上限位置点b的横坐标x
b
为420，再利用拟合后的曲线方程、方向角an和待矫正段上限位置点b点横坐标x
b
计算两条切线ac和bc的方程，根据步骤4求两切线方程的交点c，求解结果如图3所示，由图3可知，两条切线ac和bc所包络的范围，就是图像检测级配曲线与传统筛分曲线存在差异粒径范围。
[0094]
根据步骤5进行初步b样条曲线拟合，本实施例中，设b样条曲线阶数β为2，再根据步骤6获得初步矫正后的堆石料级配曲线，根据步骤7计算初步矫正后的相对误差，初步矫正后的相对误差为2.86％，再根据步骤8，利用粒子群算法寻找最优参数，本实施例中，学习率c1＝c2＝0.9，粒子数量a为10，最大迭代次数n为20，参数寻优结果如图4所示，由上图4可
以看出，粒子群算法给出的最优方向角an为68.26
°
，最优切点位置b的横坐标为481.2，相对剩余误差为1.837％。
[0095]
曲线的最终矫正效果如图5所示，由图5可知，矫正曲线与筛分级配曲线几乎完全重合，且矫正段与非矫正段节点处平滑连接，故可将矫正级配曲线视为堆石料样本的标准级配曲线。
[0096]
由以上运行实施例可知，本发明提出的一种基于b样条曲线和粒子群算法的堆石料图像检测级配曲线矫正方法，在大量历史筛分数据的基础上，能够让图像检测出的级配曲线很好的贴合筛分曲线，矫正后的两条曲线的相对剩余误差低于2％，满足设计要求，从而达到曲线矫正的效果，验证了本发明中所提算法的有效性，使图像级配检测技术在土石坝填筑过程中得以应用，从而降低施工成本并极大提高施工效率。