Smooth Interpolation Design with Consideration of Corner Tolerance Constraints for Robotics

: This paper presents a novel method for interpolation design that ensures the continuity of a velocity proﬁle and satisﬁes a speciﬁed corner tolerance constraint. The method uses an S-shaped proﬁle to generate trajectories for each line segment in the task space. The velocity proﬁles of each segment are overlapped to control the smoothness of the corners and reduce the cycle time. This study deﬁned an overlapping time parameter that is associated with the corner tolerance and the cycle time. Moreover, a corner tolerance constraint equation was derived that can allow for a given tolerance to be satisﬁed. This constraint equation enables the use of the proposed velocity proﬁle overlap (VPO) method to specify corner tolerances for each corner of the trajectory. The proposed method was compared against the conventional acceleration/deceleration after interpolation (ADAI) method. The results demonstrate that the proposed VPO method can achieve higher accuracy and lower cycle time than the ADAI method.


Introduction
Robots are widely used in various industrial applications, such as pick and place, welding, grinding, polishing, and processing. Appropriate trajectory planning is critical in these applications to prevent extremely high acceleration in the velocity profile, thereby avoiding excessive torque output, which can cause machine vibration and damage components. An acceleration/deceleration (ACC/DEC) constraint should therefore be applied when designing interpolators to ensure smooth robot movements. An interpolator generates position inputs for the servo drives in each sampling period for a given trajectory command.
Robot interpolation methods are usually divided into joint space and task space planning approaches. In the joint space planning method, inverse kinematics are used to compute moving angles for each joint. Each joint angle is then simultaneously interpolated to generate appropriate command positions to cause each joint to reach the task position. This approach is used in pick-and-place applications in which only the initial and final positions and the posture of the end-effector are relevant factors. In this approach, the continuity of the velocity and posture of the end-effector during the motion trajectory are not considered; thus, this approach may be unsuitable for grinding, welding, or polishing processes. In the task space method of planning the end-effector trajectory, both the continuity of velocity and posture rotation are considered. The trajectory of an end-effector typically comprises many line segments, and designing the velocity profile at the intersections of these line segments is the key challenge in the design of an interpolator.
Several approaches have been used to achieve smooth motion along a trajectory [1], such as designing trapezoidal or bell-shaped velocity profiles by using cubic, quartic, or quintic polynomials [2]. These profiles constrain the level of acceleration and jerk, allowing the vibrations of the robot to be used for useful work [3][4][5][6][7]. Rossi et al. [4] proposed that the VPO method. Finally, the VPO method was compared with the ADAI method [5] to demonstrate its advantages. The conclusions are presented in Section 6.

Linear Motion Planning
Point-to-point (P2P) and linear motion are the two most common types of robot trajectories. P2P is often used in pick-and-place applications. For P2P motion, only the initial and final positions are set; the trajectory is not specified. Linear motion requires designing an end-effector trajectory in the task space. The interpolated points at the endeffector are transformed into joint commands through inverse kinematics to achieve linear motion. Because velocity profiles designed for linear motion must include both translation and orientation commands, these commands should by synchronized in the interpolator (Figure 1a). Table 1 presents the toolpath information for Figure 1a, namely the position, orientation, and feed rate of the robot manipulator. The orientation is represented by A, B, and C, which indicate rotations around the X-, Y-, and Z-axes, respectively (Figure 1b).
Appl. Sci. 2023, 13, x FOR PEER REVIEW 3 of 19 was then derived, as described in Section 4, to enable recalculating the overlapping time and satisfy the corner tolerance constraint. Section 5 presents the results of simulations and experiments conducted on the HIWIN RT605 robot manipulator to validate the effectiveness of the VPO method. Finally, the VPO method was compared with the ADAI method [5] to demonstrate its advantages. The conclusions are presented in Section 6.

