A Complementary Filter Design on SE(3) to Identify Micro-Motions during 3D Motion Tracking

In 3D motion capture, multiple methods have been developed in order to optimize the quality of the captured data. While certain technologies, such as inertial measurement units (IMU), are mostly suitable for 3D orientation estimation at relatively high frequencies, other technologies, such as marker-based motion capture, are more suitable for 3D position estimations at a lower frequency range. In this work, we introduce a complementary filter that complements 3D motion capture data with high-frequency acceleration signals from an IMU. While the local optimization reduces the error of the motion tracking, the additional accelerations can help to detect micro-motions that are useful when dealing with high-frequency human motions or robotic applications. The combination of high-frequency accelerometers improves the accuracy of the data and helps to overcome limitations in motion capture when micro-motions are not traceable with 3D motion tracking system. In our experimental evaluation, we demonstrate the improvements of the motion capture results during translational, rotational, and combined movements.


Introduction
Optical motion capture technology is a valuable tool when quantifying human movements and clinical, biomechanical, and industrial applications require high system accuracy [1]. While motion capture is often tied to high costs, over the last years less expensive systems have been developed to measure human or robotic motion including goniometers [2], accelerometers [3][4][5], inertia-based and electromagnetic sensors [6][7][8]. A variety of research has been conducted to report characteristics of selected sensor systems and/or to validate new technical devices. In industrial applications such as tooling, specific motion paths We previously presented a method to filter motion capture signals and to eliminate the gravitational biasing induced by the mass of an object [20]. In this paper, we propose an experimental set-up capable of capturing micro-and macro-motions. In addition, we compare our results with an external load cell, serving as the ground truth.

Multimodal Sensor Fusion for Motion Tracking
In order to capture micro-and macro-motions of an object, we propose a complementary filter algorithm. To validate the assumptions and compare the measured motions with a ground truth, we designed a 3D printed object that includes motion capture markers, an IMU, and a loadcell (see Figure 1). As explained below, the loadcell will be used to sense apparent forces derived from non-inertial accelerations, serving as ground truth for comparison with purely kinematic information derived from IMU and Motion Tracker. A marker-based motion capture system (VZ4000, Phoenix Technologies Inc., Vancouver, BC, Canada) is used to sense the 3D location of eight infrared markers rigidly attached to the object. To track the micro-motion the IMU (will be presented in Section 3.1, 3-axis gyroscope & 3-axis accelerometer) is located inside the object to measure the change in rotation and the accelerations. A loadcell (ATI mini 40, range F X,Y : ±80 N, F Z : ±240 N, T X,YZ : ±4 Nm) is placed on top of the object to capture the dynamic vibration (i.e., forces) during the movement. The loadcell serves as the ground truth to compare the estimated acceleration based on the complementary filter algorithm. As shown in Figure 1, four major reference frames shall be considered: • task frame T : located at the center of one of the side face of the box. • loadcell frame L: attached to the center of the mounting-plate of the loadcell. • IMU frame: located at the center of the box. • global frame S: defined by the motion tracker, via external markers mounted on the workbench. The x, y, and z-axis of task frame, loadcell frame, and IMU frame are simply remains aligned together as shown in Figure 1.
Marker-based optical motion capture is a widely used technology to track 3D positions and orientations. Marker occlusions are a commonly known problem as each marker must be visible to at least two cameras in each frame to establish its position. Dorfmuller, in [21], used an extended Kalman filter (EKF) to predict the missing markers based on previous marker information, while Welch et al. in [22], used the EKF to resolve occlusions based on the skeletal model. These methods require manual intervention or become ineffective in cases where markers are missing for an extended period of time. In this section, we will propose an algorithm capable of overcoming these shortcomings.

