WO2023234279A1

WO2023234279A1 - Motion analysis method and device

Info

Publication number: WO2023234279A1
Application number: PCT/JP2023/020000
Authority: WO
Inventors: 江山本; 泰暉石垣; 星喜金; 裕太嶋根
Original assignee: 国立大学法人東京大学
Priority date: 2022-05-30
Filing date: 2023-05-30
Publication date: 2023-12-07

Abstract

The present invention achieves motion analysis with respect to an object including a tool or orthosis that deforms flexibly. In a hybrid link system composed of a rigid body link system and at least one flexible link, the rigid body link system consisting of a plurality of rigid body links including a floating base link, the flexible link is approximated by a model based on a Cosserat theory to define a motion equation of the hybrid link system, and dynamics computation is performed on the basis of the motion equation.

Description

Motion analysis method and device

The present invention relates to a motion analysis method and device.

High-speed calculation algorithms for the kinematics and dynamics of articulated rigid link systems have been developed in the field of robotics. This has been applied to muscle tension analysis based on a human musculoskeletal model. In other words, by modeling the human skeletal system as a rigid link system and approximating muscles as wires that connect the skeleton, we can use motion data measured by motion capture to reconstruct three-dimensional movement of the skeletal system and the tension generated in the wires. Through the calculation, the muscle tension occurring during exercise can be estimated (Non-Patent Document 1).

Equipment made of flexible materials is used in many sports such as golf, fencing, and ice hockey, and it is known that the elastic deformation of the equipment is used skillfully during competition. Sports prosthetics made of leaf springs are widely used not only in sports for able-bodied people but also for sports for people with disabilities. In order to analyze skills in such sports situations, a musculoskeletal calculation model that takes into account the flexibly deforming structure and the deforming forces generated therein is required. Conventional motion analysis is applied to objects that can be approximated by rigid link systems, and does not take into account the deformation of tools used in sports or the deformation of braces worn by people with disabilities.

The finite element method, which is known for calculating deformation of flexible structures, is often used under static conditions, requires a large amount of calculation, and is not suitable for calculating deformation in real time. On the other hand, under the Cosserat Rod theory (Non-Patent Document 2) that efficiently handles rod shapes, Piecewise Constant Strain (PCS) models (Non-Patent Documents 3 and 4) have been proposed. For the PCS model, a recursive dynamic calculation algorithm similar to that for rigid link systems has been proposed, and it is also well compatible with rigid link systems. If it is possible to model the flexible deformation of tools and braces in more detail, this model can be incorporated into a rigid link system to create a hybrid system consisting of a rigid link system consisting of multiple rigid links and at least one flexible link. It becomes possible to perform link system kinematic analysis.
Y. Nakamura, K. Yamane, Y. Fujita, and I. Suzuki. Somatosensory computation for man-machine interface from motion-capture data and musculoskeletal human model. Trans. Rob., 21(1):58-66, Feb 2005 . S. S. Antman. Nonlinear Problems of Elasticity, Vol. 107. 2005. F. Renda and et. al. Discrete cosserat approach for multisection soft manipulator dynamics. IEEE Transactions on Robotics, Vol. 34, No. 6, pp. 1518-1533, 2018. F. Renda and et. al. A geometric and unified approach for modeling soft-rigid multi-body systems with lumped and distributed degrees of freedom. In 2018 IEEE International Conference on Robotics and Automation (ICRA), pp. 1567-1574. IEEE , 2018. V. Sonneville and et. al. Geometrically exact beam finite element formulation on the special euclidean group se (3). Computer Methods in Applied Mechanics and Engineering, Vol. 268, pp. 451-474, 2014. A. Goswami and et. al. Humanoid robotics: A reference. Springer, 2019. MWWalker et al., "Efficient dynamic computer simulation of robotic mechanisms," 1982. T. Sugihara et al., "Whole-body cooperative balancing of humanoid robot using COG Jacobian," in Proceedings of IEEE/RSJ International Conference on Intelligent Robots and Systems, 2002, pp. 2575-2580. K. Yamamoto, "Robust walking by resolved viscoelasticity control explicitly considering structure-variability of a humanoid," in 2017 IEEE International Conference on Robotics and Automation (ICRA). IEEE, 2017, pp. 3461-3468. H. Kaminaga et al., "Mechanism and control of whole-body electrohydrostatic actuator driven humanoid robot hydra," in International Symposium on Experimental Robotics. Springer, 2016, pp. 656-665. T. Ohashi, Y. Ikegami, K. Yamamoto, W. Takano and Y. Nakamura, Video Motion Capture from the Part Confidence Maps of Multi-Camera Images by Spatiotemporal Filtering Using the Human Skeletal Model, 2018 IEEE/RSJ International Conference on Intelligent Robots and Systems (IROS), Madrid, 2018, pp. 4226-4231.

The purpose of the present invention is to realize kinematic analysis of objects including flexibly deformable tools and braces.

The present invention relates to a motion analysis method,
In a hybrid link system consisting of a rigid link system consisting of a plurality of rigid links including a floating base link and at least one flexible link,
The flexible link is approximated by a model based on Cosserat theory,
The equation of motion of the hybrid link system is

Defined by
Dynamic calculations are performed based on the equation of motion.
here,
M is the inertia matrix;
b is bias vector;
q¨ is the generalized acceleration of the mixed link;
q ₀ is the generalized coordinate of the base link;
θ is the joint vector of the rigid link with n degrees of freedom;
q _S is the generalized coordinate of the flexible link approximated by a model based on Cosserat theory;
M ₀ is the inertia matrix associated with the base link;
M _R is the inertia matrix associated with the rigid link;
M _S is the inertia matrix associated with the flexible link;
M _0R is the inertia matrix due to the connection relationship between the base link and the rigid link;
M _0S is an inertia matrix due to the connection relationship between the base link and the flexible link;
M _RS is the inertia matrix due to the connection relationship between rigid links and flexible links;
b ₀ is the bias vector associated with the base link:
b _R is the bias vector associated with the rigid link;
b _S is the bias vector associated with the flexible link;
τ _R is the generalized force of the rigid link;
τ _s is the generalized force including the viscoelastic force of the flexible link;
J _C,i is the Jacobian matrix that maps the generalized velocity to the velocity of the contact point;
f _C,i is the contact force.

In one aspect, the rigid link system is a skeletal model.
In one aspect, the rigid link system is humanoid.
In one aspect, the flexible link is a prosthetic orthosis.
In one aspect, the flexible link is a tool made of a flexible material used in sports.

In one embodiment, the flexible link is a model based on Cosserat theory in which strain, which is the amount of deformation, is assumed to be constant regardless of position, or a model that is discretized into a plurality of segments, and in each segment, the strain is This is a model that assumes that it is constant regardless of position.

In one aspect, the dynamic calculation is an inverse dynamic calculation,
Measure the movement of the target consisting of a hybrid link system and obtain movement data,
The motion data includes first motion data that is motion data of a rigid link system and second motion data that is motion data of a flexible link,
Based on the first motion data, obtain generalized coordinates, generalized velocity, generalized acceleration of the base link, generalized coordinates, generalized velocity, and generalized acceleration of the rigid link system,
obtaining generalized coordinates, generalized velocity, and generalized acceleration of the flexible link based on the second motion data;
Generalized coordinates, generalized velocity, and generalized acceleration of the hybrid link system, a measured or estimated contact force, an inertia matrix, a bias vector, and a Jacobian matrix that maps the generalized velocity to the velocity of the contact point; Using the rigidity parameter and viscosity parameter of the flexible link, the generalized force of the rigid link and the generalized force of the flexible link are calculated.

In one embodiment, in the forward calculation of the inverse dynamics calculation, at the connection between the rigid link and the flexible link, the velocity twist η _i of the i-1th rigid link and the i-th flexible link are calculated using the following formula. A transformation is performed with the velocity twist η _i-1 of the segment of the flexible link.

here,
Velocity twist is a vector that combines translational velocity and rotational velocity,
G _i-1,i is a matrix representing coordinate transformation between the i-1th rigid link and the i-th flexible link or segment of the flexible link.
In one aspect, in the forward calculation of the inverse dynamics calculation, at the connection part between the rigid link and the flexible link, the acceleration twist η ^· _i of the i-1st rigid link and the i-th flexible link or the segment of the flexible link are calculated. A conversion of acceleration twist η ^· _i-1 is performed.
Here, the acceleration twist is a vector that combines translational acceleration and rotational acceleration.

In one aspect, in the backward direction calculation of the inverse dynamics calculation, at the connection between the rigid link and the flexible link, the following formula is used to represent the translational force and rotational moment of the link connected to the link of interest. A transformation is performed between the wrench vector F _i and the wrench vector F _i-1 representing the translational force and rotational moment of the link of interest.

here,
G _i-1,i is a matrix representing coordinate transformation between wrench vector F _i and wrench vector F _i-1 .

The present invention relates to a motion analysis device that analyzes the motion of an object,
The target is a hybrid link system consisting of a rigid link system consisting of a plurality of rigid links including a floating base link and at least one flexible link, the flexible link being approximated by a model based on Cosserat theory. has been
The device includes a processing section and a storage section,
The storage unit stores the equation of motion of the hybrid link system,

