[发明专利]一种钢构件圆柱体拟合方法

权利要求书

1.一种钢构件圆柱体拟合方法，其特征在于，包括以下步骤：

步骤1：在待测的一个大型钢结构的剖面圆上设置多个测量点，并布置高精度激光仪器采集测量点的三维坐标；

步骤2：采用转站算法，把不同测量点测量的三维坐标转换到一个公共坐标系下；

步骤3：利用步骤2转站后的坐标点计算拟合平面方程，得到平面方程参数，在拟合运算过程中，算法能自动剔除误差较大的数据，获取最佳拟合方程参数，以此平面作为剖面圆度的评定平面；

步骤4：将测量点坐标投影到步骤3确定的评定平面上，并计算各点投影坐标以及到评定平面的距离；

步骤5：利用平面圆拟合算法将步骤4投影到同一平面的坐标点作为圆周上的点进行平面圆拟合，得到圆的拟合半径和圆心；

步骤6：重复步骤1-5，对多个剖面圆进行处理，在步骤5平面圆拟合的基础上，将多个剖面圆采用最小二乘法进行柱体圆度计算即进行圆柱拟合；

步骤2中所述把不同测量点测量的三维坐标转换到一个公共坐标系下的具体方法为：用两组对应的测量点数据计算出对应的转站关系参数，两组测量点数据中均有测量点编号进行对应，算法输入两组带测量点编号的空间三维数据后输出转站关系，其中包括有3个平移量和3个旋转量，采用奇异值分解算法，构造H矩阵，H矩阵为两组测量点数据去重心化后的乘积，设目标坐标为P_T(X_T,Y_T,Z_T)，待转换的坐标为P_S(X_S,Y_S,Z_S)，H为3*3的矩阵，设P₁(n*3)＝[X_S,Y_S,Z_S],S＝0,1,2...n，设P₂(n*3)＝[X_T,Y_T,Z_T],T＝0,1,2...n，则H＝P₁^T·P₂，求H矩阵奇异值，可直接求转换矩阵：R＝VU^T，U和V为左右奇异矩阵，需要判断R的行列式值是否小于0，如果小于0则为反射矩阵，需要做转换，将V[0,3]的值乘以-1，再次计算R；n表示坐标点个数；

设目标坐标点云为Center₂(x₂,y₂,z₂)，待转换坐标为Center₁(x₁,y₁,z₁)

则平移向量为

有旋转矩阵后，可求得转换角度，这里以欧拉角为例，旋转矩阵R如下所示：

可以求得旋转角为：

θ_x＝atan2(R₃₂,R₃₃)

θz＝atan2(R₂₁,R₁₁)

其中：θ_x为x的旋转角，θ_y为y的旋转角，θ_z为z的旋转角，R_ij表示旋转矩阵第i行j列的值，为自转角，θ为章动角，ψ为旋进角；

步骤3中所述平面方程参数的计算方法为：

设平面方程为

ax+by+cz+d＝0

由点云P_i(x_i,y_i,z_i)、重心Center(x₀,y₀,z₀)，i表示点云P_i中任意一点；

得到去重心化的点云P_i(X_i,Y_i,Z_i)

构造矩阵D_n*3[X_n Y_n Z_n],n＝0,1,2,3...n

求矩阵的奇异值USVT＝D.svd()

D的最小奇异值对应的奇异向量即为平面方程参数a、b、c

d＝-(ax₀+by₀+cz₀)

求坐标点到平面的距离，求中误差，设限差为3倍中误差，将剔除距离超过3倍中误差的点，然后重新拟合平面；

其中，求点(x₀,y₀,z₀)到平面的距离公式如下：

步骤5中所述平面圆拟合方法如下：

采用非线性最小二乘法，假设空间圆参数为观测值，初始值可以取上面的线性平差结果，或者圆心取重心，半径取点云到重心的距离平均值；

根据距离公式构造误差方程

泰勒展开线性化，设

展开后的法方程为

V＝BX-l

则B为[(x₀-x_i)/D,(y₀-y_i)/D,(z₀-z_i)/D,-1],i＝0,1,...,n

l为[R-D]

求点云平面拟合方程，构造约束条件方程，其中

C＝[a b c]

W_x＝[-(d+ax₀+by₀+cz₀)]

d为点(x₀,y₀,z₀)到平面的距离

得法方程

X为方程解，K_s为联系数向量，B^T、C^T分别表示B、C的转置，解方程得

解方程得圆心坐标的改正值(dx₀,dy₀,dz₀,dR)

用改正值修正初值，做迭代计算，当X的改正值小于限差结束迭代，这里限差取10^-6；

半径R为各点到圆心的距离平均值，每一个坐标点到圆心的距离，再减去半径R就是它的挠度，平均挠度为所有挠度的平均值，最大挠度为所有挠度值中的最大值；

所述步骤6中采用最小二乘法进行柱体圆度计算的方法如下：

取圆柱面上的点到轴线上的投影点的距离为圆柱半径，以此构造方程，具体方法如下：

设中心轴线方程为

x＝x₀+at

y＝y₀+bt

z＝z₀+ct

式中(a,b,c)为空间直线的单位方向矢量，(x₀,y₀,z₀)为空间直线上离原点最近的点，t为直线上任意点到(x₀,y₀,z₀)的距离，为了唯一表示直线，定义a0，若a＝0则b0；若a＝0且b＝0，则c0；a，b，c不可能同时为0；

过观测点P_i(x_i,y_i,z_i)与中心轴线垂直的平面方程

ax+by+cz+D_i＝0

其中D_i＝-(ax_i+by_i+cz_i)

将轴线方程带入得：

t＝-(ax₀+by₀+cz₀+D_i)＝a(x_i-x₀)+b(y_i-y₀)+c(z_i-z₀)

中心轴线与该垂直平面的交点坐标(x_p,y_p,z_p)：

观测点P_i(x_i,y_i,z_i)与交点坐标的距离为：

得到各点误差方程为

其中，R为圆柱半径，7个参数分别为(a,b,c,x₀,y₀,z₀,R)，在以下两个条件下求7个参数

a²+b²+c²＝1

令误差方程线性化

其中

初始坐标点可以取点云重心即所有坐标点的加权平均值，半径取坐标点到重心距离的平均值,解方程即得到7个参数的改正值，(a,b,c)取1，改正值即为初始值和迭代后值的差值，修正初值后迭代计算，迭代过程保证a、b、c不同时为0，当改正值小于限差10^-6时，迭代结束。