Damped Least Square Algorithm for the Filtering Motion Tracker's Signal
In kinematics, the set of all possible three-dimensional rotations and translations of a rigid body can be described via 4 × 4 homogeneous matrices constituting the Special Euclidean matrix Lie group, denoted by SE(3) for 3-dimensional Euclidean spaces [23,24]. In particular, the pose of a rigid object in task space T can be described in reference frame S via a homogeneous transformation matrix G: where p ∈ R 3 and R ∈ R 3×3 , respectively, are the translation vector and the rotation matrix representing the axes of the object frame T in global coordinates S. The matrix Lie group SO(3) of Special Orthogonal matrices, is the group of all rotations about the origin of three-dimensional Euclidean space R 3 under the operation of composition: The associated Lie algebra so(3), i.e., the tangent space to the group SO(3) at the identity matrix, is the set of antisymmetric matrices and is defined as A hat operator · : R 3 → so(3) maps 3D vectors (such as body angular velocities) into skew-symmetric matrices. In particular, every 3D vector ω = [ω 1 ω 2 ω 3 ] T is into one-to-one correspondence with a skew-symmetric matrix: The inverse of the "hat" operator is referred to as "vee" operator and denoted by (·) ∨ . The skew symmetric ω ∈ so(3) is the matrix logarithm of rotation matrix R: ω = log(R) and R = exp( ω). Occasionally, we shall use log ∨ as shorthand for (log(·)) ∨ , i.e., ω = log ∨ (R).
In this section, eight infrared markers are located on a rigid object with the relative position known based on computer-aided design (CAD) drawing (see Figure 1). The position of ith marker in the task space T is constant and can be defined as The position of the ith marker with respect to the global frame S changes with G and can be determined as One way to solve the marker occlusion problem, is the damped least square ( Figure 2) method which considers each marker as a spring with stiffness k i . When the position of the marker changes, there is a difference between its previous position m i (t − 1) and its current position m i (t), and consequently k i . The time interval between two samples t n ≤ t ≤ t n+1 is divided into N intervals with the same length. If the marker deviates from its correct position, due to the stiffness a force will pull each marker back to the correct position. The larger the interval N, the higher the accuracy of the algorithm.
Let G * be the estimate pose of the tool frame T , then from Equation (6), the estimated position of the ith marker is given as we can define this as where F S i is the force in the reference space, M T i is the torque in task space of marker ith, and ζ T ∈ R 3 is the center of mass in task space T measured from the mechanical design. The angular velocity of the tool generated by this change can be defined as [23] where ω T ∈ R 3 is the tools angular velocity vector. The damping of rotation b r will be identified later. To avoid inconsistencies due to numerical integration, geometric integration is preferred, following the work in [23]: The exponential exp(ω T ∆T) can be computed via the Rodrigues' formula [25], which does not drift from SO (3): where Similar to the rotation, a change in translation of the origin of T in the reference frame S for each time interval: where b p is damping of linear translation and N is the number of intervals between two samples of data.
The translation between current and previous position can be defined as For convergence with a time constant τ, set where L is a characteristic length of the object, can be defined as where ζ T is the center of mass of the object and N mrk is the number of markers. For implementation purposes the psuedocode is presented in Algorithm 1: Algorithm 1: Damped least square algorithm for the filtering motion tracker's signal.

Experimental Verification of the Filtered Motion Capture Signal
Measuring Object's position: To test and validate the filter, eight optical markers are attached to a rigid object (see Figure 1). Four additional markers from the origin of the global coordinate frame S. The position of the markers is defined via the computer-aided design (CAD) drawing. The object is moved randomly involving translations and rotations with the motion capture system tracking at 50 Hz. All markers are visible during this experiment. To verify the filter algorithm, the position data of four markers is manually manipulated i.e., the markers 3rd can be considered invisible to the camera system for a continuous time of 20s. Figure 3 shows that the complementary filter can reconstruct the missing markers. The experimental test to verify the damped least square filter algorithm is presented in [20].

Complementary Filter on SE(3)
In the previous section, a model is developed to filter the motion capture data. The motion capture can detect large rotations and translations at low frequencies, but micro-motions may stay undetected when the amplitudes are too small and the frequencies too high. This section proposes an algorithm capable to capture both the micro-and macro-motions on SE(3).
Complementary filters traditionally apply when redundant measurements are available [26] and combining all available information in order to minimize the instrumentation error. They are also used to fuse information from multiple sensors and correct the output based on the sensed inputs. Although the Kalman filters can be extended (EKF) to nonlinear cases, they fail in capturing the nonlinear structure (e.g., rotations of a rigid body) which can be detected using complementary filters [13][14][15]27]. Therefore, the question arises how to properly measure micro-motions that are characterized by small amplitudes and high frequencies. IMUs can be used to measure orientations, accelerations, and the magnetic field.
Both the IMUs and motion capture can detect the orientation of an object in space. While motion capture system have high 3D resolution at reasonably high acquisition frequency they are outperformed by IMU in terms of temporal and technical resolution. Accelerometers and gyroscopes may have high temporal resolution but are known to show signal drift over time. A combination of both systems could improve the acquisition of slow major movements (3D motion capture) and micro-motions (IMUs) at the same time.