is stored,
The processing unit is configured to perform dynamic calculations based on the equation of motion.
here,
M is the inertia matrix;
b is bias vector;
q¨ is the generalized acceleration of the mixed link;
q ₀ is the generalized coordinate of the base link;
θ is the joint vector of the rigid link with n degrees of freedom;
q _S is the generalized coordinate of the flexible link approximated by a model based on Cosserat theory;
M ₀ is the inertia matrix associated with the base link;
M _R is the inertia matrix associated with the rigid link;
M _S is the inertia matrix associated with the flexible link;
M _0R is the inertia matrix due to the connection relationship between the base link and the rigid link;
M _0S is an inertia matrix due to the connection relationship between the base link and the flexible link;
M _RS is the inertia matrix due to the connection relationship between rigid links and flexible links;
b ₀ is the bias vector associated with the base link:
b _R is the bias vector associated with the rigid link;
b _S is the bias vector associated with the flexible link;
τ _R is the generalized force of the rigid link;
τ _s is the generalized force including the viscoelastic force of the flexible link;
J _C,i is the Jacobian matrix that maps the generalized velocity to the velocity of the contact point;
f _C,i is the contact force.

In one aspect, the rigid link system is a skeletal model.
In one aspect, the rigid link system is humanoid.
In one aspect, the flexible link is a prosthetic orthosis.
In one aspect, the flexible link is a tool made of a flexible material used in sports.
In one embodiment, the flexible link is a model based on Cosserat theory in which strain, which is the amount of deformation, is assumed to be constant regardless of position, or a model that is discretized into a plurality of segments, and in each segment, the strain is This is a model that assumes that it is constant regardless of position.

In one aspect, the processing unit includes a motion measurement means, an inverse kinematics means, and an inverse dynamics means,
The motion measuring means measures the motion of an object made of a hybrid link system to obtain motion data, and the motion data includes first motion data that is motion data of a rigid link system and motion data of a flexible link. Some second movement data is included,
The inverse kinematics means are:
a first inverse kinematics means for acquiring generalized coordinates, generalized velocity, generalized acceleration of the base link, generalized coordinates, generalized velocity, and generalized acceleration of the rigid link system based on the first motion data; ,
a second inverse kinematics means for acquiring generalized coordinates, generalized velocity, and generalized acceleration of the flexible link based on the second motion data,
The inverse dynamics means converts the generalized coordinates, generalized velocity, and generalized acceleration of the hybrid link system, the measured or estimated contact force, the inertia matrix, the bias vector, and the generalized velocity into the velocity of the contact point. The generalized force of the rigid link and the generalized force of the flexible link are calculated using the Jacobian matrix mapped to , and the rigidity parameter and viscosity parameter of the flexible link.

In one embodiment, the inverse kinematics means calculates the velocity twist η _i of the i-1st rigid link and the i-th It includes twist speed conversion means for performing a conversion of the speed twist η _i-1 of the flexible link or a segment of the flexible link.

here,
Velocity twist is a vector that combines translational velocity and rotational velocity,
G _i-1,i is an adjoint matrix representing coordinate transformation between the i-1th rigid link and the i-th flexible link or segment of the flexible link.
In one embodiment, the inverse dynamics means calculates the acceleration twist η ^· _i of the i-1st rigid link and the i-th flexible link or the flexible link at the connection portion between the rigid link and the flexible link in the forward direction calculation. It includes twist acceleration conversion means for converting the segment acceleration twist η ^· _i-1 .
Here, the acceleration twist is a vector that combines translational acceleration and rotational acceleration.

In one aspect, the inverse dynamics means calculates the translational force and rotational moment of the link connected to the link of interest at the connection between the rigid link and the flexible link using the following formula in the inverse direction calculation. and a wrench vector F _i-1 _representing the translational force and rotational moment of the link of interest.

here,
G _i-1,i is an adjoint matrix representing coordinate transformation between wrench vector F _i and wrench vector F _i-1 .

In order to perform a motion analysis of a hybrid link system consisting of a rigid link system consisting of a plurality of rigid links including a floating base link and at least one flexible link, the present invention includes:
The present invention relates to a computer program configured to cause a computer to function as the storage unit and the processing unit.

In the present invention, in a hybrid link system consisting of a rigid link system consisting of a plurality of rigid links including a floating base link and at least one flexible link, the flexible link is approximated by a model based on Cosserat theory. , it is possible to realize kinematic analysis of objects including flexibly deformable tools and braces.

The left figure is a schematic diagram of a hybrid link system consisting of a rigid link system and a flexible link, and the right figure is a schematic diagram of a humanoid equipped with a flexible prosthetic leg. FIG. 2 is a schematic diagram of a PCS model. It is a figure explaining the forward direction calculation and reverse direction calculation of the inverse dynamics calculation in the connection part of a rigid body link and a flexible body. (a) is a perspective view of a sports prosthesis, (b) shows visualization of the PCS model of the sports prosthesis through simulation, and (c) is a schematic diagram of the PCS model of the sports prosthesis. The left figure is a diagram explaining the kinematics and dynamics of the PCS model of a sports prosthesis, and the right figure shows the arrangement of optical markers installed on the sports prosthesis. The inverse kinematics calculation results for the prosthetic leg are shown. The left figure shows when the Group (b) load is changed from 2.5 to 30.00 kg, and the left figure shows when the Group (e-1) load is changed from 2.5 to 30.00 kg. Indicate the case. FIG. 1 is a block diagram of a motion analysis device according to the present embodiment. FIG. 2 is a block diagram showing a dynamics calculation unit of the motion analysis device according to the present embodiment. 3 is a flowchart showing a motion analysis method according to the present embodiment. FIG. 2 is a block diagram illustrating estimation of viscoelasticity of a sports prosthetic leg. FIG. 2 is a block diagram showing the flow of motion analysis of a person wearing a prosthetic leg. It is a diagram showing the structure and generalized coordinates of a hybrid link system. The definition of marker coordinates and the relationship between marker coordinates and link or segment coordinates are shown.

[1] Overall configuration of kinematic analysis method and device according to this embodiment The object of kinematic analysis according to this embodiment is a rigid link system consisting of a plurality of rigid links including a floating base link, at least one flexible link, This is a hybrid link system consisting of (Figure 1). Examples of the rigid link system include a skeletal model or a musculoskeletal model (right figure in FIG. 1) and a humanoid (left figure in FIG. 1). In human musculoskeletal models and humanoids, every body part has a base link that is not fixed to the environment. When the rigid link system is a skeletal model, the target is, for example, an athlete with a tool made of a flexible material, or an athlete wearing a sports brace (for example, a prosthetic leg) made of a leaf spring. . Examples of tools made of flexible materials include golf clubs, ice hockey sticks, and fencing swords. When the rigid link system is a humanoid, the object may be, for example, a leaf spring attached to the legs of the humanoid.

One of the features of this embodiment is high-speed flexible deformation calculation of flexible links (tools and braces) in a hybrid link system. Until now, the finite element method has been widely used to analyze structures that undergo flexible deformation. In general, the finite element method calculates minute deformations and stresses that occur in structures, and when dealing with deformations such as sports prosthetics and flexible deformation tools, it is difficult to use even if a high-performance computer is used. It takes several hours to several days. On the other hand, structures such as leaf spring prostheses and golf clubs can be approximated as beams or rods, and by applying the Cosserat rod theory studied in the field of materials mechanics, flexible deformation can be calculated quickly. can. This makes it possible to measure three-dimensional motion using kinematic calculations of the hybrid link system.

Sports prostheses made from golf clubs, ice hockey sticks, and leaf springs can be elastically deformed in three dimensions. This embodiment discloses a method for estimating the viscoelasticity of a tool or a prosthetic leg that undergoes three-dimensional elastic deformation based on modeling of a flexible deformable body that undergoes three-dimensional elastic deformation. This makes it possible to estimate with high accuracy the elastic forces that occur in tools and prosthetic legs during exercise.

This embodiment further discloses simultaneous estimation of the elastic force and muscle tension of a tool by dynamic calculation of a flexible-rigid mixed link system including floating links. This enables real-time motion analysis that takes into account the flexible deformation of tools and braces and musculoskeletal interaction forces. For example, by displaying the muscle tension estimation results in real time, when used for training, it is possible to realize muscle tension estimation that also takes into account the flexible deformation of prosthetic legs and tools. In particular, in the movement of a multi-link system such as human walking, it is common to model the base link as a floating link that is not fixed to the environment and calculate the contact with the environment. Even in link systems, calculations including floating links are now possible. The results can also be applied to proposals for forms that effectively utilize the elastic deformation of tools, as well as proposals for optimal shapes and material properties for tools and sports prosthetics tailored to individuals.

In order to analyze the motion of an object using inverse kinematics or inverse dynamics, it is necessary to obtain motion data of the object. Motion data of the object can be acquired by a motion capture system. The motion of the target is defined by time-series data of the pose of the target acquired by motion capture. It is only necessary to obtain time-series data of poses of a rigid link system (skeletal model) of the target (hybrid link system) and time-series data of poses of a flexible link (for example, a prosthetic leg). For example, as shown in Figure 5 and Table 1, a large number of optical markers are attached to predetermined positions on the prosthetic leg, and the pose of the prosthetic leg at time t is determined from a set of positions (coordinates) of the multiple optical markers at time t. The pose of the prosthetic leg at time t+1 can be determined from a set of positions (coordinates) of multiple optical markers at time t+1, and the pose at time t and the pose at time t+1 can be determined. The difference represents the deformation of the prosthesis.

