A procedure for the evaluation and compensation ... - Semantic Scholar
Noname manuscript No. (will be inserted by the editor)
A procedure for the evaluation and compensation of form errors by means of global isometric registration with subsequent local reoptimization Laura Klein · Tobias Wagner · Christoph Buchheim · Dirk Biermann
October 11, 2013
Abstract Stresses remaining in the component after sheet metal forming processes can result in complex form errors, such as springback and torsions. In order to compensate these process-induced deformations, the local and global deformations have to be analyzed. Hence, an appropriate comparison between the actually manufactured and the target design is required. For this purpose, the surface of the actual workpiece is scanned and the so-obtained scan points have to be assigned to corresponding points of the target shape defined by the workpiece model. From these correspondences, a field of deformation vectors can be computed which represents the basis for the compensation strategy. The task of finding appropriate correspondences is called registration. It is usually solved using rigid transformations, i. e., translation and rotation. Due to the locality, strength and complexity of the deformations, rigid transformations are usually not sufficient. As a more flexible alternative, a procedure for non-rigid registration is presented in this paper. Therein, isometry, i.e., the conservation of distances between corresponding points within an appropriate neighborhood structure, is defined as the objective function. The procedure consists of three steps: definition of the neighborhood structure, global registration, and local reoptimization. The main focus of the paper is set to the latter, where an adapted gradient descent method also allowing projections into the triangles of the target shape is presented Christoph Buchheim · Laura Klein Fakult¨ at f¨ ur Mathematik, TU Dortmund, Vogelpothsweg 87, 44227 Dortmund E-mail: {christoph.buchheim,laura.klein}@tu-dortmund.de Dirk Biermann · Tobias Wagner Institut f¨ ur Spanende Fertigung, TU Dortmund, Baroper Straße 301, 44227 Dortmund E-mail: {biermann,wagner}@isf.tu-dortmund.de
and experimentally validated. With these three steps, an assignment between both shapes can be calculated, even for strong local deformations and coarse triangular meshes representing the workpiece model. Keywords Forming · Springback analysis · Non-rigid isometric registration · Local optimization · Quadratic assignment problem 1 Introduction An undesired effect in sheet metal forming is the springback of the deformed sheet after the process [14]. For controlling this springback, the local deformations have to be measured and analyzed to allow respective process and tool modifications to be derived [10, 22]. A modern approach for automatically determining a continuous deformation of the tool or the NC paths for toolmaking has already been presented [4, 5]. In this approach, the local deformations of the workpiece are directly used for adjusting the design of the tool, such that the sheet is overbent in this regions. Hence, the aim of this paper is to find appropriate correspondences to determine the local deformations as basis for the subsequent process and tool modifications. In order to find appropriate correspondences, methods to calculate assignments between the scan points measured on the actual workpiece and the surfaces of the target geometry are described. These methods define and solve an optimization problem formalizing the requirements of the practical application. They are organized within a three-step procedure. The improvements of this procedure with respect to stadard methods are based on two advantages: 1. The special characteristics of springback and torsions are considered by not minimizing the distance
2
L. Klein, T. Wagner, C. Buchheim, D. Biermann
(a)
(b)
Fig. 1: Comparison of the assignments (black lines) to the yellow target shape obtained by minimizing (a) the distance or (b) the isometry error (cf. equation 3) of the assigned point pairs. between both shapes. Instead, the distances between neighboring pairs of points on either shape are optimized in order to obtain an isometric mapping. 2. The discrete scan points representing the manufactured shape can be assigned to arbitrary points on the surface of the workpiece model. The first advantage emerges from the fact that distancebased assignments may lead to clustered mappings from the curves or kinks introduced by springback. Considering for instance the convex curve shown in Fig. 1, the mapping based on a minimization of the distance assigns a large part of the points on the knee of the workpiece to the flat part of the yellow target shape, which obviously does not represent an appropriate mapping. The second advantage is based on a higher flexibility with respect to the representation of the target geometry. As mentioned before, the presented registration approach reduces the problem to three subproblems which can be handled subsequently. First of all, a reasonable neighborhood for assessing the isometric matching is required. Secondly, an initial assignment is computed using a global approach. At this stage, a fixed set of points on the target model is considered. Thirdly, a local reoptimization based on a gradient descent method refines the initial solution based on the new objective function. In this step, the scan points can be arbitrarily assigned to the surface of the target geometry. Based on the proposed reduction to standard problem formulations, established algorithms for the neighborhood structure and the quadratic assignment problem can be adopted. Hence, the focus of this paper is set on the local reoptimization step. In the experiments, a simple heuristic for the neighborhood structure and an assignment based on a previously developed algorithm are applied. The procedure leads to significant improvements with respect to the assignment of the actual and the target shape.
2 Related work A comprehensive survey of rigid and non-rigid registration methods including those approaches has been published very recently [20]. Herein, only the most relevant methods with respect to the here given proposal are provided. The summary is categorized into methods performing either a rigid or a non-rigid transformation. As the former improve the matching only by translation and rotation of the actual shape, these methods can be mainly distinguished based on the way of generating and refining the correspondences between the points of the actual and the target shape for the distance computation. The definition of a neighborhood structure is not required. The latter methods are more complex and usually have many degrees of freedom. They therefore combine aspects of all three steps making up the proposed procedure.
2.1 Rigid registration The advantage of the restriction to rigid transformations is the reduction of degrees of freedom and the availability of a suitable analytical formulation, which significantly increases the efficiency of the methods. The rigid registration problem is mainly solved by iterative methods. In a first step, pair-wise correspondences between the points of the actual workpiece and the points of the target design are determined. In the second step, the optimal rigid transformation with respect to sum of squared distances between these pairs is computed. This is possible in an efficient way by means of a quaternion-based algorithm, see, e.g., [12]. After the transformation, the correspondences are refined and another transformation is computed. These two steps are repeated until the termination criterion is met, e.g., until the improvement falls below a predefined threshold.
Evaluation and compensation of form errors by means of isometric registration with local reoptimization
A very common approach is Iterative Closest Point (ICP) [3, 8], where each scan point of the actual shape is assigned to its closest point on the target design, either by projection on the parametrically defined workpiece surface [19] or by minimizing the distance over all points in case of discrete point data or triangulations [21]. This approach provides good solutions in cases of a well-chosen initial solution and small deviations between the actual and the target shape. Otherwise, ICP can lead to bad matchings because of local minima for inappropriate correspondences. To avoid this problem, an alternative approach was developed in the Collaborative Research Center (SFB 708) [5]. In this method, the model of the target design is discretized, i.e., replaced by a second point cloud, and a distance-minimal bipartite matching between these two point clouds is calculated. If both point clouds are similar and reasonably uniform, this approach can improve the matchings of ICP. The main weakness of this approach is its requirement for an appropriate discretization of the target shape. 2.2 Non-rigid registration
3 Procedure for non-rigid isometric registration In the following, it is assumed that the actual shape is represented by a point cloud, which can, e.g., be received by optical or tactile scanning of the workpiece. The model of the target shape is represented in terms of a triangular mesh. Even in case of a parametrical model, most CAD tools can export such kind of representation. The objective of the non-rigid registration approach is a mapping of the points within the cloud to their respective correspondences on the triangular mesh. The above mentioned disadvantages of current registration methods arise from indirect and heuristic approaches for the isometric registration. In this paper, a new approach is proposed by exactly formulating the isometry within a proper objective function. This makes the registration task accessible to methods from mathematical optimization [17]. As basis for computing the objective function, the differences between the Euclidean distances (2-norm) of two points pi , pj ∈ P of the point cloud and two corresponding points zi and zj on the target shape are computed, i.e. d(pi , pj ; zi , zj ) := |kpi − pj k2 − kzi − zj k2 |.
One possibility to relieve the restriction to rigid transformations is to allow the scan points to be independently mapped to the target shape by affine transformations. To prevent the geometry from degenerating, regulative constraints or repair steps minimizing the difference between transformations of neighboring points are added to the registration problem [1]. As a result, a quadratic minimization problem, yielding a transformation matrix for every single point measured on the surface of the manufactured workpiece, has to be solved. The number of variables increases significantly and efficient optimal algorithms are no longer available. To overcome the runtime problem, many approaches are based on local geometric features and a grouping of the scan points into clusters that undergo the same rigid transformation [5]. These approaches, however, require a merge of the clusters after the individual transformations. Other methods explicitly optimize isometry-based objective functions [6, 9, 16], but use more complex distance formulations, such as geodesic distance [9] or Gromov-Hausdorff distance [6, 16]. Based on the distance pairs, the points are embedded into a Euclidean space by multidimensional scaling (MDS) or enhanced variants thereof [6]. As an effect, the non-rigid registration problem transforms into a rigid one on the embedding space. The computation of the geodesic distance tableaux and the Gromov-Hausdorff distance, however, is extremely inefficient [6, 18].
3
(1)
With points pi fixed, the positions of the points zi on the target shape are optimized, i.e., arbitrarily moved on the target shape. The actual optimization objective is to minimize the sum of squared errors: X min d(pi , pj ; zi , zj )2 . (2) ij
The steps and the dependencies of the optimization procedure are summarized in Fig. 2. The choice of a neighborhood is worthwhile because of the new objective function which takes distance differences of two points on the shapes into account. As the target is to find a preferably isometric and non-rigid solution it is necessary to preserve distances of close-by situated points on the one hand, but not necessarily of distant points on the other hand. This idea is integrated by means of the neighborhood structure. In the second step, on a thinned set of points with the outlined neighborhood structure, an initial global assignment is calculated with the purpose of avoiding rough errors. In particular, this global assignment aims at preventing the objects to twist against each other, since later it is hard to overcome local minima and to disentangle the correspondences. The third step is a local reoptimization step improving the global registration. The assignment of the scanned points is not longer restricted to the discrete point set on the target shape as assignments are moved to arbitrary points on the respective surface.
4
L. Klein, T. Wagner, C. Buchheim, D. Biermann
Neighborhood structure defines objective function
[2-node-connected graph]
Global registration [Quadratic assignment problem]
defines objective function Local reoptimization
yields initial solution
[Gradient descent method]
More sophisticated approaches come from combinatorial optimization and graph theory. Consider a complete graph G = (P, E), whose vertices represent the points pi . The edge set F of any subgraph H = (P, F ) of G defines a selection of pairs, namely the points pi , pj being connected by an edge in F . Therefore, the neighborhood structure is formulated via an expedient subgraph H whose edge set F defines the pairs being considered in the objective function: X min d(pi , pj ; zi , zj )2 . (3) ij∈F
Fig. 2: Overview of the three steps and the related optimization problems including their dependencies.
The great benefit is the adjustment of discretization errors coming from thinning and the discrete assignment problem in the global step. In the following, the three steps are described in detail.
3.1 Neighborhood structure For a reasonable modeling of both, the non-rigidness and the isometry, it is not sensitive to take into account distances between all pairs of points on the manufactured workpiece. Different bendings naturally lead to different distances between two widely spaced points on the manufactured and their corresponding points on the target shape, therefore, their distances must not be preserved. Rather it is important to preserve the distances between two near-by situated points to preserve multiple mappings to angular points. A second technical argument for not taking into account all pairs of points is the high complexity of the quadratic assignment problem, which is proposed to be applied for finding an appropriate initial global solution and which is explained in Subsection 3.2. For both reasons, the use of a reduced neighborhood structure is proposed. As one obvious possibility for a neighborhood structure one can consider only pairs of points (pi , pj ) with a distance |pi − pj | < r, with a specified radius r > 0. Unfortunately, the choice of r is problematic. On the one hand, a too large radius value will not resolve the mentioned disadvantages, especially not the undesired rigidness. The choice of a too small value for r on the other hand runs the risk that some points might have too less neighbored points or even be isolated at all. Then again, isometry cannot be guaranteed.
To reduce the rigidness of the assignment, a minimal number or length of the edges chosen in H is demanded. In order to prevent the graph from internal twisting, however, each pair of nodes should be connected by at least two node disjoint paths, meaning that this pair is still connected if one arbitrary node is removed. This leads to the requirement of an at least 2-node-connected subgraph with minimal edge length. This problem is NP-hard [11], i.e., verifiable optimum solutions can only be found for very small graphs. Fortunately, there is no need for an exact solution of this problem for the given purpose. For a reasonable neighborhood structure, it is also possible to use graphs being denser than 2-nodeconnected, even though less and shorter edges improve computation times and the non-rigidness of the final solution. 3.2 Quadratic assignment problem For an initial global registration, only the triangle vertices Q := {q1 , . . . , qm } of the target mesh are allowed to be chosen as correspondences. As a consequence, each point pi ∈ P has to be assigned to exactly one point qj ∈ Q and vice versa. This leads to the restriction that both point sets P and Q have to be of the same cardinality, i.e., of the same number of points, which can usually be achieved by a curvature-adaptive reduction [13] of the scan point set P . Using binary variables xik ∈ {0, 1} which indicate whether point pi is assigned to point qk (xik = 1) or not (xik = 0), this global registration step can be modeled as a quadratic assignment problem (QAP) [7, 15]: XX min d(pi , pj ; qk , ql )2 xik xjl (4) ij∈F kl
subject to X xik = 1
for all i
k
X
xik = 1
for all k
xik ∈ {0, 1}
for all i, k
i
(5)
Evaluation and compensation of form errors by means of isometric registration with local reoptimization
Fig. 3: Assignment of the scan points on the manufactured shape to the vertices of the triangular mesh representing the target shape (gray). In the areas of a rough resolution of the mesh (black lines), high errors with respect to isometry (red color) can be observed. In the objective function (4) of the QAP, the respective distance costs d(pi , pk ; qj , ql )2 are used, where the arbitrary points zi of equation (3) are replaced by the triangle vertices qi . The distance costs are taken into account if two requirements are satisfied: 1. The points pi and pk have to be connected by the edges of the neighborhood defining graph. Hence, the first sum only includes the points connected by edges of the reasonable neighborhood structure F defined before. 2. The point pi on the workpiece is assigned to the point qk on the target shape, and pj is assigned to ql , respectively. This is encoded by the product of the binary variables xik and xjl . The constraints (5) express the requirement of the oneto-one assignment of a point pi ∈ P to a point qj ∈ Q. Also the QAP is NP-hard, but it is well studied in the field of discrete optimization. A comprehensive overview is given in [2], where exact and heuristic approaches are discussed. For the large point clouds faced in practice, mainly heuristic methods are applied.
3.3 Local reoptimization For the computation of a global assignment based on QAP solvers, a discretization of the target shape is necessary, as only correspondences between fixed sets can be considered – one for the scan points on the workpiece and one for the vertices of the triangular mesh representing the target shape. Moreover, both sets have to be of the same cardinality. In this context, the accuracy of the discretization of the target shape can become a problem. If the resolution is too fine, the number of variables becomes too high and the problem is not solvable in practice. A rough discretization on the other hand excludes isometric registrations from the outset, in particular if the points on the manufactured shape are spread differently compared to the triangular vertices
5
of the target mesh. Fig. 3 illustrates this effect on the basis of a rough triangular mesh of a 2.5-dimensional metal sheet, which is also used for computational tests in Section 4. In order to overcome the discretization problems and to restore isometry, a reoptimization of the discrete assignment is presented in the following. Its advantage consists in a mapping to arbitrary points on the inner surface of the triangles representing the target mesh. This allows an adaptation with respect to differences in dispersion and density of the initial point sets. To accomplish this, the problem X min d(pi , pj ; zi , zj )2 ij∈F (6) ˜ for all pi ∈ P zi ∈ Q ˜ now denotes the infinite set of all is solved, where Q positions on these triangle surfaces. For solving the problem defined in (6), an adapted iterative gradient descent method was developed. Its sequence of operations is depicted in Fig. 4. In every iteration step, each point pi ∈ P is related to exactly one triangle of the target mesh, whereby injectivity is not mandatory. In the descent step, every correspondence point zi is moved in the direction h = −∇f of steepest ascent with respect to the objective function f minimized in equation 6. If zi would move over the border of its related triangle, it is projected to the border for this iteration and it is related to the neighbored triangle in the next iteration step. In the case where the triangle border also is the border of the mesh, i.e. no neighboring triangle exists, the point remains within the current triangle. The single assignment steps of the gradient descent method are further illustrated in Fig. 5. In each of the four pictures, the four vertices of the upper two triangles represent a very small point cloud of an imaginary workpiece, whereas the two yellow triangles below represent the target mesh. In particular, the upper four points shall be assigned to the yellow triangle surfaces, whereby the isometry with respect to each pair of points shall be optimized. The coloring of the upper triangle illustrates the assignment quality, where a red color indicates a bad and a green color an isometric assignment. In case of assignments that are visually hard to relate to specific triangles, the green semicircles visualize to which of the yellow triangles the points are currently assigned. In the first two iterations, the green arrows point in the direction p of the corresponding gradient, i.e., the direction in which the corresponding points on the target shape are moved. In cases of points being situated on the border and a gradient pointing outside the related triangle, the movement is limited by the
6
L. Klein, T. Wagner, C. Buchheim, D. Biermann
For all points pi
Initial solution
Descent step: move assigned point zi in gradient direction h on the inner surface of the triangle (stop at border)
Relate pi to triangle
Yes
Max. number of iterations
No
?
No
Work-
(Change triangle)
piece border reached ?
Yes
Triangle border reached ?
Yes
Return current solution
No
Fig. 4: Sequence of operation diagram for the local reoptimization algorithm. Algorithm 1 Adaptive determination of the step length using backtracking
Fig. 5: Example of the gradient descent method for the local reoptimization.
border. The triangle assignment changes in the next iteration. This is illustrated from iteration 1 to iteration 2 in Fig. 5. In an unconstrained continuous domain, a gradient descent method guarantees an improvement of the objective value within a specific local neighborhood. In the special class of problem considered for reoptimization, however, the length of the feasible search directions are locally restricted by the triangle boundaries. Hence, the update rule x = x + αh with step length α and search direction h cannot be applied. Sophisticated gradient methods using additional information for the update of h, such as quasi-Newton or trust-region approaches [17], do not necessarily improve the convergence speed. The possible option to reduce the step length α in order to restrict the number of triangle swaps would result in slow convergence and long runtimes. In contrast, a too high value of α would result in updates far from the actual search direction. Hence, the idea of sufficient decrease [17, p. 32 ff.] was implemented using an adaptive step length algorithm based on a backtracking line search [17, p. 37]. This approach is shown in Algorithm 1.
Require: α0 /* Initial step length */ c /* Sufficient decrease parameter */ ρ /* Contraction factor */ 1: α = α0 /* Initialize step size */ 2: while f (x + αh) > f (x) + cαh∇f (x) do 3: α = ρα /* Decrease step length */ 4: end while 5: Return α
The idea behind the backtracking algorithm is to start with a large step length α0 which is reduced until a sufficient agreement with the expected gradient can be observed. This sufficient agreement is evaluated by testing in line 2 whether the improvement promised by the parameter-wise functional gradient ∇f (x) is obtained to a predefined extent. In this line, f (x + αh) is the isometry value obtained for the current test sample. By multiplying the actual step αh with the parameterwise functional gradient ∇f (x), a linear extrapolation of the gradient in search direction is performed. If at least a ratio of the improvement predicted by this extrapolation is obtained, the step length is accepted. The minimum ratio for acceptance is defined by a constant c ∈ [0, 1]. Otherwise the step length is contracted by a factor ρ ∈ ]0, 1[. This concept is particularly suited for the structure of the reoptimization problem. If the number of local projections due to triangle swaps is too high, the linear extrapolation of the gradient will fail and the step length is reduced. By accepting α the first time the sufficient decrease condition is met, too low step lengths are also avoided.
Evaluation and compensation of form errors by means of isometric registration with local reoptimization
4 Experiments In the experimental part of this paper, different aspects of the proposed procedure for evaluating form errors are analyzed. Firstly, it is shown that isometry (cf. equation (1)) is a practically relevant objective function for evaluating the error of the corresponding point pairs on both shapes. Secondly, the convergence of the proposed local reoptimization strategy is evaluated. For the neighborhood structure and the initial solution, simple baseline solutions are implemented. In order to evaluate the potential of the latter, however, the effect of the initial solution on the quality of the final assignment is also assessed.
Table 1: Comparison of the performance metrics for different assignments on the hat profile. Obj. Algo. Optimal ICP Reopt
The test workpieces for the experiments are shown in Fig. 6. They were chosen in order to cover practically relevant applications of different complexity. The hat profile is characterized by a simple shape, but a strong springback. It is thus suitable to analyze the capability of isometry for the evaluation of non-rigid deformations. The fuel tank cover represents a more complex workpiece with stress-induced deformations, in particular at the distant bulges. With regard to the research questions, two different representations of these test workpieces were used. In the first variant, the number of discrete points on both shapes is equal. Due to the maintenance of the topology during the practically motivated, but artificially computed, deformation of the target mesh, the true correspondences (pi ; zi∗ ) are known. Hence, the sum of squared errors (SSE) X kzi − zi∗ k2 (7) ˜ zi ∈Q
˜ can be computed. This of the current approximation Q SSE can be used to evaluate the suitability of different metrics to assess a given assignment of corresponding points. The knowledge about the structure of the optimal solution is used for evaluating metrics not having this knowledge. In the second variant, the number of points within the point data representing the manufactured shape is significantly higher than the one in the target mesh. Whereas the original number of points in the CAD output of the hat profile is n = 1.123, this was reduced to n = 396 using error-adaptive filtering techniques [13]. For the fuel tank cover, the reduction resulted in n = 5.041 instead of n = 11.949 points. The so-obtained workpieces were used to assess the improvements by the local reoptimization in case of a rough representation of the CAD model.
Distance
Isometry
SSE
24.392 (3.) 15.834 (1.) 16.148 (2.)
12.651 (1.) 992.028 (3.) 177.749 (2.)
0.000 (1.) 2918.169 (3.) 2843.576 (2.)
Table 2: Comparison of the performance metrics for different assignments on the fuel tank cover. Obj Algo.
4.1 Setup
7
Optimal ICP Reopt
Distance
Isometry
SSE
23.767 (3.) 12.003 (1.) 12.875 (2.)
2.606 (1.) 26513.600 (3.) 3486.250 (2.)
0.000 (1.) 25518.195 (3.) 22784.717 (2.)
Based on preliminary experiments, the parameters for the determination of the adaptive step length according to Algorithm 1 were set to α0 = 0.01, c = 0.2, and ρ = 0.9. The reoptimization algorithm was iterated until either the norm of the gradient became smaller than ε = 1e-6 or until the iteration limit nmax = 300 was reached. As mentioned before, the neighborhood structure and the global registration steps of the procedure were implemented using simple heuristics. For the former, the graph structure defined by the edges of the triangular mesh (cf. for instance Fig. 3) was considered. This resulted in an at least 2-node-connected graph. Lacking efficient methods for globally solving the quadratic assignment problems for the considered test workpieces, the initial solutions were generated using a simple projection based on the minimum distances between the shapes. In this context, point-based and triangle-based correspondences were distinguished. Both approaches were based on the iterative closest point (ICP) algorithm [3] which determines correspondences and then calculates a rigid transformation (shift and rotation) in order to minimize the distance between these correspondences in an iterative manner. Within this framework, the point-based correspondences were computed by assigning the points representing the manufactured shape to the closest points on the target shape. The triangle-based correspondences were found by directly projecting the points in normal direction to the original mesh.
8
Fig. 6: Overview of the workpieces used in the experiments. The local deformations are highlighted using red color.
4.2 Results The results for the evaluation of the first research question whether isometry is a suitable objective function for evaluating the error of the corresponding point pairs are summarized in Table 1 for the hat profile and in Table 2 for the fuel tank cover. In the columns, the values of the nearest-neighbor-distance used in ICP (Distance) are compared to the ones of isometry on the basis of three exemplary computed correspondences. These correspondences are listed in the rows and were computed using the triangle-based variant of ICP. The initial correspondences are labeled with ICP, whereas the correspondences after the refinement by means of the local reoptimization with respect to isometry are labeled with Reopt. In addition, the SSE (cf. equation (7)) with respect to the true correspondences are shown in the last column. The corresponding assignment also acts as an example in the rows of the tables (Optimal). The ranking of the different correspondences implied by the respective performance indicator in the column is provided in brackets. In both tables, the ranking resulting from isometry is identical to the one of the SSE. Moreover, the SSE value can be reduced by refining the solution of ICP using the proposed reoptimization procedure. Hence, isometry is a suitable approximation for the SSE which is actually unknown in practical applications. In contrast, the nearest-neighbor distance used in ICP resulted in extreme errors with
L. Klein, T. Wagner, C. Buchheim, D. Biermann
respect to the actual correspondences. It is no appropriate indicator for evaluating non-rigid deformations of the manufactured shape. The convergence of the adaptive step length approach within the local reoptimization algorithm is analyzed by means of convergence graphs showing the progression of the objective value, i. e., isometry in the considered case, over the iterations of the algorithm. In Figure 7 (hat profile) and Figure 8 (fuel tank cover ), these graphs are compared for the initial correspondences obtained by the point-based and triangle-based approach of ICP. In contrast to the results shown in the tables, the convergence graphs are based on the reduced target meshes. In all cases, the proposed reoptimization algorithm was able to continuously improve the correspondences. In particular in the first iterations, significant improvements were obtained. Nevertheless, the locality of the gradient-based reoptimization algorithm resulted in a strong dependence on the initial solution. The final result based on the point-based correspondences led to a premature convergence to a local optimum indicated by the end of the respective lines due to a reduction of the gradient’s norm below ε = 1e-6. These local optima were worse compared to ones found based on the triangle-based correspondences – 680 compared to 187 for the hat profile and 9540 compared to 3956 for the fuel tank cover. Starting from the solution of the point-based ICP approach, however, the reoptimization resulted in better solutions than the triangle-based ICP. Hence, it can be stated that the reoptimization algorithm successfully refines the initial solutions, but strongly depends on the quality of these correspondences.
5 Conclusion and Outlook In case of process-induced form errors, such as springback and torsions, complex local differences between the actual and the desired shape exist. As shown in the experiments, distance-based measures are not capable of finding appropriate correspondences between both shapes in these cases. These correspondences, however, are required for automatically adjusting the tool geometry or the respective NC paths [4, 5] Hence, a procedure for non-rigid registration was presented in this paper. This procedure was capable of calculating a suitable assignment between both shapes, even for nonuniform triangular meshes representing the workpiece model. To accomplish this, isometry was first defined as an objective function. Using experiments with known true correspondences, it was shown that isometry copes
Evaluation and compensation of form errors by means of isometric registration with local reoptimization
Triangle-based corresp. Point-based corresp.
9
Triangle-based corresp. Point-based corresp. 105 Isometry
Isometry
104
103 104
0
20
40
60
80
100
120
140
0
Number of iterations
Fig. 7: Convergence of the reoptimization algorithm on the reduced point set of the hat profile.
20
40
60
80
100
120
140
Number of iterations
Fig. 8: Convergence of the reoptimization algorithm on the reduced point set of the fuel tank cover. References
with non-rigid transformations. The problem of nonuniform triangular meshes was approached by a subsequent refinement step which was implemented by means of an adapted gradient method. Based on two exemplary workpieces, it was demonstrated that the proposed algorithm successfully converges to the neighboring local optima. As a result, improved correspondences for analyzing and compensating the form errors are available. The locality of the reoptimization approach resulted in a strong dependence on the initial solution. Hence, specialized algorithms for solving the global registration step of the procedure (cf. Fig. 2) have to be developed. The problem formulation as a quadratic assignment problem was already presented in this paper. Moreover, isometry depends on an appropriate neighborhood structure. The triangulation-based approach used in this paper already provided satisfying results. As alternatives, a minimum 2-node-connected graph, radiusand k-nearest-neighbor-based approaches can be considered. A comparison of these approaches is planned in the future. At the moment, the parameters of the adaptive step length algorithm were identified manually. The potential of setting these parameters will be exploited by a systematic parameter tuning. Acknowledgements This work is funded as subproject C4 of the Collaborative Research Center “3D-Surface Engineering” (SFB 708) by the German Research Foundation (DFG). Elias Kuthe and Cesaire Fondjo are acknowledged for their participation in the implementation and experimental evaluation of the presented local reoptimization algorithm.
1. Amberg, B., Romdhani, S., Vetter, T.: Optimal step nonrigid ICP algorithms for surface registration. In: IEEE Conference on Computer Vision and Pattern Recognition, pp. 1 –8 (2007). DOI 10.1109/CVPR.2007.383165 2. Anstreicher, K.: Recent advances in the solution of quadratic assignment problems. Mathematical Programming B 97, 27–42 (2003) 3. Besl, P., McKay, N.: A method for registration of 3-d shapes. IEEE Transactions on Pattern Analysis and Machine Intelligence 14(2), 239–256 (1992) 4. Biermann, D., Sacharow, A., Surmann, T., Wagner, T.: Direct free-form deformation of NC programs for surface reconstruction and form-error compensation. Production Engineering - Research and Development 4(5), 501–507 (2010) 5. Biermann, D., Surmann, T., Sacharow, A., Skutella, M., Theile, M.: Automated analysis of the form error caused by springback in metal sheet forming. In: Proceedings of the 3rd International Conference on Manufacturing Engineering, pp. 737–746 (2008) 6. Bronstein, A., Bronstein, M., Kimmel, R.: Generalized multidimensional scaling: A framework for isometryinvariant partial surface matching. Proceedings of the National Academy of Sciences 103(5), 1168–1172 (2006) 7. C ¸ ela, E.: The quadratic assignment problem - theory and algorithms. Kluwer Academic Publishers, Dordrecht, The Nederlands (1998) 8. Chen, Y., Medioni, G.: Object modeling by registration of multiple range images. In: Proceedings of the IEEE International Conference on Robotics and Automation, pp. 2724–2729 vol.3 (1991) 9. Elad, A., Kimmel, R.: On bending invariant signatures for surfaces. IEEE Transactions on Pattern Analysis and Machine Intelligence 25(10), 1285–1295 (2003) 10. Gan, W., Wagoner, R.: Die design method for sheet springback. International Journal of Mechanical Sciences 46(7), 1097–1113 (2004) 11. Garey, M.R., Johnson, D.S.: Computers and Intractability: A Guide to the Theory of NP-Completeness. W. H. Freeman & Co., New York, NY, USA (1979)
10 12. Horn, B.K.P.: Closed-form solution of absolute orientation using unit quaternions. Journal of the Optical Society of America A 4(4), 629–642 (1987) 13. Hussain, M.: Volume and normal field based simplification of polygonal models. journal of Information Science and Engineering 29(2), 267–279 (2013) 14. Kleiner, M., Tekkaya, A., Chatti, S., Hermes, M., Weinrich, A., Ben-Khalifa, N., Dirksen, U.: New incremental methods for springback compensation by stress superposition. Production Engineering - Research and Development 3(2), 137–144 (2009) 15. Lawler, E.L.: The quadratic assignment problem. Management Science 9(4), 586–599 (1963). DOI 10.1287/ mnsc.9.4.586 16. M´ emoli, F., Sapiro, G.: A theoretical and computational framework for isometry invariant recognition of point cloud data. Foundations of Computational Mathematics 5(3), 313–347 (2005) 17. Nocedal, J., Wright, S.: Numerical Optimization, 2nd edn. Springer (2006) 18. Sacharow, A., Balzer, J., Biermann, D., Surmann, T.: Non-rigid isometric icp: A practical registration method for the analysis and compensation of form errors in production engineering. Computer-Aided Design 43(12), 17581768 (2011). DOI 10.1016/j.cad.2011.07.007 19. Selimovic, I.: Improved algorithms for the projection of points on NURBS curves and surfaces. Computer Aided Geometric Design 23, 439–445 (2006) 20. Tam, G.K.L., Cheng, Z.Q., Lai, Y.K., Langbein, F.C., Liu, Y., Marshall, D., Martin, R.R., Sun, X.F., Rosin, P.L.: Registration of 3d point clouds and meshes: A survey from rigid to non-rigid. IEEE Transactions on Visualization and Computer Graphics not yet assigned (early access) (2013). DOI 10.1109/TVCG.2012.310 21. Wagner, T., Michelitsch, T., Sacharow, A.: On the design of optimisers for surface reconstruction. In: D. Thierens, et al. (eds.) Proc. 9th Annual Genetic and Evolutionary Computation Conference (GECCO 2007), London, UK, pp. 2195–2202. ACM, New York, NJ (2007). DOI 10. 1145/1276958.1277379 22. Weiher, J., Rietman, B., Kose, K., Ohnimus, S., Petzoldt, M.: Controlling springback with compensation strategies. In: AIP Conference Proceedings, pp. 1011–1015 (2004)
L. Klein, T. Wagner, C. Buchheim, D. Biermann
A procedure for the evaluation and compensation of form errors by means of global isometric registration with subsequent local reoptimization Laura Klein · Tobias Wagner · Christoph Buchheim · Dirk Biermann
October 11, 2013
Abstract Stresses remaining in the component after sheet metal forming processes can result in complex form errors, such as springback and torsions. In order to compensate these process-induced deformations, the local and global deformations have to be analyzed. Hence, an appropriate comparison between the actually manufactured and the target design is required. For this purpose, the surface of the actual workpiece is scanned and the so-obtained scan points have to be assigned to corresponding points of the target shape defined by the workpiece model. From these correspondences, a field of deformation vectors can be computed which represents the basis for the compensation strategy. The task of finding appropriate correspondences is called registration. It is usually solved using rigid transformations, i. e., translation and rotation. Due to the locality, strength and complexity of the deformations, rigid transformations are usually not sufficient. As a more flexible alternative, a procedure for non-rigid registration is presented in this paper. Therein, isometry, i.e., the conservation of distances between corresponding points within an appropriate neighborhood structure, is defined as the objective function. The procedure consists of three steps: definition of the neighborhood structure, global registration, and local reoptimization. The main focus of the paper is set to the latter, where an adapted gradient descent method also allowing projections into the triangles of the target shape is presented Christoph Buchheim · Laura Klein Fakult¨ at f¨ ur Mathematik, TU Dortmund, Vogelpothsweg 87, 44227 Dortmund E-mail: {christoph.buchheim,laura.klein}@tu-dortmund.de Dirk Biermann · Tobias Wagner Institut f¨ ur Spanende Fertigung, TU Dortmund, Baroper Straße 301, 44227 Dortmund E-mail: {biermann,wagner}@isf.tu-dortmund.de
and experimentally validated. With these three steps, an assignment between both shapes can be calculated, even for strong local deformations and coarse triangular meshes representing the workpiece model. Keywords Forming · Springback analysis · Non-rigid isometric registration · Local optimization · Quadratic assignment problem 1 Introduction An undesired effect in sheet metal forming is the springback of the deformed sheet after the process [14]. For controlling this springback, the local deformations have to be measured and analyzed to allow respective process and tool modifications to be derived [10, 22]. A modern approach for automatically determining a continuous deformation of the tool or the NC paths for toolmaking has already been presented [4, 5]. In this approach, the local deformations of the workpiece are directly used for adjusting the design of the tool, such that the sheet is overbent in this regions. Hence, the aim of this paper is to find appropriate correspondences to determine the local deformations as basis for the subsequent process and tool modifications. In order to find appropriate correspondences, methods to calculate assignments between the scan points measured on the actual workpiece and the surfaces of the target geometry are described. These methods define and solve an optimization problem formalizing the requirements of the practical application. They are organized within a three-step procedure. The improvements of this procedure with respect to stadard methods are based on two advantages: 1. The special characteristics of springback and torsions are considered by not minimizing the distance
2
L. Klein, T. Wagner, C. Buchheim, D. Biermann
(a)
(b)
Fig. 1: Comparison of the assignments (black lines) to the yellow target shape obtained by minimizing (a) the distance or (b) the isometry error (cf. equation 3) of the assigned point pairs. between both shapes. Instead, the distances between neighboring pairs of points on either shape are optimized in order to obtain an isometric mapping. 2. The discrete scan points representing the manufactured shape can be assigned to arbitrary points on the surface of the workpiece model. The first advantage emerges from the fact that distancebased assignments may lead to clustered mappings from the curves or kinks introduced by springback. Considering for instance the convex curve shown in Fig. 1, the mapping based on a minimization of the distance assigns a large part of the points on the knee of the workpiece to the flat part of the yellow target shape, which obviously does not represent an appropriate mapping. The second advantage is based on a higher flexibility with respect to the representation of the target geometry. As mentioned before, the presented registration approach reduces the problem to three subproblems which can be handled subsequently. First of all, a reasonable neighborhood for assessing the isometric matching is required. Secondly, an initial assignment is computed using a global approach. At this stage, a fixed set of points on the target model is considered. Thirdly, a local reoptimization based on a gradient descent method refines the initial solution based on the new objective function. In this step, the scan points can be arbitrarily assigned to the surface of the target geometry. Based on the proposed reduction to standard problem formulations, established algorithms for the neighborhood structure and the quadratic assignment problem can be adopted. Hence, the focus of this paper is set on the local reoptimization step. In the experiments, a simple heuristic for the neighborhood structure and an assignment based on a previously developed algorithm are applied. The procedure leads to significant improvements with respect to the assignment of the actual and the target shape.
2 Related work A comprehensive survey of rigid and non-rigid registration methods including those approaches has been published very recently [20]. Herein, only the most relevant methods with respect to the here given proposal are provided. The summary is categorized into methods performing either a rigid or a non-rigid transformation. As the former improve the matching only by translation and rotation of the actual shape, these methods can be mainly distinguished based on the way of generating and refining the correspondences between the points of the actual and the target shape for the distance computation. The definition of a neighborhood structure is not required. The latter methods are more complex and usually have many degrees of freedom. They therefore combine aspects of all three steps making up the proposed procedure.
2.1 Rigid registration The advantage of the restriction to rigid transformations is the reduction of degrees of freedom and the availability of a suitable analytical formulation, which significantly increases the efficiency of the methods. The rigid registration problem is mainly solved by iterative methods. In a first step, pair-wise correspondences between the points of the actual workpiece and the points of the target design are determined. In the second step, the optimal rigid transformation with respect to sum of squared distances between these pairs is computed. This is possible in an efficient way by means of a quaternion-based algorithm, see, e.g., [12]. After the transformation, the correspondences are refined and another transformation is computed. These two steps are repeated until the termination criterion is met, e.g., until the improvement falls below a predefined threshold.
Evaluation and compensation of form errors by means of isometric registration with local reoptimization
A very common approach is Iterative Closest Point (ICP) [3, 8], where each scan point of the actual shape is assigned to its closest point on the target design, either by projection on the parametrically defined workpiece surface [19] or by minimizing the distance over all points in case of discrete point data or triangulations [21]. This approach provides good solutions in cases of a well-chosen initial solution and small deviations between the actual and the target shape. Otherwise, ICP can lead to bad matchings because of local minima for inappropriate correspondences. To avoid this problem, an alternative approach was developed in the Collaborative Research Center (SFB 708) [5]. In this method, the model of the target design is discretized, i.e., replaced by a second point cloud, and a distance-minimal bipartite matching between these two point clouds is calculated. If both point clouds are similar and reasonably uniform, this approach can improve the matchings of ICP. The main weakness of this approach is its requirement for an appropriate discretization of the target shape. 2.2 Non-rigid registration
3 Procedure for non-rigid isometric registration In the following, it is assumed that the actual shape is represented by a point cloud, which can, e.g., be received by optical or tactile scanning of the workpiece. The model of the target shape is represented in terms of a triangular mesh. Even in case of a parametrical model, most CAD tools can export such kind of representation. The objective of the non-rigid registration approach is a mapping of the points within the cloud to their respective correspondences on the triangular mesh. The above mentioned disadvantages of current registration methods arise from indirect and heuristic approaches for the isometric registration. In this paper, a new approach is proposed by exactly formulating the isometry within a proper objective function. This makes the registration task accessible to methods from mathematical optimization [17]. As basis for computing the objective function, the differences between the Euclidean distances (2-norm) of two points pi , pj ∈ P of the point cloud and two corresponding points zi and zj on the target shape are computed, i.e. d(pi , pj ; zi , zj ) := |kpi − pj k2 − kzi − zj k2 |.
One possibility to relieve the restriction to rigid transformations is to allow the scan points to be independently mapped to the target shape by affine transformations. To prevent the geometry from degenerating, regulative constraints or repair steps minimizing the difference between transformations of neighboring points are added to the registration problem [1]. As a result, a quadratic minimization problem, yielding a transformation matrix for every single point measured on the surface of the manufactured workpiece, has to be solved. The number of variables increases significantly and efficient optimal algorithms are no longer available. To overcome the runtime problem, many approaches are based on local geometric features and a grouping of the scan points into clusters that undergo the same rigid transformation [5]. These approaches, however, require a merge of the clusters after the individual transformations. Other methods explicitly optimize isometry-based objective functions [6, 9, 16], but use more complex distance formulations, such as geodesic distance [9] or Gromov-Hausdorff distance [6, 16]. Based on the distance pairs, the points are embedded into a Euclidean space by multidimensional scaling (MDS) or enhanced variants thereof [6]. As an effect, the non-rigid registration problem transforms into a rigid one on the embedding space. The computation of the geodesic distance tableaux and the Gromov-Hausdorff distance, however, is extremely inefficient [6, 18].
3
(1)
With points pi fixed, the positions of the points zi on the target shape are optimized, i.e., arbitrarily moved on the target shape. The actual optimization objective is to minimize the sum of squared errors: X min d(pi , pj ; zi , zj )2 . (2) ij
The steps and the dependencies of the optimization procedure are summarized in Fig. 2. The choice of a neighborhood is worthwhile because of the new objective function which takes distance differences of two points on the shapes into account. As the target is to find a preferably isometric and non-rigid solution it is necessary to preserve distances of close-by situated points on the one hand, but not necessarily of distant points on the other hand. This idea is integrated by means of the neighborhood structure. In the second step, on a thinned set of points with the outlined neighborhood structure, an initial global assignment is calculated with the purpose of avoiding rough errors. In particular, this global assignment aims at preventing the objects to twist against each other, since later it is hard to overcome local minima and to disentangle the correspondences. The third step is a local reoptimization step improving the global registration. The assignment of the scanned points is not longer restricted to the discrete point set on the target shape as assignments are moved to arbitrary points on the respective surface.
4
L. Klein, T. Wagner, C. Buchheim, D. Biermann
Neighborhood structure defines objective function
[2-node-connected graph]
Global registration [Quadratic assignment problem]
defines objective function Local reoptimization
yields initial solution
[Gradient descent method]
More sophisticated approaches come from combinatorial optimization and graph theory. Consider a complete graph G = (P, E), whose vertices represent the points pi . The edge set F of any subgraph H = (P, F ) of G defines a selection of pairs, namely the points pi , pj being connected by an edge in F . Therefore, the neighborhood structure is formulated via an expedient subgraph H whose edge set F defines the pairs being considered in the objective function: X min d(pi , pj ; zi , zj )2 . (3) ij∈F
Fig. 2: Overview of the three steps and the related optimization problems including their dependencies.
The great benefit is the adjustment of discretization errors coming from thinning and the discrete assignment problem in the global step. In the following, the three steps are described in detail.
3.1 Neighborhood structure For a reasonable modeling of both, the non-rigidness and the isometry, it is not sensitive to take into account distances between all pairs of points on the manufactured workpiece. Different bendings naturally lead to different distances between two widely spaced points on the manufactured and their corresponding points on the target shape, therefore, their distances must not be preserved. Rather it is important to preserve the distances between two near-by situated points to preserve multiple mappings to angular points. A second technical argument for not taking into account all pairs of points is the high complexity of the quadratic assignment problem, which is proposed to be applied for finding an appropriate initial global solution and which is explained in Subsection 3.2. For both reasons, the use of a reduced neighborhood structure is proposed. As one obvious possibility for a neighborhood structure one can consider only pairs of points (pi , pj ) with a distance |pi − pj | < r, with a specified radius r > 0. Unfortunately, the choice of r is problematic. On the one hand, a too large radius value will not resolve the mentioned disadvantages, especially not the undesired rigidness. The choice of a too small value for r on the other hand runs the risk that some points might have too less neighbored points or even be isolated at all. Then again, isometry cannot be guaranteed.
To reduce the rigidness of the assignment, a minimal number or length of the edges chosen in H is demanded. In order to prevent the graph from internal twisting, however, each pair of nodes should be connected by at least two node disjoint paths, meaning that this pair is still connected if one arbitrary node is removed. This leads to the requirement of an at least 2-node-connected subgraph with minimal edge length. This problem is NP-hard [11], i.e., verifiable optimum solutions can only be found for very small graphs. Fortunately, there is no need for an exact solution of this problem for the given purpose. For a reasonable neighborhood structure, it is also possible to use graphs being denser than 2-nodeconnected, even though less and shorter edges improve computation times and the non-rigidness of the final solution. 3.2 Quadratic assignment problem For an initial global registration, only the triangle vertices Q := {q1 , . . . , qm } of the target mesh are allowed to be chosen as correspondences. As a consequence, each point pi ∈ P has to be assigned to exactly one point qj ∈ Q and vice versa. This leads to the restriction that both point sets P and Q have to be of the same cardinality, i.e., of the same number of points, which can usually be achieved by a curvature-adaptive reduction [13] of the scan point set P . Using binary variables xik ∈ {0, 1} which indicate whether point pi is assigned to point qk (xik = 1) or not (xik = 0), this global registration step can be modeled as a quadratic assignment problem (QAP) [7, 15]: XX min d(pi , pj ; qk , ql )2 xik xjl (4) ij∈F kl
subject to X xik = 1
for all i
k
X
xik = 1
for all k
xik ∈ {0, 1}
for all i, k
i
(5)
Evaluation and compensation of form errors by means of isometric registration with local reoptimization
Fig. 3: Assignment of the scan points on the manufactured shape to the vertices of the triangular mesh representing the target shape (gray). In the areas of a rough resolution of the mesh (black lines), high errors with respect to isometry (red color) can be observed. In the objective function (4) of the QAP, the respective distance costs d(pi , pk ; qj , ql )2 are used, where the arbitrary points zi of equation (3) are replaced by the triangle vertices qi . The distance costs are taken into account if two requirements are satisfied: 1. The points pi and pk have to be connected by the edges of the neighborhood defining graph. Hence, the first sum only includes the points connected by edges of the reasonable neighborhood structure F defined before. 2. The point pi on the workpiece is assigned to the point qk on the target shape, and pj is assigned to ql , respectively. This is encoded by the product of the binary variables xik and xjl . The constraints (5) express the requirement of the oneto-one assignment of a point pi ∈ P to a point qj ∈ Q. Also the QAP is NP-hard, but it is well studied in the field of discrete optimization. A comprehensive overview is given in [2], where exact and heuristic approaches are discussed. For the large point clouds faced in practice, mainly heuristic methods are applied.
3.3 Local reoptimization For the computation of a global assignment based on QAP solvers, a discretization of the target shape is necessary, as only correspondences between fixed sets can be considered – one for the scan points on the workpiece and one for the vertices of the triangular mesh representing the target shape. Moreover, both sets have to be of the same cardinality. In this context, the accuracy of the discretization of the target shape can become a problem. If the resolution is too fine, the number of variables becomes too high and the problem is not solvable in practice. A rough discretization on the other hand excludes isometric registrations from the outset, in particular if the points on the manufactured shape are spread differently compared to the triangular vertices
5
of the target mesh. Fig. 3 illustrates this effect on the basis of a rough triangular mesh of a 2.5-dimensional metal sheet, which is also used for computational tests in Section 4. In order to overcome the discretization problems and to restore isometry, a reoptimization of the discrete assignment is presented in the following. Its advantage consists in a mapping to arbitrary points on the inner surface of the triangles representing the target mesh. This allows an adaptation with respect to differences in dispersion and density of the initial point sets. To accomplish this, the problem X min d(pi , pj ; zi , zj )2 ij∈F (6) ˜ for all pi ∈ P zi ∈ Q ˜ now denotes the infinite set of all is solved, where Q positions on these triangle surfaces. For solving the problem defined in (6), an adapted iterative gradient descent method was developed. Its sequence of operations is depicted in Fig. 4. In every iteration step, each point pi ∈ P is related to exactly one triangle of the target mesh, whereby injectivity is not mandatory. In the descent step, every correspondence point zi is moved in the direction h = −∇f of steepest ascent with respect to the objective function f minimized in equation 6. If zi would move over the border of its related triangle, it is projected to the border for this iteration and it is related to the neighbored triangle in the next iteration step. In the case where the triangle border also is the border of the mesh, i.e. no neighboring triangle exists, the point remains within the current triangle. The single assignment steps of the gradient descent method are further illustrated in Fig. 5. In each of the four pictures, the four vertices of the upper two triangles represent a very small point cloud of an imaginary workpiece, whereas the two yellow triangles below represent the target mesh. In particular, the upper four points shall be assigned to the yellow triangle surfaces, whereby the isometry with respect to each pair of points shall be optimized. The coloring of the upper triangle illustrates the assignment quality, where a red color indicates a bad and a green color an isometric assignment. In case of assignments that are visually hard to relate to specific triangles, the green semicircles visualize to which of the yellow triangles the points are currently assigned. In the first two iterations, the green arrows point in the direction p of the corresponding gradient, i.e., the direction in which the corresponding points on the target shape are moved. In cases of points being situated on the border and a gradient pointing outside the related triangle, the movement is limited by the
6
L. Klein, T. Wagner, C. Buchheim, D. Biermann
For all points pi
Initial solution
Descent step: move assigned point zi in gradient direction h on the inner surface of the triangle (stop at border)
Relate pi to triangle
Yes
Max. number of iterations
No
?
No
Work-
(Change triangle)
piece border reached ?
Yes
Triangle border reached ?
Yes
Return current solution
No
Fig. 4: Sequence of operation diagram for the local reoptimization algorithm. Algorithm 1 Adaptive determination of the step length using backtracking
Fig. 5: Example of the gradient descent method for the local reoptimization.
border. The triangle assignment changes in the next iteration. This is illustrated from iteration 1 to iteration 2 in Fig. 5. In an unconstrained continuous domain, a gradient descent method guarantees an improvement of the objective value within a specific local neighborhood. In the special class of problem considered for reoptimization, however, the length of the feasible search directions are locally restricted by the triangle boundaries. Hence, the update rule x = x + αh with step length α and search direction h cannot be applied. Sophisticated gradient methods using additional information for the update of h, such as quasi-Newton or trust-region approaches [17], do not necessarily improve the convergence speed. The possible option to reduce the step length α in order to restrict the number of triangle swaps would result in slow convergence and long runtimes. In contrast, a too high value of α would result in updates far from the actual search direction. Hence, the idea of sufficient decrease [17, p. 32 ff.] was implemented using an adaptive step length algorithm based on a backtracking line search [17, p. 37]. This approach is shown in Algorithm 1.
Require: α0 /* Initial step length */ c /* Sufficient decrease parameter */ ρ /* Contraction factor */ 1: α = α0 /* Initialize step size */ 2: while f (x + αh) > f (x) + cαh∇f (x) do 3: α = ρα /* Decrease step length */ 4: end while 5: Return α
The idea behind the backtracking algorithm is to start with a large step length α0 which is reduced until a sufficient agreement with the expected gradient can be observed. This sufficient agreement is evaluated by testing in line 2 whether the improvement promised by the parameter-wise functional gradient ∇f (x) is obtained to a predefined extent. In this line, f (x + αh) is the isometry value obtained for the current test sample. By multiplying the actual step αh with the parameterwise functional gradient ∇f (x), a linear extrapolation of the gradient in search direction is performed. If at least a ratio of the improvement predicted by this extrapolation is obtained, the step length is accepted. The minimum ratio for acceptance is defined by a constant c ∈ [0, 1]. Otherwise the step length is contracted by a factor ρ ∈ ]0, 1[. This concept is particularly suited for the structure of the reoptimization problem. If the number of local projections due to triangle swaps is too high, the linear extrapolation of the gradient will fail and the step length is reduced. By accepting α the first time the sufficient decrease condition is met, too low step lengths are also avoided.
Evaluation and compensation of form errors by means of isometric registration with local reoptimization
4 Experiments In the experimental part of this paper, different aspects of the proposed procedure for evaluating form errors are analyzed. Firstly, it is shown that isometry (cf. equation (1)) is a practically relevant objective function for evaluating the error of the corresponding point pairs on both shapes. Secondly, the convergence of the proposed local reoptimization strategy is evaluated. For the neighborhood structure and the initial solution, simple baseline solutions are implemented. In order to evaluate the potential of the latter, however, the effect of the initial solution on the quality of the final assignment is also assessed.
Table 1: Comparison of the performance metrics for different assignments on the hat profile. Obj. Algo. Optimal ICP Reopt
The test workpieces for the experiments are shown in Fig. 6. They were chosen in order to cover practically relevant applications of different complexity. The hat profile is characterized by a simple shape, but a strong springback. It is thus suitable to analyze the capability of isometry for the evaluation of non-rigid deformations. The fuel tank cover represents a more complex workpiece with stress-induced deformations, in particular at the distant bulges. With regard to the research questions, two different representations of these test workpieces were used. In the first variant, the number of discrete points on both shapes is equal. Due to the maintenance of the topology during the practically motivated, but artificially computed, deformation of the target mesh, the true correspondences (pi ; zi∗ ) are known. Hence, the sum of squared errors (SSE) X kzi − zi∗ k2 (7) ˜ zi ∈Q
˜ can be computed. This of the current approximation Q SSE can be used to evaluate the suitability of different metrics to assess a given assignment of corresponding points. The knowledge about the structure of the optimal solution is used for evaluating metrics not having this knowledge. In the second variant, the number of points within the point data representing the manufactured shape is significantly higher than the one in the target mesh. Whereas the original number of points in the CAD output of the hat profile is n = 1.123, this was reduced to n = 396 using error-adaptive filtering techniques [13]. For the fuel tank cover, the reduction resulted in n = 5.041 instead of n = 11.949 points. The so-obtained workpieces were used to assess the improvements by the local reoptimization in case of a rough representation of the CAD model.
Distance
Isometry
SSE
24.392 (3.) 15.834 (1.) 16.148 (2.)
12.651 (1.) 992.028 (3.) 177.749 (2.)
0.000 (1.) 2918.169 (3.) 2843.576 (2.)
Table 2: Comparison of the performance metrics for different assignments on the fuel tank cover. Obj Algo.
4.1 Setup
7
Optimal ICP Reopt
Distance
Isometry
SSE
23.767 (3.) 12.003 (1.) 12.875 (2.)
2.606 (1.) 26513.600 (3.) 3486.250 (2.)
0.000 (1.) 25518.195 (3.) 22784.717 (2.)
Based on preliminary experiments, the parameters for the determination of the adaptive step length according to Algorithm 1 were set to α0 = 0.01, c = 0.2, and ρ = 0.9. The reoptimization algorithm was iterated until either the norm of the gradient became smaller than ε = 1e-6 or until the iteration limit nmax = 300 was reached. As mentioned before, the neighborhood structure and the global registration steps of the procedure were implemented using simple heuristics. For the former, the graph structure defined by the edges of the triangular mesh (cf. for instance Fig. 3) was considered. This resulted in an at least 2-node-connected graph. Lacking efficient methods for globally solving the quadratic assignment problems for the considered test workpieces, the initial solutions were generated using a simple projection based on the minimum distances between the shapes. In this context, point-based and triangle-based correspondences were distinguished. Both approaches were based on the iterative closest point (ICP) algorithm [3] which determines correspondences and then calculates a rigid transformation (shift and rotation) in order to minimize the distance between these correspondences in an iterative manner. Within this framework, the point-based correspondences were computed by assigning the points representing the manufactured shape to the closest points on the target shape. The triangle-based correspondences were found by directly projecting the points in normal direction to the original mesh.
8
Fig. 6: Overview of the workpieces used in the experiments. The local deformations are highlighted using red color.
4.2 Results The results for the evaluation of the first research question whether isometry is a suitable objective function for evaluating the error of the corresponding point pairs are summarized in Table 1 for the hat profile and in Table 2 for the fuel tank cover. In the columns, the values of the nearest-neighbor-distance used in ICP (Distance) are compared to the ones of isometry on the basis of three exemplary computed correspondences. These correspondences are listed in the rows and were computed using the triangle-based variant of ICP. The initial correspondences are labeled with ICP, whereas the correspondences after the refinement by means of the local reoptimization with respect to isometry are labeled with Reopt. In addition, the SSE (cf. equation (7)) with respect to the true correspondences are shown in the last column. The corresponding assignment also acts as an example in the rows of the tables (Optimal). The ranking of the different correspondences implied by the respective performance indicator in the column is provided in brackets. In both tables, the ranking resulting from isometry is identical to the one of the SSE. Moreover, the SSE value can be reduced by refining the solution of ICP using the proposed reoptimization procedure. Hence, isometry is a suitable approximation for the SSE which is actually unknown in practical applications. In contrast, the nearest-neighbor distance used in ICP resulted in extreme errors with
L. Klein, T. Wagner, C. Buchheim, D. Biermann
respect to the actual correspondences. It is no appropriate indicator for evaluating non-rigid deformations of the manufactured shape. The convergence of the adaptive step length approach within the local reoptimization algorithm is analyzed by means of convergence graphs showing the progression of the objective value, i. e., isometry in the considered case, over the iterations of the algorithm. In Figure 7 (hat profile) and Figure 8 (fuel tank cover ), these graphs are compared for the initial correspondences obtained by the point-based and triangle-based approach of ICP. In contrast to the results shown in the tables, the convergence graphs are based on the reduced target meshes. In all cases, the proposed reoptimization algorithm was able to continuously improve the correspondences. In particular in the first iterations, significant improvements were obtained. Nevertheless, the locality of the gradient-based reoptimization algorithm resulted in a strong dependence on the initial solution. The final result based on the point-based correspondences led to a premature convergence to a local optimum indicated by the end of the respective lines due to a reduction of the gradient’s norm below ε = 1e-6. These local optima were worse compared to ones found based on the triangle-based correspondences – 680 compared to 187 for the hat profile and 9540 compared to 3956 for the fuel tank cover. Starting from the solution of the point-based ICP approach, however, the reoptimization resulted in better solutions than the triangle-based ICP. Hence, it can be stated that the reoptimization algorithm successfully refines the initial solutions, but strongly depends on the quality of these correspondences.
5 Conclusion and Outlook In case of process-induced form errors, such as springback and torsions, complex local differences between the actual and the desired shape exist. As shown in the experiments, distance-based measures are not capable of finding appropriate correspondences between both shapes in these cases. These correspondences, however, are required for automatically adjusting the tool geometry or the respective NC paths [4, 5] Hence, a procedure for non-rigid registration was presented in this paper. This procedure was capable of calculating a suitable assignment between both shapes, even for nonuniform triangular meshes representing the workpiece model. To accomplish this, isometry was first defined as an objective function. Using experiments with known true correspondences, it was shown that isometry copes
Evaluation and compensation of form errors by means of isometric registration with local reoptimization
Triangle-based corresp. Point-based corresp.
9
Triangle-based corresp. Point-based corresp. 105 Isometry
Isometry
104
103 104
0
20
40
60
80
100
120
140
0
Number of iterations
Fig. 7: Convergence of the reoptimization algorithm on the reduced point set of the hat profile.
20
40
60
80
100
120
140
Number of iterations
Fig. 8: Convergence of the reoptimization algorithm on the reduced point set of the fuel tank cover. References
with non-rigid transformations. The problem of nonuniform triangular meshes was approached by a subsequent refinement step which was implemented by means of an adapted gradient method. Based on two exemplary workpieces, it was demonstrated that the proposed algorithm successfully converges to the neighboring local optima. As a result, improved correspondences for analyzing and compensating the form errors are available. The locality of the reoptimization approach resulted in a strong dependence on the initial solution. Hence, specialized algorithms for solving the global registration step of the procedure (cf. Fig. 2) have to be developed. The problem formulation as a quadratic assignment problem was already presented in this paper. Moreover, isometry depends on an appropriate neighborhood structure. The triangulation-based approach used in this paper already provided satisfying results. As alternatives, a minimum 2-node-connected graph, radiusand k-nearest-neighbor-based approaches can be considered. A comparison of these approaches is planned in the future. At the moment, the parameters of the adaptive step length algorithm were identified manually. The potential of setting these parameters will be exploited by a systematic parameter tuning. Acknowledgements This work is funded as subproject C4 of the Collaborative Research Center “3D-Surface Engineering” (SFB 708) by the German Research Foundation (DFG). Elias Kuthe and Cesaire Fondjo are acknowledged for their participation in the implementation and experimental evaluation of the presented local reoptimization algorithm.
1. Amberg, B., Romdhani, S., Vetter, T.: Optimal step nonrigid ICP algorithms for surface registration. In: IEEE Conference on Computer Vision and Pattern Recognition, pp. 1 –8 (2007). DOI 10.1109/CVPR.2007.383165 2. Anstreicher, K.: Recent advances in the solution of quadratic assignment problems. Mathematical Programming B 97, 27–42 (2003) 3. Besl, P., McKay, N.: A method for registration of 3-d shapes. IEEE Transactions on Pattern Analysis and Machine Intelligence 14(2), 239–256 (1992) 4. Biermann, D., Sacharow, A., Surmann, T., Wagner, T.: Direct free-form deformation of NC programs for surface reconstruction and form-error compensation. Production Engineering - Research and Development 4(5), 501–507 (2010) 5. Biermann, D., Surmann, T., Sacharow, A., Skutella, M., Theile, M.: Automated analysis of the form error caused by springback in metal sheet forming. In: Proceedings of the 3rd International Conference on Manufacturing Engineering, pp. 737–746 (2008) 6. Bronstein, A., Bronstein, M., Kimmel, R.: Generalized multidimensional scaling: A framework for isometryinvariant partial surface matching. Proceedings of the National Academy of Sciences 103(5), 1168–1172 (2006) 7. C ¸ ela, E.: The quadratic assignment problem - theory and algorithms. Kluwer Academic Publishers, Dordrecht, The Nederlands (1998) 8. Chen, Y., Medioni, G.: Object modeling by registration of multiple range images. In: Proceedings of the IEEE International Conference on Robotics and Automation, pp. 2724–2729 vol.3 (1991) 9. Elad, A., Kimmel, R.: On bending invariant signatures for surfaces. IEEE Transactions on Pattern Analysis and Machine Intelligence 25(10), 1285–1295 (2003) 10. Gan, W., Wagoner, R.: Die design method for sheet springback. International Journal of Mechanical Sciences 46(7), 1097–1113 (2004) 11. Garey, M.R., Johnson, D.S.: Computers and Intractability: A Guide to the Theory of NP-Completeness. W. H. Freeman & Co., New York, NY, USA (1979)
10 12. Horn, B.K.P.: Closed-form solution of absolute orientation using unit quaternions. Journal of the Optical Society of America A 4(4), 629–642 (1987) 13. Hussain, M.: Volume and normal field based simplification of polygonal models. journal of Information Science and Engineering 29(2), 267–279 (2013) 14. Kleiner, M., Tekkaya, A., Chatti, S., Hermes, M., Weinrich, A., Ben-Khalifa, N., Dirksen, U.: New incremental methods for springback compensation by stress superposition. Production Engineering - Research and Development 3(2), 137–144 (2009) 15. Lawler, E.L.: The quadratic assignment problem. Management Science 9(4), 586–599 (1963). DOI 10.1287/ mnsc.9.4.586 16. M´ emoli, F., Sapiro, G.: A theoretical and computational framework for isometry invariant recognition of point cloud data. Foundations of Computational Mathematics 5(3), 313–347 (2005) 17. Nocedal, J., Wright, S.: Numerical Optimization, 2nd edn. Springer (2006) 18. Sacharow, A., Balzer, J., Biermann, D., Surmann, T.: Non-rigid isometric icp: A practical registration method for the analysis and compensation of form errors in production engineering. Computer-Aided Design 43(12), 17581768 (2011). DOI 10.1016/j.cad.2011.07.007 19. Selimovic, I.: Improved algorithms for the projection of points on NURBS curves and surfaces. Computer Aided Geometric Design 23, 439–445 (2006) 20. Tam, G.K.L., Cheng, Z.Q., Lai, Y.K., Langbein, F.C., Liu, Y., Marshall, D., Martin, R.R., Sun, X.F., Rosin, P.L.: Registration of 3d point clouds and meshes: A survey from rigid to non-rigid. IEEE Transactions on Visualization and Computer Graphics not yet assigned (early access) (2013). DOI 10.1109/TVCG.2012.310 21. Wagner, T., Michelitsch, T., Sacharow, A.: On the design of optimisers for surface reconstruction. In: D. Thierens, et al. (eds.) Proc. 9th Annual Genetic and Evolutionary Computation Conference (GECCO 2007), London, UK, pp. 2195–2202. ACM, New York, NJ (2007). DOI 10. 1145/1276958.1277379 22. Weiher, J., Rietman, B., Kose, K., Ohnimus, S., Petzoldt, M.: Controlling springback with compensation strategies. In: AIP Conference Proceedings, pp. 1011–1015 (2004)
L. Klein, T. Wagner, C. Buchheim, D. Biermann