High-Frequency Motion Capture Imu Sensor Specifications
In this paper, a customized IMU described in the following is used together with motion capture to improve the acquisition of slow major movements (3D motion capture) and micro-motions simultaneously. The IMU sensor unit is mainly composed by (i) a 9-axis magneto-inertial sensor, the LSM9DS0 manufactured by STmicroelectronics, Inc., which embeds 3 axial accelerometers, gyroscopes and magnetometers. The LSM9DS0 is designed with a small outline (a 4 × 4 mm LGA24 package) and provides digital access to sensor data via I2C serial communication; (ii) a 16-bit microcontroller, the dsPIC33FJ64GP202 by Microchip, Inc., which collects the data from the sensor, and sends it to the Bluetooth module for wireless transmission. The microcontroller also processes the data with a quaternion-based complementary filter [28] in order to estimate the 3D orientation of the sensor, and transmits both raw data and computed quaternions with a sampling frequency of 200 Hz; (iii) a Bluetooth module, the SPBT2632C2A, by STmicroelectronics, Inc., guarantees the wireless communication from the sensor unit to a remote laptop via a Virtual COM Port (VCP); (iv) a rechargeable Li-Po battery with a capacity of 260 mAh, which assures the power supply with at least 4 h of continuous data streaming. The overall dimensions of the IMU sensor unit are 35 × 35 × 22 mm (see Figure 4).

Complementary Filter
Based on the complementary filter on SO(3) for dynamic attitude estimation presented by Campolo et al. [23], the measurement inputs are the rigid body orientation obtained from the motion capture system and the angular velocity measured by the IMU to estimate the kinematics of the rigid object on SE(3). An experimental validation is proposed to verify the filter algorithm will be presented in the next section. The following identities are used.

•
The rotation R MT and R I MU denotes the relative orientation of the instrumented box T with respect to the global frame S which is measured by using, respectively, the motion tracker and IMU. • Vectors p MT presents the translation of T frame in global frame S. • ω gyr denotes the angular velocity of the rigid object measured by the IMU. • a T denotes the output of the accelerometers embedded in the IMU and are expressed in moving frame coordinates. • f L v denotes the oscillation force which is caused by the plate of the loadcell.
Given a rigid object with a defined object frame T , the rotation R S MT measured by using the motion tracker presents the relative orientation of the rigid object frame T with respect to the global frame S. The gray block (i) in Figure 5 shows how raw data from motion captureR MT are geometrically preprocessed to generate a filtered output R MT to be fused with IMU angular velocity ω gyr . Such a sensor fusion takes place in the SO(3) complementary filter shown in the block (ii) in Figure 5. Specifically, the rotation R I MU represents the estimated orientation of the rigid object and is defined by the following time-continuous equations, where k is the only gain and, as in any complementary filter, it can be optimally tuned based on the noise characteristics of the signal [26].
So far, direct integration of angular velocities ω gyr led to a estimation of R I MU complementing the existing estimate R MT from motion capture. Data from accelerometers a T , however, cannot be directly integrated because they are expressed in a moving frame and, further more, contain a gravitational bias (In Figure 5, this is expressed by the additional gray-color term R −1 g S , where R is the actual (and unknown) 3D body orientation, estimated by R I MU . Because of estimation errors, R −1 R I MU is not exactly the identity matrix, this is expressed by the "tilde" variableg S which is supposed to be close to g S .). By means of the 3D orientation estimates R I MU , accelerations can be transformed in space-fixed coordinates, where the gravitational bias can be simply subtracted: as shown in the non-gray area in Figure 5. It should be noted that the gravitational bias is contact in space-fixed frame and equal to g S = [0 0 − 9.8] T m/s 2 . At this point, the space-fixed accelerations (without biased) a S could be directly integrated with respect to time to produce an estimate of linear velocities v S . A second integration, in cascade, would lead to an estimation of space-fixed positions p S . However, errors in the estimations of the gravitational bias (directly due to attitude errors in R I MU ) and noise would rapidly lead to a drift of such estimates, in typically rapid time. In order to limit the drift by fusing positions as estimated from the motion capture. The cascade of linear, complementary filter to derive velocity and position estimates is shown in block (iii) of Figure 5.

