本技术涉及岩体力学实验，具体涉及一种岩石节理面剪切破坏区识别及体积计算方法。
背景技术：
1、由于长期地质作用，工程岩体中不可避免地存在大量节理裂隙等不连续面，这些不连续面纵横交错，显著削弱了完整岩石的力学性能。节理表面形貌特征是影响节理抗剪强度的重要因素之一，只有充分掌握节理表面的几何特征信息并对其定量化表征，才能准确地揭示节理的剪切力学行为。因此，对节理面剪切破坏区的精确识别和体积计算，有助于更好地预测工程岩体的稳定性。
2、针对上述问题，已有一些专利提出了与之相关的剪切破坏面积计算方法。例如，公开号为cn102749046b的中国专利公开了一种岩体结构面直剪试验中剪切面积的测量方法，其主要步骤为：采集岩体结构面照片并预处理；将结构面照片转换为黑白二值图片；计算白色部分占整个图片的面积百分比；根据拍照时采用的刻度尺作为参照物，计算出岩体结构面实际总面积。又如公开号为cn107144243a的中国专利公开了一种测量岩体结构面剪切破坏面积的方法和系统，该方法首先获取剪切前后的结构面点云数据，然后将两者对齐并计算相应点的高度偏差，再确定高度偏差大于0的点云数及其占总数的百分比w；根据百分比w和读取的直剪试验前岩体结构面面积s计算岩体结构面剪切破坏面。
3、以上两种方法操作简单，成本较低，但均不能直接用于计算节理面剪切破坏体积。
技术实现思路
1、本技术为了解决上述技术问题，提出了如下技术方案：
2、第一方面，本技术实施例提供了一种岩石节理面剪切破坏区识别及体积计算方法，包括：
3、采集剪切前后岩石节理面点云数据，并删除具有相同坐标的重复点；
4、计算点云数据中心并移至坐标原点得到平移后的点云数据；
5、根据法线向量与三维坐标系z轴的夹角及叉积结果，旋转点云数据使其法向与z轴方向一致；
6、将点云数据在xoy平面投影后，根据边界与x轴的夹角旋转点云数据，并调整边界至坐标轴；
7、通过点云配准方法使剪切前后两组点云精确配准；
8、根据点云数据的高差标记剪切破坏区，计算剪切前后点云围成的体积，即为剪切破坏体积。
9、在一种可能的实现方式中，所述采集剪切前后岩石节理面点云数据，包括：
10、将剪切前后的节理表面清理干净，并在岩石节理面四周涂抹反射率极低的材料；
11、通过三维激光扫描设备从不同位置对节理面进行多次分区扫描，得到剪切前后节理面三维点云数据。
12、在一种可能的实现方式中，所述计算点云数据中心并移至坐标原点得到平移后的点云数据，包括：
13、计算岩石节理面点云在x轴上的平均值在y轴上的平均值在z轴上的平均值
14、然后将所有点的坐标值(xi，yi，zi)减去即把三维点云中心移动到坐标原点，得到平移后的点云数据。
15、在一种可能的实现方式中，所述根据法线向量与三维坐标系z轴的夹角及叉积结果，旋转点云数据使其法向与z轴方向一致，包括：
16、计算平移后点云数据的协方差矩阵，矩阵形式如下：
17、
18、式中，σxx为x列的方差；σyy为y列的方差；σzz为z列的方差；σxy、σyx为x和y列之间的协方差；σxz、σzx是x和z列之间的协方差；σzy、σyz为y和z列之间的协方差。其中，方差计算公式为：
19、
20、协方差计算公式为：
21、
22、式中，n为样本数；分别为样本x列和y列的平均值。
23、利用qr分解协方差矩阵，分解后对角线上的元素为特征值，正交矩阵的列向量为特征向量，其中，最小特征值对应的特征向量即为点云平面的法线向量；
24、然后计算法线向量与向量(0，0，1)的夹角得到旋转角θ，以及两向量的叉积得到旋转轴向量，并对旋转轴向量进行归一化得到旋转轴的单位向量；
25、再利用罗德里格斯旋转公式对平移后的点云进行旋转，得到法线方向与z轴方向一致的点云，即旋转后的三维点云，罗德里格斯旋转公式如下：
26、v′＝vcosθ+(u·v)u(1-cosθ)+(u×v)sinθ
27、式中，v'为旋转后的空间向量；v为要旋转的空间向量；θ为旋转角度；u为旋转轴的单位向量。
28、在一种可能的实现方式中，所述将点云数据在xoy平面投影后，根据边界与x轴的夹角旋转点云数据，并调整边界至坐标轴，包括：
29、将旋转后三维点云中的z值设置为0，获得投影到xoy平面内的二维点云数据，然后寻找投影后点云的凸包；
30、得到凸包后，按照每个点相对于凸包中心的角度对点云排序；
31、再计算排序后点云中的第一个点(xj，yj)和第二个点(xj+1，yj+1)连线的斜率((yj+1-yj)/(xj+1-xj))，继续求出arctan((yj+1-yj)/(xj+1-xj))，该值为连线与x轴的夹角α；
32、以z轴为旋转轴，将所有点云数据旋转-α角度，并将旋转后点云数据所有点的y值减去y轴上的最小值，所有点的x值减去x轴上的最小值，即将点云左边界和下边界分别移动到了y轴和x轴上；
33、对剪切前后点云数据进行以上处理后，得到粗配准的两组点云数据。
34、在一种可能的实现方式中，所述寻找投影后点云的凸包，包括：
35、找到点云中最左端点a和最右端点b，这两个点一定在凸包上；
36、以a和b为端点形成一条线段，将所有点分为在线段ab上方和下方的两组；
37、对于每一组，找到离ab线段最远的点c，则点c也在凸包上，该过程需利用点到线段的距离公式：
38、
39、式中，(x1，y1)为a点坐标，(x2，y2)为b点坐标，任一点坐标为(x，y)；
40、分别对ac和cb两条线段递归地应用上述步骤，直到没有点在线段上方。
41、在一种可能的实现方式中，所述通过点云配准方法使剪切前后两组点云精确配准，包括：
42、以剪切后粗配准的节理面点云数据为源点云，以剪切前粗配准的节理面点云数据为目标点云；
43、通过最近点迭代算法使源点云不断接近目标点云，直到收敛，得到精配准的两组点云数据。
44、在一种可能的实现方式中，所述根据点云数据的高差标记剪切破坏区，计算剪切前后点云围成的体积，即为剪切破坏体积，包括：
45、分别对精配准的两组点云数据进行插值；
46、然后计算这两组点云数据在相同坐标点(xi，yi)处zi的差值δz；
47、设定一个阈值，δz超过该阈值的被视为发生剪切破坏，并标记出来；
48、再通过数值积分计算剪切前后两组点云所围成的体积，即为剪切破坏体积。
49、在一种可能的实现方式中，插值采用三次多项式插值，公式如下：
50、
51、式中，aij为待定系数，需要通过已知数据点及其在边界点的导数来确定。
52、在本技术实施例中，对剪切前后节理面点云数据进行一系列点云数据优化处理，然后进行精配准，确保了剪切破坏区识别及体积计算的准确性和可靠性，能够精确识别剪切破坏区域并计算其体积，相比传统的基于照片或二维数据的计算方法更为精确。不仅适用于实验室的岩石力学试验，还可以推广应用于实际工程中的岩体结构面分析。
1.一种岩石节理面剪切破坏区识别及体积计算方法，其特征在于，包括：
2.根据权利要求1所述的岩石节理面剪切破坏区识别及体积计算方法，其特征在于，所述采集剪切前后岩石节理面点云数据，包括：
3.根据权利要求1所述的岩石节理面剪切破坏区识别及体积计算方法，其特征在于，所述计算点云数据中心并移至坐标原点得到平移后的点云数据，包括：
4.根据权利要求3所述的岩石节理面剪切破坏区识别及体积计算方法，其特征在于，所述根据法线向量与三维坐标系z轴的夹角及叉积结果，旋转点云数据使其法向与z轴方向一致，包括：
5.根据权利要求4所述的岩石节理面剪切破坏区识别及体积计算方法，其特征在于，所述将点云数据在xoy平面投影后，根据边界与x轴的夹角旋转点云数据，并调整边界至坐标轴，包括：
6.根据权利要求5所述的岩石节理面剪切破坏区识别及体积计算方法，其特征在于，所述寻找投影后点云的凸包，包括：
7.根据权利要求5或6所述的岩石节理面剪切破坏区识别及体积计算方法，其特征在于，所述通过点云配准方法使剪切前后两组点云精确配准，包括：
8.根据权利要求7所述的岩石节理面剪切破坏区识别及体积计算方法，其特征在于，所述根据点云数据的高差标记剪切破坏区，计算剪切前后点云围成的体积，即为剪切破坏体积，包括：
9.根据权利要求8所述的岩石节理面剪切破坏区识别及体积计算方法，其特征在于，插值采用三次多项式插值，公式如下：
