一种林木应力波断层成像的速度误差校正方法与流程

本发明属于木材无损检测技术领域，具体是涉及一种林木应力波断层成像的速度误差校正方法。

背景技术：

面对木材资源匮乏的现状，如何提高木材利用率，促进我国林业现代化的可持续发展，已成为林业工作者关注的焦点。木材无损检测对于发现林木内部缺陷，防止腐朽扩散，提高木材利用率具有重要意义。应力波无损检测技术具有低成本、操作方便、无辐射，不受被测木材的大小形状影响等优点，已经在林木无损检测中得到了广泛的应用。虽然现在国际上主流的林木应力波无损检测产品如德国的ARBOTOM和PICUS以及匈牙利的FAKOPP都提供了二维或三维断层成像功能，但成像效果和检测精度仍然有待提高。

技术实现要素：

发明目的：为了克服现有技术中存在的不足，本发明提供一种林木应力波断层成像的速度误差校正方法，利用速度误差校正机制对速度反演结果进行优化，利用该方法可以准确检测出木材内部是否存在缺陷以及识别缺陷的面积及形状，与其他方法相比，本发明方法准确率高。

技术方案：为实现上述目的，本发明的一种林木应力波断层成像的速度误差校正方法，包括以下步骤：

S1生成木材横截面模拟图，对该木材横截面模拟图进行网格划分；

S2计算所述木材横截面模拟图中网格单元的传播速度；

S3计算所述木材横截面模拟图中网格单元的传播速度相对误差值；

S4根据投票机制对所述木材横截面模拟图中非边缘网格单元的传播速度相对误差值进行校正；

S5根据校正后的所述木材横截面模拟图中网格单元的传播速度相对误差值对该木材横截面模拟图中网格单元进行填色，形成该木材横截面的断层成像。

进一步地，对所述木材横截面模拟图进行网格划分基准是以该木材横截面模拟图轮廓为内切圆，对该内切圆的外切正四边形进行网格划分，设该外切正四边形的边长为a，则划分成a*a个网格单元，设所述木材横截面中网格单元的个数为m。

进一步地，所述步骤S2包括以下步骤：

S21利用应力波时间测定仪测量应力波在木材内部的传播时间tj，tj表示第j条应力波传播路径的传播时间，j＝(1,2,…,n)；

S22计算应力波传播路径在所述木材横截面模拟图中网格单元上的分割距离dij， i＝(1,2,…,m)，dij表示第j条应力波传播路径在木材横截面模拟图中第i个网格单元上的分割距离，j＝(1,2,…,n)；

S23通过方程组(一)计算所述木材横截面模拟图中网格单元的慢度pi，pi表示木材横截面模拟图中第i个网格单元的慢度，i＝(1,2,…,m)：

S24通过公式(一)计算所述木材横截面模拟图中网格单元的传播速度：

其中，vi表示木材横截面模拟图中第i个网格单元的传播速度，i＝(1,2,…,m)。

进一步地，所述步骤S22中，基于零点法求交点坐标的方法确定应力波传播路径与所述木材横截面模拟图中网格单元的交点，具体步骤如下：确定应力波传播路径直线的方程式；确定所述木材横截面模拟图中网格线的方程式；将应力波传播路径直线的方程式与所述木材横截面模拟图中网格线的方程式作差构成新的方程式；求解该新的方程式的零点，该零点即为应力波传播路径与所述木材横截面模拟图中网格线的交点，该应力波传播路径对应的发射传感器和接收传感器的位置坐标即为该应力波传播路径与所述木材横截面模拟图边界的两个交点；以该应力波传播路径对应的发射传感器的位置为起点，以该应力波传播路径对应的接收传感器的位置为终点，将该应力波传播路径与所述木材横截面模拟图中网格线的交点按照该交点与起点之间的距离排序，计算排序好的相邻两个交点之间的欧式距离，该欧式距离即为该应力波传播路径在所述木材横截面模拟图中网格单元上的分割距离。

进一步地，所述步骤S3中，根据公式(二)计算所述木材横截面模拟图中网格单元的传播速度相对误差值：

其中，δi表示所述木材横截面模拟图中第i个网格单元的传播速度相对误差值，vi表示木材横截面模拟图中第i个网格单元的传播速度，vref表示预设参考速度， i＝(1,2,…,m)。

进一步地，所述步骤S4包括以下步骤：

S41遍历校正所述木材横截面模拟图中非边缘网格单元的传播速度相对误差值；

S42判断非边缘网格单元的上、下、左、右四个邻居网格单元的传播速度相对误差值，如果四个邻居网格单元中有三个或三个以上的传播速度相对误差值大于或等于10％，则将该非边缘网格单元的传播速度相对误差值校正为10％；如果四个邻居网格单元中有三个或三个以上的传播速度相对误差值小于10％，则将该非边缘网格单元的传播速度相对误差值校正为[0,10％)的随机数；