The posture of the target rigid link system is specified by multiple feature points (typically joints) on the target body (skeletal model), and the three-dimensional coordinate values of the multiple feature points are acquired in each frame. Then, the operation of the target rigid link system is defined from time-series data of three-dimensional coordinate values of a plurality of feature points. For the parts of the target skeletal model (rigid link system), time-series data of joint positions and joint angles are acquired by motion capture, and the velocity and acceleration of the rigid links are calculated by inverse kinematics.

For tools that flexibly deform (flexible links, in this embodiment, the PCS model), time series data of the amount of expansion and contraction of the beam and the curvature of bending are acquired. A flexible link has a total of 6 dimensions of deformation (strain), with the amount of expansion/contraction and curvature each being 3 dimensions. For example, a PCS model in which a leaf spring prosthesis is divided into 6 segments has a deformation amount of 6x6 = 36. However, for example, in an actual prosthetic leg, the amount of expansion and contraction can be considered to be very small compared to the displacement of curvature, so during motion capture, the approximation is made so that only the three degrees of freedom of curvature are handled. If the strain of the PCS model can be obtained from motion capture, it is possible to calculate the amount equivalent to the speed of deformation from it. From there, the viscoelastic force of the PCS model can be estimated using equation (4) in [3][E]. At this time, the stiffness matrix K and viscosity matrix D must be estimated in advance.

As shown in FIG. 7, the motion analysis device according to this embodiment includes a processing section and a storage section, and the processing section and storage section are composed of one or more computers. A computer includes an input section, a processing section (processor), a storage section (memory including RAM and ROM), and an output section. The motion analysis device may include a display unit. The motion analysis device analyzes the motion of the target based on the motion data of the target acquired by the motion capture system.

The storage unit stores hybrid link system models (skeletal model or musculoskeletal model, and PCS model that defines flexible links), and motion capture is executed corresponding to each model, and the motion capture measurement Data is stored. The processing unit executes kinematic calculation based on the measurement data of the motion capture, and the calculation data of the inverse kinematic calculation based on the measurement data is stored in the storage unit.

The storage unit stores the equation of motion of the hybrid link system, and the processing unit executes dynamic calculations (inverse dynamic calculations and forward dynamic calculations) according to the equation of motion, and stores the calculated data. . The storage unit stores stiffness parameters and viscosity parameters of the flexible link necessary for dynamic calculation based on the equation of motion in the processing unit.

The storage unit stores the contact force when the object moves. The contact force is measurement data (ground reaction force) acquired by, for example, a force plate during the movement of the object. Alternatively, it may be data (ground reaction force) estimated by the processing unit during the movement of the object. The storage unit stores computer programs for the processing unit to execute various predetermined calculations. Specifically, the program includes a program to perform motion capture and obtain measurement data, a program to perform inverse kinematic calculations based on measurement data, and a program to perform dynamic calculations based on processed data and stored data. (inverse dynamics calculation, forward dynamics calculation), program to estimate stiffness parameters and viscosity parameters of flexible links, and to perform contact force estimation calculation when estimating contact force. programs, etc. are stored.

The processing section of the motion analysis device is configured to execute predetermined processing according to a predetermined program, and the processing results are stored in the storage section. In this embodiment, the processing section typically includes an inverse kinematics calculation section, an inverse dynamics calculation section, and a forward dynamics calculation section. It should be noted that the processing section in FIG. 7 clearly shows the inverse kinematics calculation section and the inverse dynamics calculation section, while FIG. 8 clearly shows the inverse dynamics calculation section and the forward dynamics calculation section. Further, the processing section may include a muscle tension calculation section and a muscle activity calculation section.

As shown in Fig. 8, the inverse dynamics calculation unit calculates the inertia matrix based on the equation of motion (22) of the hybrid link system, according to the input of generalized coordinates q, generalized velocity q ^, and generalized acceleration q. Generalized forces τ _R and τ _S are output using M, bias vector b, contact force, and contact Jacobian matrix. The generalized force (joint torque) τ _R is used to estimate muscle tension in the rigid link system (musculoskeletal model) of the hybrid link system.

As shown in FIG. 8, the forward dynamics calculation unit calculates the inertia matrix M, the bias vector b, and the contact force according to the input of the generalized forces τ _R and τ _S based on the equation of motion (22) of the hybrid link system. , using the contact Jacobian matrix, it is possible to output the generalized acceleration q¨, and also calculate the generalized velocity q ^· and the generalized coordinate q. The generalized acceleration q¨, the generalized velocity q ^· , and the generalized coordinate q are used, for example, to control the walking motion of a humanoid.

The type of motion capture method used in this embodiment is not limited, and includes an optical motion capture method that uses optical markers to identify feature points, and an optical motion capture method that uses so-called inertial sensors such as acceleration sensors, gyroscopes, and geomagnetic sensors on the target body. Examples include a method in which the sensor is attached to the camera to acquire motion data of the target, and a so-called markerless motion capture method in which no optical marker or sensor is attached. The markerless motion capture method is advantageous from the viewpoint of not interfering with the natural movement of the subject. Examples of the markerless motion capture method include a motion capture method using a system equipped with a camera and a depth sensor, or a video motion capture method that acquires motion data by analyzing RGB images (Non-Patent Document 11). I can do it. In this embodiment, one type of motion capture method is used, but different motion capture methods may be combined. For example, video motion capture may be employed for the skeletal model, optical motion capture using a plurality of optical markers may be employed for the flexible link, and measurement data may be integrated and used. In one aspect, a motion capture system includes one or more cameras that acquire video data of a subject's motion, and one or more computers that acquire time-series data representing the subject's motion based on the video data. configured. This computer may also be used as a computer for motion analysis.

Contact force is required for dynamic calculations. In one aspect, a force plate is used to obtain ground reaction force during motion measurement. Alternatively, by detecting the contact state between the foot and the floor based on the measurement results of motion capture (for example, from marker information on the foot), that information can be used to simultaneously estimate the ground reaction force during inverse dynamics calculation. You may.

FIG. 9 shows an example of a processing process for motion analysis of an object (hybrid link system) using motion capture. The motion of a musculoskeletal model (target of a hybrid link system), including flexibly deformable tools or braces, is measured using motion capture. The motion of the target skeleton (rigid link system) is acquired by motion capture according to this embodiment. The time series data of the acquired joint angles and joint positions is acquired, and the position, angle, velocity, angular velocity, acceleration, and angular acceleration of each rigid link are acquired by inverse kinematics calculations and input into the inverse dynamics engine. At the same time, the motion data of the target tool or brace (flexible link) is acquired. Based on the motion data, the viscoelastic force of the flexible link is obtained by inverse kinematics calculation and input into the inverse dynamics engine. The inertia tensor, centrifugal force, Coriolis force, and gravity vector are also input to the inverse dynamics engine. As an external force (contact force), for example, a floor reaction force acquired by a force plate is input to the inverse dynamics engine.

In a rigid link system, more specifically, for example, by interpolating the time-series inter-frame displacement of all degrees of freedom of the skeleton obtained by motion capture using a continuous function, the displacement of all degrees of freedom of the skeleton at each frame time is calculated. , calculate its time derivative, velocity, and its time derivative, acceleration. The position, angle, velocity, angular velocity, acceleration, and angular acceleration of each link calculated from these are sent to the inverse dynamics engine, which calculates the mechanical information associated with the movement of the skeleton assuming mass, thereby matching the movement. Calculate the joint torque. Each segment of the skeleton is a rigid body, and its mass, center of gravity position, and inertia tensor can be estimated from statistical measurement information of each part of the person using physique information. Alternatively, these parameters can also be estimated by identification from motion information of the object. The physique information of the target used for estimation is acquired in advance.

Once you know the viscoelastic force of the flexible link (PCS model), you can know the force that acts on the rigid link to which it is connected from the PCS model side (for example, the force that acts on the attachment site from a leaf spring prosthesis), so you can create a form that includes this force. You can perform inverse dynamics calculations using Although the stiffness matrix K and the viscosity matrix D are assumed to be obtained before motion capture is executed, the viscoelastic force can be calculated in real time during motion capture. In the inverse dynamics calculation, the value of the input "viscoelastic force" is used as a reference value, and the "generalized force of the flexible deformation element" is calculated again in consideration of the acceleration term in the inverse dynamics calculation. Muscle tension is estimated by estimating the viscoelastic force that occurs when the prosthetic leg deforms and combining it with inverse dynamics calculations.

Contact force is required for dynamic calculations. In one embodiment, a force plate is used to obtain ground reaction force during motion measurement. Alternatively, by detecting the contact state between the foot and the floor based on the measurement results of motion capture (for example, from marker information on the foot), that information can be used to simultaneously estimate the ground reaction force during inverse dynamics calculation. You may.

Joint torque is obtained by inverse dynamics calculation, and using the joint torque, wire tension in a musculoskeletal model including wires imitating muscles is obtained by optimization calculation (quadratic programming or linear programming). For example, the tension of a wire distributed over a whole body modeled on muscles is calculated using bias/weighted quadratic programming. Regarding the calculation of this wire tension, reference can be made to Non-Patent Document 1. By obtaining the measured values of the force distribution when antagonist muscles are used according to the classified movements, and using the bias and weights based on those values, we obtain a solution that better approximates the actual muscle tension. be able to. In acquiring muscle tension, measurement data from an electromyograph may be taken into consideration. When optimizing muscle tension in a hybrid link system, the same evaluation function as in the optimization calculation in a conventional rigid link system can be used. However, since the types of muscles considered in the evaluation function may differ depending on the defective region, it is necessary to omit the corresponding muscle term from the evaluation function depending on the defective region of the subject.