Linear Motion Planning
Point-to-point (P2P) and linear motion are the two most common types of robot trajectories. P2P is often used in pick-and-place applications. For P2P motion, only the initial and final positions are set; the trajectory is not specified. Linear motion requires designing an end-effector trajectory in the task space. The interpolated points at the end-effector are transformed into joint commands through inverse kinematics to achieve linear motion. Because velocity profiles designed for linear motion must include both translation and orientation commands, these commands should by synchronized in the interpolator (Figure 1a). Table 1 presents the toolpath information for Figure 1a, namely the position, orientation, and feed rate of the robot manipulator. The orientation is represented by A, B, and C, which indicate rotations around the X-, Y-, and Z-axes, respectively (Figure 1b).
A robot's orientation is determined by a combination of rotations in the X, Y, and Z directions. Various methods can be used to describe the orientation of a rigid body, such as the Euler angle, roll-pitch-yaw angle, angle-axis representation, unit quaternion, and Cayley-Rodrigues parameter methods. In this study, the equivalent angle-axis representation was adopted; this method is similar to the quaternion method but involves fewer parameters. The rotation matrix in the equivalent angle-axis method can be represented by a unit vector (u) and an angle (ϕ) of revolution about the u vector.
The parameters u and ϕ are used to rotate the orientation of the end-effector from XYZ coordinates to X′Y′Z′ coordinates (Figure 2a). To obtain u and ϕ, rotation matrices for A robot's orientation is determined by a combination of rotations in the X, Y, and Z directions. Various methods can be used to describe the orientation of a rigid body, such as the Euler angle, roll-pitch-yaw angle, angle-axis representation, unit quaternion, and Cayley-Rodrigues parameter methods. In this study, the equivalent angle-axis representation was adopted; this method is similar to the quaternion method but involves fewer parameters. The rotation matrix in the equivalent angle-axis method can be represented by a unit vector (u) and an angle (φ) of revolution about the u vector.
The parameters u and φ are used to rotate the orientation of the end-effector from XYZ coordinates to X Y Z coordinates (Figure 2a). To obtain u and φ, rotation matrices for the XYZ and X Y Z coordinates are defined as R 1 and R 2 . The matrices R 1 and R 2 can be Appl. Sci. 2023, 13, 8789 4 of 18 obtained from two arbitrary orientation commands (A, B, C) on the toolpath. A general rotation matrix R can be calculated using Equation (1).
To convert a rotation matrix R back to an orientation command, the rotation matrix can be inverted by using Equations (2)-(4).
To transform R 1 into R 2 , R 1 must be multiplied by a rotation matrix R 12 . The relationship between R 1 and R 2 can be expressed as follows, where R 12 represents the rotation matrix from R 1 to R 2 : The matrices R 1 and R 2 can be obtained from the orientation commands P i and P i+1 . Matrix R 12 can then be calculated from Equation (6). However, the parameters u and φ are still unknown. The equivalent angle-axis representation can be employed to determine matrix R 12 by applying R φ,u .  11 To convert a rotation matrix R back to an orientation command, the rotation matrix can be inverted by using Equations (2)-(4). (2) ( ) ( ) 23 To transform R1 into R2, R1 must be multiplied by a rotation matrix R12. The relationship between R1 and R2 can be expressed as follows, where R12 represents the rotation matrix from R1 to R2: The matrices R1 and R2 can be obtained from the orientation commands Pi and Pi+1. Matrix R12 can then be calculated from Equation (6). However, the parameters u and ϕ are still unknown. The equivalent angle-axis representation can be employed to determine matrix R12 by applying R ϕ ,u. As shown in Figure 2b, the angle between the u vector and the i-k plane is represented by γ, and β represents the angle between the u vector and the k-axis. The projection of Appl. Sci. 2023, 13, 8789 5 of 18 the u vector onto the i-, j-, and k-axes is represented by u i , u j , and u k , respectively. The following relationships can be determined from Figure 2b: The rotation matrix R φ,u , which describes the rotation φ around vector u, can be obtained as follows. First, vector u is rotated by angle γ on the i-axis to vector u located on the ik plane; the corresponding rotation matrix is represented by R γ,i . Vector u is then rotated by angle −β on the j-axis to vector u , which is located on the k-axis. This rotation matrix is represented by R − β,j . Vector u is then rotated by φ on the k-axis; the rotation matrix is R φ,k . The reverse sequence of rotations and their respective opposite angles can be expressed as R β,j and R −γ,i . After this sequence of rotations, the matrix R φ,u can be represented as where sγ = sinγ, cγ = cosγ, sβ = sinβ, cβ = cosβ, sφ = sinφ, and cφ = cosφ.
By substituting Equations (7)-(10) into Equation (12), the following equation can be derived, where the abbreviation vφ = 1 − cosφ is used: 2 vφ + cφ u x u y vφ − u z sφ u x u z vφ + u y sφ u x u y vφ + u z sφ u y 2 vφ + cφ u y u z vφ − u x sφ u x u z vφ − u y sφ u y u z vφ + u x sφ u z 2 vφ + cφ The rotational angle and the unit vector can be obtained from Equations (14) and (15). After u and φ are obtained, φ can be used as the variable for interpolation.
As shown in Figure 1, the toolpath can be used to compute the translation distance L and the rotational angle φ determined using Equations (1)- (15). The S-shaped design method can be employed to generate a velocity profile for L or φ. The velocity profile can be divided into five phases ( Figure 3) by Equation (16), where S i represents the linear motion distance L or rotational angle φ of the end-effector.
where J m represents the maximum jerk, V m represents the maximum feed rate, and V s and V e represent initial and final velocity (usually 0), respectively.
where Jm represents the maximum jerk, Vm represents the maximum feed rate, and Vs and Ve represent initial and final velocity (usually 0), respectively. For example, consider the toolpath from Pi to Pi+1 for the position command for translation from (368, 0, 293.5) to (368, 200, 100) and the orientation command for rotation from (180, 0, 90) to (150, 0, 80). The interpolation should simultaneously plan both the translation and orientation commands. Table 2 presents the constraints on the interpolation parameters, where Vlf and Vof represent the maximum feed rate for translation and orientation, respectively, and Jlm and Jom represent the maximum jerk for translation and orientation, respectively. The velocity profile for the translation is planned first; Equation (17) is used to determine if the maximum jerk limit Jlm is exceeded. If so, the acceleration time can be adjusted using Equation (18), and the orientation can then be planned in accordance with the translation results. If the angular jerk still does not exceed the limit, planning is complete. However, if the angular jerk exceeds the limit, orientation must be planned first, and translation must be allocated afterwards in the same manner. In this example, the translation distance is 278.28 mm and the rotational angle is 31.58°; the results indicate that translation can be planned first.  For example, consider the toolpath from P i to P i+1 for the position command for translation from (368, 0, 293.5) to (368, 200, 100) and the orientation command for rotation from (180, 0, 90) to (150, 0, 80). The interpolation should simultaneously plan both the translation and orientation commands. Table 2 presents the constraints on the interpolation parameters, where V lf and V of represent the maximum feed rate for translation and orientation, respectively, and J lm and J om represent the maximum jerk for translation and orientation, respectively. The velocity profile for the translation is planned first; Equation (17) is used to determine if the maximum jerk limit J lm is exceeded. If so, the acceleration time can be adjusted using Equation (18), and the orientation can then be planned in accordance with the translation results. If the angular jerk still does not exceed the limit, planning is complete. However, if the angular jerk exceeds the limit, orientation must be planned first, and translation must be allocated afterwards in the same manner. In this example, the translation distance is 278.28 mm and the rotational angle is 31.58 • ; the results indicate that translation can be planned first.  Figure 4 reveals that the resulting single program command achieves both translation and orientation. After the interpolation points (P x , P y , P z ) and the rotation matrix (R φ,u R 1 ) for each sampling period are obtained, they can be substituted into Equation (19) to obtain the rotation matrix ( 0 6 R). The joint commands can then be solved through inverse kinematics and used to command the motor drives, achieving synchronized translation and orientation for linear motion.  Figure 4 reveals that the resulting single program command achieves both translation and orientation. After the interpolation points (Px, Py, Pz) and the rotation matrix (R ϕ ,uR1) for each sampling period are obtained, they can be substituted into Equation (19) to obtain the rotation matrix ( ). The joint commands can then be solved through inverse kinematics and used to command the motor drives, achieving synchronized translation and orientation for linear motion.