S43经过一次遍历校正后，判断所述木材横截面模拟图中网格单元的传播速度相对误差值是否发生变化，若发生变化则进入S41，否则进入S5。

进一步地，所述步骤S5中，根据遍历校正后的所述木材横截面模拟图中网格单元的传播速度相对误差值对所述木材横截面模拟图中网格单元进行填色，若所述木材横截面模拟图中网格单元的传播速度相对误差值小于10％，则将对应的网格单元填充为颜色一；若所述木材横截面模拟图中网格单元的传播速度相对误差值大于或等于10％，则将对应的网格单元填充为颜色二，形成该木材横截面的断层成像。

有益效果：本发明与现有技术比较，具有的优点是：本发明方法首先具体是通过 LSQR(Least Square QR-factorization)最小平方QR分解法求解速度反演方程组，通过反演计算得到木材横截面上各网格单元的应力波传播速度，然后求出各网格单元的传播速度相对误差值，接着通过速度误差校正机制对传播速度相对误差值进行误差校正，最后根据校正后的网格单元的传播速度相对误差值对木材横截面进行填色，形成该木材横截面的断层成像；本发明方法在传统的应力波速度反演成像算法基础上，利用速度误差校正机制对反演结果进行优化，形成的木材横截面断层成像准确率高，可以准确检测出木材内部是否存在缺陷以及识别缺陷的面积及形状。

附图说明

图1是本发明方法流程图。

图2是应力波木材无损检测实验平台。

图3是采用本发明方法对被测木材进行试验得到的木材横截面的断层成像。

图4是图3中采用的被测木材的同一横截面实际缺陷情况图。

具体实施方式

下面结合附图对本发明作更进一步的说明。

本发明提出了一种林木应力波断层成像的速度误差校正方法，参照图1，包括以下步骤，S1生成木材横截面模拟图，对该木材横截面模拟图进行网格划分；

S2计算所述木材横截面模拟图中网格单元的传播速度；

S3计算所述木材横截面模拟图中网格单元的传播速度相对误差值；

S4根据投票机制对所述木材横截面模拟图中非边缘网格单元的传播速度相对误差值进行校正；

S5根据校正后的所述木材横截面模拟图中网格单元的传播速度相对误差值对该木材横截面模拟图中网格单元进行填色，形成该木材横截面的断层成像。

实施例：

在2016年度江苏省普通高效专业学位研究生实践创新计划项目-林木内部应力波传播仿真及成像算法研究(项目编号为SJLX16_0499)中，通过实验来实现本实施例，本实施例中，采用12个传感器，首先生成木材横截面模拟图，对该模拟图进行网格划分，对木材横截面模拟图进行网格划分基准是以该木材横截面模拟图轮廓为内切圆，对该内切圆的外切正四边形进行网格划分，设该外切正四边形的边长为a，则划分成a*a个网格单元，设所述木材横截面中网格单元的个数为m；如图3所示，图3中在该木材横截面轮廓上等分12等分，将该木材横截面轮廓模拟成正圆(模拟后的正圆上等分点显示不均匀)，将该正圆的外切正四边形进行网格划分，设该外切正四边形的边长为17，则划分成17*17＝289个网格单元，设所述木材横截面中网格单元的个数m小于289个；

参照图2，采用应力波木材无损检测实验平台测试得到应力波在木材内部的传播时间，图2中包括被测木材试样1、固定在木材试样1某个横截面上的一组传感器、用于敲击传感器的脉冲锤2，与传感器连接的应力波时间测试仪3以及与应力波测试仪连接的PC电脑4，PC电脑上设有分析软件系统；首先选择木材试样1的某个横截面，测量其半径及周长，根据其半径及周长将横截面轮廓12等分，分别在12等分点处(等分点数量根据具体实验需要而定)固定传感器，规定1号传感器为横截面水平正方向的传感器，逆时针依次顺序编号1～12号传感器，分别将12个传感器均匀固定在木材试样1横截面上，使用脉冲锤2从1号传感器开始敲击，每个传感器被敲击数次，记录PC电脑显示的度数并取其平均值作为最终数据，得到应力波在木材内部的传播时间，如表1所示，表1中显示了各发送端传感器和各接收端传感器之间的应力波传播时间，得到144 条应力波传播路径的传播时间，其中存在12个时间是发送端传感器和接收端传感器为同一个传感器时的时间，时间为0，将其删除，还剩132条应力波传播路径的传播时间，可以标记为t1～t132，得到时间列向量为n＝132；

表1应力波传播时间矩阵

