-
本文中的研究对象为Leica At-930激光跟踪仪。以测量位置靶球中心为原点可以建立不确定度椭球坐标系(uncertainty ellipsoid coordinate system,UCS),在此坐标系(w, u, v)下不确定度椭球包含的范围可以描述为方程:
$ E(w, u, v)=\frac{w^{2}}{a^{2}}+\frac{u^{2}}{b^{2}}+\frac{v^{2}}{c^{2}}-1 \leqslant 0 $
(1) 式中,w, u, v是测量点在UCS三坐标轴w, u, v上的坐标值,a,b,c分别是该不确定椭球的三轴半轴长。图 1为UCS下不确定度椭球示意图。
在实际测量中,靶球的位置会随着测量目标的变化而变化,即各测量点分别存在其独立的UCS,因此, 测量数据往往采用统一的基于跟踪仪的测量系统坐标系(measurement coordinate system,MCS)。对于任意UCS其坐标轴单位向量在MCS上可以描述为:
$ \left\{\begin{array}{l} \boldsymbol{u}=\left[\begin{array}{lll} \sin \alpha & -\cos \alpha & 0 \end{array}\right]^{\mathrm{T}} \\ \boldsymbol{v}=\left[\begin{array}{lll} \cos \beta \cos \alpha & \cos \beta \sin \alpha & -\sin \beta \end{array}\right]^{\mathrm{T}} \\ \boldsymbol{w}=\left[\begin{array}{lll} \sin \beta \cos \alpha & \sin \beta \sin \alpha & \cos \beta \end{array}\right]^{\mathrm{T}} \end{array}\right. $
(2) 式中,α为方位角,β为天顶角,w为测量的激光方向,v为w与z平面上垂直于w的方向,u为遵循右手定则垂直于w与v的方向。由UCS到MCS的旋转矩阵则可以描述为:
$ {}_{\operatorname{UCS}}^{\operatorname{MCS}} \boldsymbol{R}=[\boldsymbol{w}, \boldsymbol{u}, \boldsymbol{v}] $
(3) 图 2为激光跟踪仪MCS与坐标系UCS的转换示意图。
对于测量点P,在其位置构建的UCS在MCS上的位置矢量可以描述为:
$ ^{\operatorname{MCS}}_{\operatorname{UCS}} \boldsymbol{P}_{O}=\left[x_{0}, y_{0}, z_{0}\right] $
(4) 式中,x0,y0,z0为在该位置累积测量采样分布的期望,以其近似为理论靶球中心。据此两坐标系变换的变换矩阵可以描述为:
$ ^{\text {MCS}}_{\text {UCS}} \boldsymbol{T}=\left[ {\begin{array}{*{20}{c}} {_{{\rm{UCS}}}^{{\rm{MCS}}}\mathit{\boldsymbol{R}}}&{_{{\rm{UCS}}}^{{\rm{MCS}}}{\mathit{\boldsymbol{P}}_O}}\\ {\begin{array}{*{20}{c}} 0&0&0 \end{array}}&1 \end{array}} \right] $
(5) $ { }_{\text {MCS}}^{\text {UCS}} \boldsymbol{T}=\left[ {\begin{array}{*{20}{c}} {_{{\rm{UCS}}}^{{\rm{MCS}}}{\mathit{\boldsymbol{R}}^{\rm{T}}}}&{ - _{{\rm{UCS}}}^{{\rm{MCS}}}{\mathit{\boldsymbol{R}}^{\rm{T}}}{\kern 1pt} {\kern 1pt} _{{\rm{UCS}}}^{{\rm{MCS}}}{\mathit{\boldsymbol{P}}_O}}\\ {\begin{array}{*{20}{c}} 0&0&0 \end{array}}&1 \end{array}} \right] $
(6) 据此依据通用坐标系转换公式:
$ \left[\begin{array}{c} { }^{A} \boldsymbol{P} \\ 1 \end{array}\right]={ }_{B}^{A} \boldsymbol{T} \cdot\left[\begin{array}{c} { }^{B} \boldsymbol{P} \\ 1 \end{array}\right] $
(7) 式中,A, B作为通用坐标系表达均可以替换为MCS或UCS,据此任意测量位置在MCS下的累积采样位置数据均可以通过坐标变换转换为UCS下进行标准不确定度椭球的分布分析;而分析获得的不确定度椭球同样可以转换为MCS下测量场不同位置的椭球分布。
-
基于固定设备与采样点位的单点重复位置测量与基于蒙特卡洛等算法的单点位随机仿真测点均是获取单点位置不确定度分布的重要方法。这些方法可以获得对于单点的大量测量/仿真位置数据,实测过程中由于场地与设备等因素的扰动会获得极少数较大的误差点,而在算法仿真中同样会由于必要的扰动函数产生少量较大误差点。在进行不确定度包络椭球计算时,这些较大误差点会对包络计算形成干扰,使得最终获得的包络椭球过大。本研究中,引入了孤立森林算法来对测量/仿真数据进行采样筛选。
基于本课题问题研究,孤立森林算法模型可以通用化描述为{D, F, H, fs, rs, Nt, Nsub, nDFR},其中D={P1, P2, …, PN}为数据空间,包含所有的待分析筛选的数据,其数据个数为N,对于每个数据点P,其维度为d,在本研究问题中,D即为所有测量/仿真位置数据的集合,维度d=3;F={T1, T2, …, TNt}为孤立森林算法划分出的森林,其中T为森林中的独立数据树,Nt为每座森林中树的个数;对于每个数据树T包含集合Dsub={Pr, 1, Pr, 2, …, Pr, i…Pr, Nsub}的全部数据,并通过孤立函数fs划分为二叉树结构,Nsub为每棵树的数据个数,Pr, i为数据空间中随机出的第i个数据;Dsub为每棵树的数据集合,其为D的随机样本量是Nsub的子集;孤立函数fs可以描述为:
$ \left\{\begin{aligned} \boldsymbol{D}_{\mathrm{sub}, \mathrm{l}, h+1}=& f_{\mathrm{s}, \mathrm{l}}\left(\boldsymbol{D}_{\mathrm{sub}, h}, S_{\mathrm{r}, h}\right)=\left\{\boldsymbol{P}_{\mathrm{r}, i} \mid \boldsymbol{P}_{\mathrm{r}, i}\left(r_{d}\right) \leqslant\right.\\ &\left.S_{\mathrm{r}, h}, \boldsymbol{P}_{\mathrm{r}, i}\left(r_{d}\right) \in \boldsymbol{D}_{\mathrm{sub}, h}\right\} \\ \boldsymbol{D}_{\mathrm{sub}, \mathrm{r}, h+1}=& ∁\ \ \ _{{D_{\mathrm{sub}, h}}} \boldsymbol{D}_{\mathrm{sub}, \mathrm{l}, h+1} \end{aligned}\right. $
(8) 式中,Dsub, l, h与Dsub, r, h分别代表该数据树在距离根部深度为h(即h层)的左右数据子集;Sr, h为该数据树在h层数据孤立分割的标准,其取值为[min(Dsub, h), max(Dsub, h)]范围内的随机值;因为D为d维的数据集,在孤立分割时,取随机整数rd∈[1, d]作为分割判定时的比较维度,Pr, i(rd)即为数据Pr, i在比较维度rd下的数据分量。在筛选计算中, 将后续的Dsub, l, h+1, 或者Dsub, l, h+1分别作为新的子数据集Dsub, h+1, 采用孤立函数fs进行递归迭代划分,直到最终层次的Dsub样本量为1,即完成数据树的构建。H为每个数据点在森林中每个树上距离根部的深度的平均值的集合;rs为算法的筛选次数,即森林的个数;nDFR则为每次筛选时判定标准,即数据距离根节点的深度(deep from root, DFR),即当H(Pi) < nDFR时,Pi为较大误差点, 当被舍弃。最终经过孤立森林算法筛选后的数据空间D0则为后续不确定度椭球模型的数据基础。
-
通过坐标转换与数据筛选后的测量/仿真位置数据可以通过PSO基于(1)式进行最小包络椭球模型的建立。在UCS下,椭球中心为原点,因此最小包络的控制因素则为椭球三轴的半轴长。定义粒子X={a, b, c},解的搜索空间为3维。在解空间随机生成n个粒子,每个粒子Xi均为该系列数据的潜在一个包络椭球解。随着迭代时刻t的推进,每个粒子均在解空间以一定速度Vi运动探索以获取最优的包络解。每个粒子的运动规则遵循:
$ \begin{gathered} \boldsymbol{V}_{i, t+1}=\omega \cdot \boldsymbol{V}_{i, t}+c_{1} r_{1}\left(L_{i, t}-\boldsymbol{X}_{i, t}\right)+ \\ c_{2} r_{2}\left(\boldsymbol{G}_{t}-\boldsymbol{X}_{i, t}\right) \end{gathered} $
(9) $ \boldsymbol{X}_{i, t+1}=\boldsymbol{X}_{i, t}+\boldsymbol{V}_{i, t+1} $
(10) 式中,Xi, t与Vi, t为第i个粒子在t时刻的位置与运动速度;Li, t为第i个粒子到t时刻其个体运动轨迹中的最优解,代表着粒子的个体认知;Gt为t时刻所有粒子运动轨迹中的最优解,代表着粒子的群体交流与共识;ω为惯性因子,c1和c2为学习因子,控制着粒子向自身经验和群体经验学习的倾向,r1和r2为随机因子,为[0, 1]之间的随机数。在每次运动迭代完成后,粒子Xi, t均以适应值函数fe来进行评价,通常适应值越小,该粒子距离理论最优解越近。本研究中以包络椭球解的体积作为评价标准,因此对于粒子Xi={ai, bi, ci}:
$ f_{\mathrm{e}}\left(\boldsymbol{X}_{i}, r_{\mathrm{e}}\right)=\left\{\begin{array}{l} \frac{4}{3} {\rm{ \mathit{ π} }} a_{i} b_{i} c_{i},\left(n_{\mathrm{in}} \geqslant N_{0} \cdot r_{\mathrm{e}}\right) \\ -1,(\text {others}) \end{array}\right. $
(11) 式中,re为椭球模型对数据D0的包络比例,N0为D0中数据的个数,nin为D0中处于包络范围内的数据个数,可以描述为:
$ n_{\mathrm{in}}=\overset{N_0}{\underset{i}{\operatorname{argnum}}}\left(\boldsymbol{P}_{i}, E_{\boldsymbol{X}_{i}}\left(\boldsymbol{P}_{i}\right) \leqslant 0\right),\left(\boldsymbol{P}_{i} \in \boldsymbol{D}_{0}\right) $
(12) 式中,argnum函数的输出结果为满足后续条件的数据Pi个数,EXi(Pi)为点Pi代入粒子Xi条件下的包络椭球模型E(见(1)式)计算其相对于椭球包络面的位置。通过适应值函数,在每次粒子群运动迭代后对个体历史最优解Pi与全局最优解G进行更新:
$ \boldsymbol{L}_{i, t+1}=\left\{\begin{array}{l} \boldsymbol{X}_{i, t+1},\left(f_{\mathrm{e}}\left(\boldsymbol{L}_{i, t}\right)>f_{\mathrm{e}}\left(\boldsymbol{X}_{i, t+1}\right),\right. \\ \ \ \ \ \ \ \ \ \ \ \ \ \left.f_{\mathrm{e}}\left(\boldsymbol{X}_{i, t+1}\right) \neq-1\right) \\ \boldsymbol{L}_{i, t},(\text {others}) \end{array}\right. $
(13) $ \boldsymbol{G}_{t}=\overset{n}{\underset{i=1}{\operatorname{argmin}}}\left(\boldsymbol{L}_{i, t}, f_{\mathrm{e}}\left(\boldsymbol{L}_{i, t}\right)\right) $
(14) 式中,argmin函数的输出结果为使得适应值函数fe最小的粒子Li。
在本研究问题算法执行流程为对所有粒子在粒子空间进行随机位置初始化,通过(9)式~(14)式对每个粒子的运动位置与速度的迭代计算完成最优包络的探索。区别于传统优化问题,最小包络的探索在(11)式~ (13)式中存在对包络数据比例的判定,因此为实现优化探索的准确性与效率,需要选择合适的初始化粒子位置范围与粒子陷入过小包络模型局部位置后的扰动逃逸模型。粒子初始化位置为:
$ \boldsymbol{X}_{i}=\left(w_{\max }, u_{\max }, v_{\max }\right) \cdot r_{x}\left(l_{0}\right) $
(15) 式中,wmax,umax,vmax为UCS下所有点在三坐标轴下的绝对值极值,rx为[1, l0]范围的随机值,l0为空间放大系数。扰动逃逸函数在当前粒子无法满足包络要求时可以替代其位置更新函数,描述为:
$ \boldsymbol{X}_{i, t+1}=\boldsymbol{X}_{i, t} \cdot e_{\mathrm{s}}+\boldsymbol{V}_{i, t+1} \mid f_{\mathrm{e}}\left(\boldsymbol{X}_{i, t}\right)=-1 $
(16) 式中,es为逃逸系数。
激光测量系统不确定度最小包络椭球模型研究
Research on uncertainty minimum ellipsoid envelope model of laser measurement system
-
摘要: 为了对激光测量系统的测量误差3维空间分布进行有效评估, 以特定点测量或仿真的大量位置数据为基础, 采用孤立森林算法对点云进行了异常数据筛除。基于误差椭球理论, 引入粒子群优化算法, 针对有效数据建立了最小包络椭球的不确定度模型; 采用测量场与单点不确定度的坐标系变换, 将不确定度最小包络椭球模型应用于测量场景内不确定度场的空间分布分析; 通过单点以及10m量级范围空间场景实测数据的测试, 该模型可以高效地筛选有效采样数据, 并依据需求进行不同程度的最小包络椭球计算, 得到相应的不确定度。结果表明, 基于测量位置数据, 该模型可以高效准确地描述单点位置的3维不确定度范围, 并能够有效地再现测量空间内的不确定度分布, 在4.7m的测量距离、94.2%的筛选后有效数据、97.5%的包络比例下, 计算获得不确定度范围为三轴长4.95μm, 18.39μm和30.53μm的椭球。该最小包络椭球不确定度模型在基于实测的理论模型验证、设备状态与测量场景环境分析, 以及测量布局设计等方面具有着重要的价值。Abstract: In order to effectively evaluate the three-dimensional spatial distribution of measurement errors of the laser measurement system, a new uncertainty model based on the calculation of the minimum envelope ellipsoid was proposed. Based on the measured or simulated location data, the isolated forest algorithm was introduced to filter out abnormal data in the point cloud. With the valid data, the minimum envelope ellipsoid uncertainty model was established based on the particle swarm optimization and the error ellipsoid theory. By the coordinate system transformation between the measurement field and the single point uncertainty, the minimum envelope ellipsoid model was applied to the spatial uncertainty distribution analysis. Through the test of the measured data of a single point and a 10m-level space scene, the model can efficiently screen the sampling data and perform different levels of minimum envelope ellipsoid calculations according to the requirements. And then the corresponding uncertainty can be obtained. The results show that based on the measurement position data, the model can efficiently and accurately describe the three-dimensional uncertainty range of a single point position, and can effectively reproduce the uncertainty distribution in the measurement space. With the experimental conditions of a measurement distance 4.7m, the effective data after screening 94.2%, and an envelope ratio 97.5%, an ellipsoid with an uncertainty range of three-axis length 4.95μm, 18.39μm, 30.53μm is obtained by the model calculation. The minimum envelope ellipsoid uncertainty model has important value in many aspects, such as theoretical model verification based on actual measurement, equipment state and measurement environment analysis, and measurement layout design.
-
-
[1] SWYT D A. Length and dimensional measurements at NIST[J]. Journal of research of the National Institute of Standards and Techno-logy, 2001, 106(1): 1-23. doi: 10.6028/jres.106.002 [2] SHI H Y, GUO T, WANG D, et al. Power line suspension point location method based on laser point cloud[J]. Laser Technology, 2020, 44(3): 364-370(in Chinese). [3] SUI Sh Ch, ZHU X Sh. Digital measurement technique for evaluating aircraft final assembly quality[J]. Scientia Sinica Technologica, 2020, 50: 1449-1460(in Chinese). doi: 10.1360/SST-2020-0049 [4] DENG Zh P, LI S G, HUANG X. Coordinate transformation uncertainty analysis and reduction using hybrid reference system for aircraft assembly[J]. Assembly Automation, 2018, 38(4): 487-496. doi: 10.1108/AA-08-2017-097 [5] MEI Z, MAROPOULOS P G. Review of the application of flexible, measurement-assisted assembly technology in aircraft manufacturing[J]. Proceedings of the Institution of Mechanical Engineers, 2014, B228(10): 1185-1197. [6] CHEN Z H, DU F Zh, TANG X Q. Research on uncertainty in mea-surement assisted alignment in aircraft assembly[J]. Chinese Journal of Aeronautics, 2013, 26(6): 1568-1576. doi: 10.1016/j.cja.2013.07.037 [7] DENG Zh Ch, WU Zh Y, YANG J G. Point cloud uncertainty analysis for laser radar measurement system based on error ellipsoid model[J]. Optics and Lasers in Engineering, 2016, 79: 78-84. doi: 10.1016/j.optlaseng.2015.11.010 [8] COX M G, HARRIS P M. Measurement uncertainty and traceability[J]. Measurement Science and Technology, 2006, 17(3): 533-540. doi: 10.1088/0957-0233/17/3/S13 [9] ZHANG F M, QU X H. Fusion estimation of point sets from multiple stations of spherical coordinate instruments utilizing uncertainty estimation based on Monte Carlo[J]. Measurement Science Review, 2012, 12(2): 40-45. [10] REN Y, LIN J R, ZHU J G, et al. Coordinate transformation uncertainty analysis in large-scale metrology[J]. IEEE Transactions on Instrumentation and Measurement, 2015, 64(9): 2380-2388. doi: 10.1109/TIM.2015.2403151 [11] ZHU J K, LI L J, LIN X Zh. Research on the measurement field planning of lidar measurement system[J]. Laser Technology, 2021, 45(1): 99-104(in Chinese). [12] BERGSTRÖM P, EDLUND O. Robust registration of point sets using iteratively reweighted least squares[J]. Computational Optimization and Applications, 2014, 58(3): 543-561. doi: 10.1007/s10589-014-9643-2 [13] WANG Q, HUANG P, LI J X, et al. Uncertainty evaluation and optimization of INS installation measurement using Monte Carlo method[J]. Assembly Automation, 2015, 35(3): 221-233. doi: 10.1108/AA-08-2014-070 [14] JIN Zh J, YU C J, LI J X, et al. Configuration analysis of the ERS points in large-volume metrology system[J]. Sensors, 2015, 15(9): 24397-24408. doi: 10.3390/s150924397 [15] PREDMORE C R. Bundle adjustment of multi-position measurements using the Mahalanobis distance[J]. Precision Engineering, 2010, 34(1): 113-123. doi: 10.1016/j.precisioneng.2009.05.003 [16] CALKINS J M. Quantifying coordinate uncertainty fields in coupled spatial measurement systems[D]. Virginia, USA: Virginia Polytechnic Institute and State University, 2002, 1: 48. [17] LIU Y, SUN Sh Y. Laser point cloud denoising based on principal component analysis and surface fitting[J]. Laser Technology, 2020, 44(4): 497-502(in Chinese). [18] CHEN H S, MA H Zh, CHU X N, et al. Anomaly detection and critical attributes identification for products with multiple operating conditions based on isolation forest[J]. Advanced Engineering Informatics, 2020, 46: 101139. doi: 10.1016/j.aei.2020.101139 [19] SUSTO G A, BEGHI A, MCLOONE S. Anomaly detection through on-line isolation Forest: An application to plasma etching[C] // Proceeding of the 28th Annual SEMI Advanced Semiconductor Ma-nufacturing Conference (ASMC). New York, USA: IEEE, 2017: 89-94. [20] LIU F T, TING K M, ZHOU Zh H. Isolation-based anomaly detection[J]. ACM Transactions on Knowledge Discovery from Data, 2012, 6(1): 1-39. [21] CHEN W R, YUN Y H, WEN M, et al. Representative subset selection and outlier detection via isolation forest[J]. Analytical Methods, 2016, 8(39): 7225-7231. doi: 10.1039/C6AY01574C [22] KENNEDY J, EBERHART R. Particle swarm optimization[C] // Proceeding of IEEE International Conference on Neural Networks. New York, USA: IEEE, 1995: 1942-1948. [23] POLI R, KENNEDY J, BLACKWELL T. Particle swarm optimization[J]. Swarm Intelligence, 2007, 1(1): 33-57. doi: 10.1007/s11721-007-0002-0 [24] ALAM S, DOBBIE G, KOH Y S, et al. Research on particle swarm optimization based clustering: A systematic review of literature and techniques[J]. Swarm and Evolutionary Computation, 2014, 17: 1-13. doi: 10.1016/j.swevo.2014.02.001 [25] LI Y, XING Y, FANG C, et al. An experiment-based method for focused ion beam milling profile calculation and process design[J]. Sensors and Actuators, 2019, A286: 78-90. [26] LI Y, GOSÁLVEZ M A, PAL P, et al. Particle swarm optimization-based continuous cellular automaton for the simulation of deep reactive ion etching[J]. Journal of Micromechanics and Microengineering, 2015, 25(5): 055023. doi: 10.1088/0960-1317/25/5/055023 [27] CHEN S Q, ZHANG H Y, ZHAO Ch M, et al. Point cloud registration method based on particle swarm optimization algorithm improved by beetle antennae algorithm[J]. Laser Technology, 2020, 44(6): 678-683(in Chinese).