Corner Smoothing
In addition to single-line segment motion in positioning applications, corner interpolation must also be considered for tracking multi-line segment motion in polishing and grinding applications. Corner smoothing techniques can be employed to improve machining efficiency and reduce vibration at junctions. The motion of end-effector translation and orientation should be designed to ensure smoothness and continuity throughout the entire trajectory. Conventional methods for generating a smooth trajectory without excessive ACC/DEC can be categorized as ADBI and ADAI methods [5]. The ADBI method constrains the acceleration and the jerk of the trajectory; this prevents specifying the corner tolerance. However, the ADBI method can cause velocity discontinuities at block junctions, which can result in vibration. As shown in Figure 5, the feed rate (tangential velocity) is distributed across x-and y-axes in accordance with the direction of movement. The axis velocities Vx and Vy at the junction are not continuous.

Corner Smoothing
In addition to single-line segment motion in positioning applications, corner interpolation must also be considered for tracking multi-line segment motion in polishing and grinding applications. Corner smoothing techniques can be employed to improve machining efficiency and reduce vibration at junctions. The motion of end-effector translation and orientation should be designed to ensure smoothness and continuity throughout the entire trajectory. Conventional methods for generating a smooth trajectory without excessive ACC/DEC can be categorized as ADBI and ADAI methods [5]. The ADBI method constrains the acceleration and the jerk of the trajectory; this prevents specifying the corner tolerance. However, the ADBI method can cause velocity discontinuities at block junctions, which can result in vibration. As shown in Figure 5, the feed rate (tangential velocity) is distributed across xand y-axes in accordance with the direction of movement. The axis velocities V x and V y at the junction are not continuous. The ADAI method, in which digital convolution is used eliminate velocity discontinuities on each axis, can resolve this problem ( Figure 6). Corner tolerance occurs because the velocities of blocks N1 and N2 overlap. Digital convolution is performed using Equation (20), where = / , represents the sample time and represents the time constant of the ACC/DEC. The input signal represents a velocity command. The ADAI method, in which digital convolution is used eliminate velocity discontinuities on each axis, can resolve this problem ( Figure 6). Corner tolerance occurs because the velocities of blocks N 1 and N 2 overlap. Digital convolution is performed using Equation (20), where N = τ/T s , T s represents the sample time and τ represents the time constant of the ACC/DEC. The input signal x(kT s ) represents a velocity command. Figure 5. ADBI method.
The ADAI method, in which digital convolution is used eliminate velocity discontinuities on each axis, can resolve this problem ( Figure 6). Corner tolerance occurs because the velocities of blocks N1 and N2 overlap. Digital convolution is performed using Equation (20), where = / , represents the sample time and represents the time constant of the ACC/DEC. The input signal represents a velocity command.
Although the ADAI method is commonly used for interpolators, its applicability is limited because it can only be applied to the entire trajectory; the tolerance for each corner cannot be specified. To achieve both corner smoothing and a specified tolerance function, this paper proposes a VPO algorithm that calculates the overlapping time of line segments to produce blended S-shaped ACC/DEC profiles. This approach improves cycle times and enables the trajectory to be systematically designed to meet different accuracy requirements.
An example of the overlapping time of line segments is presented in Figure 7, which depicts two identical blocks with velocity blending. First, blocks 1 and 2 are planned using S-shaped ACC/DEC. Block 2 is then blended into block 1, and Tol is the overlapped area of the two blocks; this represents as overlapping time. To calculate the overlapping time, the OVLP parameter is defined as a percentage. If OVLP is 100%, the overlapping time is 2Trr, and the junction of the two blocks occurs exactly at half of Vm. Therefore, 2Vc is equal to Although the ADAI method is commonly used for interpolators, its applicability is limited because it can only be applied to the entire trajectory; the tolerance for each corner cannot be specified. To achieve both corner smoothing and a specified tolerance function, this paper proposes a VPO algorithm that calculates the overlapping time of line segments to produce blended S-shaped ACC/DEC profiles. This approach improves cycle times and enables the trajectory to be systematically designed to meet different accuracy requirements.
An example of the overlapping time of line segments is presented in Figure 7, which depicts two identical blocks with velocity blending. First, blocks 1 and 2 are planned using S-shaped ACC/DEC. Block 2 is then blended into block 1, and T ol is the overlapped area of the two blocks; this represents as overlapping time. To calculate the overlapping time, the OVLP parameter is defined as a percentage. If OVLP is 100%, the overlapping time is 2T rr , and the junction of the two blocks occurs exactly at half of V m . Therefore, 2V c is equal to V m . The relationship between the OVLP parameter and the velocity can be expressed as follows: where V c is the velocity junction between the two blocks. The relationship between V c and the overlapping time (T ol ) can be calculated using the A-T diagram in Figure 7. For an S-shaped velocity profile, A m = V m /T rr , where A m represents the maximum acceleration and T rr is equal to the time of acceleration from 0 to A m within the ACC/DEC period. A c represents V c at the junction and can be calculated using Equation (22).
The parameter V c can be calculated by integrating the area of the A-T diagram and is given as follows: Vm. The relationship between the OVLP parameter and the velocity can be expressed as follows: where Vc is the velocity junction between the two blocks. The relationship between Vc and the overlapping time (Tol) can be calculated using the A-T diagram in Figure 7. For an S-shaped velocity profile, Am = Vm/Trr, where Am represents the maximum acceleration and Trr is equal to the time of acceleration from 0 to Am within the ACC/DEC period. Ac represents Vc at the junction and can be calculated using Equation (22).
The parameter Vc can be calculated by integrating the area of the A-T diagram and is given as follows: Tol can be calculated from OVLP by substituting Equation (21) into Equation (24) as follows: T ol can be calculated from OVLP by substituting Equation (21) into Equation (24) as follows: Equation (25) indicates that T ol is a function of OVLP and T rr . Because T rr is already known, OVLP can be adjusted to quickly control the smoothing of the corner, improving the cycle time. The velocity profiles on the xand y-axes for different values of OVLP are shown in Figure 8c,d, in which OVLP is equal to 0% and 80%, respectively. When OVLP is 80%, the starting time for the N 2 block moves forward by T ol on both the xand y-axes. The VPO approach ensures that the velocity is continuous at the junction for each axis, and the cycle time decreases from 1.6 to 1.38 s. However, a corner tolerance can occur at the junction if OVLP increases (Figure 8a). The next section describes the derivation of the CTC equation, which can be used to systematically evaluate the corner tolerance on the basis of the overlapping time.
is 80%, the starting time for the N2 block moves forward by Tol on both the x-and y-axes. The VPO approach ensures that the velocity is continuous at the junction for each axis, and the cycle time decreases from 1.6 to 1.38 s. However, a corner tolerance can occur at the junction if OVLP increases (Figure 8a). The next section describes the derivation of the CTC equation, which can be used to systematically evaluate the corner tolerance on the basis of the overlapping time.