计算应力波传播路径在所述木材横截面模拟图中网格单元上的分割距离dij， i＝(1,2,…,m)，dij表示第j条应力波传播路径在木材横截面模拟图中第i个网格单元上的分割距离，j＝(1,2,…,n)，基于零点法求交点坐标的方法确定应力波传播路径与所述木材横截面模拟图中网格单元的交点，具体步骤如下：确定应力波传播路径直线的方程式；确定所述木材横截面模拟图中网格线的方程式；将应力波传播路径直线的方程式与所述木材横截面模拟图中网格线的方程式作差构成新的方程式；求解该新的方程式的零点，该零点即为应力波传播路径与所述木材横截面模拟图中网格线的交点，该应力波传播路径对应的发射传感器和接收传感器的位置坐标即为该应力波传播路径与木材横截面模拟图边界的两个交点；以该应力波传播路径对应的发射传感器的位置为起点，以该应力波传播路径对应的接收传感器的位置为终点，将该应力波传播路径与木材横截面模拟图中网格线的交点按照该交点与起点之间的距离排序，计算排序好的相邻两个交点之间的欧式距离，该欧式距离即为该应力波传播路径在木材横截面模拟图中网格单元上的分割距离；上述132条应力波传播路径经过木材横截面模拟图中各网格(有部分网格不经过)，最终计算得到应力波传播路径在木材横截面模拟图中网格单元上的分割距离dij， dij表示第j条应力波传播路径在木材横截面模拟图中第i个网格单元上的分割距离， i＝(1,2,…,m)，j＝(1,2,…,n)，得到距离矩阵本实施例中取n＝132；

通过列方程组(一)计算木材横截面模拟图中网格单元的慢度pi，pi表示木材横截面模拟图中第i个网格单元的慢度，i＝(1,2,…,m)：

通过LSQR(Least Square QR-factorization)最小平方QR分解法求解速度反演方程组，得到木材横截面模拟图中m个网格单元的慢度p1～pm；设方程组(一)为T＝DP，n＝132，m为木材横截面模拟图中网格单元的个数， LSQR(Least Square QR-factorization)最小平方QR分解法中对于方程组AX＝B，调用函数 X＝lsqr(A，B)，就可以求得X，因此对于T＝DP方程组(一)通过最小平方QR分解法可以求得木材横截面模拟图中m个网格单元的慢度p1～pm；

通过公式(一)计算木材横截面模拟图中网格单元的传播速度：

其中，vi表示木材横截面模拟图中第i个网格单元的传播速度，i＝(1,2,…,m)，得到木材横截面模拟图中m个网格单元的慢度v1～vm；

计算木材横截面模拟图中网格单元的传播速度相对误差值：根据公式(二)计算木材横截面模拟图中网格单元的传播速度相对误差值：

其中，δi表示所述木材横截面模拟图中第i个网格单元的传播速度相对误差值，vi表示木材横截面模拟图中第i个网格单元的传播速度，vref表示预设参考速度， i＝(1,2,…,m)，最终得到木材横截面模拟图中m个网格单元的传播速度相对误差值δ1～δm；

根据投票机制对木材横截面模拟图中非边缘网格单元的传播速度相对误差值进行校正：遍历校正木材横截面模拟图中所有非边缘网格单元的传播速度相对误差值，判断非边缘网格单元的上、下、左、右四个邻居网格单元的传播速度相对误差值，如果四个邻居网格单元中有三个或三个以上的传播速度相对误差值大于或等于10％，则将该非边缘网格单元的传播速度相对误差值校正为10％；如果四个邻居网格单元中有三个或三个以上的传播速度相对误差值小于10％，则将该非边缘网格单元的传播速度相对误差值校正为[0,10％)的随机数；经过一次遍历校正后，判断木材横截面模拟图中网格单元的传播速度相对误差值是否发生变化，若发生变化则再次进行遍历校正，直至相对误差值不再发生变化；

最后根据校正后的木材横截面模拟图中网格单元的传播速度相对误差值对该木材横截面模拟图中网格单元进行填色，形成该木材横截面的断层成像，若木材横截面模拟图中网格单元的传播速度相对误差值小于10％，则将对应的网格单元填充为颜色一；若木材横截面模拟图中网格单元的传播速度相对误差值大于或等于10％，则将对应的网格单元填充为颜色二，形成该木材横截面的断层成像，本实施例中将颜色一设置为成淡灰色，将颜色二设置为深灰色，得到如图3所示的木材横截面的断层成像，在图3中可以清楚显示木材横截面上的缺陷位置和大小，将通过本发明方法得到的木材横截面的断层成像(图3所示)与采用的被测木材的同一横截面实际缺陷情况图作对比，被测木材的同一横截面实际缺陷情况图如图4所示，可以发现本发明方法得到的缺陷位置和大小和实际缺陷位置和大小基本一致，通过本发明方法可以准确检测出木材内部是否存在缺陷以及识别缺陷的面积及形状。

以上所述仅是本发明的优选实施方式，应当指出：对于本技术领域的技术人员来说，在不脱离本发明原理的前提下，还可以做出若干改进和润饰，这些改进和润饰也应视为本发明的保护范围。