Experimental Verification
In this section, the ability of the proposed complementary filter to detect both micro-motions and macro-motions will be demonstrated experimentally. On the one hand, micro-motions entail small amplitude but high-frequency content which is typical of vibratory phenomena. To this end, a 6-degrees of freedom (dof) loadcell will be used as ground-truth to detect accelerations induced by the vibratory phenomena. On the other hand, macro-motions entail large amplitude but low frequency kinematics which is easily detected via motion trackers. We shall show that by fusing both motion-tracking information as well as accelerometers, a complementary filter is able to reliably detect both micro-and macro-motions.
Although the filter structure in Section 2 is general and application-independent, the specific hardware (motion tracker and IMU) along with specific sampling rates, reflects applications involving human motion. In particular, targeting human activities of daily living, the instrumented object was mounted on a leaf-spring and held by one of the authors while performing arbitrary, oscillatory hand-motions. Such motions are deemed representative of activities of daily leaving. The leaf spring mechanism was used to amplify natural and comfortable hand motions, at least in one direction.
Throughout the experimental condition the IMU, motion capture and loadcell recorded the movement (for more information and animations readers are referred to the supplementary data online). The instrumented object is mounted on a leaf spring allowing only rotations as shown in Figure 6. In this set-up, the instrumented object is attached to a leaf spring fixed to a base. The leaf spring is activated and oscillates until it gets back to its neutral position. Then, the whole setup is moved in space to allow both translations and rotations. In order to validate the complementary filter algorithm, the sampling frequency was set to: 200 Hz, 50 Hz, and 1000 Hz for the IMU, motion capture, and load cell, respectively.
The motion capture acquisition frequency was purposely set to 50 Hz to prove the concept of the complementary filter. The filter algorithm is able to detect the movements which the motion capture cannot such as the vibration with high frequencies and small amplitudes. Motion capture can detect large rotations and translations at lower frequencies but micro-motions may stay undetected when the amplitudes are too small and the frequencies too high. The measured force of the loadcell is the result of the acceleration and its moments of inertia. They are given by the hardware specifications and a multiplication with the estimated accelerations from the complementary filter will allow us to compare the estimated force with the measured force from the loadcell (see Equation (20)).

Equation (17) is time-continuous and any solution starting on SO(3) is guaranteed to remain on SO(3).
However, this would no longer be the case if (17) were to be directly discretized, e.g., by a numerical simulator, and the solution will likely drift away from SO(3).
Campolo et al. [23] proposed a technique that avoids this phenomenon and the following pseudocode in Algorithm 2 describes a possible numerical implementation, Algorithm 2: Complementary filter algorithm on SO(3).
Where f I MU = 200 Hz, f MT = 50 Hz is the sampling frequency of the IMU and motion capture, respectively; the notation · denotes the ceiling function maps a real number to the smallest following integer.
The complementary filter uses the motion capture data as reference for the IMU ( An estimated force from the complementary filter is defined as where L G T is an homogeneous transformation matrix computed from the CAD of the instrumented box. Therefore, the error between the measured force f L v by loadcell and the estimated force f L est by the complementary filter in the loadcell frame L can be determined as Figure 7 shows how well the estimated forces f L est match the force detected by the loadcell f L v . A linear regression of measured forces vs. estimated forces resulted in a goodness of fit with R-squared of 93%, 70%, and 96% for the x-, y-, and z-components, respectively. Bland-Altman plots [29] were computed to analyze the agreement between the estimated and the measured forces, as shown in Figure 8. Furthermore, the determination coefficient (r2) and the squared sum errors (SSE) were computed for each experiment and each force direction (Fx, Fy, Fz) to measure the discrepancy between the estimated and measured forces. The experimental results show strong coherence and small discrepancy between the estimated and measured forces. The direct comparison between the estimated and measured force shows quasi-linear relationship across all conditions. The complementary filter increases the accuracy of the estimated values of acceleration.   Finally, Figure 9 shows a comparison of z-axis position, velocity, and accelerations as detected by the motion tracker and as estimated by the SE(3) filter in Figure 5. While position estimates are almost indistinguishable (R-squared > 99% when MT positions data are linearly regressed with IMU position estimates), a decay in similarity can be seen when comparing velocities (R-squared < 65%) and even more when it comes to accelerations (R-squared < 10%).

Conclusions
In this paper, we have discussed the integration of two measurement systems to detect micro-and macro-motions. Multimodal sensor fusion is conducted using a complementary filter specialized to the nonlinear structure of the SE(3) group of rigid body rotations. This work aimed to improve the estimation of acceleration, free from gravitational bias using a complementary filter. Both the numerical implementation and the experimental validation have been provided as a verification and to complement the theoretical approach. The proposed filter was experimentally tested and validated based on the force measurements from a loadcell.
Optical and IMU technologies can be either used as standalone systems for specific applications but the strength lies in the their combination to deepen our understanding of movement assessment. Optical technology has traditionally been used to the laboratory setting; inertial technology takes movement analysis out to the field. As such, the combination of optical and inertial technologies will contribute to the performance assessments both in sports and industrial tasks.
Potential applications of the proposed approach include clinical movement analysis to identify micro-motions that are clinically relevant [30] as well as tooling task, i.e., capturing kinematics of instrumented tools during manual operations then transferring human skills data to robot [9,31,32].