CTC Algorithm
A corner tolerance occurs when the velocity profiles of line segments overlap. The maximum corner tolerance occurs at half of the overlapping time (Tol/2). Figure 9 presents the overlapping area of Figure 7. Figure 9a shows the trajectory of block 1 and block 2 at the corner after overlapping; P1, P4, and P3 and are produced from P1, P2, and P3. Here, the distance from P2 to P4 is the maximum corner tolerance. Figure 9b presents the velocity profiles of block 1 and block 2 and the overlapped velocity profile, where Vx1 and Vx2 represent the velocity profiles on the x-axis of block 1 and block 2, respectively, and the blue solid line formed by P1, P4, and P3 represents the velocity profile of the two blocks after overlapping. The maximum corner tolerance occurs at point P4.

CTC Algorithm
A corner tolerance occurs when the velocity profiles of line segments overlap. The maximum corner tolerance occurs at half of the overlapping time (T ol /2). Figure 9 presents the overlapping area of Figure 7. Figure 9a shows the trajectory of block 1 and block 2 at the corner after overlapping; P 1 , P 4 , and P 3 and are produced from P 1 , P 2 , and P 3 . Here, the distance from P 2 to P 4 is the maximum corner tolerance. Figure 9b presents the velocity profiles of block 1 and block 2 and the overlapped velocity profile, where V x1 and V x2 represent the velocity profiles on the x-axis of block 1 and block 2, respectively, and the blue solid line formed by P 1 , P 4 , and P 3 represents the velocity profile of the two blocks after overlapping. The maximum corner tolerance occurs at point P 4 . The three-dimensional corner tolerance (ε) can be calculated from the x, y, and z parameters as follows: The three-dimensional corner tolerance (ε) can be calculated from the x, y, and z parameters as follows: In Figure 9a, ε x1 is the distance from P 1 to P 2 for block 1 in the x direction and can be calculated by integrating the area of BDE in Figure 9b; this area represents the distance from P 5 to P 2 in the X direction. To obtain ε x1 , the S 1 term in Equation (16) can be replaced by Equation (27), where t is substituted with T ol /2 and J is substituted with J x1 , obtaining Equation (28). S = Jt 3 /6 (27) When considering the overlap effect, P 5 should move to a distance projected on the x-axis relative to P 4 , represented as ε x2 . The value of ε x2 represents the area under ACD and is given by Equation (29), where J x2 represents the jerk of N 2 along the x direction.
The distance between P 4 and P 2 on the X-axis is represented by ε x and can be calculated using the following equation.
Here, ε x represents the corner tolerance ε projected on the x-axis. Similarly, the corner tolerance components in the y and z directions can be calculated using Equations (31) and (32).
Summing the squares of Equation (33) produces Equation (34) After jerk and corner tolerance have been determined, T ol can be obtained using Equation (35), and the CTC equation can be expressed as follows: In Section 2, the S-shaped design method was adopted to plan the translation and orientation. In Section 3, the VPO method was applied to plan a line segment, solving the discontinuity problem. By contrast with the conventional ADAI approach, the VPO method can explicitly specify the corner tolerance for each corner. Finally, the designed velocity profile can be integrated to obtain positions and substituted into Equation (19) to generate the command points for each axis, completing the interpolation process.

