CN113524235B - 一种面向动态交互的对象肢体时变刚度辨识方法及装置 - Google Patents

一种面向动态交互的对象肢体时变刚度辨识方法及装置 Download PDF

Info

Publication number
CN113524235B
CN113524235B CN202110798706.7A CN202110798706A CN113524235B CN 113524235 B CN113524235 B CN 113524235B CN 202110798706 A CN202110798706 A CN 202110798706A CN 113524235 B CN113524235 B CN 113524235B
Authority
CN
China
Prior art keywords
moment
limb
force
time
varying
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
CN202110798706.7A
Other languages
English (en)
Other versions
CN113524235A (zh
Inventor
张琴
陈豫生
刘济汉
熊蔡华
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Huazhong University of Science and Technology
Original Assignee
Huazhong University of Science and Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Huazhong University of Science and Technology filed Critical Huazhong University of Science and Technology
Priority to CN202110798706.7A priority Critical patent/CN113524235B/zh
Publication of CN113524235A publication Critical patent/CN113524235A/zh
Application granted granted Critical
Publication of CN113524235B publication Critical patent/CN113524235B/zh
Priority to US17/858,028 priority patent/US11696687B2/en
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/0048—Detecting, measuring or recording by applying mechanical forces or stimuli
    • A61B5/0057—Detecting, measuring or recording by applying mechanical forces or stimuli by applying motion other than vibrations, e.g. rolling, rubbing, applying a torque, tribometry
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/103—Detecting, measuring or recording devices for testing the shape, pattern, colour, size or movement of the body or parts thereof, for diagnostic purposes
    • A61B5/107—Measuring physical dimensions, e.g. size of the entire body or parts thereof
    • A61B5/1071—Measuring physical dimensions, e.g. size of the entire body or parts thereof measuring angles, e.g. using goniometers
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/103—Detecting, measuring or recording devices for testing the shape, pattern, colour, size or movement of the body or parts thereof, for diagnostic purposes
    • A61B5/11—Measuring movement of the entire body or parts thereof, e.g. head or hand tremor, mobility of a limb
    • A61B5/1121—Determining geometric values, e.g. centre of rotation or angular range of movement
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/103—Detecting, measuring or recording devices for testing the shape, pattern, colour, size or movement of the body or parts thereof, for diagnostic purposes
    • A61B5/11—Measuring movement of the entire body or parts thereof, e.g. head or hand tremor, mobility of a limb
    • A61B5/1126—Measuring movement of the entire body or parts thereof, e.g. head or hand tremor, mobility of a limb using a particular sensing technique
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B25—HAND TOOLS; PORTABLE POWER-DRIVEN TOOLS; MANIPULATORS
    • B25J—MANIPULATORS; CHAMBERS PROVIDED WITH MANIPULATION DEVICES
    • B25J13/00—Controls for manipulators
    • B25J13/08—Controls for manipulators by means of sensing devices, e.g. viewing or touching devices
    • B25J13/085—Force or torque sensors
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B25—HAND TOOLS; PORTABLE POWER-DRIVEN TOOLS; MANIPULATORS
    • B25J—MANIPULATORS; CHAMBERS PROVIDED WITH MANIPULATION DEVICES
    • B25J13/00—Controls for manipulators
    • B25J13/08—Controls for manipulators by means of sensing devices, e.g. viewing or touching devices
    • B25J13/087—Controls for manipulators by means of sensing devices, e.g. viewing or touching devices for sensing other physical parameters, e.g. electrical or chemical properties
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B25—HAND TOOLS; PORTABLE POWER-DRIVEN TOOLS; MANIPULATORS
    • B25J—MANIPULATORS; CHAMBERS PROVIDED WITH MANIPULATION DEVICES
    • B25J18/00—Arms
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B25—HAND TOOLS; PORTABLE POWER-DRIVEN TOOLS; MANIPULATORS
    • B25J—MANIPULATORS; CHAMBERS PROVIDED WITH MANIPULATION DEVICES
    • B25J9/00—Programme-controlled manipulators
    • B25J9/16—Programme controls
    • B25J9/1628—Programme controls characterised by the control loop
    • B25J9/1633—Programme controls characterised by the control loop compliant, force, torque control, e.g. combined with position control
    • B—PERFORMING OPERATIONS; TRANSPORTING
    • B25—HAND TOOLS; PORTABLE POWER-DRIVEN TOOLS; MANIPULATORS
    • B25J—MANIPULATORS; CHAMBERS PROVIDED WITH MANIPULATION DEVICES
    • B25J9/00—Programme-controlled manipulators
    • B25J9/16—Programme controls
    • B25J9/1628—Programme controls characterised by the control loop
    • B25J9/1653—Programme controls characterised by the control loop parameters identification, estimation, stiffness, accuracy, error analysis
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/20—Analysis of motion
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/70—Determining position or orientation of objects or cameras
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06V—IMAGE OR VIDEO RECOGNITION OR UNDERSTANDING
    • G06V40/00—Recognition of biometric, human-related or animal-related patterns in image or video data
    • G06V40/20—Movements or behaviour, e.g. gesture recognition
    • G06V40/23—Recognition of whole body movements, e.g. for sport training
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B90/00—Instruments, implements or accessories specially adapted for surgery or diagnosis and not covered by any of the groups A61B1/00 - A61B50/00, e.g. for luxation treatment or for protecting wound edges
    • A61B90/06—Measuring instruments not otherwise provided for
    • A61B2090/064—Measuring instruments not otherwise provided for for measuring force, pressure or mechanical tension
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B90/00—Instruments, implements or accessories specially adapted for surgery or diagnosis and not covered by any of the groups A61B1/00 - A61B50/00, e.g. for luxation treatment or for protecting wound edges
    • A61B90/06—Measuring instruments not otherwise provided for
    • A61B2090/064—Measuring instruments not otherwise provided for for measuring force, pressure or mechanical tension
    • A61B2090/066—Measuring instruments not otherwise provided for for measuring force, pressure or mechanical tension for measuring torque
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B90/00—Instruments, implements or accessories specially adapted for surgery or diagnosis and not covered by any of the groups A61B1/00 - A61B50/00, e.g. for luxation treatment or for protecting wound edges
    • A61B90/06—Measuring instruments not otherwise provided for
    • A61B2090/067—Measuring instruments not otherwise provided for for measuring angles
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B90/00—Instruments, implements or accessories specially adapted for surgery or diagnosis and not covered by any of the groups A61B1/00 - A61B50/00, e.g. for luxation treatment or for protecting wound edges
    • A61B90/39—Markers, e.g. radio-opaque or breast lesions markers
    • A61B2090/3904—Markers, e.g. radio-opaque or breast lesions markers specially adapted for marking specified tissue
    • A61B2090/3916—Bone tissue
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B34/00—Computer-aided surgery; Manipulators or robots specially adapted for use in surgery
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B90/00—Instruments, implements or accessories specially adapted for surgery or diagnosis and not covered by any of the groups A61B1/00 - A61B50/00, e.g. for luxation treatment or for protecting wound edges
    • A61B90/50—Supports for surgical instruments, e.g. articulated arms
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30196—Human being; Person
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30204—Marker

Abstract

本发明提供了一种面向动态交互的对象肢体时变刚度辨识方法及装置，属于肢体刚度测量领域，方法包括：采集对象肢体位移与实测力数据组合或角度与实测力矩数据组合；基于二阶阻抗模型构建的时变动力学系统，利用线性变参数法，替代时变阻抗参数并重建回复力/力矩表达；根据基函数展开的动态交互力/力矩表达式使用时变动力学参数对变量权重、动态交互力/力矩和回复力/力矩进行迭代辨识；依照替代时变阻抗参数的表达式，使用变量权重及动态交互力/力矩解算时变刚度。本发明不仅提高了时变刚度辨识技术的精度而且扩展了时变刚度辨识技术的应用场景。

Description

一种面向动态交互的对象肢体时变刚度辨识方法及装置
技术领域
本发明属于肢体刚度测量领域，更具体地，涉及一种面向动态交互的对象肢体时变刚度辨识方法及装置。
背景技术
在日常生活中面对丰富多样的运动任务，如行走、跑步、上下楼梯等下肢任务，抑或是接球、举重物这样的上肢任务，生物个体总能表现出卓越的运动稳定性，这得益于生物个体丰富的组织构成和快速高效的神经调控机制。当生物个体的肢体与环境或物体进行动态交互时，中枢神经系统会根据不同的任务需求以调节肢体的刚度特性，从而保持动态交互过程中肢体末端或关节的稳定。
现有刚度研究分为两类，即时不变刚度研究与时变刚度研究。时不变刚度研究假设生物个体的刚度不随时间发生变化，显然这一过程不符合生物个体在与环境或物体动态交互过程中刚度实时变化的特征。因此，研究生物个体时变刚度特征不仅有助于探明生物个体时变刚度的调节机理，将生物个体时变刚度调节机制引入到机器人系统的人机交互控制、自适应变阻抗控制、变阻抗驱动器的研究中还能够显著提升机器人系统的性能表现。生物个体时变刚度的辨识是这类研究的前提，提出一种快速精准的面向动态交互的时变刚度辨识算法十分重要。
目前时变刚度辨识技术有短数据段法、集总数据法、线性变参数法和基函数展开法等。其中，短数据段法是将非平稳的数据集分割成大量被认定为系统参数时不变的准平稳短数据段，然后将具有相似准平稳特征的短数据段分组并用时不变模型辨识从而得到时变刚度。由于短数据段法是基于短段数据准平稳的假设，这一前提使其只能适用于辨识慢速变化的肢体刚度特征，且无法得到肢体刚度的全局时变特征。集总数据法采用多达上百个具有相同时变特征的试验数据构建集总数据，再使用从同一时刻的集总数据中估计出的时不变模型描述时变系统在该时刻的特征。由此可以看出基于集总数据的辨识方法的关键在于进行大量具有相同时变特征的实验。不难发现，短数据段法和集总数据法具有相同的特征，即都需要基于大量且相同时变行为的数据集进行时变刚度参数的辨识，而这一特征导致了实验过程冗长、所需数据量大且多次时变行为的一致性难以保证。
不同于上述两种方法，线性变参数法和基函数展开法能够通过单次实验数据直接辨识出肢体的时变刚度。其中，线性变参数法将时变系统参数建模为以一个或多个随时间变化的调度变量为自变量的函数，这种做法的前提是确定该调度变量直接影响系统时变特征，且调度变量必须能通过直接方法或间接方法得到。基函数展开法则是直接将时变系统参数建模为一系列以时间自变量的基函数的加权和，在基函数形式确定后，该方法将时变系统方程转换为以各阶基函数时不变权重值为未知参数的时不变系统方法，以时不变系统辨识方法辨识出这些未知参数后即可合成得到原时变系统参数。上述两种辨识方法的局限性在于只能用于肢体不与环境交互的情景。
发明内容
现有的短数据段法和集总数据法需要大量时变行为的数据集进行时变刚度参数的辨识，导致多次时变行为的一致性难以保证；同时线性变参数法和基函数展开法只适用于肢体不与环境交互的情景，因此，针对现有技术的缺陷，本发明的目的在于提供一种面向动态交互的对象肢体时变刚度辨识方法及装置，旨在解决现有的时变刚度辨识技术精度不够和无法用于动态交互场景的问题。
为实现上述目的，一方面，本发明提供了一种面向动态交互的对象肢体时变刚度辨识方法，包括以下步骤：
在随机微扰动运动环境下，对象保持肢体姿态不变，通过改变肌肉状态，采集对象肢体位移与实测力数据组合或角度与实测力矩数据组合；
基于替代时变阻抗参数的回复力/力矩的表达式以及采用基函数展开的动态交互力/力矩的表达式，根据对象肢体位移与实测力数据组合或角度与实测力矩数据组合，采用时变动力学的参数对动态交互力/力矩和回复力/力矩进行迭代辨识；
利用参数迭代辨识过程中得到的线性系数矩阵及动态交互力/力矩，依照替代时变阻抗参数的回复力/力矩表达式解算对象肢体时变刚度；
具体地，采集对象肢体位移与实测力数据组合或角度与实测力矩数据组合的方法为：
在协作型机械臂构造的三维空间随机微扰动运动环境，将对象肢体与协作型机械臂工具端固连在一起；
协作型机械臂的工具端固定了一定数量的反光标记点，光学运动捕捉系统会收集反光标记点的三维位置，以得到对象肢体的位移/角度信息；
对象肢体末端或关节相邻的肢体与躯干会固定一定数量的反光标记点用于记录肢体姿态；
机械臂工具端连接六维力/力矩传感器用于采集对象肢体与机械臂工具端之间的力/力矩信息；
对象正前方放置一显示器，显示机械臂工具端的六维力/力矩传感器在实验过程中的力/力矩的变化，为对象肢体提供视觉反馈；
对象肢体时变刚度数据采集过程中，在协作型机械臂构造的随机微扰动运动环境下，对象的肢体末端或关节与协作型机械臂固定连接在一起，保持肢体姿态不变，通过改变肌肉的肌肉状态改变对象肢体的时变刚度；
时变刚度通过协作型机械臂工具端的六维力/力矩传感器所记录的交互力/力矩反映，在X、Y、Z三个方向(或旋转轴)上，先增加正向的交互力/力矩后减小，再增加负向的交互力/力矩再减小，以此规律完成三个方向(或旋转轴)的力/力矩循环即为一个动态力/力矩交互周期；
在规定的扰动时间内，对象肢体需要完成至少一个动态力/力矩变换周期，每个周期的完成时间不受限制；
按照上述方法采集到的对象肢体位移/角度信息及力/力矩信息能够用于辨识对象肢体时变刚度。
具体地，替代时变参数的回复力/力矩的表达式获取方法，包括以下步骤：
根据对象肢体时变阻抗特性，利用二阶阻抗模型，构建对象肢体时变动力学系统，其表达式为：
Figure BDA0003163790340000031
其中，I为对象肢体的惯性参数；B(t)为对象肢体的时变阻尼参数；K(t)为对象肢体的时变刚度参数；x(t)为由光学运动捕捉系统采集的对象肢体位移/角度；fr(t)为对象肢体的回复力/力矩数据；
Figure BDA0003163790340000032
为对象肢体位移/角度的一阶导；
Figure BDA0003163790340000033
为对象肢体位移/角度的二阶导；t为时间。
利用线性变参数法，将动力学模型中对象肢体的时变阻尼参数以及对象肢体的时变刚度参数构建为以调度变量动态交互力/力矩为自变量的函数加权和，其表达式为：
Figure BDA0003163790340000034
其中，B(t)为对象肢体的时变阻尼参数，K(t)是对象肢体的时变刚度参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，ai、bi为第i阶函数的权重，P为所用函数最大阶数；t为时间；
基于替代时变阻抗参数的表达式的对象肢体时变动力学模型的表达式为：
Figure BDA0003163790340000041
其中，I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，ai、bi为第i阶函数的权重，P为所用函数最大阶数，fr(t)为对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度，
Figure BDA0003163790340000042
是对象肢体位移/角度的一阶导，
Figure BDA0003163790340000043
是对象肢体位移/角度的二阶导；t为时间。
构建中间变量xg(i,t)和
Figure BDA0003163790340000044
表达式分别为：
xg(i,t)＝gi(fv(t))x(t)
Figure BDA0003163790340000045
其中，xg(i,t)是由函数gi(fv(t))和对象肢体位移/角度x(t)构建的中间变量，gi(fv(t))为以动态交互力/力矩；fv(t)为自变量的函数，
Figure BDA0003163790340000046
是由函数gi(fv(t))和对象肢体位移/角度的一阶导
Figure BDA0003163790340000047
构建的中间变量，i是gi函数的阶数，t为时间。
基于上述替代时变阻抗参数的表达式及变量，重建回复力/力矩的表达式为：
Figure BDA0003163790340000048
其中，I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，a0…P、b0…P为第0…P阶函数的权重，P为所用函数最大阶数，fr(t)是对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度,β是由对象肢体的惯性参数I、第0…P阶函数的权重a0…P、b0…P构成的线性系数矩阵，
Figure BDA0003163790340000049
是由变量xg(i,t)、
Figure BDA00031637903400000410
对象肢体位移/角度的二阶导
Figure BDA00031637903400000411
构成的位移/角度参数矩阵，xg(i,t)是由函数gi(fv(t))和对象肢体位移/角度x(t)构建的中间变量，
Figure BDA00031637903400000412
是由函数gi(fv(t))和对象肢体位移/角度的一阶导
Figure BDA00031637903400000413
构建的中间变量，i是gi函数的阶数，t为时间。
具体地，将动态交互力/力矩进行基函数展开，采用扰动持续区间内的基函数序列以及基函数的权重表示动态交互力/力矩，其表达式为：
Figure BDA00031637903400000414
其中，fv(t)为动态交互力，hj(t)为扰动持续区间内的基函数序列，cj代表第j阶基函数的权重，Q为所使用基函数的最大阶数；t为时间。
优选地，采用时变动力学参数对动态交互力/力矩和回复力/力矩进行迭代辨识的方法，包括以下步骤：
(1)采集的对象肢体位移与实测力数据组合或角度与实测力矩数据组合，并初始化回复力/力矩的估计值为0；
(2)基于对象肢体位移与实测力数据组合或角度与实测力矩数据组合，结合当前迭代下的回复力/力矩的估计值，根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取动态交互力/力矩的中间变量；
(3)根据基函数展开的动态交互力/力矩的表达式，采用动态交互力/力矩的中间变量辨识基函数的权重；
(4)利用基函数的权重以及扰动持续区间内的基函数序列，获取当前迭代下的动态交互力/力矩的估计值；
(5)根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取回复力/力矩的中间变量；
(6)基于替代时变阻抗参数的回复力/力矩的表达式，采用回复力/力矩的中间变量，获取线性系数矩阵；
(7)根据位移参数矩阵以及当前迭代下获取的线性系数矩阵，获取当前迭代下的回复力/力矩的估计值；
(8)转至步骤(2)，直至实际迭代次数达到预定的最大迭代次数，输出最后一次迭代下获取的线性系数矩阵、回复力/力矩的估计值和动态交互力/力矩的估计值。
具体地，采用时变动力学参数对动态交互力/力矩和回复力/力矩进行迭代辨识的方法的各步骤表达式为：
步骤:(1)：初始化，fm(t)、x(t)、
Figure BDA0003163790340000051
和
Figure BDA0003163790340000052
已知，令fr(t)的初始估计值为
Figure BDA0003163790340000053
步骤(2)(3)(4)：对第n次迭代，计算中间变量
Figure BDA0003163790340000054
根据公式(6)，使
Figure BDA0003163790340000055
从而辨识出权重cj，并重建得到fv(t)的估计值
Figure BDA0003163790340000056
步骤(5)(6)(7)：计算中间变量
Figure BDA0003163790340000057
根据公式(4)、公式(5)，使
Figure BDA0003163790340000058
从而辨识出参数β，并重建得到fr(t)的估计值
Figure BDA0003163790340000059
步骤(8)：令n＝n+1，重新回到步骤2，直至n超过设定的最大迭代次数时迭代终止。
其中，fm(t)为实测力/力矩，fv(t)为动态交互力/力矩，fr(t)为回复力/力矩，
Figure BDA0003163790340000061
为回复力/力矩的初值，
Figure BDA0003163790340000062
为中间变量，
Figure BDA0003163790340000063
为第n-1次迭代的回复力/力矩的估计值，
Figure BDA0003163790340000064
为中间变量，β是由对象肢体的惯性参数Iend、第0…P阶函数的权重a0…P、b0…P构成的线性系数矩阵，
Figure BDA0003163790340000065
是由变量xg(i,t)、
Figure BDA0003163790340000066
对象肢体位移/角度的二阶导
Figure BDA0003163790340000067
构成的位移/角度参数矩阵，
Figure BDA0003163790340000068
为第n次迭代的实测力/力矩的估计值，
Figure BDA0003163790340000069
为第n次迭代的动态交互力/力矩的估计值，
Figure BDA00031637903400000610
为第n次迭代回复力/力矩的估计值；i,j为函数的阶数，t为时间。
利用参数迭代辨识过程中得到的线性系数矩阵及动态交互力/力矩，依照替代时变阻抗参数的回复力/力矩表达式解算对象肢体时变刚度；
另一方面，本发明提供了一种面向动态交互的对象肢体时变刚度辨识装置，包括：数据采集系统、回复力/力矩表达式获取系统、参数迭代辨识系统和动态交互力/力矩表达式获取系统；
数据采集系统用于在随机微扰动运动环境下，对象保持肢体姿态不变，通过改变肌肉状态，采集对象肢体位移与实测力数据组合或角度与实测力矩数据组合；
回复力/力矩表达式获取系统用于根据对象肢体的时变动力学系统，使用线性变参数法替代对象肢体时变阻抗参数重建回复力/力矩的表达；
参数迭代辨识系统用于基于替代时变阻抗参数的回复力/力矩的表达式以及采用基函数展开的动态交互力/力矩的表达式，根据对象肢体位移与实测力数据组合或角度与实测力矩数据组合，采用时变动力学的参数对动态交互力/力矩和回复力/力矩进行迭代辨识；
动态交互力/力矩表达式获取系统用于将动态交互力/力矩进行基函数展开，采用扰动持续区间内的基函数序列以及基函数的权重表示动态交互力/力矩。
优选地，数据采集系统包括协作型机械臂、反光标记点、光学运动捕捉系统、六维力/力矩传感器和显示器；
协作型机械臂的工具端固定设置有反光标记点、且连接有六维力/力矩传感器，使用时通过连接装置与肢体末端/关节进行固定连接；
光学运动捕捉系统用于收集反光标记点的三维位置，获取对象肢体的位移/角度信息；
反光标记点用于记录肢体姿态；
协作型机械臂用于构建三维空间随机微扰动运动环境；
六维力/力矩传感器用于采集对象肢体与机械臂工具端之间的实测力/力矩信息；
显示器为对象提供动态交互过程中六维力/力矩传感器实时数据的视觉反馈。
优选地，回复力/力矩表达式获取系统替代时变阻抗参数的回复力/力矩的表达式获取方法，包括以下步骤：
根据对象肢体时变阻抗特性，利用二阶阻抗模型，构建对象肢体的时变动力学系统；
利用线性变参数法，将动力学模型中对象肢体的时变阻尼参数以及对象肢体的时变刚度参数构建为以调度变量为自变量的函数加权和，对回复力/力矩进行变形转换，替代时变阻抗参数；其中，动态交互力/力矩为调度变量；
具体地，回复力/力矩表达式获取系统替代时变参数的回复力/力矩的表达式获取方法，包括以下步骤：
根据对象肢体时变阻抗特性，利用二阶阻抗模型，构建对象肢体时变动力学系统，其表达式为：
Figure BDA0003163790340000071
其中，I为对象肢体的惯性参数；B(t)为对象肢体的时变阻尼参数；K(t)为对象肢体的时变刚度参数；x(t)为由光学运动捕捉系统采集的对象肢体位移/角度；fr(t)为对象肢体的回复力/力矩数据；
Figure BDA0003163790340000072
为对象肢体位移/角度的一阶导；
Figure BDA0003163790340000073
为对象肢体位移/角度的二阶导；t为时间。
利用线性变参数法，将动力学模型中对象肢体的时变阻尼参数以及对象肢体的时变刚度参数构建为以调度变量动态交互力/力矩为自变量的函数加权和，其表达式为：
Figure BDA0003163790340000074
其中，B(t)为对象肢体的时变阻尼参数，K(t)是对象肢体的时变刚度参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，ai、bi为第i阶函数的权重，P为所用函数最大阶数，t为时间；
基于替代时变阻抗参数的表达式的对象肢体时变动力学模型的表达式为：
Figure BDA0003163790340000075
其中，I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，ai、bi为第i阶函数的权重，P为所用函数最大阶数，fr(t)为对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度，
Figure BDA0003163790340000076
是对象肢体位移/角度的一阶导，
Figure BDA0003163790340000077
是对象肢体位移/角度的二阶导，t为时间。
构建中间变量xg(i,t)和
Figure BDA0003163790340000081
表达式分别为：
xg(i,t)＝gi(fv(t))x(t)
Figure BDA0003163790340000082
其中，xg(i,t)是由函数gi(fv(t))和对象肢体位移/角度x(t)构建的中间变量；gi(fv(t))为以动态交互力/力矩；fv(t)为自变量的函数；
Figure BDA0003163790340000083
是由函数gi(fv(t))和对象肢体位移/角度的一阶导
Figure BDA0003163790340000084
构建的中间变量；i是gi函数的阶数；t为时间。
基于上述替代时变阻抗参数的表达式及变量，重建回复力/力矩的表达式为：
Figure BDA0003163790340000085
Figure BDA0003163790340000086
β＝[I b0 b1 … bP a0 a1 … aP]T
其中，I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，a0…P、b0…P为第0…P阶函数的权重，P为所用函数最大阶数，fr(t)是对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度,β是由对象肢体的惯性参数I、第0…P阶函数的权重a0…P、b0…P构成的线性系数矩阵，
Figure BDA0003163790340000087
是由变量xg(i,t)、
Figure BDA0003163790340000088
对象肢体位移/角度的二阶导
Figure BDA0003163790340000089
构成的位移/角度参数矩阵，xg(i,t)是由函数gi(fv(t))和对象肢体位移/角度x(t)构建的中间变量，
Figure BDA00031637903400000810
是由函数gi(fv(t))和对象肢体位移/角度的一阶导
Figure BDA00031637903400000811
构建的中间变量，i是gi函数的阶数，t为时间。
优选地，参数迭代辨识系统辨识动态交互力/力矩和回复力/力矩的方法，包括以下步骤：
(1)采集的对象肢体位移与实测力数据组合或角度与实测力矩数据组合，并初始化回复力/力矩的估计值为0；
(2)基于对象肢体位移与实测力数据组合或角度与实测力矩数据组合，结合当前迭代下的回复力/力矩的估计值，根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取动态交互力/力矩的中间变量；
(3)根据基函数展开的动态交互力/力矩的表达式，采用动态交互力/力矩的中间变量辨识基函数的权重；
(4)利用基函数的权重以及扰动持续区间内的基函数序列，获取当前迭代下的动态交互力/力矩的估计值；
(5)根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取回复力/力矩的中间变量；
(6)基于替代时变阻抗参数的回复力/力矩的表达式，采用回复力/力矩的中间变量，获取线性系数矩阵；
(7)根据位移参数矩阵以及当前迭代下获取的线性系数矩阵，获取当前迭代下的回复力/力矩的估计值；
(8)转至步骤(2)，直至实际迭代次数达到预定的最大迭代次数，输出最后一次迭代下获取的线性系数矩阵、回复力/力矩的估计值和动态交互力/力矩的估计值。
具体地，参数迭代辨识系统采用时变动力学参数对动态交互力/力矩和回复力/力矩进行迭代辨识的方法的各步骤表达式为：
步骤(1)：初始化，fm(t)、x(t)、
Figure BDA0003163790340000091
和
Figure BDA0003163790340000092
已知，令fr(t)的初始估计值为
Figure BDA0003163790340000093
步骤(2)(3)(4)：对第n次迭代，计算中间变量
Figure BDA0003163790340000094
使
Figure BDA0003163790340000095
从而辨识出权重cj，并重建得到fv(t)的估计值
Figure BDA0003163790340000096
步骤(5)(6)(7)：计算中间变量
Figure BDA0003163790340000097
使
Figure BDA0003163790340000098
从而辨识出参数β，并重建得到fr(t)的估计值
Figure BDA0003163790340000099
步骤(8)：令n＝n+1，重新回到步骤2，直至n超过设定的最大迭代次数时迭代终止。
其中，fm(t)为实测力/力矩，fv(t)为动态交互力/力矩，fr(t)为回复力/力矩，
Figure BDA00031637903400000910
为回复力/力矩的初值，
Figure BDA00031637903400000911
为中间变量，
Figure BDA00031637903400000912
为第n-1次迭代的回复力/力矩的估计值，
Figure BDA00031637903400000913
为中间变量，β是由对象肢体的惯性参数Iend、第0…P阶函数的权重a0…P、b0...P构成的线性系数矩阵，
Figure BDA00031637903400000914
是由变量xg(i,t)、
Figure BDA00031637903400000915
对象肢体位移/角度的二阶导
Figure BDA00031637903400000916
构成的位移参数矩阵，
Figure BDA00031637903400000917
为第n次迭代的实测力/力矩的估计值，
Figure BDA00031637903400000918
为第n次迭代的动态交互力/力矩的估计值，
Figure BDA00031637903400000920
为第n次迭代回复力/力矩的估计值，i,j为函数的阶数，t为时间。
优选地，动态交互力/力矩表达式获取系统对动态交互力/力矩进行基函数展开的步骤为：采用扰动持续区间内的基函数序列以及基函数的权重表示动态交互力/力矩。
具体地，采用基函数展开的动态交互力/力矩的表达式为：
Figure BDA00031637903400000919
其中，fv(t)为动态交互力，hj(t)为扰动持续区间内的基函数序列，cj代表第j阶基函数的权重，Q为所使用基函数的最大阶数，t为时间。
总体而言，通过本发明所构思的以上技术方案与现有技术相比，具有以下有益效果：
本发明提供的面向动态交互的对象肢体时变刚度辨识方法，是基于线性变参数法、基函数展开法和参数迭代辨识方法构建的，通过线性变参数法将对象肢体时变刚度研究问题转换为时不变研究问题，为时变刚度的辨识提供了新的思路；提供的参数迭代辨识能够依靠单次扰动数据完整准确地辨识出对象肢体时变刚度，极大地减小了时变刚度研究的复杂性，更为重要地突破了现有依靠单次实验数据的时变刚度辨识技术无法应用于与环境动态交互情景的局限性。总体而言，本发明提供的对象肢体时变刚度辨识方法提高了时变刚度辨识技术的精度。
本发明提供的对象肢体时变刚度辨识的运动数据采集方法能够统一地指导对象肢体的时变刚度研究的实验数据采集，实验范式简单、操作性高，能够完整地改变对象肢体的刚度特性，且记录的数据能够精准地解算出对象肢体的时变刚度。
本发明利用多自由度协作型机械臂构建微扰动环境以进行对象肢体时变刚度测量的方法，省去了开发专用型测量装置的复杂过程，降低了对象肢体时变刚度测量的难度，并且拓展了利用协作型机械臂进行刚度测量的应用场景，使对象上、下肢的时变刚度测量能够借助协作型机械臂进行研究。
附图说明
图1是本发明实施例提供的面向动态交互的对象肢体时变刚度辨识方法示意图；
图2是本发明实施例提供的扰动试验示意图；
图3是本发明实施例提供的实验过程中力循环变化图；
图4是本发明实施例提供的上肢末端时变刚度椭球示意图；
图5是本发明实施例提供的上肢末端时变刚度曲线示意图。
具体实施方式
为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。
一方面，本发明提供了一种面向动态交互的对象肢体时变刚度的辨识方法，包括对象肢体时变刚度辨识的运动数据采集方法和对象肢体时变刚度辨识方法；
对象肢体时变刚度辨识的运动数据采集方法的过程为：在协作机械臂构造的随机微扰动运动情景下通过调整肌肉激活状态改变对象肢体阻抗特性，并利用物理传感器获取对象肢体运动数据；
对象肢体时变刚度辨识方法包括线性变参数法、基函数展开法和迭代辨识方法；
其中，线性变参数法将对象肢体时变阻抗参数转换成时不变参数形式，再用基函数展开法描述对象肢体的动态交互力，然后使用迭代辨识方法实现对象肢体回复力和交互力的解耦，最终实现对对象肢体时变刚度的精确辨识。
本发明提出的对象肢体时变刚度辨识方法不仅能够获取肢体末端的时变刚度，如上肢末端，还能够解算对象肢体各个关节的时变刚度。此外，本发明不限于使用对象肢体的动态交互力与位移数据组合对对象肢体的刚度进行辨识，当对象肢体用于刚度辨识的数据库为力矩与角度的数据组合时，本发明所提出的方法同样适用，只需将力与力矩进行替换，位移与角度进行替换即可。
具体内容如下：
(1)对象肢体时变刚度辨识的运动数据采集
本发明使用协作型机械臂以构造三维空间随机微扰动运动环境，将对象肢体与协作型机械臂工具端固连在一起。协作型机械臂的工具端固定了一定数量的反光标记点，光学运动捕捉系统会收集反光标记点的三维位置，以得到对象肢体的位移/角度信息。对象肢体与躯干会固定一定数量的反光标记点用于记录肢体姿态；
机械臂工具端连接六维力/力矩传感器用于采集对象肢体与机械臂工具端之间的力/力矩信息；
对象正前方放置一显示器，显示机械臂工具端的六维力/力矩传感器在实验过程中的力/力矩的变化，为对象肢体提供视觉反馈；
对象肢体时变刚度数据采集过程中，在协作型机械臂构造的随机微扰动运动环境下，对象的肢体末端或关节与协作型机械臂固定连接在一起，保持肢体姿态不变，通过改变肌肉的肌肉状态改变对象肢体的时变刚度；
时变刚度通过协作型机械臂工具端的六维力/力矩传感器所记录的交互力/力矩反映，在X、Y、Z三个方向(或旋转轴)上，先增加正向的交互力/力矩后减小，再增加负向的交互力/力矩再减小，以此规律完成三个方向(或旋转轴)的力/力矩循环即为一个动态力/力矩交互周期；
在规定的扰动时间内，对象肢体需要完成至少一个动态力/力矩变换周期，每个周期的完成时间不受限制；
按照上述方法采集到的对象肢体位移/角度信息及力/力矩信息能够用于辨识对象肢体的时变刚度。
(2)对象肢体时变动力学模型的建立
为了表征对象肢体的时变刚度特性，具有结构简单且能较全面地反映对象肢体时变阻抗特征的二阶阻抗模型被应用于对象肢体刚度研究，其动力学模型通常建模为：
Figure BDA0003163790340000121
其中，t为时间，I为对象肢体的惯性参数；B(t)为对象肢体的时变阻尼参数；K(t)为对象肢体的时变刚度参数；x(t)为由光学运动捕捉系统采集的对象肢体位移/角度；fr(t)为对象肢体的回复力/力矩数据；
Figure BDA0003163790340000122
为对象肢体位移/角度的一阶导；
Figure BDA0003163790340000123
为对象肢体位移/角度的二阶导；由上述动力学模型可知，微扰动环境下的对象肢体位移/角度x(t)和回复力/力矩fr(t)是辨识末端刚度的关键数据；其中，对象肢体位移/角度x(t)一般通过光学运动捕捉系统进行获取；回复力/力矩fr(t)则为机械臂工具端力/力矩传感器所记录数据的一部分；动态交互过程中机械臂工具端的力/力矩传感器所记录的数据为实测力/力矩fm(t)，其与动态交互力/力矩fv(t)、回复力/力矩fr(t)之间的关系为描述为：
fm(t)＝fr(t)+fv(t)(2)
因此，对实测力/力矩fm(t)进行解耦获取动态交互过程中对象肢体的回复力/力矩fr(t)和动态交互力/力矩fv(t)是解决对象肢体时变刚度辨识的关键。
(3)对象肢体时变刚度的辨识
(3.1)对象肢体时变阻抗参数的线性表达：
线性变参数法的思想是：通过先验知识在被辨识的系统确定出某些直接影响系统时变特性的时变量，作为调度变量，然后将未知的时变参数建模为调度变量的一系列函数的加权和。在动态交互过程中，对象肢体所有肌肉的活动状态是难以直接获得的，但它们的共同作用使得对象肢体与外部环境之间表现出状态的交互力/力矩fv(t)，因此动态交互力/力矩fv(t)直接影响系统刚度特征。根据上面的分析，可以将动态交互力/力矩fv(t)设定为调度变量，确定调度变量后，将对象肢体时变刚度K(t)和对象肢体时变阻尼B(t)构建为以调度变量为自变量的一系列函数加权和：
Figure BDA0003163790340000124
其中，B(t)为对象肢体的时变阻尼参数，K(t)是对象肢体的时变刚度参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，ai、bi为第i阶函数的权重，P为所用函数最大阶数；t为时间；
将公式(3)代入公式(1)中，可得到：
Figure BDA0003163790340000131
其中，I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，ai、bi为第i阶函数的权重，P为所用函数最大阶数，fr(t)为对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度，
Figure BDA0003163790340000132
是对象肢体位移/角度的一阶导，
Figure BDA0003163790340000133
是对象肢体位移/角度的二阶导，t为时间。
构建变量xg(i,t)和
Figure BDA0003163790340000134
令
Figure BDA0003163790340000135
其中，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度，
Figure BDA00031637903400001312
是对象肢体位移/角度的一阶导，i是gi函数的阶数，t为时间。
则式(4)可转换为如下形式：
Figure BDA0003163790340000136
其中，I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，a0…P、b0…P为第0…P阶函数的权重，P为所用函数最大阶数，fr(t)是对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移,β是由对象肢体的惯性参数I、第0…P阶函数的权重a0…P、b0…P构成的线性系数矩阵，
Figure BDA0003163790340000137
是由变量xg(i,t)、
Figure BDA0003163790340000138
对象肢体位移/角度的二阶导
Figure BDA0003163790340000139
构成的位移/角度参数矩阵，xg(i,t)是由gi(fv(t))和x(t)构建的中间变量，
Figure BDA00031637903400001310
是由gi(fv(t))和
Figure BDA00031637903400001311
构建的中间变量，i是gi函数的阶数，t为时间。
由此，线性变参数法将时变系统参数的辨识问题转化为时不变未知参数的辨识问题。
(3.2)动态交互力/力矩fv(t)的基函数展开：
为了能够构建动态交互力/力矩fv(t)的时不变表达形式，基函数展开法被用于动态交互力的表达，对动态交互力/力矩fv(t)进行基函数展开如下式所示：
Figure BDA0003163790340000141
其中，fv(t)为动态交互力/力矩，hj(t)为扰动持续区间内的基函数序列，cj代表第j阶基函数的权重，Q为所使用基函数的最大阶数，t为时间。
由于基函数hj(t)只与时间有关，当基函数形式与阶数确定时，借助线性回归方法即可得到权重值c0,c1,…,cQ。
(3.3)时变动力学系统的参数迭代辨识
基于以上对动态交互力/力矩fv(t)与对象回复力/力矩fr(t)的描述，用逐步迭代的辨识方法从实测力/力矩数据fm(t)中对这两部分力进行解耦，迭代辨识方法步骤如下：
1.初始化，fm(t)、x(t)、
Figure BDA0003163790340000142
和
Figure BDA0003163790340000143
已知，令fr(t)的初始估计值为
Figure BDA0003163790340000144
2.对第n次迭代，计算中间变量
Figure BDA0003163790340000145
根据公式(7)，使
Figure BDA0003163790340000146
从而辨识出权重cj，并重建得到fv(t)的估计值
Figure BDA0003163790340000147
3.计算中间变量
Figure BDA0003163790340000148
根据公式(5)、公式(6)，使
Figure BDA0003163790340000149
从而辨识出参数β，并重建得到fr(t)的估计值
Figure BDA00031637903400001410
4.令n＝n+1，重新回到步骤2，直至n超过设定的最大迭代次数时迭代终止。
其中，fm(t)为实测力/力矩，fv(t)为动态交互力/力矩，fr(t)为回复力/力矩，
Figure BDA00031637903400001411
为回复力/力矩的初值，
Figure BDA00031637903400001412
为中间变量，
Figure BDA00031637903400001413
为第n-1次迭代的回复力/力矩的估计值，
Figure BDA00031637903400001414
为中间变量，β是由对象肢体的惯性参数I、第0…P阶函数的权重a0…P、b0…P构成的矩阵，
Figure BDA00031637903400001415
是由变量xg(i,t)、
Figure BDA00031637903400001416
对象肢体位移/角度的二阶导
Figure BDA00031637903400001417
构成的矩阵，
Figure BDA00031637903400001418
为第n次迭代的实测力/力矩的估计值，
Figure BDA00031637903400001419
为第n次迭代的动态交互力/力矩的估计值，
Figure BDA00031637903400001420
为第n次迭代回复力/力矩的估计值，i,j为函数的阶数，t为时间。
经过上述的迭代辨识算法，就能从参数β中提取出式(3)中的权重值，回代权重值及动态交互力的值于式(3)从而得到动态交互过程中的对象肢体时变刚度K(t)。
另一方面，本发明提供了一种面向动态交互的对象肢体时变刚度辨识装置，包括：数据采集系统、回复力/力矩表达式获取系统、参数迭代辨识系统和动态交互力/力矩表达式获取系统；
(1)数据采集系统
数据采集系统用于在随机微扰动运动环境下，对象保持肢体姿态不变，通过改变肌肉状态，采集对象肢体位移与实测力数据组合或角度与实测力矩数据组合，构建对象肢体时变刚度辨识的数据库；
数据采集系统包括协作型机械臂、反光标记点、光学运动捕捉系统、六维力/力矩传感器和显示器；
协作型机械臂的工具端固定设置有反光标记点、且连接有六维力/力矩传感器，使用时通过连接装置与肢体末端/关节进行固定连接；
光学运动捕捉系统用于收集反光标记点的三维位置，获取对象肢体的位移/角度信息；
反光标记点用于记录肢体姿态；
协作型机械臂用于构建三维空间随机微扰动运动环境；
六维力/力矩传感器用于采集对象肢体与机械臂工具端之间的实测力/力矩信息。
显示器为对象提供动态交互过程中六维力/力矩传感器实时数据的视觉反馈。
(2)回复力/力矩表达式获取系统
回复力/力矩表达式获取系统用于根据对象肢体的时变动力学系统，使用线性变参数法替代对象肢体时变阻抗参数重建回复力/力矩的表达；
回复力/力矩表达式获取系统替代时变阻抗参数的回复力/力矩的表达式获取方法，包括以下步骤：
根据对象肢体时变阻抗特性，利用二阶阻抗模型，构建对象肢体的时变动力学系统；
利用线性变参数法，将动力学模型中对象肢体的时变阻尼参数以及对象肢体的时变刚度参数构建为以调度变量为自变量的函数加权和，对回复力/力矩进行变形转换，替代时变阻抗参数；其中，动态交互力/力矩为调度变量；
具体地，回复力/力矩表达式获取系统替代时变参数的回复力/力矩的表达式获取方法，包括以下步骤：
根据对象肢体时变阻抗特性，利用二阶阻抗模型，构建对象肢体时变动力学系统，其表达式为：
Figure BDA0003163790340000151
其中，t为时间，I为对象肢体的惯性参数；B(t)为对象肢体的时变阻尼参数；K(t)为对象肢体的时变刚度参数；x(t)为由光学运动捕捉系统采集的对象肢体位移/角度；fr(t)为对象肢体的回复力/力矩数据；
Figure BDA0003163790340000161
为对象肢体位移/角度的一阶导；
Figure BDA0003163790340000162
为对象肢体位移/角度的二阶导。
利用线性变参数法，将动力学模型中对象肢体的时变阻尼参数以及对象肢体的时变刚度参数构建为以调度变量动态交互力/力矩为自变量的函数加权和，其表达式为：
Figure BDA0003163790340000163
其中，B(t)为对象肢体的时变阻尼参数，K(t)是对象肢体的时变刚度参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，ai、bi为第i阶函数的权重，P为所用函数最大阶数，t为时间；
基于替代时变阻抗参数的表达式的对象肢体时变动力学模型的表达式为：
Figure BDA0003163790340000164
其中，t为时间；I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，ai、bi为第i阶函数的权重，P为所用函数最大阶数，fr(t)为对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度，
Figure BDA0003163790340000165
是对象肢体位移/角度的一阶导，
Figure BDA0003163790340000166
是对象肢体位移/角度的二阶导。
构建中间变量xg(i,t)和
Figure BDA0003163790340000167
表达式分别为：
xg(i,t)＝gi(fv(t))x(t)
Figure BDA0003163790340000168
其中，xg(i,t)是由函数gi(fv(t))和对象肢体位移/角度x(t)构建的中间变量，
Figure BDA0003163790340000169
是由函数gi(fv(t))和对象肢体位移/角度的一阶导
Figure BDA00031637903400001610
构建的中间变量，i是gi函数的阶数，t为时间。
基于上述替代时变阻抗参数的表达式及变量，重建回复力/力矩的表达式为：
Figure BDA00031637903400001611
Figure BDA00031637903400001612
β＝[I b0 b1 … bP a0 a1 … aP]T
其中，I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，a0…P、b0…P为第0…P阶函数的权重，P为所用函数最大阶数，fr(t)是对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度,β是由对象肢体的惯性参数I、第0…P阶函数的权重a0…P、b0…P构成的线性系数矩阵，
Figure BDA0003163790340000171
是由变量xg(i,t)、
Figure BDA0003163790340000172
对象肢体位移/角度的二阶导
Figure BDA0003163790340000173
构成的位移/角度参数矩阵，xg(i,t)是由函数gi(fv(t))和对象肢体位移/角度x(t)构建的中间变量，
Figure BDA0003163790340000174
是由函数gi(fv(t))和对象肢体位移/角度的一阶导
Figure BDA0003163790340000175
构建的中间变量，i是gi函数的阶数，t为时间。
(3)参数迭代辨识系统
参数迭代辨识系统用于基于替代时变阻抗参数的回复力/力矩的表达式以及采用基函数展开的动态交互力/力矩的表达式，根据对象肢体位移与实测力数据组合或角度与实测力矩数据组合，采用时变动力学的参数对动态交互力/力矩和回复力/力矩进行迭代辨识；
参数迭代辨识系统辨识动态交互力/力矩和回复力/力矩的方法，包括以下步骤：
(1)采集的对象肢体位移与实测力数据组合或角度与实测力矩数据组合，并初始化回复力/力矩的估计值为0；
(2)基于对象肢体位移与实测力数据组合或角度与实测力矩数据组合，结合当前迭代下的回复力/力矩的估计值，根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取动态交互力/力矩的中间变量；
(3)根据基函数展开的动态交互力/力矩的表达式，采用动态交互力/力矩的中间变量辨识基函数的权重；
(4)利用基函数的权重以及扰动持续区间内的基函数序列，获取当前迭代下的动态交互力/力矩的估计值；
(5)根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取回复力/力矩的中间变量；
(6)基于替代时变阻抗参数的回复力/力矩的表达式，采用回复力/力矩的中间变量，获取线性系数矩阵；
(7)根据位移参数矩阵以及当前迭代下获取的线性系数矩阵，获取当前迭代下的回复力/力矩的估计值；
(8)转至步骤(2)，直至实际迭代次数达到预定的最大迭代次数，输出最后一次迭代下获取的线性系数矩阵、回复力/力矩的估计值和动态交互力/力矩的估计值。
具体地，参数迭代辨识系统采用时变动力学参数对动态交互力/力矩和回复力/力矩进行迭代辨识的方法的各步骤表达式为：
步骤(1):初始化，fm(t)、x(t)、
Figure BDA0003163790340000181
和
Figure BDA0003163790340000182
已知，令fr(t)的初始估计值为
Figure BDA0003163790340000183
步骤(2)(3)(4):对第n次迭代，计算中间变量
Figure BDA0003163790340000184
使
Figure BDA0003163790340000185
从而辨识出权重cj，并重建得到fv(t)的估计值
Figure BDA0003163790340000186
步骤(5)(6)(7):计算中间变量
Figure BDA0003163790340000187
使
Figure BDA0003163790340000188
从而辨识出参数β，并重建得到fr(t)的估计值
Figure BDA0003163790340000189
步骤(8)：令n＝n+1，重新回到步骤2，直至n超过设定的最大迭代次数时迭代终止。
其中，fm(t)为实测力/力矩，fv(t)为动态交互力/力矩，fr(t)为回复力/力矩，
Figure BDA00031637903400001810
为回复力/力矩的初值，
Figure BDA00031637903400001811
为中间变量，
Figure BDA00031637903400001812
为第n-1次迭代的回复力/力矩的估计值，
Figure BDA00031637903400001813
为中间变量，β是由对象肢体的惯性参数Iend、第0…P阶函数的权重a0…P、b0…P构成的线性系数矩阵，
Figure BDA00031637903400001814
是由变量xg(i,t)、
Figure BDA00031637903400001815
对象肢体位移/角度的二阶导
Figure BDA00031637903400001816
构成的位移参数矩阵，
Figure BDA00031637903400001817
为第n次迭代的实测力/力矩的估计值，
Figure BDA00031637903400001818
为第n次迭代的动态交互力/力矩的估计值，
Figure BDA00031637903400001819
为第n次迭代回复力/力矩的估计值；i,j为函数的阶数，t为时间。
(4)动态交互力/力矩表达式获取系统
动态交互力/力矩表达式获取系统用于将动态交互力/力矩进行基函数展开，采用扰动持续区间内的基函数序列以及基函数的权重表示动态交互力/力矩。
动态交互力/力矩表达式获取系统对动态交互力/力矩进行基函数展开的步骤为：采用扰动持续区间内的基函数序列以及基函数的权重表示动态交互力/力矩。
具体地，采用基函数展开的动态交互力/力矩的表达式为：
Figure BDA00031637903400001820
其中，fv(t)为动态交互力，hj(t)为扰动持续区间内的基函数序列，cj代表第j阶基函数的权重，Q为所使用基函数的最大阶数，t为时间。
实施例
本实施例是使用一种六自由度的协作型机械臂对人体上肢末端施加三维随机扰动，从而构造随机微扰动的运动情景。通过改变肌肉的激活状态改变上肢末端的刚度特性，利用物理传感器收集扰动过程中上肢末端的位移与回复力数据，根据提出的时变刚度辨识方法计算得到三维的上肢末端的时变刚度。整体测量流程如图1所示，具体包含以下步骤：
(1)上肢末端运动数据采集
受试者坐在一张硬质椅子上，六自由度的协作型机械臂的握把位于受试者右肩旋转中线前方390mm、下方100mm的位置，如图2所示；在受试者握住机械臂工具端所固定的握把后，启动机械器并使机械臂执行扰动，扰动的最大幅值设定为10mm，频带范围设定为0～5Hz，参数的设定能确保扰动过程为后续的末端动力学辨识提供的数据具有准确性。受试者保持固定的上肢姿态，通过调节肌肉的激活状态以改变肢体末端的刚度特性，机械臂工具端固定的力传感器收集上肢因受到扰动而产生的回复力和动态交互力数据，外部的光学三维运动捕捉系统收集肢体末端的位移数据。在运动数据采集的过程中，受试者正前方放置一显示器，显示机械臂工具端的六维力/力矩传感器在实验过程中的力/力矩的变化，为受试者提供视觉反馈。受试者需要在X、Y、Z三个方向上分别进行范围在±30N的力循环，根据显示器所提供的力反馈信息，先增加正向的交互力后减小，再增加负向的交互力再减小，以此规律完成三个方向的力循环即为一个动态力交互周期。在20s的扰动时间内，受试者需要完成至少一个动态力变换周期，每个周期的完成时间不受限制，实测力的力循环如图3所示。
(2)上肢末端时变动力学模型的建立
为了表征人体上肢末端的时变刚度特性，具有结构简单且能较全面地反映人体上肢末端时变阻抗特征的二阶阻抗模型被应用于上肢末端刚度辨识的研究，其动力学模型通常建模为：
Figure BDA0003163790340000191
其中，t为时间，I是人体上肢末端的惯性参数，B(t)是人体上肢末端的时变阻尼参数，K(t)是人体上肢末端的时变刚度参数，x(t)是由光学运动捕捉系统采集的人体上肢末端的位移，fr(t)是人体上肢末端的回复力数据，
Figure BDA0003163790340000192
是人体上肢末端位移的一阶导，
Figure BDA0003163790340000193
是人体上肢末端位移的二阶导。由上述动力学模型可知，微扰动环境下的人体上肢末端的位移x(t)和上肢末端的回复力fr(t)是辨识末端刚度的关键数据。其中，人体上肢末端的位移x(t)一般通过光学运动捕捉系统进行获取，上肢末端的回复力fr(t)则是机械臂工具端力/力矩传感器所记录数据的一部分。动态交互过程中机械臂工具端的力/力矩传感器所记录的数据为实测力fm(t)，其与动态交互力fv(t)、回复力fr(t)之间的关系被描述为：
fm(t)＝fr(t)+fv(t) (2)
因此，对实测力fm(t)进行解耦获取动态交互过程中人体上肢末端的回复力fr(t)和动态交互力fv(t)是解决人体上肢末端时变刚度辨识的关键。
(3)人体上肢末端时变刚度的解算
(3.1)人体上肢末端时变阻抗参数的线性表达：
线性变参数法的思想是，通过先验知识在被辨识的系统中确定出某些直接影响系统时变特性的时变量，作为调度变量，然后将未知的时变参数建模为调度变量的一系列函数的加权和。在动态交互过程中，上肢所有肌肉的活动状态是难以直接获得的，但它们的共同作用使得人体上肢末端与外部环境之间表现出动态的交互力fv(t)，因此动态交互力fv(t)直接影响系统刚度特征。根据上面的分析，可以将动态交互力fv(t)设定为调度变量，确定调度变量后，将上肢末端时变刚度Kend(t)和人体上肢末端时变阻尼Bend(t)构建为以调度变量为自变量的一系列函数加权和：
Figure BDA0003163790340000201
其中，B(t)为人体上肢末端的时变阻尼参数；K(t)是人体上肢末端的时变刚度参数；gi(fv(t))为以动态交互力fv(t)为自变量的函数；ai、bi为第i阶函数的权重；P为所用函数最大阶数；t为时间；
将式(3)代入式(1)，可得到：
Figure BDA0003163790340000202
其中，Iend是人体上肢末端的惯性参数，gi(fv(t))为以动态交互力fv(t)为自变量的函数，ai、bi为第i阶函数的权重，P为所用函数最大阶数，fr(t)为人体上肢末端的回复力，x(t)是由光学运动捕捉系统采集的人体上肢末端的位移，
Figure BDA0003163790340000203
是人体上肢末端位移的一阶导，
Figure BDA0003163790340000204
是人体上肢末端位移的二阶导，t为时间。
构建变量xg(i,t)和
Figure BDA0003163790340000205
令
Figure BDA0003163790340000206
则式(4)可转换为如下形式：
Figure BDA0003163790340000211
其中，Iend是人体的惯性参数，gi(fv(t))为以动态交互力fv(t)为自变量的函数，a0…P、b0…P为第0…P阶函数的权重，P为所用函数最大阶数，fr(t)是人体的回复力，x(t)是由光学运动捕捉系统采集的人体位移,β是由人体的惯性参数I、第0…P阶函数的权重a0…P、b0…P构成的线性系数矩阵，
Figure BDA0003163790340000212
是由变量xg(i,t)、
Figure BDA0003163790340000213
人体位移的二阶导
Figure BDA0003163790340000214
构成的位移参数矩阵，xg(i,t)是由gi(fv(t))和x(t)构建的中间变量，
Figure BDA0003163790340000215
是由gi(fv(t))和
Figure BDA0003163790340000216
构建的中间变量，i是gi函数的阶数，t为时间。
由此，线性变参数法将时变系统参数的辨识问题转化为时不变未知参数的辨识问题。本实施例选择使用0～3阶切比雪夫多项式作为描述上肢末端时变阻抗参数与动态交互力之间关系的函数gi(fv(t))，但函数gi(fv(t))的形式不限于切比雪夫多项式。
(3.2)动态交互力fv(t)的基函数展开：
为了能够构建动态交互力fv(t)的时不变表达形式，基函数展开法被用于动态交互力的表达，对动态交互力fv(t)进行基函数展开如下式所示：
Figure BDA0003163790340000217
其中，fv(t)为动态交互力，hj(t)为扰动持续区间内的基函数序列，cj代表第j阶基函数的权重，Q为所使用基函数的最大阶数，t为时间。由于动态交互力具有平滑的特征，B样条(B-spline)基函数由于能良好地描述平滑变化量，被选作基函数hj(t)。基于对动态交互力辨识精度的需求，需要使用较高阶数的基函数，但过高的阶数会导致未知参数数量的增加反而降低整体辨识结果。本实施例将B样条基函数的最大阶数Q设为30，此阶数设置能大体保证被描述变量的整体变化规律以及部分时间区间内的细微变化能被良好的估计。
(3.3)上肢末端时变动力学系统的参数迭代辨识
基于以上对动态交互力fv(t)与人体回复力fr(t)的描述，用逐步迭代的辨识方法从上肢末端实测力数据fm(t)中对这两部分力进行解耦，迭代辨识算法步骤如下：
1.初始化，fm(t)、x(t)、
Figure BDA0003163790340000218
和
Figure BDA0003163790340000219
已知，令fr(t)的初始估计值为
Figure BDA00031637903400002110
2.对第n次迭代，计算中间变量
Figure BDA0003163790340000221
根据公式(7)，使
Figure BDA0003163790340000222
从而辨识出权重cj，并重建得到fv(t)的估计值
Figure BDA0003163790340000223
3.计算中间变量
Figure BDA0003163790340000224
根据公式(5)、公式(6)，使
Figure BDA0003163790340000225
从而辨识出参数β，并重建得到fr(t)的估计值
Figure BDA0003163790340000226
4.令n＝n+1，重新回到步骤2，直至n超过设定的最大迭代次数时迭代终止。
其中，fm(t)为实测力，fv(t)为动态交互力，fr(t)为回复力，
Figure BDA0003163790340000227
为回复力的初值，
Figure BDA0003163790340000228
为中间变量，
Figure BDA0003163790340000229
为第n-1次迭代的回复力的估计值，
Figure BDA00031637903400002210
为中间变量，β是由人体上肢末端的惯性参数Iend、第0…P阶函数的权重a0…P、b0…P构成的矩阵，
Figure BDA00031637903400002211
是由变量xg(i,t)、
Figure BDA00031637903400002212
人体上肢末端位移的二阶导
Figure BDA00031637903400002213
构成的矩阵，
Figure BDA00031637903400002214
为第n次迭代的实测力的估计值，
Figure BDA00031637903400002215
为第n次迭代的动态交互力的估计值，
Figure BDA00031637903400002216
为第n次迭代回复力的估计值，i,j为函数的阶数，t为时间。
经过上述的迭代辨识算法，得到参数β中的各项参数、各项回复力fr(t)和动态交互力fv(t)，将上述各项数据回传至公式(3)，便能得到动态交互过程中的上肢末端时变刚度Kend(t)；通过对所获取的刚度矩阵进行奇异值分解得到如图4所示的上肢末端的刚度椭球示意图和如图5所示的上肢末端刚度的曲线图K(t)。
与现有技术相比，本发明存在以下优势：
本发明提供的面向动态交互的对象肢体时变刚度辨识方法，是基于线性变参数法、基函数展开法和参数迭代辨识的对象肢体时变刚度辨识方法构建的，通过线性变参数法将对象肢体时变刚度研究问题转换为时不变研究问题，为时变刚度的辨识提供了新的思路；提供的参数迭代辨识能够依靠单次扰动数据完整准确地辨识出对象肢体时变刚度，极大地减小了时变刚度研究的复杂性，更为重要地突破了现有依靠单次实验数据的时变刚度辨识技术无法应用于与环境动态交互情景的局限性。总体而言，本发明提供的对象肢体时变刚度辨识方法提高了时变刚度辨识技术的精度。
本发明提供的对象肢体时变刚度辨识的运动数据采集方法能够统一地指导对象肢体的肢体末端或关节的时变刚度研究的实验数据采集，实验范式简单、操作性高，能够完整地改变对象肢体的肢体末端或关节的刚度特性，且记录的数据能够精准地解算出对象肢体的时变刚度。
本发明利用多自由度协作型机械臂构建微扰动环境以进行对象肢体时变刚度测量的方法，省去了开发专用型测量装置的复杂过程，降低了对象肢体时变刚度测量的难度，并且拓展了利用协作型机械臂进行刚度测量的应用场景，使对象上、下肢的时变刚度测量能够借助协作型机械臂进行刚度研究。
本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。

Claims (10)

1.一种面向动态交互的对象肢体时变刚度辨识方法，其特征在于，包括以下步骤：
在随机微扰动运动环境下，对象保持肢体姿态不变，通过改变肌肉状态，采集对象肢体位移与实测力数据组合或角度与实测力矩数据组合；
基于替代时变阻抗参数的回复力/力矩表达式以及基函数展开的动态交互力/力矩的表达式，根据对象肢体位移与实测力数据组合或角度与实测力矩数据组合，采用时变动力学参数对动态交互力/力矩和回复力/力矩进行迭代辨识；
利用参数迭代辨识过程中得到的线性系数矩阵及动态交互力/力矩，依照替代时变阻抗参数的回复力/力矩表达式解算对象肢体时变刚度；
其中，替代时变阻抗参数的回复力/力矩的表达式获取方法，包括以下步骤：
根据对象肢体时变阻抗特性，利用二阶阻抗模型，构建对象肢体的时变动力学模型；
利用线性变参数法，将时变动力学模型中对象肢体的时变阻尼参数以及对象肢体的时变刚度参数构建为以调度变量为自变量的函数加权和，对回复力/力矩进行变形转换，获取替代时变阻抗参数的回复力/力矩的表达式；其中，动态交互力/力矩为调度变量；
采用基函数展开的动态交互力/力矩的表达式的获取方法为：
将动态交互力/力矩进行基函数展开，采用扰动持续区间内的基函数序列以及基函数的权重表示动态交互力/力矩。
2.根据权利要求1所述的对象肢体时变刚度辨识方法，其特征在于，采用时变动力学参数对动态交互力/力矩和回复力/力矩进行迭代辨识的方法，包括以下步骤：
(1)采集的对象肢体位移与实测力数据组合或角度与实测力矩数据组合，并初始化回复力/力矩的估计值为0；
(2)基于对象肢体位移与实测力数据组合或角度与实测力矩数据组合，结合当前迭代下的回复力/力矩的估计值，根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取动态交互力/力矩的中间变量；
(3)根据基函数展开的动态交互力/力矩的表达式，采用动态交互力/力矩的中间变量辨识基函数的权重；
(4)利用基函数的权重以及扰动持续区间内的基函数序列，获取当前迭代下的动态交互力/力矩的估计值；
(5)根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取回复力/力矩的中间变量；
(6)基于替代时变阻抗参数的回复力/力矩的表达式，采用回复力/力矩的中间变量，获取线性系数矩阵；
(7)根据位移参数矩阵以及当前迭代下获取的线性系数矩阵，获取当前迭代下的回复力/力矩的估计值；
(8)转至步骤(2)，直至实际迭代次数达到预定的最大迭代次数，输出最后一次迭代下获取的线性系数矩阵、回复力/力矩的估计值和动态交互力/力矩的估计值。
3.根据权利要求1或2所述的对象肢体时变刚度辨识方法，其特征在于，采集对象肢体位移与实测力数据组合或角度与实测力矩数据组合的方法为：
在随机微扰动运动环境下，对象肢体保持姿态不变，将对象肢体的肢体末端或关节与协作型机械臂连接；
在扰动时间内，对象肢体完成至少一个动态力/力矩交互周期，采集对象肢体位移与力数据组合或角度与力矩数据组合；
其中，一个动态力/力矩为在协作型机械臂的X、Y和Z三个方向/旋转轴均按先增加正向的交互力/力矩后减小为0，再增加负向的交互力/力矩后再减小为0的规律施加交互力/力矩。
4.根据权利要求1所述的对象肢体时变刚度辨识方法，其特征在于，替代时变阻抗参数的回复力/力矩的表达式为：
Figure FDA0003163790330000021
Figure FDA0003163790330000022
β＝[I b0 b1 … bP a0 a1 … aP]T
其中，I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，a0…P、b0…P为第0…P阶函数的权重，P为所用函数最大阶数，fr(t)是对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度,β是由对象肢体的惯性参数I、第0…P阶函数的权重a0…P、b0…P构成的线性系数矩阵，
Figure FDA0003163790330000023
是由变量xg(i,t)、
Figure FDA0003163790330000024
对象肢体位移/角度的二阶导
Figure FDA0003163790330000025
构成的位移参数矩阵，xg(i,t)是由函数gi(fv(t))和对象肢体位移/角度x(t)构建的中间变量，
Figure FDA0003163790330000026
是由函数gi(fv(t))和对象肢体位移/角度的一阶导
Figure FDA0003163790330000027
构建的中间变量，i是gi函数的阶数，t为时间。
5.根据权利要求1所述的对象肢体时变刚度辨识方法，其特征在于，采用基函数展开的动态交互力/力矩的表达式为：
Figure FDA0003163790330000031
其中，fv(t)为动态交互力，hj(t)为扰动持续区间内的基函数序列，cj代表第j阶基函数的权重，Q为所使用基函数的最大阶数，t为时间。
6.一种面向动态交互的对象肢体时变刚度辨识装置，其特征在于，包括：数据采集系统、回复力/力矩表达式获取系统、参数迭代辨识系统和动态交互力/力矩表达式获取系统；
数据采集系统用于在随机微扰动运动环境下，对象保持肢体姿态不变，通过改变肌肉状态，采集对象肢体位移与实测力数据组合或角度与实测力矩数据组合；
回复力/力矩表达式获取系统用于根据对象肢体的时变动力学系统，使用线性变参数法替代对象肢体时变阻抗参数重建回复力/力矩的表达；
参数迭代辨识系统用于基于替代时变参数的回复力/力矩的表达式以及基函数展开的动态交互力/力矩的表达式，根据对象肢体位移/角度信息以及实测力/力矩信息，采用时变动力学的参数对动态交互力/力矩和回复力/力矩进行迭代辨识；
动态交互力/力矩表达式获取系统用于将动态交互力/力矩进行基函数展开，采用扰动持续区间内的基函数序列以及基函数的权重表示动态交互力/力矩。
7.根据权利要求6所述的对象肢体时变刚度辨识装置，其特征在于，所述数据采集系统包括协作型机械臂、反光标记点、光学运动捕捉系统、六维力/力矩传感器和显示器；
协作型机械臂的工具端固定设置有反光标记点、且连接有六维力/力矩传感器，使用时通过连接装置与肢体末端/关节固定连接；
光学运动捕捉系统用于收集反光标记点的三维位置，获取对象肢体的位移/角度信息；
反光标记点用于记录肢体姿态；
协作型机械臂用于构建三维空间随机微扰动运动环境；
六维力/力矩传感器用于采集对象肢体与机械臂工具端之间的实测力/力矩信息；
显示器为对象提供动态交互过程中六维力/力矩传感器实时数据的视觉反馈。
8.根据权利要求6所述的对象肢体时变刚度辨识装置，其特征在于，回复力/力矩表达式获取系统获取替代时变阻抗参数的回复力/力矩表达式的方法，包括以下步骤：
根据对象肢体时变阻抗特性，利用二阶阻抗模型，构建对象肢体的时变动力学系统；
利用线性变参数法，将动力学模型中对象肢体的时变阻尼参数以及对象肢体的时变刚度参数构建为以调度变量为自变量的函数加权和，对回复力/力矩进行变形转换，替代时变阻抗参数；其中，动态交互力/力矩为调度变量。
9.根据权利要求6或8所述的对象肢体时变刚度辨识装置，其特征在于，所述替代时变阻抗参数的回复力/力矩的表达式为：
Figure FDA0003163790330000041
Figure FDA0003163790330000042
β＝[I b0 b1 … bP a0 a1 … aP]T
其中，I是对象肢体的惯性参数，gi(fv(t))为以动态交互力/力矩fv(t)为自变量的函数，a0…P、b0…P为第0…P阶函数的权重，P为所用函数最大阶数，fr(t)是对象肢体的回复力/力矩，x(t)是由光学运动捕捉系统采集的对象肢体位移/角度,β是由对象肢体的惯性参数I、第0…P阶函数的权重a0…P、b0…P构成的线性系数矩阵，
Figure FDA0003163790330000043
是由变量xg(i,t)、
Figure FDA0003163790330000044
对象肢体位移/角度的二阶导
Figure FDA0003163790330000045
构成的位移参数矩阵，xg(i,t)是由函数gi(fv(t))和对象肢体位移/角度x(t)构建的中间变量，
Figure FDA0003163790330000046
是由函数gi(fv(t))和对象肢体位移/角度的一阶导
Figure FDA0003163790330000047
构建的中间变量，i是gi函数的阶数，t为时间；
基函数展开的动态交互力/力矩的表达式为：
Figure FDA0003163790330000048
其中，fv(t)为动态交互力，hj(t)为扰动持续区间内的基函数序列，cj代表第j阶基函数的权重，Q为所使用基函数的最大阶数，t为时间。
10.根据权利要求6所述的对象肢体时变刚度辨识装置，其特征在于，参数迭代辨识系统辨识动态交互力/力矩和回复力/力矩的方法，包括以下步骤：
(1)采集的对象肢体位移与实测力数据组合或角度与实测力矩数据组合，并初始化回复力/力矩的估计值为0；
(2)基于对象肢体位移与实测力数据组合或角度与实测力矩数据组合，结合当前迭代下的回复力/力矩的估计值，根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取动态交互力/力矩的中间变量；
(3)根据基函数展开的动态交互力/力矩的表达式，采用动态交互力/力矩的中间变量辨识基函数的权重；
(4)利用基函数的权重以及扰动持续区间内的基函数序列，获取当前迭代下的动态交互力/力矩的估计值；
(5)根据实测力/力矩等于回复力/力矩与动态交互力/力矩之和，获取回复力/力矩的中间变量；
(6)基于替代时变阻抗参数的回复力/力矩的表达式，采用回复力/力矩的中间变量，获取线性系数矩阵；
(7)根据位移参数矩阵以及当前迭代下获取的线性系数矩阵，获取当前迭代下的回复力/力矩的估计值；
(8)转至步骤(2)，直至实际迭代次数达到预定的最大迭代次数，输出最后一次迭代下获取的线性系数矩阵、回复力/力矩的估计值和动态交互力/力矩的估计值。
CN202110798706.7A 2021-07-15 2021-07-15 一种面向动态交互的对象肢体时变刚度辨识方法及装置 Active CN113524235B (zh)

Priority Applications (2)

Application Number Priority Date Filing Date Title
CN202110798706.7A CN113524235B (zh) 2021-07-15 2021-07-15 一种面向动态交互的对象肢体时变刚度辨识方法及装置
US17/858,028 US11696687B2 (en) 2021-07-15 2022-07-05 Dynamic interaction-oriented subject's limb time-varying stiffness identification method and device

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202110798706.7A CN113524235B (zh) 2021-07-15 2021-07-15 一种面向动态交互的对象肢体时变刚度辨识方法及装置

Publications (2)

Publication Number Publication Date
CN113524235A CN113524235A (zh) 2021-10-22
CN113524235B true CN113524235B (zh) 2022-03-11

Family

ID=78099305

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202110798706.7A Active CN113524235B (zh) 2021-07-15 2021-07-15 一种面向动态交互的对象肢体时变刚度辨识方法及装置

Country Status (2)

Country Link
US (1) US11696687B2 (zh)
CN (1) CN113524235B (zh)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN117398085B (zh) * 2023-12-14 2024-02-27 合肥合滨智能机器人有限公司 基于力扰动的手臂阻抗计算方法及超声扫查技能采集系统

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102252836A (zh) * 2011-04-14 2011-11-23 唐德尧 旋转机器支承状态、动平衡及不对中的测试方法和装置
CN104698837A (zh) * 2014-12-11 2015-06-10 华侨大学 一种时变线性结构工作模态参数识别方法、装置及应用
CN106726341A (zh) * 2016-12-06 2017-05-31 西安交通大学 基于人机阻抗匹配模型的变刚度肘关节康复机器人及其控制方法
CN107671860A (zh) * 2017-10-27 2018-02-09 安徽工程大学 一种多运动模式变刚度约束机械臂振动控制系统和方法
CN112244833A (zh) * 2020-09-27 2021-01-22 华中科技大学 一种基于协作机械臂的人体上肢多维末端刚度测量方法

Family Cites Families (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US8928161B2 (en) * 2012-05-25 2015-01-06 Kcf Technologies, Inc. Apparatuses and methods for harvesting energy from prosthetic limbs

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN102252836A (zh) * 2011-04-14 2011-11-23 唐德尧 旋转机器支承状态、动平衡及不对中的测试方法和装置
CN104698837A (zh) * 2014-12-11 2015-06-10 华侨大学 一种时变线性结构工作模态参数识别方法、装置及应用
CN106726341A (zh) * 2016-12-06 2017-05-31 西安交通大学 基于人机阻抗匹配模型的变刚度肘关节康复机器人及其控制方法
CN107671860A (zh) * 2017-10-27 2018-02-09 安徽工程大学 一种多运动模式变刚度约束机械臂振动控制系统和方法
CN112244833A (zh) * 2020-09-27 2021-01-22 华中科技大学 一种基于协作机械臂的人体上肢多维末端刚度测量方法

Also Published As

Publication number Publication date
US20230028793A1 (en) 2023-01-26
CN113524235A (zh) 2021-10-22
US11696687B2 (en) 2023-07-11

Similar Documents

Publication Publication Date Title
CN107378944B (zh) 一种基于主成分分析法的多维表面肌电信号假手控制方法
Lapham et al. The use of artificial intelligence in the analysis of sports performance: A review of applications in human gait analysis and future directions for sports biomechanics
CN110710970B (zh) 肢体动作的识别方法、装置、计算机设备及存储介质
Tucker et al. Human preference-based learning for high-dimensional optimization of exoskeleton walking gaits
Wu et al. Using reinforcement learning to estimate human joint moments from electromyography or joint kinematics: An alternative solution to musculoskeletal-based biomechanics
CN116507276A (zh) 用于机器学习以从图像分析肌肉骨骼康复的方法和设备
Mattioli et al. Classification of EMG signals using artificial neural networks for virtual hand prosthesis control
CN113524235B (zh) 一种面向动态交互的对象肢体时变刚度辨识方法及装置
Tatinati et al. Multidimensional modeling of physiological tremor for active compensation in handheld surgical robotics
Liu et al. A new IMMU-based data glove for hand motion capture with optimized sensor layout
CN113143257A (zh) 基于个体运动行为层次模型的泛化应用系统及方法
JP4078419B2 (ja) ポインティングデバイス
Rezzoug et al. Prediction of fingers posture using artificial neural networks
KR100994408B1 (ko) 손가락 힘 추정 방법 및 추정 장치, 손가락 힘 추정을 위한근육 판별 방법 및 근육 판별 장치
Li et al. Continuous motion and time-varying stiffness estimation of the human elbow joint based on sEMG
Chalfoun et al. Forces, activation and displacement prediction during free movement in the hand and forearm
Allouch et al. Proposition, identification, and experimental evaluation of an inverse dynamic neuromusculoskeletal model for the human finger
Sitole et al. Continuous Prediction of Human Joint Mechanics using EMG Signals: A Review of Model-based and Model-free Approaches
Artal-Sevil et al. Control of a bionic hand based on neural networks and improved gesture recognition techniques using multiple EMG sensors
Saputra et al. A real-time control system of upper-limb human musculoskeletal model with environmental integration
Lueken et al. Using synthesized imu data to train a long-short term memory-based neural network for unobtrusive gait analysis with a sparse sensor setup
Li et al. A novel method for estimating continuous motion and time-varying stiffness of human elbow joint
Guo et al. Patient’s Healthy-Limb Motion Characteristic-Based Assist-As-Needed Control Strategy for Upper-Limb Rehabilitation Robots
CN111452022B (zh) 基于贝叶斯优化的上肢康复机器人主动训练参考轨迹复杂度调节方法
Starke et al. Temporal force synergies in human grasping

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant