A tolerance constrained robot path circular interpolation method for industrial SCARA robots

Abstract

The linear robot path is tangential and curvature discontinuity, which will lead to vibration and unnecessary hesitation during execution. Local corner transition method and local spline interpolation method are used in state-of-art industrial robot controller to reduce vibration, while local corner transition method cannot interpolate target points and local spline interpolation method cannot constrain chord errors. This research proposes a robot path local interpolation method that eliminates deficiencies of each method. The smoothing method satisfies all of the following requirements: G¹ continuity, target point interpolation and chord tolerance confined, shape-preserving (free of self-intersection), and unified parameterization. The generated smooth path consists of linear path and circular arc path with G¹ continuity. A geometric iterative method cooperating with local corner transition method is used to generate local interpolation path. Simulations and actual experiments verify the generated smooth path is G¹ continuous, tolerance constrained, shape-preserving, and have high computational efficiency.

Introduction

Industrial robots are widely used in automotive and electrical/electronics industry. The continuity of robot path is important to robots’ execution accuracy and efficiency. In complex robotic applications, robot path is often generated by robot offline programming software, ¹ formed of a series of short linear segments in Cartesian space. The linear path only has G⁰ continuity, which will cause unwanted vibrations, slow-downs, and low quality.^2,3

Path smoothing algorithm can be used to generate continuous and tolerance constrained motions, thus to improve the efficiency and accuracy of robot application. A target point of linear robot path can be represented by a position in R³ (3-dimensional Euclidean space) and orientation in SO(3) (the group of all rotations in R³). ⁴ Both positions and orientations in Cartesian space are needed to be smoothed.

Position path smoothing methods

The smoothing methods in CNC machining or robotics for 3-dimensional paths can be divided into four classes: global interpolation method, global approximation method (fitting), local interpolation method, and local approximation methods (corner transition). ⁴

Global methods usually use splines to interpolate^5,6 or approximate^7,8 the whole set of data points. Global methods are often used offline because it is difficult to satisfy the requirements of interpolation, chord error constrained, and no self-intersection in real-time performance. PIA (Progressive iterative approximation methods),^7,9 or the name of geometric iterative methods, are effective iterative methods used in global B-spline approximation or interpolation to avoid numerical instability. Deng and Lin ⁹ reported a 3D data fitting method of LSPIA (PIA for Least Squares). He et al. ⁷ used the ELSPIA (LSPIA method with energy term) to smooth three-axis tool path by cubic B-splines. Chord errors are always unconstrained because it is difficult to represent the chord error constraints analytically.

Local corner transition methods utilize a smooth curve to approximate linear path at corners. The method is popular in CNC tool path smoothing as the advantage of analytical computation. Zhao et al. ¹⁰ adopted a curvature-continuous B-spline to blend 3-axis linear tool path. Du et al. ¹¹ used a cubic Bezier to smooth G01 paths under constraints of curvature variation energy. Li et al. ¹² used the local transition scheme for a real-time and look-ahead interpolation of tool paths. Zhang et al. ¹³ proposed a corner smoothing method using feedrate blending under geometric and kinematic constraints. Similarly, the local corner transition methods were extended to five-axis tool path smoothing.^14–18 In robotics, the parabolic path was commonly used for corner transition of position path.^19,20 The reason is that trapezoidal speed planning can be easily realized by adopting the parabolic path. Also, quadratic Bezier path by Wang et al. ²¹ were designed in a SCARA robot control system. The local corner transition methods can generate smooth paths by simple calculation, but the smooth path cannot interpolate the target points.

Local spline interpolation methods utilize splines to interpolate target points considering adjacent path continuity. Fan et al. ²² used iterative methods to interpolate a three-axis tool path by Bezier curves. Fan et al. ²³ proposed a chord error constrained and shape-preserving local interpolation method by adjusting undesired segments. Fleisig and Spence ²⁴ used quintic polynomial splines to interpolate tool path with G² continuity. The polynomial interpolation methods are also used in robotics^3,25 to get target point interpolated smooth path. Zhao et al. ²⁶ presented a NURBS curve interpolation method for a SCARA robot. Compared with local corner transition methods, local interpolation methods can generate smooth paths that interpolating target points, but the local interpolation methods are difficult to be solved analytically. High computation efficiency is an important requirement of local interpolation methods.

Orientation path smoothing methods

The orientation path of robots is often represented by quaternions, Euler angles or axis-angle method. ²⁷ Linear orientation path is usually represented by SLERP (Spherical Linear Interpolation) curve. ²⁰ The orientation path can be smoothed by the quaternion method. Kim and Nam ²⁸ used circular blending quaternion curves to interpolate solid orientations. Then Kim and Nam ²⁹ proposed a Hermite interpolation method using two circular quaternion curves on SO(3). Park and Ravani ³⁰ generalized the concept of Bezier curves to Riemannian manifolds, which can be used for orientation interpolation. Furthermore, Kim et al. ³¹ presented a construction method of a C²-continuous B-spline quaternion curve which interpolates a sequence of orientations in SO(3).