Results of Simulations and Experiments
Experiments were conducted to evaluate the effectiveness of the proposed algorithm for the HIWIN RT605 robot manipulator. The robotic equipment was equipped with a PCbased controller with a real-time operating system enabling multitasking; the sampling time was set to 1 ms. The software included a human-machine interface, a numerical control (NC) interpreter, a kinematics module, an interpolator, and a multiaxis contour trajectory motion system. Figure 10 presents the hardware used in the experiment, which comprised Sanyo Denki servo motors and drivers. The six-axis servo drives were commanded through the EtherCAT protocol to ensure that the robotic equipment moved synchronously.
In Section 2, the S-shaped design method was adopted to plan the translation and orientation. In Section 3, the VPO method was applied to plan a line segment, solving the discontinuity problem. By contrast with the conventional ADAI approach, the VPO method can explicitly specify the corner tolerance for each corner. Finally, the designed velocity profile can be integrated to obtain positions and substituted into Equation (19) to generate the command points for each axis, completing the interpolation process.

Results of Simulations and Experiments
Experiments were conducted to evaluate the effectiveness of the proposed algorithm for the HIWIN RT605 robot manipulator. The robotic equipment was equipped with a PCbased controller with a real-time operating system enabling multitasking; the sampling time was set to 1 ms. The software included a human-machine interface, a numerical control (NC) interpreter, a kinematics module, an interpolator, and a multiaxis contour trajectory motion system. Figure 10 presents the hardware used in the experiment, which comprised Sanyo Denki servo motors and drivers. The six-axis servo drives were commanded through the EtherCAT protocol to ensure that the robotic equipment moved synchronously.