The obtained muscle tension is divided by the assumed maximum muscle tension of the muscle, and the muscle activity level is determined, and an image of the whole body musculoskeletal system is generated and visualized with the muscle color changed according to the muscle activity level. Musculoskeletal images with muscle activity are output at a predetermined frame rate and displayed on a display as a moving image. Furthermore, changes in the values of each variable (for example, joint angle, velocity, muscle tension, ground reaction force, center of gravity position, etc.) are graphed and output. These outputs are presented as analysis results in the form of images and graphs, and are used as a record of muscle and body activities during exercise, or the movements of each part of the body. In this way, it is possible to automatically and efficiently perform everything from photographing the movement of the target, obtaining the three-dimensional pose of the target during movement, and estimating and visualizing the muscle activity necessary for the movement.

Kinematic calculations and dynamic calculations of rigid link systems are well known to those skilled in the art, so detailed explanations will be omitted. It will be understood by those skilled in the art that knowledge of system kinematics calculations and dynamical calculations can be used. Also, in the following explanation, please note that the formulas are numbered independently for each chapter.

[2] Dynamic calculation and control of hybrid link system We propose a dynamic model of a hybrid link system that integrates a rigid body and a soft body. In order to apply it to humans and humanoids, we consider a floating link system and realize its forward and inverse dynamics calculations. In addition, we derive a centroid Jacobian matrix that takes into account the entire hybrid link system. Assuming a humanoid wearing a sports prosthetic leg, we will verify the dynamic calculations of the proposed hybrid link system through dynamic simulation of balance control. By modeling flexible deformation in more detail using the PCS model, we can expect application to the kinematic analysis of people with disabilities who wear sports prosthetics, as well as sports movements that utilize flexible deformation of tools, such as golf swings. In addition, if it becomes possible to make contact with the environment and calculate the center-of-gravity Jacobian matrix, it will become possible to calculate and control the dynamics of humanoids with flexible structures.

[A] PCS model (Piecewise Constant Strain Model)
[A-1] Kinematics of the PCS model Figure 2 shows a schematic diagram of the PCS model (Non-Patent Document 2). The PCS model discretizes the continuous Cosserat model (Non-Patent Document 2) into N segments, and in each segment, the 6-dimensional strain ξi∈R ⁶ (i = 1, 2,...,N), which is the amount of deformation, is constant. This is a model that assumes that Let the coordinate system of the central axis of the continuum be s∈R, and define the i-th segment to be between s=L _i-1 and s=L _i of the continuum. The constituent curve G(s) of the continuous Cosserat model at a certain time is expressed by the following equation.

However, R∈SO(3) is a rotation matrix and p∈R ³ is a position vector. Here, the time differential and spatial differential of G(s) are defined as follows.

The deformation of the continuous Cosserat model is defined as the infinitesimal displacement and strain ξ(s) of the constituent curves as shown in the following equation.

In this case, u∈R ³ and k∈R ³ represent translational and rotational distortions, respectively. Also, [・×] is a matrix representation of the vector representation of the Lie algebras se(3) and so(3).

The fact that the strain is constant within each segment is expressed as follows.

Based on this assumption and from equation (3), we obtain the following equation.

At this time, Gi(s) is defined as follows using an exponential mapping.

Here, R _i is a rotation matrix and p _i is a position vector. By using equation (7), the constituent curve frame G(s) in the i-th segment (s∈[L _i-1 ,L _i ]) is

It can be calculated as follows.

In this specification, E and O are an identity matrix and a zero matrix, respectively. The adjoint representation of the homogeneous transformation matrix G is

is defined as The adjoint expression of ξ∈se(3) is

is defined as

[A-2] The time change of the differential kinematics constitutive curve of the PCS model is given by the velocity twist η of the following equation.

Here, v∈R ³ and w∈R ³ represent translational speed and rotational speed, respectively. The velocity twist η(s) in the i segment can be calculated using η( _Li-1 ) as shown in the following equation.

Here, ^Ad G is an adjoint expression of SE(3), and T _i (s) is defined as a tangent operator of an exponential map as shown in the following equation (Non-Patent Document 5).

Further, the acceleration twist η ^· (s) can be calculated as shown in the following equation.

[A-3] Jacobian matrix of PCS model
The distortion of N segments is defined as the generalized coordinate q∈R ^6N of the PCS model as shown in the following equation.

If we rewrite equation (14) as a relational expression between speed twist η(s) and generalized speed ^q , we get

becomes. Here, J(s)=[J ₁ . . . J _N ]∈R ^6×6N is the Jacobian matrix of the PCS model (segment), and J _j ∈R ^6×6 can be calculated as follows.

Assuming that the base of the continuum is fixed, η(0) = 0, and by using the recursion in equation (12), η(s) is a polynomial of the strain ξ _j of each segment, as follows: It is expressed as follows.

here,

It is.

[A-4] Dynamics of PCS model
The equation of motion of the PCS model is expressed as the following equation (Non-Patent Document 2).

here,
M∈R ^6N×6N is the inertia matrix;
b∈R ^6N is a bias vector including Coriolis force and gravity;
τ _S (q ^· ,q)=[τ ^T _S,1 ...τ ^T _S,N ]∈R ^6N represents the internal force due to the stiffness and viscosity of the continuum.

Considering the influence of external force, equation (18) becomes as follows.

here,
M∈R ^6N×6N is the inertia matrix;
b∈R ^6N is a bias vector including Coriolis force and gravity;
τ _int (q, q^・)=[τ ^T _int,1 ...τ ^T _int,N ] ^T ∈R ^6N is the internal force due to the viscoelasticity of the continuum;
τ _ext = [τ ^T _ext,1 ...τ ^T _ext,N ] ^T ∈R ^6N represents external force;

In equation (18), the internal force τ _S is

It is expressed as Here, q ^eq is an equivalent point of the generalized coordinate vector, and K∈R ^6N×6N and D∈R ^6N×6N are a stiffness matrix and a viscosity matrix, respectively.

The equation of motion of the PCS model in Equation (18) has a form similar to the equation of motion of a rigid link system, and the recursive algorithm used for rigid link systems can be applied to the dynamics calculation of the PCS model. The recursive Newton-Euler method can be applied to inverse dynamics calculations in the PCS model (Non-Patent Documents 3 and 4). In the inverse dynamics calculation, the generalized force τ is obtained from the generalized coordinate q, its velocity q ^· , and acceleration q¨. The recursive algorithm uses the following equation as the equation of motion in the i-segment.

Here, A∈R ^6×6 is the inertia matrix of the strain acceleration ^... ξ _i , B∈R ^6×6 is the inertia matrix of the acceleration twist η^・(L _i-1 ), and c∈R ⁶ is the Coriolis force. and the external force, τ _int,i is the internal force of the i segment, and F _i+1 is the transmitted force from the i+1 segment received at the cross section of s=L _i .

