CN111027010A - 一种钢构件圆柱体拟合算法 - Google Patents

一种钢构件圆柱体拟合算法 Download PDF

Info

Publication number
CN111027010A
CN111027010A CN201911112195.8A CN201911112195A CN111027010A CN 111027010 A CN111027010 A CN 111027010A CN 201911112195 A CN201911112195 A CN 201911112195A CN 111027010 A CN111027010 A CN 111027010A
Authority
CN
China
Prior art keywords
plane
equation
point
fitting
circle
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201911112195.8A
Other languages
English (en)
Other versions
CN111027010B (zh
Inventor
唐坚强
周慧明
吴颖
游勇
李翔
华程铭
韩玉超
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Wuhan Tianheng Information Technology Co ltd
Original Assignee
Wuhan Tianheng Information Technology Co ltd
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Wuhan Tianheng Information Technology Co ltd filed Critical Wuhan Tianheng Information Technology Co ltd
Priority to CN201911112195.8A priority Critical patent/CN111027010B/zh
Publication of CN111027010A publication Critical patent/CN111027010A/zh
Application granted granted Critical
Publication of CN111027010B publication Critical patent/CN111027010B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G06F17/13—Differential equations
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Mathematical Physics (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Theoretical Computer Science (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Algebra (AREA)
  • General Engineering & Computer Science (AREA)
  • Operations Research (AREA)
  • Computing Systems (AREA)
  • Length Measuring Devices With Unspecified Measuring Means (AREA)

Abstract

本发明公开了一种钢构件圆柱体拟合算法，包括以下步骤：在待测剖面圆上设置多个测量点，并布置高精度激光仪器采集测量点的三维坐标采用转站算法，把不同测量位置测量的三维坐标，转换到一个公共坐标系下；利用转站后的坐标点计算拟合平面方程，得到平面方程参数，以此平面作为剖面圆度的评定平面；将测量点坐标投影到评定平面上，并计算各点投影坐标以及到评定平面的距离；利用线性最小二乘法将投影到同一平面的坐标点进行平面圆拟合，得到圆的拟合半径和圆心；在平面圆拟合的基础上，将多个剖面圆采用最小二乘法进行柱体圆度计算即进行圆柱拟合，通过本发明将多种算法整合为一体，对大型圆柱形钢构件尺寸测量具有很强的针对性。

Description

一种钢构件圆柱体拟合算法
技术领域
本发明涉及工业测量、钢结构检测领域，具体涉及一种钢构件圆柱体拟合算法。
背景技术
在现代建造业中，人们对大型钢结构(尺寸在十米至数十米)的加工精度要求越来越高，对建造技术提出了更高要求的同时，对建造质量的准确评估和检验也有着迫切的需求。
大型圆柱形刚钢构件广泛应用于大型桥梁，涵洞，隧道，船舶，航天，航空等领域，圆柱形钢构件的主要指标主要包括端面的平面度，柱体的圆度、挠度、直径等。相对于小型构件可用测量仪器直接测量，大型圆柱形钢构件通常都是固定的不可移动的物体，其直径往往非常巨大，且结构复杂，周围伴随着诸多遮挡物，测量空间也可能比较受限，往往很难在一个位置完成测量，而多位置测量的数据如何进行同坐标系转换又是一大难题。
目前的大型圆柱形刚构件的测量一般采用撑杆法等几何测量的方式，这种方式除需要大量人力外，检测慢，效率低，操作繁琐，主观性大，测量精度不高，自动化程度低，此外该方法需要整个构件处在开阔的场地，须保证其在某个位置能被整体测量，大大限制了测量条件。同时该方法存在的不确定因素多，即影响了建造效率，也无法对建造质量作出精确的评估。
对于大型圆柱形刚构件传统方法也只能测量圆柱端面的形状，对于端面的平整度及整个柱体的圆度、挠度无法做出测量。因此提出一套完整的测量、分析、计算方法来满足大尺寸圆柱形钢构件高精度测量的要求。
发明内容
有鉴于此，本发明提供了一种钢构件圆柱体拟合算法
本发明提供了一种钢构件圆柱体拟合算法，包括以下步骤：
步骤1：在待测一个大型钢结构的剖面圆上设置多个测量点，并布置高精度激光仪器采集测量点的三维坐标；
步骤2：采用转站算法，把不同测量点测量的三维坐标转换到一个公共坐标系下；
步骤3：利用步骤2转站后的坐标点计算拟合平面方程，得到平面方程参数，在拟合运算过程中，算法能自动剔除误差较大的数据，获取最佳拟合方程参数，以此平面作为剖面圆度的评定平面；
步骤4：将测量点坐标投影到步骤3确定的评定平面上，并计算各点投影坐标以及到评定平面的距离；
步骤5：利用平面圆拟合算法将步骤4投影到同一平面的坐标点作为圆周上的点进行平面圆拟合，得到圆的拟合半径和圆心；
步骤6：重复步骤1-5，对多个剖面圆进行处理，在步骤5平面圆拟合的基础上，将多个剖面圆采用最小二乘法进行柱体圆度计算即进行圆柱拟合。
进一步地，步骤1中所述高精度激光仪器与计算机连接，实时传输所测数据给计算机进行处理计算。
进一步地，步骤2中所述把不同测量位置测量的三维坐标转换到一个公共坐标系下的具体方法为：用两组对应的测量点数据计算出对应的转站关系参数，两组测量点数据中均有测量点编号进行对应，算法输入两组带测量点编号的空间三维数据后输出转站关系，其中包括有3个平移量和3个旋转量，采用奇异值分解算法，构造H矩阵，H矩阵为两组测量点数据去重心化后的乘积，设目标坐标为PT(XT,YT,ZT)，待转换的坐标为PS(XS,YS,ZS)，H为3*3的矩阵，设P1(n*3)＝[XS,YS,ZS],(S＝0,1,2L n)，设P2(n*3)＝[XT,YT,ZT],(T＝0,1,2Ln)，则H＝P1 T·P2，求H矩阵奇异值，可直接求转换矩阵：R＝VUT，U和V为左右奇异矩阵，需要判断R的行列式值是否小于0，如果小于0则为反射矩阵，需要做转换，将V[0,3]的值乘以-1，再次计算R；n表示坐标点个数。
设目标坐标点云为Center2(x2,y2,z2)，待转换坐标为Center1(x1,y1,z1)
则平移向量为
Figure BDA0002273053680000031
有旋转矩阵后，可求得转换角度，这里以欧拉角为例，旋转矩阵R如下所示：
Figure BDA0002273053680000032
可以求得旋转角为：
θx＝atan2(R32,R33)
Figure BDA0002273053680000033
θz＝atan2(R21,R11)
其中：θx为x的旋转角，θy为y的旋转角，θz为z的旋转角，Rij表示旋转矩阵第i行j列的值，
Figure BDA0002273053680000034
为自转角，θ为章动角，ψ为旋进角。
进一步地，步骤3中所述平面方程参数的计算方法为：
设平面方程为
ax+by+cz+d＝0
由点云Pi(xi,yi,zi)、重心Center(x0,y0,z0)i表示点云pi中任意一点；
得到去重心化的点云Pi(Xi,Yi,Zi)
构造矩阵Dn*3[Xn Yn Zn],n＝0,1,2,3Kn
求矩阵的奇异值USVT＝D.svd()
D的最小奇异值对应的奇异向量即为平面方程参数a、b、c
d＝-(ax0+by0+cz0)
求坐标点到平面的距离，求中误差，设限差为3倍中误差，将剔除距离超过3倍中误差的点，然后重新拟合平面；
其中，求点(x0,y0,z0)到平面的距离公式如下：
Figure BDA0002273053680000041
进一步地，步骤4中所述空间点投影到指定的空间平面的方法为：
由平行关系有下面投影方程：
Figure BDA0002273053680000042
(x′，y′，z′)为各坐标点在投影平面上的投影坐标；
由垂直关系求得测量点到投影平面的距离t
Figure BDA0002273053680000043
进一步地，步骤5中所述平面圆拟合方法如下：
采用非线性最小二乘法，假设空间圆参数为观测值，初始值可以取上面的线性平差结果，或者圆心取重心，半径取点云到重心的距离平均值；
根据距离公式构造误差方程
Figure BDA0002273053680000044
泰勒展开线性化，设
Figure BDA0002273053680000051
展开后的法方程为V＝BX-l
则B为[(x0-xi)/D,(y0-yi)/D,(z0-zi)/D,-1],i＝0,1,K,n
l为[R-D]
求点云平面拟合方程，构造约束条件方程，其中
C＝[a b c]
Wx＝[-(d+ax0+by0+cz0)]
d为点(x0,y0,z0)到平面的距离
得法方程
Figure BDA0002273053680000052
X为方程解，Ks为联系数向量，BT、CT分别表示B、C的转置，解方程得
Figure BDA0002273053680000053
解方程得圆心坐标的改正值(dx0,dy0,dz0,dR)
用改正值修正初值，做迭代计算，当X的改正值小于限差结束迭代，这里限差取10-6；
半径R为各点到圆心的距离平均值，每一个坐标点到圆心的距离，再减去半径R就是它的挠度，平均挠度为所有挠度的平均值，最大挠度为所有挠度值中的最大值。
进一步地，所述步骤6中采用最小二乘法进行柱体圆度计算的方法如下：
取圆柱面上的点到轴线上的投影点的距离为圆柱半径，以此构造方程，具体方法如下：
设中心轴线方程为
x＝x0+at
y＝y0+bt
z＝z0+ct
式中(a,b,c)为空间直线的单位方向矢量，(x0,y0,z0)为空间直线上离原点最近的点，t为直线上任意点到(x0,y0,z0)的距离，为了唯一表示直线，定义a>0，若a＝0则b>0；若a＝0且b＝0，则c>0；a，b，c不可能同时为0。
过观测点Pi(xi,yi,zi)与中心轴线垂直的平面方程
ax+by+ch+Di＝0
其中Di＝-(axi+byi+chi)
将轴线方程带入得：
t＝-(ax0+by0+ch0+Di)＝a(xi-x0)+b(yi-y0)+c(hi-h0)
中心轴线与该垂直平面的交点坐标：
Figure BDA0002273053680000061
观测点Pi(xi,yi,zi)与交点坐标的距离为：
Figure BDA0002273053680000062
得到各点误差方程为
Figure BDA0002273053680000063
其中，R圆柱为半径，7个参数分别为(a,b,c,x0,y0,z0,R)，在以下两个条件下求7个参数
a2+b2+c2＝1
Figure BDA0002273053680000064
令
Figure BDA0002273053680000071
误差方程线性化
Figure BDA0002273053680000072
Figure BDA0002273053680000073
Figure BDA0002273053680000074
Figure BDA0002273053680000075
Figure BDA0002273053680000076
Figure BDA0002273053680000077
Figure BDA0002273053680000078
其中
Figure BDA0002273053680000081
Figure BDA0002273053680000082
Figure BDA0002273053680000083
Figure BDA0002273053680000084
初始坐标点可以取点云重心即所有坐标点的加权平均值，半径取坐标点到重心距离的平均值,解方程即得到7个参数的改正值，(a,b,c)取1，改正值即为初始值和迭代后值的差值，修正初值后迭代计算，迭代过程保证a、b、c不同时为0，当改正值小于限差10-6时，迭代结束。
本发明提供的技术方案带来的有益效果是：通过本发明将多种算法整合为一体，对大型圆柱形钢构件尺寸测量具有很强的针对性，且适应各种实地建造环境，整套系统集成化高，便于携带，配套高精度激光测量手段，测量效率和测量精度都大大提高，解决大型圆柱形钢构件尺寸参数快速高精度测量的问题。
附图说明
图1是本发明一种一种钢构件圆柱体拟合算法的功能结构图；
图2是本发明一种钢构件圆柱体拟合算法的转站关系计算流程图；
图3是本发明一种钢构件圆柱体拟合算法的平面方程计算流程图；
图4是本发明一种钢构件圆柱体拟合算法的投影计算流程图；
图5是本发明一种钢构件圆柱体拟合算法的平面圆拟合流程图；
图6是本发明一种钢构件圆柱体拟合算法的圆柱拟合流程图；
图7是本发明一种钢构件圆柱体拟合算法的软件主界面图；
图8是本发明一种钢构件圆柱体拟合算法的计算转站关系界面图；
图9是本发明一种钢构件圆柱体拟合算法的计算转站关系结果图；
图10是本发明一种钢构件圆柱体拟合算法的平面圆度计算界面；
图11是本发明一种钢构件圆柱体拟合算法的平面圆度计算结果图；
图12是本发明一种钢构件圆柱体拟合算法的计算平面圆度限定半径的操作界面图；
图13是本发明一种钢构件圆柱体拟合算法的圆柱拟合操作软件界面图；
图14是本发明一种钢构件圆柱体拟合算法的软件圆柱拟合结果图。
具体实施方式
为使本发明的目的、技术方案和优点更加清楚，下面将结合附图对本发明实施方式作进一步地描述。
请参考图1，本发明提供了一种钢构件圆柱体拟合算法，包括以下步骤：
步骤1：在待测剖面圆上设置多个测量点，并布置高精度激光仪器采集测量点的三维坐标；高精度激光仪器与计算机连接，实时传输所测数据给计算机进行处理计算。
步骤2：采用转站算法，包括正交矩阵法，奇异值分解法，单元四元数法，对偶四元数法4种转站算法，把不同测量位置测量的三维坐标，转换到一个公共坐标系下，如图2所示，具体方法为：用两组对应的测量点数据计算出对应的转站关系参数，两组测量点数据中均有测量点编号进行对应，算法输入两组带测量点编号的空间三维数据后输出转站关系，其中包括有3个平移量和3个旋转量，本实施例采用奇异值分解算法，构造H矩阵，H矩阵为两组测量点数据去重心化后的乘积，设目标坐标为PT(XT,YT,ZT)，待转换的坐标为PS(XS,YS,ZS)，H为3*3的矩阵，设P1(n*3)＝[XS,YS,ZS],(S＝0,1,2L n)，设P2(n*3)＝[XT,YT,ZT],(T＝0,1,2L n)，则H＝P1 T·P2，求H矩阵奇异值，可直接求转换矩阵：R＝VUT，U和V为左右奇异矩阵，需要判断R的行列式值是否小于0，如果小于0则为反射矩阵，需要做转换，将V[0,3]的值乘以-1，再次计算R；
设目标坐标点云为Center2(x2,y2,z2)，待转换坐标为Center1(x1,y1,z1)
则平移向量为
Figure BDA0002273053680000101
有旋转矩阵后，可求得转换角度，这里以欧拉角为例，旋转矩阵R如下所示：
Figure BDA0002273053680000102
可以求得旋转角为：
θx＝atan2(R32,R33)
Figure BDA0002273053680000103
θz＝atan2(R21,R11)
其中：θx为x的旋转角，θy为y的旋转角，θz为z的旋转角，Rij表示旋转矩阵第i行j列的值，
Figure BDA0002273053680000111
为自转角，θ为章动角，ψ为旋进角。
步骤3：利用步骤2转站后的坐标点计算拟合平面方程，得到平面方程参数，如图3所示，在拟合运算过程中，算法能自动剔除误差较大的数据，获取最佳拟合方程参数，以此平面作为剖面圆度的评定平面；
设平面方程为
ax+by+cz+d＝0
由点云Pi(xi,yi,zi)、重心Center(x0,y0,z0)，i表示点云pi中任意一点；
得到去重心化的点云Pi(Xi,Yi,Zi)
构造矩阵Dn*3[Xn Yn Zn],n＝0,1,2,3Kn
求矩阵的奇异值USVT＝D.svd()
D的最小奇异值对应的奇异向量即为平面方程参数a、b、c
d＝-(ax0+by0+cz0)
求坐标点到平面的距离，中误差取10-6，设限差为3倍中误差，将剔除距离超过3倍中误差的点，然后重新拟合平面；
其中，求点(x0,y0,z0)到平面的距离公式如下：
Figure BDA0002273053680000112
步骤4：将测量点坐标投影到步骤3确定的评定平面上，并计算各点投影坐标以及到评定平面的距离，如图4所示，具体方法为：
由平行关系有下面投影方程
Figure BDA0002273053680000113
(x′，y′，z′)为各坐标点在投影平面上的投影坐标；
由垂直关系求得测量点到投影平面的距离t
Figure BDA0002273053680000121
步骤5：采用非线性最小二乘法将步骤4投影到同一平面的坐标点进行平面圆拟合，得到圆的拟合半径和圆心，如图5所示，具体方法如下：
假设空间圆参数为观测值，初始值可以取上面的线性平差结果，或者圆心取重心，半径取点云到重心的距离平均值；
根据距离公式构造误差方程
Figure BDA0002273053680000122
泰勒展开线性化，设
Figure BDA0002273053680000123
展开后的法方程为V＝BX-l
则B为[(x0-xi)/D,(y0-yi)/D,(z0-zi)/D,-1],i＝0,1,K,n
l为[R-D]
求点云平面拟合方程，构造约束条件方程，其中
C＝[a b c]
Wx＝[-(d+ax0+by0+cz0)]
d为点(x0,y0,z0)到平面的距离
得法方程
Figure BDA0002273053680000124
X为方程解，Ks为联系数向量，BT、CT分别表示B、C的转置，解方程得
Figure BDA0002273053680000125
解方程得圆心坐标的改正值(dx0,dy0,dz0,dR)
用改正值修正初值，做迭代计算，当X的改正值小于限差结束迭代，这里限差取10-6；
半径R为各点到圆心的距离平均值，每一个坐标点到圆心的距离，再减去半径R就是它的挠度，平均挠度为所有挠度的平均值，最大挠度为所有挠度值中的最大值。
步骤6：重复步骤1-5，对多个剖面圆进行处理，在步骤5平面圆拟合的基础上，将多个剖面圆采用最小二乘法进行柱体圆度计算即进行圆柱拟合，如图6所示，方法如下：
取圆柱面上的点到轴线上的投影点的距离为圆柱半径，以此构造方程，具体方法如下：
设中心轴线方程为
x＝x0+at
y＝y0+bt
z＝z0+ct
式中(a,b,c)为空间直线的单位方向矢量，(x0,y0,z0)为空间直线上离原点最近的点，t为直线上任意点到(x0,y0,z0)的距离，为了唯一表示直线，定义a>0，若a＝0则b>0；若a＝0且b＝0，则c>0；a，b，c不可能同时为0；
过观测点Pi(xi,yi,zi)与中心轴线垂直的平面方程
ax+by+cz+Di＝0
其中Di＝-(axi+byi+czi)
将轴线方程带入得：
t＝-(ax0+by0+cz0+Di)＝a(xi-x0)+b(yi-y0)+c(zi-z0)
中心轴线与该垂直平面的交点坐标：
Figure BDA0002273053680000141
观测点Pi(xi,yi,zi)与交点坐标的距离为：
Figure BDA0002273053680000142
得到各点误差方程为
Figure BDA0002273053680000143
其中，R为圆柱半径，7个参数分别为(a,b,c,x0,y0,z0,R)，在以下两个条件下求7个参数
a2+b2+c2＝1
Figure BDA0002273053680000144
令
Figure BDA0002273053680000145
误差方程线性化
Figure BDA0002273053680000146
Figure BDA0002273053680000147
Figure BDA0002273053680000148
Figure BDA0002273053680000149
Figure BDA00022730536800001410
Figure BDA00022730536800001411
Figure BDA00022730536800001412
其中
Figure BDA0002273053680000151
Figure BDA0002273053680000152
Figure BDA0002273053680000153
Figure BDA0002273053680000154
初始坐标点可以取点云重心即所有坐标点的加权平均值，半径取坐标点到重心距离的平均值,解方程即得到7个参数的改正值，(a,b,c)取1，改正值即为初始值和迭代后值的差值，修正初值后迭代计算，迭代过程保证a、b、c不同时为0，当改正值小于限差10-6时，迭代结束。
考虑系统运行于微软Windows及Microsoft.NET Framework 4.0环境中，根据平台统一性和可重用性原则，使用C#语言进行开发，采用微软动态链接库(Dynamic LinkLibrary)技术，实现在Windows操作系统下共享函数库。同时实现模块化，方便二次开发。最终以DLL(动态链接库)文件的方式形成成果，将各种算法封装于动态链接库中，供测量时使用。
该程序用于调用测量辅助计算模块程序，并演示各种计算方法，同时输出各类计算结果，以便于使用者检验或者比较计算结果，主界面如图7所示，程序的文件数据格式为xyz坐标，以空格隔开每行一个点坐标。
(1)转站关系计算界面如图8所示
读取数据，一次选中两个文件，打开，进行计算，会计算第一个文件的数据转到第二个文件的坐标系中参数。点击交换数据，再次计算可以反过来计算，计算结果如图9所示。
(2)平面圆度计算
输入限差，会自动以该限差作为距离中误差的倍数进行计算，勾选求投影点和挠度，会自动调用计算方法。单击“读取数据”，打开文件，读取数据，操作界面如图10所示。
勾选了显示绘图，会自动把结果绘制出来，勾选清空输出，每次计算会自动将上次输出的结果清除，关闭上次的绘图窗口。
绘图中以挠度大小点的按从小到大颜色以蓝色渐变至红色，如图11所示。
限定半径的约束平差，需要输入约束半径如图12所示。
(3)圆柱拟合
直接打开数据，进行计算，计算完后程序会自动绘制出图形，如图13所示，结果如图14所示。
在不冲突的情况下，本文中上述实施例及实施例中的特征可以相互结合，以上所述仅为本发明的较佳实施例，并不用以限制本发明，凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (7)

1.一种钢构件圆柱体拟合算法，其特征在于，包括以下步骤：
步骤1：在待测一个大型钢结构的剖面圆上设置多个测量点，并布置高精度激光仪器采集测量点的三维坐标；
步骤2：采用转站算法，把不同测量点测量的三维坐标转换到一个公共坐标系下；
步骤3：利用步骤2转站后的坐标点计算拟合平面方程，得到平面方程参数，在拟合运算过程中，算法能自动剔除误差较大的数据，获取最佳拟合方程参数，以此平面作为剖面圆度的评定平面；
步骤4：将测量点坐标投影到步骤3确定的评定平面上，并计算各点投影坐标以及到评定平面的距离；
步骤5：利用平面圆拟合算法将步骤4投影到同一平面的坐标点作为圆周上的点进行平面圆拟合，得到圆的拟合半径和圆心；
步骤6：重复步骤1-5，对多个剖面圆进行处理，在步骤5平面圆拟合的基础上，将多个剖面圆采用最小二乘法进行柱体圆度计算即进行圆柱拟合。
2.根据权利要求1所述的一种钢构件圆柱体拟合算法，其特征在于，步骤1中所述高精度激光仪器与计算机连接，实时传输所测数据给计算机进行处理计算。
3.根据权利要求1所述的一种钢构件圆柱体拟合算法，其特征在于，步骤2中所述把不同测量位置测量的三维坐标转换到一个公共坐标系下的具体方法为：用两组对应的测量点数据计算出对应的转站关系参数，两组测量点数据中均有测量点编号进行对应，算法输入两组带测量点编号的空间三维数据后输出转站关系，其中包括有3个平移量和3个旋转量，采用奇异值分解算法，构造H矩阵，H矩阵为两组测量点数据去重心化后的乘积，设目标坐标为PT(XT,YT,ZT)，待转换的坐标为PS(XS,YS,ZS)，H为3*3的矩阵，设P1(n*3)＝[XS,YS,ZS],(S＝0,1,2L n)，设P2(n*3)＝[XT,YT,ZT],(T＝0,1,2L n)，则H＝P1 T·P2，求H矩阵奇异值，可直接求转换矩阵：R＝VUT，U和V为左右奇异矩阵，需要判断R的行列式值是否小于0，如果小于0则为反射矩阵，需要做转换，将V[0,3]的值乘以-1，再次计算R；n表示坐标点个数。
设目标坐标点云为Center2(x2,y2,z2)，待转换坐标为Center1(x1,y1,z1)
则平移向量为
Figure FDA0002273053670000021
有旋转矩阵后，可求得转换角度，这里以欧拉角为例，旋转矩阵R如下所示：
Figure FDA0002273053670000022
可以求得旋转角为：
θx＝a tan 2(R32,R33)
Figure FDA0002273053670000023
θz＝a tan 2(R21,R11)
其中：θx为x的旋转角，θy为y的旋转角，θz为z的旋转角，Rij表示旋转矩阵第i行j列的值，
Figure FDA0002273053670000024
为自转角，θ为章动角，ψ为旋进角。
4.根据权利要求1所述的一种钢构件圆柱体拟合算法，其特征在于，步骤3中所述平面方程参数的计算方法为：
设平面方程为
ax+by+cz+d＝0
由点云Pi(xi,yi,zi)、重心Center(x0,y0,z0)，i表示点云pi中任意一点；
得到去重心化的点云Pi(Xi,Yi,Zi)
构造矩阵Dn*3[Xn Yn Zn],n＝0,1,2,3Kn
求矩阵的奇异值USVT＝D.svd()
D的最小奇异值对应的奇异向量即为平面方程参数a、b、c
d＝-(ax0+by0+cz0)
求坐标点到平面的距离，求中误差，设限差为3倍中误差，将剔除距离超过3倍中误差的点，然后重新拟合平面；
其中，求点(x0,y0,z0)到平面的距离公式如下：
Figure FDA0002273053670000031
5.根据权利要求1所述的一种钢构件圆柱体拟合算法，其特征在于，步骤4中所述空间点投影到指定的空间平面的方法为：
由平行关系有下面投影方程：
Figure FDA0002273053670000032
(x′，y′，z′)为各坐标点在投影平面上的投影坐标；
由垂直关系求得测量点到投影平面的距离t
Figure FDA0002273053670000033
6.根据权利要求1所述的一种钢构件圆柱体拟合算法，其特征在于，步骤5中所述平面圆拟合方法如下：
采用非线性最小二乘法，假设空间圆参数为观测值，初始值可以取上面的线性平差结果，或者圆心取重心，半径取点云到重心的距离平均值；
根据距离公式构造误差方程
Figure FDA0002273053670000041
泰勒展开线性化，设
Figure FDA0002273053670000042
展开后的法方程为
V＝BX-l
则B为[(x0-xi)/D,(y0-yi)/D,(z0-zi)/D,-1],i＝0,1,K,n
l为[R-D]
求点云平面拟合方程，构造约束条件方程，其中
C＝[a b c]
Wx＝[-(d+ax0+by0+cz0)]
d为点(x0,y0,z0)到平面的距离
得法方程
Figure FDA0002273053670000043
X为方程解，Ks为联系数向量，BT、CT分别表示B、C的转置，解方程得
Figure FDA0002273053670000044
解方程得圆心坐标的改正值(dx0,dy0,dz0,dR)
用改正值修正初值，做迭代计算，当X的改正值小于限差结束迭代，这里限差取10-6；
半径R为各点到圆心的距离平均值，每一个坐标点到圆心的距离，再减去半径R就是它的挠度，平均挠度为所有挠度的平均值，最大挠度为所有挠度值中的最大值。
7.根据权利要求1所述的一种钢构件圆柱体拟合算法，其特征在于，所述步骤6中采用最小二乘法进行柱体圆度计算的方法如下：
取圆柱面上的点到轴线上的投影点的距离为圆柱半径，以此构造方程，具体方法如下：
设中心轴线方程为
x＝x0+at
y＝y0+bt
z＝z0+ct
式中(a,b,c)为空间直线的单位方向矢量，(x0,y0,z0)为空间直线上离原点最近的点，t为直线上任意点到(x0,y0,z0)的距离，为了唯一表示直线，定义a>0，若a＝0则b>0；若a＝0且b＝0，则c>0；a，b，c不可能同时为0；
过观测点Pi(xi,yi,zi)与中心轴线垂直的平面方程
ax+by+cz+Di＝0
其中Di＝-(axi+byi+czi)
将轴线方程带入得：
t＝-(ax0+by0+cz0+Di)＝a(xi-x0)+b(yi-y0)+c(zi-z0)
中心轴线与该垂直平面的交点坐标(xp,yp,zp)：
Figure FDA0002273053670000051
观测点Pi(xi,yi,zi)与交点坐标的距离为：
Figure FDA0002273053670000052
得到各点误差方程为
Figure FDA0002273053670000053
其中，R为圆柱半径，7个参数分别为(a,b,c,x0,y0,z0,R)，在以下两个条件下求7个参数
a2+b2+c2＝1
Figure FDA0002273053670000061
令
Figure FDA0002273053670000062
误差方程线性化
Figure FDA0002273053670000063
Figure FDA0002273053670000064
Figure FDA0002273053670000065
Figure FDA0002273053670000066
Figure FDA0002273053670000067
Figure FDA0002273053670000068
Figure FDA0002273053670000069
其中
Figure FDA0002273053670000071
Figure FDA0002273053670000072
Figure FDA0002273053670000073
Figure FDA0002273053670000074
初始坐标点可以取点云重心即所有坐标点的加权平均值，半径取坐标点到重心距离的平均值,解方程即得到7个参数的改正值，(a,b,c)取1，改正值即为初始值和迭代后值的差值，修正初值后迭代计算，迭代过程保证a、b、c不同时为0，当改正值小于限差10-6时，迭代结束。
CN201911112195.8A 2019-11-14 2019-11-14 一种钢构件圆柱体拟合方法 Active CN111027010B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201911112195.8A CN111027010B (zh) 2019-11-14 2019-11-14 一种钢构件圆柱体拟合方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201911112195.8A CN111027010B (zh) 2019-11-14 2019-11-14 一种钢构件圆柱体拟合方法

Publications (2)

Publication Number Publication Date
CN111027010A true CN111027010A (zh) 2020-04-17
CN111027010B CN111027010B (zh) 2023-09-22

Family

ID=70200155

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201911112195.8A Active CN111027010B (zh) 2019-11-14 2019-11-14 一种钢构件圆柱体拟合方法

Country Status (1)

Country Link
CN (1) CN111027010B (zh)

Cited By (17)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111795651A (zh) * 2020-07-06 2020-10-20 安徽工程大学 一种运用机械手臂测量大型回转体参数的方法及设备
CN112066893A (zh) * 2020-08-14 2020-12-11 苏州杰锐思智能科技股份有限公司 测量键盘键帽高度的方法及装置
CN112465979A (zh) * 2020-11-10 2021-03-09 国网内蒙古东部电力有限公司检修分公司 一种基于剖切面几何特征的圆柱面点云拟合方法
CN112652068A (zh) * 2021-01-30 2021-04-13 上海汇像信息技术有限公司 旋转柱体3d模型的圆台拟合方法
CN112815849A (zh) * 2021-02-19 2021-05-18 三门核电有限公司 一种基于激光跟踪的核电管道建模方法
CN112902844A (zh) * 2021-02-24 2021-06-04 三门核电有限公司 一种基于激光跟踪的设备圆形端口建模方法
CN113219426A (zh) * 2021-05-21 2021-08-06 上海航天电子通讯设备研究所 一种大尺寸sar天线装配测量方法
CN113251919A (zh) * 2021-04-28 2021-08-13 中铁第四勘察设计院集团有限公司 一种基于坐标变换的圆柱体点云数据拟合方法及装置
CN114297753A (zh) * 2021-12-21 2022-04-08 中铁七局集团有限公司 基于bim的压力钢管智能安装施工方法、装置及系统
CN114322802A (zh) * 2021-12-30 2022-04-12 苏州中科行智智能科技有限公司 一种基于3d线激光相机的线径测量方法
CN114740801A (zh) * 2022-03-21 2022-07-12 成都飞机工业（集团）有限责任公司 一种用于数控设备群协同生产线安装的基坐标系创建方法
CN114782315A (zh) * 2022-03-17 2022-07-22 清华大学 轴孔装配位姿精度的检测方法、装置、设备及存储介质
CN114877849A (zh) * 2022-05-27 2022-08-09 包头钢铁（集团）有限责任公司 一种全站仪测量环形加热炉上部钢结构圆度的方法
CN115235376A (zh) * 2022-09-23 2022-10-25 国网天津市电力公司电力科学研究院 一种非接触式电缆敷设质量检测方法及检测装置
CN115310028A (zh) * 2022-07-21 2022-11-08 成都飞机工业（集团）有限责任公司 一种盲孔轴线与平面夹角测量不确定度的计算方法
CN115994403A (zh) * 2023-03-22 2023-04-21 中国水利水电第七工程局有限公司 一种基于三维圆心拟合的护筒校核方法、装置和设备
CN117928680A (zh) * 2024-03-21 2024-04-26 青岛清万水技术有限公司 换能器自动定位方法、系统、电子设备及存储介质

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101387494A (zh) * 2008-10-06 2009-03-18 天津大学 大型隧道管片构件几何量测量装置和方法
JP2009097985A (ja) * 2007-10-17 2009-05-07 Enzan Kobo:Kk 位置座標の測量方法
CN103900528A (zh) * 2012-12-28 2014-07-02 朱志洁 一种大型构件三维空间剖面圆度检测方法
CN105067011A (zh) * 2015-09-15 2015-11-18 沈阳飞机工业（集团）有限公司 一种基于视觉标定及坐标转换的测量系统整体校准方法
CN107833249A (zh) * 2017-09-29 2018-03-23 南京航空航天大学 一种基于视觉引导的舰载机着陆过程姿态预估方法
CN110095060A (zh) * 2019-03-12 2019-08-06 中建三局第一建设工程有限责任公司 基于三维扫描技术的钢结构快速质量检测方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2009097985A (ja) * 2007-10-17 2009-05-07 Enzan Kobo:Kk 位置座標の測量方法
CN101387494A (zh) * 2008-10-06 2009-03-18 天津大学 大型隧道管片构件几何量测量装置和方法
CN103900528A (zh) * 2012-12-28 2014-07-02 朱志洁 一种大型构件三维空间剖面圆度检测方法
CN105067011A (zh) * 2015-09-15 2015-11-18 沈阳飞机工业（集团）有限公司 一种基于视觉标定及坐标转换的测量系统整体校准方法
CN107833249A (zh) * 2017-09-29 2018-03-23 南京航空航天大学 一种基于视觉引导的舰载机着陆过程姿态预估方法
CN110095060A (zh) * 2019-03-12 2019-08-06 中建三局第一建设工程有限责任公司 基于三维扫描技术的钢结构快速质量检测方法

Non-Patent Citations (4)

* Cited by examiner, † Cited by third party
Title
于志明: "论刚体对坐标轴转动惯量的平移和旋转变换" *
王君刚,王解先: "圆柱体的拟合与质量检测" *
褚建春,张泽峰: "空间任意方向圆柱面拟合方法" *
郝雪,毕超,刘京亮等: "精密测量中的圆柱拟合算法研究" *

Cited By (26)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111795651A (zh) * 2020-07-06 2020-10-20 安徽工程大学 一种运用机械手臂测量大型回转体参数的方法及设备
CN112066893B (zh) * 2020-08-14 2022-06-07 苏州杰锐思智能科技股份有限公司 测量键盘键帽高度的方法及装置
CN112066893A (zh) * 2020-08-14 2020-12-11 苏州杰锐思智能科技股份有限公司 测量键盘键帽高度的方法及装置
CN112465979A (zh) * 2020-11-10 2021-03-09 国网内蒙古东部电力有限公司检修分公司 一种基于剖切面几何特征的圆柱面点云拟合方法
CN112652068A (zh) * 2021-01-30 2021-04-13 上海汇像信息技术有限公司 旋转柱体3d模型的圆台拟合方法
CN112652068B (zh) * 2021-01-30 2022-03-11 上海汇像信息技术有限公司 旋转柱体3d模型的圆台拟合方法
CN112815849A (zh) * 2021-02-19 2021-05-18 三门核电有限公司 一种基于激光跟踪的核电管道建模方法
CN112815849B (zh) * 2021-02-19 2022-07-12 三门核电有限公司 一种基于激光跟踪的核电管道建模方法
CN112902844A (zh) * 2021-02-24 2021-06-04 三门核电有限公司 一种基于激光跟踪的设备圆形端口建模方法
CN113251919A (zh) * 2021-04-28 2021-08-13 中铁第四勘察设计院集团有限公司 一种基于坐标变换的圆柱体点云数据拟合方法及装置
CN113219426A (zh) * 2021-05-21 2021-08-06 上海航天电子通讯设备研究所 一种大尺寸sar天线装配测量方法
CN113219426B (zh) * 2021-05-21 2022-11-25 上海航天电子通讯设备研究所 一种大尺寸sar天线装配测量方法
CN114297753A (zh) * 2021-12-21 2022-04-08 中铁七局集团有限公司 基于bim的压力钢管智能安装施工方法、装置及系统
CN114322802A (zh) * 2021-12-30 2022-04-12 苏州中科行智智能科技有限公司 一种基于3d线激光相机的线径测量方法
CN114782315B (zh) * 2022-03-17 2024-07-09 清华大学 轴孔装配位姿精度的检测方法、装置、设备及存储介质
CN114782315A (zh) * 2022-03-17 2022-07-22 清华大学 轴孔装配位姿精度的检测方法、装置、设备及存储介质
CN114740801B (zh) * 2022-03-21 2023-09-29 成都飞机工业（集团）有限责任公司 一种用于数控设备群协同生产线安装的基坐标系创建方法
CN114740801A (zh) * 2022-03-21 2022-07-12 成都飞机工业（集团）有限责任公司 一种用于数控设备群协同生产线安装的基坐标系创建方法
CN114877849A (zh) * 2022-05-27 2022-08-09 包头钢铁（集团）有限责任公司 一种全站仪测量环形加热炉上部钢结构圆度的方法
CN114877849B (zh) * 2022-05-27 2024-04-30 包头钢铁（集团）有限责任公司 一种全站仪测量环形加热炉上部钢结构圆度的方法
CN115310028A (zh) * 2022-07-21 2022-11-08 成都飞机工业（集团）有限责任公司 一种盲孔轴线与平面夹角测量不确定度的计算方法
CN115310028B (zh) * 2022-07-21 2023-11-10 成都飞机工业（集团）有限责任公司 一种盲孔轴线与平面夹角测量不确定度的计算方法
CN115235376A (zh) * 2022-09-23 2022-10-25 国网天津市电力公司电力科学研究院 一种非接触式电缆敷设质量检测方法及检测装置
CN115994403A (zh) * 2023-03-22 2023-04-21 中国水利水电第七工程局有限公司 一种基于三维圆心拟合的护筒校核方法、装置和设备
CN117928680A (zh) * 2024-03-21 2024-04-26 青岛清万水技术有限公司 换能器自动定位方法、系统、电子设备及存储介质
CN117928680B (zh) * 2024-03-21 2024-06-07 青岛清万水技术有限公司 换能器自动定位方法、系统、电子设备及存储介质

Also Published As

Publication number Publication date
CN111027010B (zh) 2023-09-22

Similar Documents

Publication Publication Date Title
CN111027010B (zh) 一种钢构件圆柱体拟合方法
CN108759665B (zh) 一种基于坐标转换的空间目标三维重建精度分析方法
CN104019799B (zh) 一种利用局部参数优化计算基础矩阵的相对定向方法
CN109099883A (zh) 高精度大视场机器视觉测量与标定装置及方法
Xia et al. An accurate and robust method for the measurement of circular holes based on binocular vision
CN112697041B (zh) 一种基于蒙特卡洛法的装配位姿测量精度预评估方法
CN102778224B (zh) 一种基于极坐标参数化的航空摄影测量光束法平差的方法
US20070052974A1 (en) Three-dimensional modeling from arbitrary three-dimensional curves
CN106248014A (zh) 一种基于单相片的三维坐标测量方法及装置
CN109141266B (zh) 一种钢结构测量方法及系统
CN105354850B (zh) 基于电场性质的复杂曲面零件尺寸三维匹配检测方法
CN112775935A (zh) 一种基于末端误差检测信息子集的并联机器人标定方法
CN107480356A (zh) 基于catia和激光跟踪仪的部件设计检验一体化方法
CN113702994B (zh) 一种基于刚性约束的激光跟踪仪测量精度提升方法
CN103344252B (zh) 一种航空高光谱成像系统定位误差分析方法
Wang et al. A binocular vision method for precise hole recognition in satellite assembly systems
CN109883381A (zh) 一种关节式坐标测量机的三维空间大尺寸测量方法
Stępień et al. Dimensioning method of floating offshore objects by means of quasi-similarity transformation with reduced tolerance errors
CN115493616B (zh) 一种激光跟踪姿态角现场精度的评定方法
Acero et al. Evaluation of a metrology platform for an articulated arm coordinate measuring machine verification under the ASME B89. 4.22-2004 and VDI 2617_9-2009 standards
CN115493617B (zh) 一种激光跟踪姿态角现场精度评定系统
Cui et al. Novel method of rocket nozzle motion parameters non-contact consistency measurement based on stereo vision
Wang et al. Research on improving the measurement accuracy of the AACMM based on indexing joint
Wang et al. An accurate and stable pose estimation method based on geometry for port hoisting machinery
CN109059761A (zh) 一种基于eiv模型的手持靶标测头标定方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
CB02 Change of applicant information

Address after: 430000 Wuhan Donghu New Technology Development Zone, Wuhan, Hubei Province

Applicant after: WUHAN TIANHENG INFORMATION TECHNOLOGY Co.,Ltd.

Address before: 430000 building B26, financial port, Wuhan City, Hubei Province

Applicant before: WUHAN TIANHENG INFORMATION TECHNOLOGY Co.,Ltd.

CB02 Change of applicant information
GR01 Patent grant
GR01 Patent grant