Simulation Results for Difference OVLP Parameters
To evaluate the VPO methodology, a simple polygon trajectory comprising six line segments was tested (Figure 11a). The corresponding NC code is listed in Table 3, and the interpolation parameters are provided in Table 4. The trajectory planning process must

Simulation Results for Difference OVLP Parameters
To evaluate the VPO methodology, a simple polygon trajectory comprising six line segments was tested (Figure 11a). The corresponding NC code is listed in Table 3, and the interpolation parameters are provided in Table 4. The trajectory planning process must satisfy both kinematic constraints and tolerance constraints, and simultaneous translation and orientation were required. The velocity and acceleration profiles of the tool center point (TCP) at a full stop in Case 1 and in Case 2 are shown in Figure 12. In the full stop scenario, movement stops completely at each corner. Although this case has no corner tolerance, it has the longest cycle time of approximately 6.37 s. In Case 1, all of the OVLP parameters were set to 100%, reducing the cycle time from 6.37 to 5.12 s (19.6%). However, the corner tolerance was larger than in other cases. In Case 2, the OVLP parameters were adjusted to control the smoothness of each corner. This method is more flexible than the ADAI method. The corner smoothing results for each case are presented in Figure 11b-f. In the VPO algorithm, the parameter OVLP controlled the level of overlap in the velocity profiles. Higher OVLP values result in smoother corner trajectories and shorter cycle times. However, increasing the OVLP parameter leads to a larger corner tolerance.  satisfy both kinematic constraints and tolerance constraints, and simultaneous translation and orientation were required. The velocity and acceleration profiles of the tool center point (TCP) at a full stop in Case 1 and in Case 2 are shown in Figure 12. In the full stop scenario, movement stops completely at each corner. Although this case has no corner tolerance, it has the longest cycle time of approximately 6.37 s. In Case 1, all of the OVLP parameters were set to 100%, reducing the cycle time from 6.37 to 5.12 s (19.6%). However, the corner tolerance was larger than in other cases. In Case 2, the OVLP parameters were adjusted to control the smoothness of each corner. This method is more flexible than the ADAI method. The corner smoothing results for each case are presented in Figure 11b-f. In the VPO algorithm, the parameter OVLP controlled the level of overlap in the velocity profiles. Higher OVLP values result in smoother corner trajectories and shorter cycle times. However, increasing the OVLP parameter leads to a larger corner tolerance.   satisfy both kinematic constraints and tolerance constraints, and simultaneous translation and orientation were required. The velocity and acceleration profiles of the tool center point (TCP) at a full stop in Case 1 and in Case 2 are shown in Figure 12. In the full stop scenario, movement stops completely at each corner. Although this case has no corner tolerance, it has the longest cycle time of approximately 6.37 s. In Case 1, all of the OVLP parameters were set to 100%, reducing the cycle time from 6.37 to 5.12 s (19.6%). However, the corner tolerance was larger than in other cases. In Case 2, the OVLP parameters were adjusted to control the smoothness of each corner. This method is more flexible than the ADAI method. The corner smoothing results for each case are presented in Figure 11b-f. In the VPO algorithm, the parameter OVLP controlled the level of overlap in the velocity profiles. Higher OVLP values result in smoother corner trajectories and shorter cycle times. However, increasing the OVLP parameter leads to a larger corner tolerance.