The steps of the recursive algorithm for inverse dynamics calculation are as follows.
(a) At the boundary of each segment from the forward calculation base (s=0) to the tip (s=L _N ), the constituent frame G(s) from (8), the velocity twist η, (12') from (12) ), calculate the acceleration twist η^・.
(b) Inverse calculation generalized force (18”) calculates τ _int,i from the tip to the base segment.
In the forward dynamics calculation, the acceleration q¨ of the generalized coordinate is determined from the generalized coordinate q, its velocity q ^· , and the generalized force τ. Similar to the rigid link system, the inertia matrix M and bias vector b are determined by the unit vector method (non-patent document 7), and (18) is solved.

[B] Hybrid link system with floating links [B-1] Dynamics of hybrid link system Consider a hybrid link system as shown in Figure 1. In order to apply it to humanoid and human skeletal models, we assume a floating link system in which the base link is not fixed. In floating link systems, the base link is not fixed to the environment and is actuated indirectly by contact forces.

The equation of motion of the hybrid link system is expressed as follows.

M is the inertia matrix;
b is bias vector;
q¨ is the generalized acceleration of the hybrid link;
q ₀ ∈R ⁶ is the 6-dimensional generalized coordinate of the base link;
θ∈R ⁿ is the joint vector of the rigid link with n degrees of freedom;
q _S ∈R ^6N is the generalized coordinate of the flexible link approximated by the PCS model;
M ₀ is the inertia matrix associated with the base link;
M _R is the inertia matrix associated with the rigid link;
M _S is the inertia matrix associated with the flexible link;
M _0R is the inertia matrix due to the connection relationship between the base link and the rigid link;
M _0S is an inertia matrix due to the connection relationship between the base link and the flexible link;
M _RS is the inertia matrix due to the connection relationship between rigid links and flexible links;
b ₀ is the bias vector associated with the base link:
b _R is the bias vector associated with the rigid link;
b _S is the bias vector associated with the flexible link;
τ _R is the generalized force of the rigid link;
τ _s is the generalized force including the viscoelastic force of the flexible link;
J _C,i ∈R3 ^×(6+n+6N) is the Jacobian matrix that maps the generalized velocity to the velocity of the contact point;
f _C,i ∈R ³ is the contact force.

The dynamic calculation of equation (20) requires an inertia matrix M and a bias vector b. Several methods can be used to calculate the inertia matrix M and bias vector b. It is possible to calculate the Jacobian matrix J _i and its time derivative for each link or segment of the hybrid link system based on the hybrid link structure, and use them to calculate the matrix M and vector b. For specific calculations, for example, refer to Non-Patent Document 3. Non-Patent Document 3 discloses that the result of integrating the density of a PCS model (flexible link) is used in calculating the inertia matrix M and the bias vector b. In this embodiment, a recursive algorithm for inverse dynamics calculation is considered in a hybrid link system, and this is used in a unit vector method to calculate an inertia matrix M and a bias vector b. In the calculation of the recursive algorithm of inverse dynamics, the result of integrating the density of the PCS model (flexible link) is used. When using the unit vector method, there is no need to calculate the Jacobian matrix or its time derivative; if the connection relationship between links (or segments) is known, the inertia matrix M can be calculated by repeating the inverse dynamics calculation. and bias vector b can be calculated.

[B-2] Inverse dynamics calculation As described above, the inverse dynamics calculation of the PCS model is performed using the recursive formula. We implement inverse dynamics calculations for hybrid link systems while considering appropriate transformations at the connections between rigid links and flexible links (Table 1).

[B-2-1] Forward calculation As shown in FIG. 3, consider the connection between the i-1th rigid link and the i-th flexible link segment (ie, base). Here, the subscript i is a link ID that does not distinguish between rigid links and PCS model segments. The forward calculation of the inverse dynamics calculation requires conversion from the twisting speed of rigid link i-1 to the twisting speed of base i of the PCS model. In the forward direction calculation, the transformation between the velocity twist η _i-1 of the rigid link and the base velocity twist η _i of the PCS model is obtained as follows.

Here, G _i-1,i ∈SE(3) represents the coordinate transformation between the rigid link i-1 and the base (segment) i of the PCS model, as shown in FIG. Furthermore, the conversion of acceleration twist η ^· _i−1 and η ^· _i can be obtained by first-order differentiating equation (21) with respect to time.

[B-2-1] Reverse calculation In the backward calculation, it is necessary to transform the wrench vectors F _i and F _i-1 , which is obtained as shown in the following equation.

Here, F _i and F _i-1 are six-dimensional wrench vectors (consisting of a three-dimensional translational force vector and a three-dimensional rotational moment vector).

A transformation matrix used for velocity twist transformation and wrench vector transformation in inverse kinematics calculations is obtained in advance. For example, when the flexible link is a prosthetic leg, a transformation matrix can be obtained based on the dimensions of a joint that fixes a leaf spring and geometric information about the attachment position of the leaf spring. Alternatively, for example, when initializing optical motion capture measurements, for example, have the subject assume a T-shaped pose and measure while standing still, and obtain a transformation matrix from the position information of the multiple optical markers obtained. You may.

[B-3] Forward dynamics calculation of forward dynamics calculation formula (20) is,
1) Calculation of inertia matrix M and bias vector b: As with the rigid link system, inertia matrix M and bias vector b are calculated by applying the unit vector method in inverse dynamics calculation.
2) Calculation of contact Jacobian matrix J _C,i : Calculate contact Jacobian matrix J _C,i using the unit vector method. M and b were calculated in the previous step. The column vector of the Jacobian matrix is calculated by setting the unit vector to be the contact force f _C,i , and the generalized force is calculated by inverse dynamics calculation.
3) Modeling of contact force f _C,i ∈R ³ : Assuming that the contact between the robot and the environment is a non-perfect elastic collision, the contact force f _C,i is calculated (Non-Patent Document 6).
4) Forward dynamics calculation: Calculate the generalized acceleration of the hybrid link system using the calculations in the three steps above. From equation (20), the generalized acceleration is obtained as follows.

[C] Centroid Jacobian Matrix (COG Jacobian Matrix)
[C-1] Centroid Jacobian Matrix of Hybrid Link System In the control of humanoids, the center of gravity is an important feature. When controlling a humanoid (bipedal walking robot) to maintain balance to prevent it from falling over, it is necessary to control the center of gravity to keep it as constant as possible. For example, when calculating the kinematics and dynamics of a hybrid link of a bipedal walking robot equipped with a leaf spring, the center-of-gravity Jacobian matrix of the hybrid link system is used for control to maintain such balance. Note that the use of the center-of-gravity Jacobian matrix is not limited to humanoids, and can be used for human body motion analysis, for example, when optimally calculating the shape of a leaf spring that facilitates balance.

This section describes the centroid Jacobian matrix of a hybrid link system. Let m _i and p _G,i ∈R ³ be the mass and center of gravity of each link, respectively. Here, as in the previous section, the subscript i does not distinguish between rigid links and PCS models. The center of gravity velocity p ^· _G ∈R ³ of the entire hybrid link system can be calculated as follows.

Here, J _G,i is the centroid Jacobian matrix of the i-th link; p ^· _G,i = J _G,i q ^· ; M is the total mass (M = Σ _i m _i );

The centroid Jacobian matrix J _G of the entire hybrid link system can be calculated as shown below.

In a local coordinate system in a link, the center of gravity of a rigid link is constant, and a method for calculating the center of gravity Jacobian matrix has already been established (Non-Patent Document 8). On the other hand, the center of gravity of a PCS segment moves as the volume or shape of the segment changes. Therefore, below, we will discuss the centroid Jacobian matrix of the PCS model.

[C-2] Centroid Jacobian matrix of PCS model
The center of gravity velocity of the i-segment of the PCS model can be calculated as follows.

Here, ρ(s) and A(s) are the density and cross-sectional area of the i-th segment of the PCS model. Assuming that ρ(s) and A(s) are constant in the segment, equation (26) can be simplified as follows.

Considering the time derivative of equation (27), the following relationship is established between p^・(s) and the translational velocity v of the twist vector.

Using this fact, the translational velocity v can be calculated from the translational velocity part of equation (12) as follows.

In equation (29), η(L _i-1 ) and ξ ^· i are as follows.

Here, S _S,i ∈R ^6×6 is a selection matrix that selects the i-th PCS segment ξ _i ^from the generalized velocity q.

From equations (27) to (31), the center of gravity velocity p ^· _G,i and the center of gravity Jacobian matrix can be calculated as follows.

Here, R(L _i-1 )=R(s)R ^T _i , s∈[L _i-1 ,L _i ], and equation (15) was used.

[D] Compliance optimization of a humanoid in a hybrid link system [D-1] Compliance distribution control We will examine the standing stabilization control of a humanoid modeled as a hybrid link system. Compliance optimization control (Non-Patent Document 9) is applied to a hybrid link system. Here, for simplicity, we assume that a flexible link such as a prosthesis generates active compliance by controlling τ _s . Note that this assumption does not apply in real cases where the prosthesis has passive compliance. The purpose here is to verify the forward dynamics calculation of a hybrid link system. Therefore, assuming that τ _s can be actively controlled, the driving torque τ of the hybrid link system is given as follows.

Here, consider obtaining the center of gravity by the following joint compliance control.

Here, K _q and D _q ∈R ^{(n+6N)×(n+6N)} are an elastic matrix and a viscous matrix, respectively. Design COG compliance K _G , D _G ∈R ^3×3 so that stance balance is maximized. Based on the principle of virtual work, the relationship between K _q and K _G is expressed as follows.

Here, J _G ∈ ^3×(n+6N) is a Jacobian matrix representing the relationship between COG and the generalization speed related to τ. The general solution for K _q that satisfies equation (36) is given as follows.

Here, Y is any positive definite matrix and can be regarded as the target value of the inverse matrix of the stiffness matrix. J ^♯ _G is a pseudo inverse matrix of J _G , and the viscosity D _q can be calculated in the same way.

[D-2] Stance balance control when supporting one leg In order to verify the dynamic calculations of the hybrid link system, we simulated the dynamics of balance control using compliance optimization control assuming a humanoid wearing a sports prosthesis. Here, we will consider the standing stabilization control during the single-leg support phase of the left leg in a model (right figure in FIG. 1) in which the legs of the hydrostatically driven humanoid Hydra (Non-Patent Document 10) are replaced with prosthetic legs (PCS model). In the simulation, the target value of joint compliance was set as follows.

Here, Y _k and Y _d are target values of the stiffness matrix and the viscosity matrix, respectively. k _d and d _d are positive constants.
The target viscosity matrix was set as d _d =0.2, and the target elastic matrix was verified using a positive constant k _d by changing the value of k _d to 100, 200, and 400. At k _d = 100, he lost his balance after about 8 seconds, but at k _d = 200 and 400, he was able to maintain his standing position for 12 seconds. It was also confirmed that the larger the value of k _d , the faster the center of gravity converges.