The position path and orientation path can be smoothed separately using different parameters or smoothed synchronously using unified parameters. In five-axis machining, the tool tips and tool axes are synchronized using a third curve to fit the relationship between two parameters,^6,14,24 or using the adjusting method in the smoothing process.^5,8,32 Chen and Chen ³³ proposed a synchronized fitting method for a six-axis machining robot manipulators, the orientation uses a rotation matrix and a vector to plan the tool orientations. However, the orientation errors cannot be constrained. He et al. ² proposed a G² continuous path interpolation method for SCARA robots, while the commercial robot control system does not accept user-defined splines.

Current commercial robot control systems provide local corner transition methods and local interpolation methods. ABB robotics ³⁴ and KUKA robotics ³⁵ use a distance to control the local corner transition shape. FANUC robotics ³⁶ and YASKAWA Motoman robotics ³⁷ use a scale parameter to specify the transition curve. YASKAWA Motoman robotics also provides a CR command to specify the radius of a transition arc. The two kinds of smoothing methods cannot satisfy the precision of target points and chord errors simultaneously. Moreover, the smoothing methods of orientation path are not specified in these methods.

About this research

The program path of SCARA robots can be represented by a position and a rotation angle of orientation, as the rotation axis of orientation is a constant value. This research uses a sequence of 4D points (a 3D position and a rotation angle) to represent SCARA robot path. Circular arc will be used as the smooth blending path, which has good properties of analytical arc length and curvature. The designed computational frame can be extended to smooth six-axis robot paths.

This research focuses on geometrical smoothing of robot path. This research aims to design a path smoothing method for SCARA robots. The smoothing methods can fulfill the requirements of G¹ continuity, target points interpolation, chord error constrained, shape-preserving, and high computational efficiency. Moreover, the position path and orientation path are smoothed using unified parameters, no complicated parameter synchronization operation is required. A synchronized speeding planning method based on unified parameters of positions and orientations will be presented in next research.

This research is organized as follows. Section 2 presents the local arc transition method for SCARA robot path. The transit path is formulated by a convex combination formulation. Section 3 presents SCARA robot path interpolation method based on local arc transition. A geometric iterative method is used to calculate transition parameters. Section 4 presents several simulation experiments of the proposed method and comparisons with other methods. Section 5 is the actual experiments by drawing with a TURIN SCARA robot. At last, Section 6 concludes this research.

SCARA robot path local arc transition

The representation of SCARA linear robot path

In Cartesian space, a target point of linear robot path represents the position and orientation of Tool Coordinate System (TCS) relative to Robot base Coordinate System (RCS), as shown in Figure 1. The target points of SCARA robot path are represented by 4D point $\{{\mathit{P}}_i(x_i,y_i,z_i,{\theta }_i){\}}_{i=0}^N$ . Positions of SCARA robot can be represented by the Tool Center Point (TCP) ( $x_i,y_i,z_i)$ in RCS. The rotation axis of the three rotation joints are constant, only rotation angles ( ${\theta }_i$ ) around z-axis of TCS are variables.

Figure 1.

CAD model of a SCARA robot.

The distance between two 4D points can be defined by position distance and orientation distance. Suppose a 4D point is ${\mathit{P}}_i=(x_i,y_i,z_i,{\theta }_i)$ . Denote ${\mathit{P}}_{i,\left[3\right]}=(x_i,y_i,z_i)$ , ${\mathit{P}}_{i,\left[\theta \right]}={\theta }_i$ . The position distance and orientation distance are defined by equations (1) and (2).

$$D_3({\mathit{P}}_i,{\mathit{P}}_{i-1})=|{\mathit{P}}_{i,\left[3\right]}-{\mathit{P}}_{i-1,\left[3\right]}|=\sqrt{{(x_i-x_{i-1})}^2+{(y_i-y_{i-1})}^2+{(z_i-z_{i-1})}^2}$$

(1)

$$D_{\theta }({\mathit{P}}_i,{\mathit{P}}_{i-1})=|{\mathit{P}}_{i,\left[\theta \right]}-{\mathit{P}}_{i-1,\left[\theta \right]}|=|{\theta }_i-{\theta }_{i-1}|$$

(2)

Based on the definition of 4D representation and distance, a 4D arc path for SCARA robot path will be defined in Section 2.2.

Definition of a 4D circular arc path

A circular arc in R³ can be defined by the arc center, start point, end point and arc normal. Extend this definition to 4D space: suppose $\mathit{O}$ is the center point of a 4D arc $\mathit{C}\left(u\right),$ ${\mathit{P}}_s$ is the start point, ${\mathit{P}}_e$ is the end point, and $\mathit{n}$ is the arc normal in R³. The 4D arc $\mathit{C}\left(u\right)$ can be represented by equation (3).

$$\mathit{C}\left(u\right)=\mathit{O}+\frac{\sin ((1-u)\alpha )\left(\mathit{O}{\mathit{P}}_s\right)+\sin \left(u\alpha \right)\left(\mathit{O}{\mathit{P}}_e\right)}{\sin \alpha },u\in [0,1]$$

(3)

In equation (3), $\alpha$ is the rotation angle from ${\mathit{P}}_s$ to ${\mathit{P}}_e$ . As shown in Figure 2, the rotation from ${\mathit{P}}_s$ to ${\mathit{P}}_e$ may follow minor arc or major arc. Assume the angle between 3D vector ${\left(\mathit{O}{\mathit{P}}_s\right)}_{\left[3\right]}$ and ${\left(\mathit{O}{\mathit{P}}_e\right)}_{\left[3\right]}$ is $\phi$ . The rotation angle $\alpha$ in equation (3) can be calculated by equation (4).