Simulation Results for the CTC Method
To confirm that the VPO algorithm not only satisfies the kinematic constraints but can also be used to specify the tolerance at the corner of the corner, another experiment, Case 3, was performed. The NC code for Case 3 is shown in Table 3. The velocity and acceleration of the TCP for Case 3 are shown in Figure 13a. The feed rate reached 150 mm/s with a maximum acceleration of −1430 mm/s 2 and a cycle time of 5.45 s. Figure 13b presents the profiles of axial velocity, axial acceleration, and axial jerk. The maximum axial acceleration and jerk occurred on the y-axis and were −1684 mm/s 2 and 20,000 mm/s 3 , respectively. For orientation planning, the rotational angles were calculated from Equations (1), (6), (14) and (15) Figure 14. The maximum velocity was 48.48 deg/s, the maximum acceleration was 509.7 deg/s 2 , and the maximum jerk was 5309 deg/s 3 . Hence, the position and orientation of the end-effectors of the robot manipulator satisfy the kinematic constraints. Finally, the rotation matrix and position were substituted back into Equation (19), and the joint commands were obtained through inverse kinematics. Figure 15 presents the joint commands, angular velocity, angular acceleration, and angular jerk. The maximum values of velocity, acceleration, and jerk are −69.7 deg/s, −932 deg/s 2 , and −9701 deg/s 3 , respectively.

Simulation Results for the CTC Method
To confirm that the VPO algorithm not only satisfies the kinematic constraints but can also be used to specify the tolerance at the corner of the corner, another experiment, Case 3, was performed. The NC code for Case 3 is shown in Table 3. The velocity and acceleration of the TCP for Case 3 are shown in Figure 13a. The feed rate reached 150 mm/s with a maximum acceleration of −1430 mm/s 2 and a cycle time of 5.45 s. Figure 13b presents the profiles of axial velocity, axial acceleration, and axial jerk. The maximum axial acceleration and jerk occurred on the y-axis and were −1684 mm/s 2 and 20,000 mm/s 3 , respectively. For orientation planning, the rotational angles were calculated from Equations (1), (6), (14), and (15), obtaining 17.8°, 33.2°, 49.7°, 25.1°, 22.2°, and 38.6°. The results for the x, y, and z rotations are presented in Figure 14. The maximum velocity was 48.48 deg/s, the maximum acceleration was 509.7 deg/s 2 , and the maximum jerk was 5309 deg/s 3 . Hence, the position and orientation of the end-effectors of the robot manipulator satisfy the kinematic constraints. Finally, the rotation matrix and position were substituted back into Equation (19), and the joint commands were obtained through inverse kinematics. Figure 15 presents the joint commands, angular velocity, angular acceleration, and angular jerk. The maximum values of velocity, acceleration, and jerk are −69.7 deg/s, −932 deg/s 2 , and −9701 deg/s 3 , respectively.  To validate the effectiveness of the CTC equation in the VPO algorithm, the results for Case 3 (Table 3) were analyzed. The corner tolerance is the shortest distance from the sharp corner to the ultimate trajectory. By inputting the overlapping time calculated with the CTC equation to the VPO algorithm, the corner tolerance at each corner can be correctly constrained within the specified tolerance. When the corner tolerances of A, B, C, D, and E on the trajectory in Figure 16a were set to 4.2, 3.4, 2.6, 0, and 4.0 mm, respectively, the resulting trajectories at each corner in Figure 16b-f satisfied these constraints with overlapping times of 0.237, 0.226, 0.181, 0, and 0.22 s, respectively. Hence, the simulation results confirm that the proposed VPO algorithm, unlike the conventional ADAI approach, enables separately setting corner tolerances at each corner to satisfy the specified constraints. To validate the effectiveness of the CTC equation in the VPO algorithm, the results for Case 3 (Table 3) were analyzed. The corner tolerance is the shortest distance from the sharp corner to the ultimate trajectory. By inputting the overlapping time calculated with the CTC equation to the VPO algorithm, the corner tolerance at each corner can be correctly constrained within the specified tolerance. When the corner tolerances of A, B, C, D, and E on the trajectory in Figure 16a were set to 4.2, 3.4, 2.6, 0, and 4.0 mm, respectively, the resulting trajectories at each corner in Figure 16b-f satisfied these constraints with overlapping times of 0.237, 0.226, 0.181, 0, and 0.22 s, respectively. Hence, the simulation results confirm that the proposed VPO algorithm, unlike the conventional ADAI approach, enables separately setting corner tolerances at each corner to satisfy the specified constraints.