[3] Calculation of flexible deformation of a competition prosthesis based on the model and estimation of stiffness and viscosity parameters
[A] Overview/Background Because there is a strong interaction between the mechanical properties of a prosthetic leg and body motion, prediction of body motion is important in order to design a prosthetic leg with appropriate functionality. Research is being conducted on motion prediction using simulation models that take into account human whole body motion. However, many methods do not take into account the flexible deformation of the prosthetic leg, such as those that treat the entire lower limb as a simple spring mass point model for reconstructing the dynamics and motion of the prosthetic leg, or those that model the prosthetic leg using a rigid link system. . Furthermore, motion is often limited to two dimensions. On the other hand, in the field of soft robotics, a Piece-wise Constant Strain (PCS) model (Non-Patent Documents 3 and 4) has been proposed to calculate the deformation of flexible beams and rods. This model has high extensibility and compatibility with dynamic models of multi-joint rigid link systems, and allows calculation of flexible deformation at relatively low cost. Here, for the purpose of kinematic analysis of prosthetic runners including competitive prosthetic legs, we will consider applying the PCS model to the calculation of flexible deformation of competitive prosthetic legs. First, the shape of the competition prosthesis is modeled by determining the strain parameters of the PCS model for each segment based on the actual prosthesis dimensions and image data. In addition, the amount of deformation when a load is applied to the prosthetic leg is measured using optical motion capture, and based on the measurement data, inverse kinematics and statics are calculated, and numerical optimization is performed to improve the stiffness parameters of the prosthetic leg for competitions. Estimate the viscosity parameters. Please note that in this chapter, formula numbers are assigned independently.

[B] Flexible deformation model of prosthetic leg
[B-1] Piece wise constant strain (PCS) model Figure 4 (a) shows the competition prosthesis (Sprinter 1E90, Ottobock, Germany) that is the subject of this study. This flexible deformation is modeled using a PCS model (non-patent document 3) as shown in FIG. 4(b). Figure 4(c) shows an outline of the PCS model. The PCS model is based on the discrete Cosserat model (Non-Patent Document 2), and calculates the flexible deformation of a beam or rod by dividing it into a finite number of segments. The constitutive curve of the continuous Cosserat model is expressed as follows.

Here, s is the central axis coordinate of the continuum, and R∈SO(3) and p∈R ³ represent a rotation matrix and a position vector, respectively.

In the PCS model, the displacement of the constituent curve G(s) caused by flexible deformation is defined as a six-dimensional strain vector ξ(s) as shown in the following equation.

However, k and u represent strain in the rotational direction and translational direction, respectively. In the PCS model, a complex beam structure such as an athletic prosthesis is divided into a finite number of sections, and the strain ξ is assumed to be constant in each segment. That is, the constant strain ξ _i of a certain segment i is defined as follows.

However, L _i indicates the value of the center coordinate s of segment i. Using this strain vector, the generalized coordinates q∈R ^6N of the PCS model can be defined as follows.

However, N is the number of segment divisions.

The constituent curve G(s) in each segment can be calculated from equation (2) as follows.

However, G _i (s) can be calculated as follows using an exponential mapping.

The time change of the constituent curve G(s) is defined as a twist vector η(s) as shown in the following equation.

However, w and v represent rotational and translational speeds, respectively. For each segment i, the relationship between velocity twist η(s) and generalized velocity ^q can be obtained as follows.

However, T _i (s) is an exponential map shown below,

Here, ^Ad G is an adjoint expression of SE(3), ξ is an adjoint expression of se(3), and is defined as follows.

Here, J(s)=[J ₁ ,J ₂ ,...,J _N ]∈R ^6×6N is the Jacobian matrix of the PCS model and is calculated as follows.

In this way, similar to the kinematics of an articulated rigid link system, the generalized velocity ^q can be mapped to the twist vector η using the Jacobian matrix.

[B-2] Modeling the prosthetic leg using the PCS model Calculate the strain parameters of the PCS model that reproduces the shape of the prosthetic leg for competition. First, in segment i, by substituting s=L _i and ξ=ξ _i into equation (6) and transforming it, the following equation is obtained.

however,

The following relationship was used. That is, if the homogeneous transformation matrices G(L _i-1 ) and G(L _i ) of the end faces of each segment are known, the strain can be calculated. In this embodiment, a central curve on a two-dimensional plane corresponding to the side surface is extracted from a CAD model of a competition prosthetic leg, the shape is divided into six segments that can be approximated by a constant curvature, and each segment is I asked for G _Li . Figure 4(b) shows the shape of a competition prosthetic leg drawn as a PCS model.

[C] Motion capture measurement of sports prosthesis
[C-1] Experiment to measure the deformation of a prosthetic leg In order to estimate the viscoelasticity of a competitive prosthetic leg, an experiment was conducted to measure the strain Δqj when a load F ^j ∈R ⁶ was applied to the prosthetic leg. First, they designed a jig to apply force and moment to the upper end of the prosthesis. In the experiment, a load was applied to the prosthesis, and the amount of deformation of the prosthesis and the force applied to the end points were measured using optical motion capture and a ground reaction force meter. In the experiment, 22 optical markers using retroreflective materials were installed. The right figure in FIG. 5 shows the arrangement of markers used in the experiment. In the experiment, loads were applied to the prosthetic leg in 5 to 10 steps in eight different directions, including the translational and rotational directions, and the static deformation states were measured. The coordinates of a plurality of optical markers provided at predetermined parts of the prosthetic leg are measured using an optical motion capture method, and the amount of deformation of each segment can be calculated by optimization calculation using the measured values. The marker positions and the dimensions of the prosthesis are listed in Table 1.

Table 2 lists the experimental items conducted this time.

[C-2] Inverse kinematics calculation Obtain the value of the generalized coordinate q of the PCS model by inverse kinematics calculation from the measured marker position data. In this embodiment, one aspect is a simulation of a hybrid link system consisting of a musculoskeletal model and a prosthetic leg. Considering that the prosthetic leg is connected to the limb, the upper end side is modeled as a base link, and inverse kinematics calculations are performed as a floating link system as shown in the left diagram of FIG. In the left diagram of Figure 5, we are considering a situation in which the base link (segment 1) is not fixed anywhere in the environment, and in order to represent the 6 degrees of freedom of the base link's position and orientation as variables, we have created a virtual environment and base link. It is assumed that the joint angle of the virtual joint represents the six degrees of freedom of the position and posture of the base link. The motion relationship between the twist speed η(s) and the generalized coordinate q is expressed as follows.

Here, J^(s) is the Jacobian matrix of the floating-based PCS model. The calculation results are shown in FIG.

Obtain the value of the generalized coordinate q of the PCS model by inverse kinematics calculation from the measured marker position data. Here, the subscript representing the number of the acquired experimental data is j. For example, when a load is applied in the x direction using a 20kg weight, the generalized coordinate and generalized force data number of the PCS model is set to j=1, and when the weight is increased by 10kg and a 30kg weight is applied in the x direction. The generalized coordinates of the PCS model when a load is applied to Let the data number of generalization power be j=...etc. The index of the obtained data is set as j, and the stiffness of the competition prosthetic leg is estimated by least squares from the generalized coordinate data qj of the PCS model and the ground reaction force data fj, τj.

In the next section, we will explain the estimation of stiffness and viscosity parameters of flexible links.
a) Generalized coordinates, generalized velocity, and generalized acceleration of the PCS model;
b) Information on the magnitude and direction of the load when externally applied, and at what point on the PCS model the load was applied;
c) Density of PCS model;
is used.
By multiplying the ground reaction force by the transposed matrix of the Jacobian matrix and converting it into the generalized force of the PCS model (Equation (21) described later), the relationship between the generalized coordinates and the generalized force can be obtained.

c) Regarding the density of the PCS model, the value of the inertia matrix M is used when estimating the viscosity parameter, but the value of the density of the PCS model is obtained in advance when calculating the inertia matrix M. There is a need. Relatedly, the density of the PCS model also needs to be known in inverse dynamics calculations and forward dynamics calculations. In objects in which the flexible links are made of a single material, such as a leaf spring, the density is constant throughout. If the density catalog value can be referenced from the material model number information, that value can be used. Alternatively, after measuring the mass of the flexible link with a measuring device, the volume can be calculated from 3D CAD data and the density can be calculated. The PCS model itself does not have a constant density, and can handle objects with different density values for each segment.

[D] Estimation of stiffness parameters of prosthetic leg
[D-1] Estimation of stiffness matrix of PCS model
Let τ be the generalized force of the PCS model, and consider that viscoelasticity exists between the generalized coordinate q and the generalized velocity ^q as shown in the following equation.

However, q ₀ represents the value of the generalized coordinates of the PCS model under no load, and K and D represent the stiffness matrix and viscosity matrix, respectively.

Here, we consider static balance and estimate only the stiffness. In addition, in a highly rigid structure such as a prosthetic leg for competitions, it is assumed that the strain u in the translational direction is negligibly small compared to the strain k in the rotational direction, and unless otherwise specified, the rotation Only directional strain is considered. That is, q∈R ^3N , τ∈R ^3N , K∈R ^3N×3N .

Based on equation (18), the stiffness matrix is estimated by the following optimization using experimental data.

However, W is a weight matrix, and Δqj was defined as in the following equation.

Also, since we are considering a state of static equilibrium here, τj can be calculated as follows using the principle of virtual work.

[D-2] Estimation of stiffness by semidefinite programming problem The original PCS model (Non-Patent Document 3) assumes a soft robot made of flexible material such as silicone rubber, and only a diagonal matrix is used as the stiffness matrix K. was expected. In the competition prosthetic leg used in this embodiment, considering that it is made of a material with greater rigidity and that directional interference may occur between force and deformation even in one segment, K is calculated as follows: We decided to consider it as a block diagonal matrix.

However, K _i is a stiffness matrix existing inside segment i, and is considered as a more general positive definite symmetric matrix. In order to solve the optimization of equation (19) under the condition of Ki≧0, consider a semidefinite programming problem. In a semidefinite programming problem, it is necessary to consider a linear objective function for variables, so we introduce an auxiliary variable t and rewrite equation (19) as follows.

The inequality in equation (23) can be converted to the following linear matrix inequality using Schur's complement.

[E] Estimation of viscosity parameters Please note that in this section, equation numbers are assigned independently.
[E-1]PCS model
The dynamic equation in the PCS model is as follows.

Here, τ is the internal elastic load, defined as:

Here, F(q) is the external concentrated load and is defined as follows.

[E-2] Estimation of viscosity matrix The equation of motion derived from equation (1) is as follows.

Here, the left side is calculated by inverse dynamics calculation using the inertia matrix M, and q, q ^· , q¨ are defined as follows.

Each value is calculated by numerical differentiation. p and α in Equation 8 represent the position and orientation of the base link. α is a vector that summarizes three variables, for example, Euler angles. The angular velocity ω ₀ of the base link can be calculated as follows.
Here, the differential of the rotation matrix is calculated as follows.

here,

Therefore, the viscosity matrix D is estimated by the following optimization calculation.

Here, τ _p is defined as follows.

This inequality is transformed into the following linear matrix inequality using the Schur complement.

[F] Estimating the viscoelasticity of the prosthetic leg The kinematic analysis of the hybrid link system requires the viscoelastic force of the flexible device. The viscoelastic force is unique to the flexible device (material, shape, etc.), and basically it is necessary to use an estimated value for each device. Viscoelastic force is calculated from stiffness parameters (stiffness matrix) and viscosity parameters (viscosity matrix). In order to estimate the viscoelastic force, the amount of deformation occurring in the flexible device, its rate of change, and the 6-axis force (3-dimensional translational force and 3-dimensional rotational moment) acting on a specific position on the flexible device are used. Is required.

Estimation of the viscoelastic force of a prosthetic leg, which is exemplified as a flexible device, will be explained. As shown in Figure 10, an optical motion capture system measures the movement (deformation) of the prosthetic leg (on which multiple markers are attached) when force is applied from different directions to the prosthetic leg on the force plate. Obtain the three-dimensional position p of the marker and at the same time obtain the force f using the force plate. For details on the prosthetic leg deformation measurement experiment, please refer to [3][C-1]. Please note that in this section, formula numbers are given independently.

As shown in FIG. 10, generalized coordinates q are obtained by performing inverse kinematic calculation using the three-dimensional position p of the marker. The stiffness matrix K is calculated by solving a semi-definite programming problem using the generalized coordinates q obtained by inverse kinematic calculation and the force f obtained by the force plate.

When estimating the stiffness parameters, in the case of a leaf spring prosthesis, we repeat several trials by applying a load, measuring it, and then removing the load. When compared, it was found that there were subtle variations between trials. In the equation below, it is assumed that plastic (subtle) deformation occurs in the prosthetic leg when a load is applied or removed, and the parameters are estimated including the influence of this variation.

When there is q _err due to inverse kinematics calculation or modeling error,

So, the evaluation function when including all trials is

becomes.

As shown in Figure 10, the generalized force τ is calculated by performing inverse dynamics calculation using the generalized coordinate q (and its numerical derivative) obtained by inverse kinematics calculation and the force f obtained by the force plate. calculate. The viscosity matrix D is calculated by solving a semi-definite programming problem using the estimated stiffness matrix K and the generalized force τ obtained by inverse dynamics calculation.

The evaluation function for estimating the viscosity matrix D is as follows.

Although the estimation of the viscoelasticity of the prosthetic leg has been described, the above estimation means is merely an example. When measuring the amount of deformation of a flexible device while applying a load, the measurement method is not necessarily limited to motion capture; for example, in the case of a prosthetic leg, strain sensors are attached to the surface of the structure to measure deformation. The amount may also be measured. Further, in order to measure changes in the force applied to the flexible device during deformation, a six-axis force sensor may be used instead of the force plate. In the motion analysis of a hybrid link system including a flexible device, if the viscoelasticity of the flexible device is known, that value can be used. If the flexible device is not a one-piece product but is provided as multiple identical devices of the same product/model number, the viscoelasticity obtained as a representative of one device is the common viscoelasticity (minor individual differences). can be treated as a tolerance). Furthermore, in the motion analysis of a hybrid link system, the motion analysis may be performed by giving appropriate values to the stiffness matrix and the viscosity matrix. For example, when it is desired to perform a simulation assuming a structure that is softer than the actual device, virtual values of the elastic matrix K and the viscous matrix D may be set within the simulation. The effect of using a prosthesis that is softer/harder than the one actually used can be seen by running a simulation by setting the values of the virtual elastic matrix K and viscosity matrix D. good.

[G] Motion analysis of a hybrid link system using motion capture Motion analysis of a hybrid link system using motion capture will be explained. A hybrid link system consists of a rigid part and a soft part, and the rigid part can be treated as a rigid link system, and the soft part can be treated as a PCS model. In this embodiment, a prosthetic leg is exemplified as a soft part, and for a prosthetic runner, the prosthetic leg part is modeled as a 6-segment PCS model, the human musculoskeletal part is modeled as a rigid link system, and these are integrated to form a hybrid link system. Perform motion analysis using Figure 11 shows the flow of motion analysis for a human wearing a prosthetic leg. Please note that in this chapter, formula numbers are assigned independently.

[G-1] Kinematics of a hybrid link system In a hybrid link system, the model is divided into N parts, and each part is composed of a rigid link system or a PCS model. In FIG. 12, the relationships among parts, links, and segments are defined as follows.
(a) The subscript j (j=1,...N) is the subscript of the part, and N is the number of parts.
(a) q _j is the generalized coordinate of part j, and does not distinguish whether it is a rigid part or a soft part.
(C) In the rigid part, the generalized coordinates q _j are defined by joint angle vectors, and in the soft part, the generalized coordinates q _j are defined by the strain vectors.
(D) The subscript i (i=1,...n _j ) is the index of the rigid link or PCS segment in part j, and n _j is the number of links or segments in part j.
(e) G _j,i is the position and orientation of link i in rigid part j.
(f) G _j,i (s) is the position and orientation of s (center axis) of segment i in PCS part j.
(g) G _j is the position and orientation of the base link or base segment in part j.
(h) ^j-1 G _j (q _j-1 ) is a homogeneous transformation matrix from parent part j-1 to child part j, and is a function of q _j-1 .

The recursive formula for the forward kinematics of the parent part and child part is

Furthermore, G _j,i , G _j,i (s) are expressed as

It can be calculated from

Given the homogeneous transformation matrix G ₀ representing the position and orientation of the base link and the generalized coordinates q _j of all links and segments, by recursively applying the above formula, any link in the hybrid link system can be or the position and orientation of the segment can be calculated. Assuming a floating link system in which the base link is not fixed, the generalized coordinates of the hybrid link system including G ₀ are expressed as follows.

[G-2] Jacobian matrix of hybrid link system The Jacobian matrix ^w J _j,i of G j _,i in the hybrid link system will be explained. The velocity twist ^w v _j,i of link or segment i of part j is expressed as follows.

This formula can be rewritten as follows.

Here, ^w v _j is the velocity of the base link or segment of part j, the first term on the right side represents the influence of velocity twist (p is the position vector) of the base itself of part j, and the second term is the velocity of the base link or segment of part j. The term represents the effect of velocity twist on segment i (R is a 6-dimensional square matrix with rotation matrix as block diagonal element) considering the base of part j is fixed.

The velocity twist ^j v _j,i of the base link of part j is as follows using the Jacobian matrix J _i (q _j ) of part j,

The following formula is obtained.

Here, ^w J _i (q _j ) is a Jacobian matrix obtained by converting the Jacobian matrix J _i (q _j ) to the world coordinate system. The velocity twist ^w v _j is the same as the velocity twist of the end link of the parent part j-1 ^w v _j-1,nj-1 ,

becomes.

By solving the above recursive equation, the Jacobian matrix ^w J _j,i can be obtained as follows.

^w J _k (k= 0, ..., N-1) is the Jacobian matrix of the kth part, and each matrix contained in [ ^w J ₀ ... ^w J _N-1 ] is represent. When k>i, it becomes a zero matrix, and when k≦i, it is expressed by recursive calculation.

[G-3] Motion measurement using motion capture The type of motion capture used for motion analysis of the hybrid link system is not limited, but in one embodiment, optical motion capture using markers is used. The hybrid link system, which models the human body with a prosthetic leg, consists of two parts (N = 2). One part is a human skeletal model with a missing right lower limb, and in one embodiment, it consists of 50 rigid links and corresponding joints, but does not include the three links corresponding to the right lower limb. The link that connects to the right lower limb is replaced by a cylinder that is connected to the sports prosthesis. The other part is a sports prosthetic leg, which is a PCS model consisting of six segments. In the experiment, 54 optical markers were used: 33 were attached to the human body link, and the remaining 21 were attached to the prosthetic leg.

Optical motion capture systems output three-dimensional position data of a large number of markers. When p _k ∈R ³ (i=k,..,M) is the marker position and M is the number of markers, the vector that stores all p _k is

It is defined as If p^ is the marker position measured by motion capture, the inverse kinematics problem can be formulated as follows.

When solving this optimization problem using the gradient method, a Jacobian matrix J _k at each marker position is required.

To calculate the Jacobian matrix of a marker, consider the local coordinate system G _k of the kth marker. As shown in FIG. 13, the origin of G _k is different from ⁱ p _k . ⁱ p _k is the relative position from the i-th link or segment to the k-th marker. If v _k is a twist vector in the marker coordinate system, the relationship between v _k and the generalized speed of the hybrid link system can be written as follows.