$$\{\begin{array}{c}\alpha =\phi ,\mathrm{if}({\left(\mathit{O}{\mathit{P}}_s\right)}_{\left[3\right]}\times {\left(\mathit{O}{\mathit{P}}_e\right)}_{\left[3\right]}).\mathit{n}>0\hfill \\ \multicolumn{1}{c}{\alpha =2\pi -\phi ,\mathrm{else}}\end{array}$$

(4)

Figure 2.

4D arc path: (a) minor arc and (b) major arc.

Using the unified parametrization, the position portion of 4D circular arc is a 3D arc, and the rotation portion is a continuous quadratic curve of rotation angles relative to rotation arc length.

Local arc transition and convex combination formulation

Circular arc in R³ is often used for corner smoothing, because its arc length and curvature can be analytically computed. This research will use 4D circular arcs to smooth linear SCARA robot path. Preliminary, local circular arc transition method is introduced.

As shown in Figure 3, considering three 4D points ${\mathit{Q}}_{\mathit{i}-1},{\mathit{Q}}_{\mathit{i}},{\mathit{Q}}_{\mathit{i}+1}$ . A 4D circular arc ${\mathit{C}}_{\mathit{i}}\left(u\right)$ is used to transit ${\mathit{Q}}_{\mathit{i}-1},{\mathit{Q}}_{\mathit{i}},{\mathit{Q}}_{\mathit{i}+1}$ with G¹ continuity. By deduction, the arc center ${\mathit{O}}_{\mathit{i}}$ , start point ${\mathit{E}}_{\mathit{i},0}$ , end point ${\mathit{E}}_{\mathit{i},1}$ of the 4D arc $\mathit{C}\left(u\right)$ can be determined by the continuous condition and transition ratio as in equation (5). In which, ${\varphi }_i$ is the position angle between vector ${\mathit{Q}}_i{\mathit{Q}}_{i-1}$ and ${\mathit{Q}}_i{\mathit{Q}}_{i+1}$ , Transition ratio $r_i=\frac{D_3({\mathit{Q}}_i,{\mathit{E}}_{i,0})}{D_3({\mathit{Q}}_i,{\mathit{Q}}_{i-1})}\in (0,1)$ , $a_i=\frac{D_3({\mathit{Q}}_i,{\mathit{Q}}_{i-1})}{D_3({\mathit{Q}}_i,{\mathit{Q}}_{i+1})}$ and $b_i=\frac{1}{1+\cos {\varphi }_i}$ . $|.|$ is the modular operation for 3-dimetional points.

$$\{\begin{array}{c}{\mathit{E}}_{\mathit{i},0}={\mathit{Q}}_{\mathit{i}}+{\mathit{r}}_{\mathit{i}}({\mathit{Q}}_{\mathit{i}-1}-{\mathit{Q}}_{\mathit{i}})\hfill \\ \multicolumn{1}{c}{{\mathit{E}}_{\mathit{i},1}={\mathit{Q}}_{\mathit{i}}+{\mathit{a}}_{\mathit{i}}{\mathit{r}}_{\mathit{i}}({\mathit{Q}}_{\mathit{i}+1}-{\mathit{Q}}_{\mathit{i}})}\\ \multicolumn{1}{c}{{\mathit{O}}_{\mathit{i}}={\mathit{Q}}_{\mathit{i}}+{\mathit{b}}_{\mathit{i}}({\mathit{Q}}_{\mathit{i}}{\mathit{E}}_{\mathit{i},0}+{\mathit{Q}}_{\mathit{i}}{\mathit{E}}_{\mathit{i},0})}\\ \multicolumn{1}{c}{{\mathit{n}}_{\mathit{i}}=\frac{{({\mathit{Q}}_{\mathit{i}+1}-{\mathit{Q}}_{\mathit{i}})}_{\left[3\right]}\times {({\mathit{Q}}_{\mathit{i}-1}-{\mathit{Q}}_{\mathit{i}})}_{\left[3\right]}}{|{({\mathit{Q}}_{\mathit{i}+1}-{\mathit{Q}}_{\mathit{i}})}_{\left[3\right]}\times {({\mathit{Q}}_{\mathit{i}-1}-{\mathit{Q}}_{\mathit{i}})}_{\left[3\right]}|}}\end{array}$$

(5)

Figure 3.

Local arc transition.

The arc transition curve can be built by equations (3)–(5) once the transition ratio $r_i$ is determined. Note that the rotation angle in equation (4) is computed by ${\phi }_i=\pi -{\varphi }_i$ .

The following will use a simplified representation for the transition circular arc. By the theory of convex combination, the transition arc ${\mathit{C}}_i\left(u\right)$ can be represented by equation (6). ²

$${\mathit{C}}_i\left(u\right)=F_i\left(u\right){\mathit{Q}}_{i-1}+J_i\left(u\right){\mathit{Q}}_i+K_i\left(u\right){\mathit{Q}}_{i+1},u\in [0,1]$$

(6)

In which, $F_i\left(u\right)$ , ${\mathit{J}}_i\left(u\right)$ and $K_i\left(u\right)$ are basis functions that satisfy $F_i\left(u\right)+J_i\left(u\right)+K_i\left(u\right)=1$ . By equations (3)–(5), the basic functions can be represented by equation (7).

$$\{\begin{array}{c}{F}_i\left(u\right)=\frac{\sin {\phi }_i-\sin ((1-u){\phi }_i)\cos {\phi }_i-\sin \left(u{\phi }_i\right)}{\sin {\phi }_i(1-\cos {\phi }_i)}{r}_i\hfill \\ \multicolumn{1}{c}{J_i\left(u\right)=1-F_i\left(u\right)-K_i\left(u\right)}\\ \multicolumn{1}{c}{K_i\left(u\right)=a_i{F}_i\left(u\right)}\end{array}$$

(7)

The convex combination based representation can simply computation and internal storage in this research, and the smoothing method can be easily extended to other transition curves.

Denoting ${\beta }_i=F_i(0.5)$ . If ${\beta }_i$ is known first, the transition ratio $r_i$ can be determined by $r_i=\frac{(1-\cos {\phi }_i)}{1-\cos \left(\frac{{\phi }_i}{2}\right)}{\beta }_i$ . Substituting $u=0.5$ into equation (6), parameter middle point ${\mathit{R}}_i$ can be represented by equation (8).

$${\mathit{R}}_i={\mathit{C}}_i(0.5)={\mathit{Q}}_i+{\beta }_i({\mathit{Q}}_i{\mathit{Q}}_{i-1}+a_i{\mathit{Q}}_i{\mathit{Q}}_{i+1})$$

(8)

The parameter ${\beta }_i$ is named as transition parameter and will be determined in Section 3.

The local arc transition method can be used to generate G¹ continuous smooth path. To generate target point interpolation smooth path, this research will use an iterative method to construct a virtual linear path, and then use local arc transition method to get target point interpolated smooth path. The details will be introduced in Section 3.

SCARA robot path interpolation based on local arc transition

In this research, a geometric iterative method is used to generate virtual linear path (named pseudo path) and transition parameters, and then smooth path is constructed by the pseudo path and transition parameters using local arc transition. The smooth path should interpolate target points while satisfying chord tolerance constraints and free of self-intersection.

Algorithm overview

The illustration of the linear path, pseudo path, and smoothing path are shown in Figure 4(a). The first pseudo path $\{{\mathit{Q}}_i^0{\}}_{i=0}^N$ is the same with input linear path $\{{\mathit{P}}_i{\}}_{i=0}^N$ . The linear path $\{{\mathit{Q}}_i^k{\}}_{i=0}^N$ is the pseudo path constructed in the iteration. The linear path $\{{\mathit{Q}}_i^M{\}}_{i=0}^N$ is the final pseudo path that satisfies all the requirements. The smooth path is constructed by the convex combination of the final pseudo path and transition parameters, using local arc transition method in Section 2. The smooth path is consisted of linear paths and arcs with G¹ continuity.

Figure 4.

Illustration and flowchart of SCARA robot path smoothing method: (a) illustration of pseudo path and (b) flowchart of robot path smoothing.

The flowchart of the smoothing method is illustrating in Figure 4(b). First input linear path $\{{\mathit{P}}_i{\}}_{i=0}^N$ and smoothing parameters (data point tolerance, chord error tolerance, and shape-preserving parameters). Then a geometric iteration method is used to compute iterative pseudo path $\{{\mathit{Q}}_i^k{\}}_{i=0}^N$ and transition parameters $\{{\beta }_i^k{\}}_{i=1}^{N-1}$ (step1∼step 5). The iteration is stopped if the data point tolerance of positions and orientations are both satisfied. If the iteration is failed or exceeds the maximum iteration number, a refinement method (step 6) that the smooth result of last iteration is used. At last (step 7 and step 8), the smooth path is constructed by local arc transition of final pseudo path $\{{\mathit{Q}}_i^M{\}}_{i=0}^N$ and transition parameters $\{{\beta }_i^M{\}}_{i=1}^{N-1}$ . This research will focus on describing the geometric iterative method in Section 3.2.

Geometric iterative method for local arc interpolation

Except the G¹ continuity requirement, the smoothing path should satisfy three other requirements: (1) Target points interpolation within given data point tolerance (include position data point tolerance and orientation data point tolerance). (2) Chord error confined within given chord tolerance. (3) Free of self-intersection by given shape-preserving parameters.

The requirements and objective function will be introduced in Section 3.2.1. The algorithm of computing pseudo path and transition parameters will be introduced in Section 3.2.2, and at last two adjusting vector selection methods are introduced in Section 3.2.3.

Requirements and objective function

Suppose the input linear path is $\{{\mathit{P}}_i{\}}_{i=0}^N$ , the position data point tolerance is ${\epsilon }_{\max }$ and the orientation data point tolerance is $o_{\max }$ , the position chord tolerance is ${\delta }_{\max }$ , and a blank ratio $\rho$ is defined to avoid self-intersection. The above mentioned three requirements can be represented by a mathematical description.

Suppose the constructed final pseudo path and transition parameters are $\{{\mathit{Q}}_i^M{\}}_{i=0}^N$ and $\{{\beta }_i^M{\}}_{i=0}^N$ . For index $i=1,2,\dots N-1$ , ${\mathit{C}}_i^M\left(u\right)$ is the arc that interpolates ${\mathit{P}}_i$ at ${\mathit{R}}_i^M$ . By equation (8), ${\mathit{R}}_i^M={\mathit{C}}_i^M(0.5)={\mathit{Q}}_i^M+{\beta }_i^M({\mathit{Q}}_i^M{\mathit{Q}}_{i-1}^M+a_i^M{\mathit{Q}}_i^M{\mathit{Q}}_{i+1}^M)$ , $a_i^M=\frac{D_3({\mathit{Q}}_i^M,{\mathit{Q}}_{i-1}^M)}{D_3({\mathit{Q}}_i^M,{\mathit{Q}}_{i+1}^M)}$ .

The interpolation requirements can be represented using equations (9) and (10).

$$D_3({\mathit{R}}_i^M,{\mathit{P}}_i)\le {\epsilon }_{\max }$$

(9)

$$D_{\theta }({\mathit{R}}_i^M,{\mathit{P}}_i)\le o_{\max }$$

(10)

The position chord error around the i th corner can be developed by equation (11) as mentioned in He et al. ² and Fan et al. ²²

$${\beta }_i^M\le \frac{{\delta }_{\max }\left|{\mathit{Q}}_{i,\left[3\right]}^M{\mathit{Q}}_{i+1,\left[3\right]}^M\right|}{|{\mathit{Q}}_{i,\left[3\right]}^M{\mathit{Q}}_{i+1,\left[3\right]}^M\times {\mathit{Q}}_{i,\left[3\right]}^M{\mathit{Q}}_{i-1,\left[3\right]}^M|}$$

(11)

As for the shape-preserving constraint, as introduced in Huang et al. ¹⁷ a blank ratio $\rho \in (0,1)$ defines blank region between two adjacent transition arcs. The transition ratios $r_i^k$ and $a_i^k{r}_i^k$ should be smaller than half of $(1-\rho )$ . After derivation, the shape-preserving constraint can be represented by equation (12).

$${\beta }_i^M\le \min \{\frac{\rho (1-\cos \left(\frac{{\phi }_i^M}{2}\right))}{2(1-\cos {\phi }_i^M)},\frac{\rho (1-\cos \left(\frac{{\phi }_i^M}{2}\right))}{2a_i^M(1-\cos {\phi }_i^M)}\}$$

(12)

To reach the requirements of equations (9)–(12), a geometric iteration method is designed to update iterative pseudo path and transition parameter. In each iteration, a constrained optimization problem is built to compute iterative pseudo path and transition parameters.

Firstly, the objective function at the k th ( $k=0,1,\dots M$ ) iteration and i th ( $i=1,2,\dots N-1$ ) target point is the minimum of data errors. Suppose ${\mathit{R}}_i^k={\mathit{C}}_i^k(0.5)={\mathit{Q}}_i^k+{\beta }_i^k({\mathit{Q}}_i^k{\mathit{Q}}_{i-1}^k+a_i^k{\mathit{Q}}_i^k{\mathit{Q}}_{i+1}^k).$ The position data error is $f_{i,p}^k=D_3({\mathit{R}}_i^k,{\mathit{P}}_i)$ , and the orientation data error is $f_{i,o}^k=D_{\theta }({\mathit{R}}_i^k,{\mathit{P}}_i)$ . The objection function can be written by equation (13). The coefficient s is a weight between position and orientation.

Then, the chord tolerance constraint and shape-preserving constraint are used to constrain ${\beta }_i^k$ in the k th iteration by equations (14) and (15).

Furthermore, the iteration must be convergent for maximum position data error and maximum orientation data error, which means the k th data error at each target point should be smaller than the maximum data error at ( $k-1)$ th iteration. These two constraints can be represented by equations (16) and (17).

Based on the above analysis, a constrained optimization problem can be built for the i th target point in the k th iteration. The objection function is equation (13), under constraints of equations (14)–(17).

$$\mathrm{Minimize}(f_{i,p}^k+{\mathrm{sf}}_{i,o}^k)$$

(13)

$${\beta }_i^k\le \frac{{\delta }_{\max }\left|{\mathit{Q}}_{i,\left[3\right]}^k{\mathit{Q}}_{i+1,\left[3\right]}^k\right|}{|{\mathit{Q}}_{i,\left[3\right]}^k{\mathit{Q}}_{i+1,\left[3\right]}^k\times {\mathit{Q}}_{i,\left[3\right]}^k{\mathit{Q}}_{i-1,\left[3\right]}^k|}$$

(14)

$${\beta }_i^k\le min\{\frac{\rho (1-\cos \left(\frac{{\phi }_i^k}{2}\right))}{2(1-\cos {\phi }_i^k)},\frac{\rho (1-\cos \left(\frac{{\phi }_i^k}{2}\right))}{2(1-\cos {\phi }_i^k)}.\frac{1}{a_i^k}\},a_i^k=\frac{D_3({\mathit{Q}}_i^k,{\mathit{Q}}_{i-1}^k)}{D_3({\mathit{Q}}_i^k,{\mathit{Q}}_{i+1}^k)}$$

(15)

$$D_3({\mathit{R}}_i^k,{\mathit{P}}_i)<\mathrm{ma}{x}_{0<i<N}\left(D_3({\mathit{R}}_i^{k-1},{\mathit{P}}_i)\right)$$

(16)

$$D_{\theta }({\mathit{R}}_i^k,{\mathit{P}}_i)<\mathrm{ma}{x}_{0<i<N}\left(D_{\theta }({\mathit{R}}_i^{k-1},{\mathit{P}}_i)\right)$$

(17)

By derivation, the constrained optimization problem is a quadratic optimization problem of ${\beta }_i^k$ . As the transition parameters $\{{\beta }_i^k{\}}_{i=1}^{N-1}$ and iterative pseudo path $\{{\mathit{Q}}_i^k{\}}_{i=0}^N$ are both variables in each iteration. The general iterative method, such as Newton’s method, cannot be used directly. In the following, the pseudo path will be computed firstly, and then the transition parameters are computed by solving a general quadratic optimization problem.

Details of computing pseudo path and transition parameters

To solve the constrained optimization problem, the pseudo path in the k th iteration is updated by equation (18).

$${\mathit{Q}}_i^k={\mathit{Q}}_i^{k-1}+{\mu }^k{\mathrm{\Delta }}_i^k,i=1,2,\dots N-1$$

(18)

In which, ${\mu }^k$ is a step, ${\mathrm{\Delta }}_i^k$ is an adjusting vector. Two methods of calculating adjusting vector and the range of step will be introduced in Section 3.2.3. An inner iteration which describes details of step 1.2 is designed by updating ${\mu }^k$ .

The range of ${\mu }^k$ can be estimated by adjusting vectors. The details will be introduced in Section 3.3.3. Then the pseudo path $\{{\mathit{Q}}_i^k{\}}_{i=0}^N$ can be computed by equation (18). Furthermore, transition parameters $\{{\beta }_i^k{\}}_{i=1}^{N-1}$ can be computed by solving the constrained optimization problem (equations (13)–(17)). If there are no solutions and current iteration number dose not exceed maximum iteration number, the step is updated by one-dimension search method.

Adjusting vector selection and step range estimation

Based on the methods of LSPIA ⁹ and WPIA ³⁸ with constant coefficient. This section proposes two adjusting vector selection methods.

Equation (19) can be used to represent the target point $\{{\mathit{P}}_i{\}}_{i=0}^N$ , iterative pseudo path $\{{\mathit{Q}}_i^k{\}}_{i=0}^N$ , coefficient matrix ${\mathit{A}}^k$ , and adjusting matrix ${\mathbf{\Delta }}^k$ . In which, ${\gamma }_i^k=a_i^k{\beta }_i^k$ , ${\alpha }_i^k=1-{\beta }_i^k-{\gamma }_i^k$ .

$$\mathit{P}=\left[\begin{array}{c}{\mathit{P}}_0\hfill \\ \multicolumn{1}{c}{{\mathit{P}}_{\mathit{N}}}\end{array}\right],{\mathit{Q}}^{\mathit{k}}=\left[\begin{array}{c}{\mathit{Q}}_0^{\mathit{k}}\hfill \\ \multicolumn{1}{c}{{\mathit{Q}}_{\mathit{N}}^{\mathit{k}}}\end{array}\right],{\mathit{A}}^{\mathit{k}}=\left(\begin{array}{cccccc}1\hfill & 0\hfill & 0\hfill & 0\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{{\beta }_1^k}& {\alpha }_1^k\hfill & {\gamma }_1^k\hfill & 0\hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{0}& {\beta }_2^k\hfill & {\alpha }_2^k\hfill & \dots \hfill & \dots \hfill & 0\hfill \\ \multicolumn{1}{c}{⋮}& \ddots \hfill & \ddots \hfill & \ddots \hfill & \ddots \hfill & ⋮\hfill \\ \multicolumn{1}{c}{}& \hfill & \hfill & {\beta }_{N-1}^k\hfill & {\alpha }_{N-1}^k\hfill & {\gamma }_{N-1}^k\hfill \\ \multicolumn{1}{c}{0}& 0\hfill & 0\hfill & \dots \hfill & 0\hfill & 1\hfill \end{array}\right),{\mathit{\Delta }}^{\mathit{k}}=\left[\begin{array}{c}{\mathrm{\Delta }}_0^{\mathit{k}}\hfill \\ \multicolumn{1}{c}{{\mathrm{\Delta }}_{\mathit{N}}^{\mathit{k}}}\end{array}\right]$$

(19)

The first method is VCWPIA (Variable Coefficient Weighted Progressive Iterative Approximation).The adjusting vector can be represented by equation (20). The adjusting vector of VCWPIA considers the error vector of one target point, as shown in Figure 5(a).

$${\mathbf{\Delta }}^k=\mathit{P}-{\mathit{A}}^{k-1}{\mathit{Q}}^{k-1}$$

(20)

Figure 5.

Adjusting vectors computation of VCWPIA and VCLSPIA: (a) VCWPIA method and (b) VCLSPIA method.

The second method is VCLSPIA (Variable Coefficient Least-Squares Progressive Iterative Approximation). The adjusting vector of VCLSPIA considers the error vector of adjacent three target points, which can be represented by equation (21). The diagram can be seen in Figure 5(b).

$${\mathbf{\Delta }}^k={{\mathit{A}}^{k-1}}^T(\mathit{P}-{\mathit{A}}^{k-1}{\mathit{Q}}^{k-1}),$$

(21)

The range of step ${\mu }^k$ are computed according to the eigenvalue of ${\mathit{A}}^{k-1}$ or ${{\mathit{A}}^{k-1}}^T{\mathit{A}}^{k-1}$ . ⁹

For the two methods, the VCWPIA is prefer to be used as an interpolation method while VCLSPIA as an approximation method. An example is used to compare smoothing path generated by the two methods using same tolerance constraints.

As shown in Figure 6, the input linear path has five target points. The position path using VCWPIA and VCLSPIA are illustrated in Figure 6(a), the rotation angle path relative to the parameter of position arc length are plotted in Figure 6(b). It demonstrates that the smooth path generated by VCLSPIA is closer to linear path than by VCWPIA.

Figure 6.

Smoothing results comparison between VCWPIA and VCLSPIA: (a) position point path before and after robot path smoothing and (b) rotation angle path before and after robot path smoothing.

In the simulations of Section 4, the method of VCWPIA will be used, because under the same chord tolerance, the arc radius of position path by VCWPIA are bigger than by VCLSPIA, which can run faster with a smaller curvature.

Simulations

The proposed robot path smoothing method for SCARA robots is implemented in a self-developed software, using Visual Studio 2013, with a Window 8.1 Pro OS, 3.60 GHZ CPU and 16.0 GB memory. The figures are plotted in MATLAB R2014a. The linear robot paths are generated by a robot offline programming and simulation software “HiperMOS”. ³⁹

In this section, the path before and after smoothing are illustrated and analyzed in Section 4.1. Then the smoothing method is compared with G¹ smoothing results of Yaskawa Motoman robotics and ABB robotics.

SCARA robot path before and robot path smoothing

The path of a planar “dolphin” shaped model is generated by HiperMOS with 251 target points. The orientations are changed with the position path. In the smoothing algorithm, the position data point tolerance is set to 0.001 mm, the orientation data point tolerance is set to 1°, the chord tolerance is set to 0.1 mm. The shape-preserving parameter is set to $\rho =0.1$ .

Figure 7 is the position paths before and after robot path smoothing. The linear robot path is plotted with black segments. The smooth robot path is consisting of segments (green lines) and arcs (red curves). It can be seen that the smooth position path has good smoothness and is shape-preserving.

Figure 7.

Position robot path before and after smoothing: (a) overview, (b) marker A and (c) marker B.

The rotation angle about position arc length after smoothing is G1 continuity. Figure 8 verifies the smoothness of rotation angles. Moreover, the first derivation of robot path about position arc length ( $\frac{\mathrm{dx}}{\mathrm{ds}},\frac{\mathrm{dy}}{\mathrm{ds}},\frac{d\theta }{\mathrm{ds}}$ ) is computed and illustrating in Figure 9. $\frac{\mathrm{dz}}{\mathrm{ds}}$ is not illustrated in this case because this example is a planar path ( $\frac{\mathrm{dz}}{\mathrm{ds}}=0$ for all points). It can be seen that the first derivation before smoothing is G⁰ discontinuous. After smoothing, the first derivation for $x,y,\theta$ is G⁰ continuous. It demonstrates that the path after smoothing is G¹ continuous, both for positions and for orientations.

Figure 8.

Rotation angles before and after smoothing: (a) overview and (b) amplified view around marker C.

Figure 9.

First derivation about position path: (a) overview and (b) amplified view around marker D.

The data point tolerance is satisfied both for positions and rotation angles, which can be seen in Table 1. The maximum position data error and maximum orientation data error are both reduced in each iteration and reach the data point tolerance requirements at the fourth iteration.

Table 1.

Data errors in the iterations.

	Data point tolerance	Maximum data errors
		Iteration 1	Iteration 2	Iteration 3	Iteration 4
Position [mm]	0.00100	0.33987	0.00946	0.00169	0.00031
Rotation angle [degree]	1.00000	14.87591	1.42286	0.59554	0.59554

The actual position chord errors are plotted in Figure 10. The maximum chord error is 0.09998 mm, and the average chord error is 0.05057 mm. It can be seen that the smooth path is chord tolerance constrained.

Figure 10.

Chord tolerance and real chord errors.

The computational efficiency is estimated by Visual Studio 2013. The time of smoothing the 251 target points is 0.00162 s (average time for one point is 6.732e-6 s). It demonstrates that the proposed smoothing method has high computational efficiency, thus can be used for real-time motion control.

Comparison with smoothing method of Motoman

Yaskawa Motoman robotics use a code “CR” to approach linear path by local arc transition. The arc radius can be specified by the end-user. The simulation will compare the smoothing path generated by CR and by the proposed method.

A planar path of a character “Gong” is used. The position data point tolerance is set to 0.1 mm. The rotation angle tolerance is set to 1°. The chord tolerance is set to 3 mm. The CR is simulated run by the software MotoSimEG-VRC 5.2.

Figure 11(a) and (b) are the position path and orientation path generated by CR and by this research. The smooth path generated by the proposed method interpolates target points, while the smooth path by CR has target point data errors. Moreover, the rotation angles about position arc length of this research are smoother than of CR path.

Figure 11.

Smoothing results comparison between CR path of Yaskawa Motoman and this research: (a) position path comparison, (b) rotation angle path comparison, (c) curvature profiles, and (d) chord error profiles.

Figure 11(c) demonstrates the curvature profiles generated by this research are smoother by CR path. Moreover, by comparison of chord errors in Figure 11(d), the proposed smooth path is chord error confined, while the CR path is chord error uncontrolled by the end-user.

The comparison demonstrates that smooth path generated in this research is easy to control errors, and the rotation angle path is smoother than by CR code.

Comparison with smoothing method of ABB

ABB robotics uses a parabola curve to transit corners. The transition path is generated by specifying a transition level. ³⁴ The real target point error cannot be constrained by the end-user directly. This section compares the smoothing results of ABB Z200 with smoothing results in this research.

A bigger character “Gong” is used. The smoothing parameters are set the same as the experiment in Section 4.2. The Z200 path is generated by simulation software RobotStudio 6.08 and a robot IRB910SC.

Figure 12(a) and (b) are the position path and rotation path generated by Z200 and this research. The smooth path generated by this research interpolates target points while Z200 path transits linear path at corners. Moreover, smooth path by this research is closer to the linear path than by ABB Z200 path.

Figure 12.

Smoothing results comparison between Z200 path of ABB and this research: (a) position path comparison, (b) rotation angle path comparison, (c) curvature profiles, and (d) chord error profiles.

Figure 12(c) is the comparison of curvature profile. It can be seen that curvature of smooth path by this research is more evenly and smaller than by Z200 path. The actual chord errors are computed and compared in Figure 12(d). It can be seen that chord errors of this research are confined by the tolerance 3 mm, while the maximum chord error of Z200 is 9.77634 mm. It demonstrates that the generated smooth path by this research has better accuracy and can be easier controlled by the end-user.

To summarize the above simulations, all the mentioned requirements are satisfied. Comparisons with smoothing methods of Motoman and ABB, the generated smooth paths by this research are smoother and tolerance constrained than other methods.

Experimental results

The generated G¹ smooth path is actually verified by drawing through a SCARA robot. Three types of robot path (the linear path, the corner transition path generated by the robot controller, and the G¹ smooth path generated by this research) are drawing together in the experiments. The drawing results and actual errors are compared.

Experiment platform

A “Heart” shaped planar linear path with nine target point is designed in the offline programming software “HiperMOS”. To smooth the linear path, the position data point tolerance is set to 0.01 mm, and the chord tolerance is set to 2 mm. The orientation is not changed in the experiment.

The linear path and G¹ smooth path are simulated in software HiperMOS with a TURIN SCARA robot (STH030S/E-3KG-600 mm). Figure 13(a) is the illustration of the SCARA robot and G¹ smooth path in HiperMOS. Figure 13(b) is the actual SCARA robot and the setup of drawing. A pen is installed on the end effector of robot by a clamp. A piece of A4 paper is fixed in the work table.

Figure 13.

Drawing by a SCARA robot: (a) simulation in HiperMOS and (b) actual running.

Experiment results

The three types of path are drawn by three different color of pens. Firstly, the linear path is drawn by a black pen, then the corner transition path generated by TURIN robot controller is drawn by a blue pen, at last the G¹ smooth path is drawn by a red pen. Two drawings are presented in Figure 14, corresponding to experiments performed for two different speed settings of robot.

Figure 14.

Three types of path by robot drawing: (a) first experiment and (b) second experiment.

A first experiment is designed with speed of 40 mm/s, and the acceleration of 200 mm/s². The corner transition path is generated by robot controller with maximum transition. The drawing results of three types of path are shown in Figure 14(a). It can be seen that the blue transition path do not interpolate the target points. By measurement, the position error is about 1 mm. while the red G¹ smooth path interpolates the targets points under the chord tolerance of 2 mm. It demonstrates that the resulting path has advantage of target point interpolation and chord errors constrained.

A second experiment is designed with speed of 80 mm/s, and the acceleration of 400 mm/s². The corner transition path is generated by the same method with the first experiment. Figure 14(b) presents the drawing results of three types of path. It can be seen that the transition path changes transition shape with the change of speed and acceleration. By measurement, the position error is about 3 mm. However, the red G¹ smooth path has same shape as in the first experiments.

The two experiments show the advantage of the smooth path relative to corner transition path of TURIN robot controller, i.e. easy-controlled tolerance and unchanged shape with different speed. It indicates that the G¹ smooth path is fit for precision-demanding application.

Conclusion

This research optimizes any given SCARA linear robot path to achieve $G^1$ continuity for both position and orientation, to increase precision and efficiency. A convex combination formulation is used to represent local arc transition path. A geometric iterative method is proposed to generate local arc transition smooth path. Compared with the existing works, the proposed method has the following advantages: (1) The generated smooth path is G¹ continuous, both for position path and for orientation path with unified parameterization. (2) The resulting smooth path satisfies high precision and shape-preserving requirements simultaneously, which is suitable for high-precision-demanding robotic applications. (3) The smoothing method has high computational efficiency, thus can be used for general industrial implementation and real-time motion control.

Simulations and experiments demonstrate the effectiveness of the proposed SCARA robot path arc interpolation. Comparison with similar methods demonstrates the proposed method can generate smoother path while satisfying all the requirements.

At the consideration of future work, the proposed method will be used for six-axis robot path smoothing.