Experimental Comparison of VPO and ADAI
The simulation results discussed in the previous sections clearly indicate that the OVLP parameter affects the cycle time and corner tolerance. In this section, experiments were conducted for the HIWIN robot RT605 using real trajectories; the results are shown in Figure 17. Trajectories were obtained using the VPO method and the conventional interpolation method (ADAI). The feed rate, tolerance, and maximum acceleration of the TCP were the same for both algorithms. For the ADAI method with a time constant of 500 ms and the VPO method with the acceleration time set to 0.25 s and OVLP to 74%, the contour and tracking errors in the x, y, and z directions were almost identical ( Figure  17c,d). However, the cycle time of the VPO method was only 8.98 s, 10.6% less than that of the ADAI method. However, if the acceleration time and OVLP parameter in the VPO method were adjusted to 0.35 s and 58% (to maintain a similar cycle time), the contour error [25] at point A decreased from 0.1841 mm to 0.1566 mm, a reduction of approximately 15.22% (Figure 17e,f). This result indicates that the trajectory planning of the VPO

Experimental Comparison of VPO and ADAI
The simulation results discussed in the previous sections clearly indicate that the OVLP parameter affects the cycle time and corner tolerance. In this section, experiments were conducted for the HIWIN robot RT605 using real trajectories; the results are shown in Figure 17. Trajectories were obtained using the VPO method and the conventional interpolation method (ADAI). The feed rate, tolerance, and maximum acceleration of the TCP were the same for both algorithms. For the ADAI method with a time constant of 500 ms and the VPO method with the acceleration time set to 0.25 s and OVLP to 74%, the contour and tracking errors in the x, y, and z directions were almost identical (Figure 17c,d). However, the cycle time of the VPO method was only 8.98 s, 10.6% less than that of the ADAI method. However, if the acceleration time and OVLP parameter in the VPO method were adjusted to 0.35 s and 58% (to maintain a similar cycle time), the contour error [25] at point A decreased from 0.1841 mm to 0.1566 mm, a reduction of approximately 15.22% (Figure 17e,f). This result indicates that the trajectory planning of the VPO algorithm is excellent. Table 5 reveals that the cycle time of the VPO method was 9.19-12.30% lower than that of the ADAI method. For Cases A, B, and C, the cycle times of the VPO method were approximately 8.98, 9.18, and 6.42 s, respectively. Table 6 reveals that the contour error in the VPO method was 6-15% lower than that of the ADAI method. For Cases D, E, and F, the maximum contour errors of the VPO method were approximately 0.1566, 0.1693, and 0.1883 mm, respectively. The statistical results presented in Tables 5 and 6 demonstrate that the VPO method outperforms ADAI in terms of both feasibility and performance.

Conclusions
This study proposed a VPO algorithm for designing robotic interpolators that uses the velocity overlap method on each axis to remove the discontinuities caused by projecting the velocity profile of the TCP onto each axis. A superior balance between smoothness and cycle time can be achieved by adjusting the overlapping time. This paper described the derivation of the CTC equation, which enables using the VPO method to individually control the tolerance for each corner. Other interpolation design methods, such as the conventional ADAI approach, do not incorporate these features. The experimental results indicate that the VPO algorithm outperformed the ADAI method by reducing the cycle time by 6.08-15.22% if the parameters of the methods were adjusted to obtain approximately equal contour error. Experiments were also conducted to confirm that the contour errors could be reduced by 9.19-12.30% if the cycle times of both algorithms were adjusted to be approximately equal. The results confirm the effectiveness of the VPO method.

Data Availability Statement:
The data presented in this study are available upon request from the corresponding author.