Calculate ^w J _k from the above formula.

The Jacobian matrix for all marker positions is obtained as follows.

[G-4] Motion analysis of hybrid link system As shown in Figure 11, by performing inverse kinematics calculations based on marker position information obtained from motion capture, generalized coordinates, generalized velocity, and generalized acceleration can be calculated. Calculate. By performing inverse dynamics using generalized coordinates, generalized velocity, and generalized acceleration, we can calculate the generalized force τ _S of the PCS model, the joint torque (generalized force of a rigid body) τ _R of the human skeleton model, and the contact force f Get _C.

The dynamics of the hybrid link system can be expressed by the following equation.

The optimization function for inverse dynamics calculations is as follows.
By solving a quadratic programming problem, the generalized force τ _s of the PCS model, the joint torque (generalized force of a rigid body) τ _R of the human skeleton model, and the contact force f _C are estimated. (See Figure 11). ^ref τ _s is calculated using the elasticity matrix K and viscosity matrix D estimated in advance. The muscle tension f _m can be estimated using the joint torque (generalized force of a rigid body) τ _R of the human skeleton model.

Claims

In a hybrid link system consisting of a rigid link system consisting of a plurality of rigid links including a floating base link and at least one flexible link,
The flexible link is approximated by a model based on Cosserat theory,
The equation of motion of the hybrid link system is

Defined by
performing a dynamic calculation based on the equation of motion;
Motion analysis method.
here,
M is the inertia matrix;
b is bias vector;
q¨ is the generalized acceleration of the mixed link;
q 0 is the generalized coordinate of the base link;
θ is the joint vector of the rigid link with n degrees of freedom;
q S is the generalized coordinate of the flexible link approximated by a model based on Cosserat theory;
M 0 is the inertia matrix associated with the base link;
M R is the inertia matrix associated with the rigid link;
M S is the inertia matrix associated with the flexible link;
M 0R is the inertia matrix due to the connection relationship between the base link and the rigid link;
M 0S is an inertia matrix due to the connection relationship between the base link and the flexible link;
M RS is the inertia matrix due to the connection relationship between rigid links and flexible links;
b 0 is the bias vector associated with the base link:
b R is the bias vector associated with the rigid link;
b S is the bias vector associated with the flexible link;
τ R is the generalized force of the rigid link;
τ s is the generalized force including the viscoelastic force of the flexible link;
J C,i is the Jacobian matrix that maps the generalized velocity to the velocity of the contact point;
f C,i is the contact force.
The flexible link is a model based on Cosserat theory that assumes that strain, which is the amount of deformation, is constant regardless of position, or a model that is discretized into multiple segments, and in each segment, strain is constant regardless of position. This is a model that assumes that
The motion analysis method according to claim 1.
The dynamic calculation is an inverse dynamic calculation,
Measure the movement of the target consisting of a hybrid link system and obtain movement data,
The motion data includes first motion data that is motion data of a rigid link system and second motion data that is motion data of a flexible link,
Based on the first motion data, obtain generalized coordinates, generalized velocity, generalized acceleration of the base link, generalized coordinates, generalized velocity, and generalized acceleration of the rigid link system,
obtaining generalized coordinates, generalized velocity, and generalized acceleration of the flexible link based on the second motion data;
Generalized coordinates, generalized velocity, and generalized acceleration of the hybrid link system, a measured or estimated contact force, an inertia matrix, a bias vector, and a Jacobian matrix that maps the generalized velocity to the velocity of the contact point; Calculating the generalized force of the rigid link and the generalized force of the flexible link using the rigidity parameter and viscosity parameter of the flexible link,
The motion analysis method according to claim 1.
In the forward calculation of the inverse dynamics calculation, at the connection between the rigid link and the flexible link, the velocity twist η i of the i-1st rigid link and the segment of the i-th flexible link or flexible link are calculated using the following formula. With the velocity twist η i-1 , the transformation of is performed,

The motion analysis method according to claim 3.
here,
Velocity twist is a vector that combines translational velocity and rotational velocity,
G i-1,i is a matrix representing coordinate transformation between the i-1th rigid link and the i-th flexible link or segment of the flexible link.
In the forward calculation of the inverse dynamics calculation, at the connection between the rigid link and the flexible link, the acceleration twist η・ i of the i-1st rigid link and the acceleration twist η・ i of the i-th flexible link or segment of the flexible link are calculated. -1 and the conversion is performed,
The motion analysis method according to claim 4.
Here, the acceleration twist is a vector that combines translational acceleration and rotational acceleration.
In the backward calculation of the inverse dynamics calculation, at the connection between the rigid link and the flexible link, use the following formula to calculate the wrench vector F i that represents the translational force and rotational moment of the link connected to the link of interest. , wrench vector F i-1 representing the translational force and rotational moment of the link of interest, is transformed into

The motion analysis method according to claim 3.
here,
G i-1,i is a matrix representing coordinate transformation between wrench vector F i and wrench vector F i-1 .
A motion analysis device that analyzes the motion of a target,
The target is a hybrid link system consisting of a rigid link system consisting of a plurality of rigid links including a floating base link and at least one flexible link, the flexible link being approximated by a model based on Cosserat theory. has been
The device includes a processing section and a storage section,
The storage unit stores the equation of motion of the hybrid link system,

is stored,
The processing unit is configured to perform a dynamic calculation based on the equation of motion.
Motion analysis device.
here,
M is the inertia matrix;
b is bias vector;
q¨ is the generalized acceleration of the mixed link;
q 0 is the generalized coordinate of the base link;
θ is the joint vector of the rigid link with n degrees of freedom;
q S is the generalized coordinate of the flexible link approximated by a model based on Cosserat theory;
M 0 is the inertia matrix associated with the base link;
M R is the inertia matrix associated with the rigid link;
M S is the inertia matrix associated with the flexible link;
M 0R is the inertia matrix due to the connection relationship between the base link and the rigid link;
M 0S is an inertia matrix due to the connection relationship between the base link and the flexible link;
M RS is the inertia matrix due to the connection relationship between rigid links and flexible links;
b 0 is the bias vector associated with the base link:
b R is the bias vector associated with the rigid link;
b S is the bias vector associated with the flexible link;
τ R is the generalized force of the rigid link;
τ s is the generalized force including the viscoelastic force of the flexible link;
J C,i is the Jacobian matrix that maps the generalized velocity to the velocity of the contact point;
f C,i is the contact force.
The flexible link is a model based on Cosserat theory that assumes that strain, which is the amount of deformation, is constant regardless of position, or a model that is discretized into multiple segments, and in each segment, strain is constant regardless of position. This is a model that assumes that
The motion analysis device according to claim 7.
The processing unit includes a motion measurement means, an inverse kinematics means, and an inverse dynamics means,
The motion measuring means measures the motion of an object made of a hybrid link system to obtain motion data, and the motion data includes first motion data that is motion data of a rigid link system and motion data of a flexible link. Some second movement data is included,
The inverse kinematics means are:
a first inverse kinematics means for acquiring generalized coordinates, generalized velocity, generalized acceleration of the base link, generalized coordinates, generalized velocity, and generalized acceleration of the rigid link system based on the first motion data; ,
a second inverse kinematics means for acquiring generalized coordinates, generalized velocity, and generalized acceleration of the flexible link based on the second motion data,
The inverse dynamics means converts the generalized coordinates, generalized velocity, and generalized acceleration of the hybrid link system, the measured or estimated contact force, the inertia matrix, the bias vector, and the generalized velocity into the velocity of the contact point. Calculate the generalization force of the rigid link and the generalization force of the flexible link using the Jacobian matrix mapped to the , and the stiffness parameter and viscosity parameter of the flexible link.
The motion analysis device according to claim 7.
In the forward direction calculation, the inverse kinematics means calculates the velocity twist η i of the i-1st rigid link and the i-th flexible link or the flexible link at the connection part between the rigid link and the flexible link using the following formula. the velocity twist η i-1 of the segment of and twist velocity conversion means for performing the transformation of

The motion analysis device according to claim 9.
here,
Velocity twist is a vector that combines translational velocity and rotational velocity,
G i-1,i is a matrix representing coordinate transformation between the i-1th rigid link and the i-th flexible link or segment of the flexible link.
In the forward calculation, the inverse dynamics means calculates the acceleration twist η · i of the i-1st rigid link and the acceleration twist η of the i-th flexible link or segment of the flexible link at the connection portion between the rigid link and the flexible link. - includes twist acceleration conversion means for performing the conversion of i-1 ;
The motion analysis device according to claim 10.
Here, the acceleration twist is a vector that combines translational acceleration and rotational acceleration.
In the reverse direction calculation, the inverse dynamics means calculates a wrench vector F representing the translational force and rotational moment of the link connected to the link of interest using the following formula at the connection between the rigid link and the flexible link. i and a wrench vector F i-1 representing the translational force and rotational moment of the link of interest;

The motion analysis device according to claim 9.
here,
G i-1,i is a matrix representing coordinate transformation between wrench vector F i and wrench vector F i-1 .
In order to perform motion analysis of a target of a hybrid link system consisting of a rigid link system consisting of a plurality of rigid links including a floating base link and at least one flexible link,
A computer program configured to cause a computer to function as the storage unit and processing unit according to any one of claims 7 to 12.