WO2017031718A1 - 弹性物体变形运动的建模方法 - Google Patents

弹性物体变形运动的建模方法 Download PDF

Info

Publication number
WO2017031718A1
WO2017031718A1 PCT/CN2015/088126 CN2015088126W WO2017031718A1 WO 2017031718 A1 WO2017031718 A1 WO 2017031718A1 CN 2015088126 W CN2015088126 W CN 2015088126W WO 2017031718 A1 WO2017031718 A1 WO 2017031718A1
Authority
WO
WIPO (PCT)
Prior art keywords
simulated
mesh model
elastic object
deformation
tetrahedral mesh
Prior art date
Application number
PCT/CN2015/088126
Other languages
English (en)
French (fr)
Inventor
王滨
黄惠
伍龙华
Original Assignee
中国科学院深圳先进技术研究院
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 中国科学院深圳先进技术研究院 filed Critical 中国科学院深圳先进技术研究院
Priority to PCT/CN2015/088126 priority Critical patent/WO2017031718A1/zh
Priority to US15/392,688 priority patent/US10394979B2/en
Publication of WO2017031718A1 publication Critical patent/WO2017031718A1/zh

Links

Images

Classifications

    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T13/00—Animation
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/11—Complex mathematical operations for solving equations, e.g. nonlinear equations, general mathematical optimization problems
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10—Complex mathematical operations
    • G06F17/16—Matrix or vector computation, e.g. matrix-matrix or matrix-vector multiplication, matrix factorization
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00—Computer-aided design [CAD]
    • G06F30/20—Design optimisation, verification or simulation
    • G06F30/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06F—ELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00—Details relating to CAD techniques
    • G06F2111/10—Numerical modelling

Definitions

  • the invention relates to the field of simulation modeling technology, and in particular to a modeling method for deformation motion of an elastic object.
  • the traditional modeling method is mainly to apply external force with known strength at different positions of the object by means of the force sensing device, and then obtain the shape change of the object under stable conditions under different external force conditions through the position tracking device, and finally use machine learning, probability statistics
  • the theoretical method is used to establish the relationship between stress and strain, so as to obtain the mathematical physical model of the object to be measured for simulation.
  • the motion model built by applying the existing elastic object motion modeling method is not realistic enough.
  • Embodiments of the present invention provide a modeling method for deformation motion of an elastic object, which is used to establish a realistic deformation model of an elastic object, and the method includes:
  • a simulated tetrahedral mesh model for simulating the deformation motion of the elastic object is established
  • each iteration cycle performs the following operations: obtaining the reference shape of the simulated tetrahedral mesh model corresponding to the current material property coefficient of the elastic object; According to the current elastic object material property coefficient and the reference shape of the corresponding simulated tetrahedral mesh model, the simulated tetrahedral mesh model is driven to simulate the elastic object deformation motion from the same initial deformation, and the simulated deformation of the simulated tetrahedral mesh model is obtained.
  • Motion sequence calculating the position of the simulated deformation motion sequence and the tracking deformation motion sequence Deviating; updating the material property coefficient of the elastic object along a direction that reduces the positional deviation; until the material property coefficient at the minimum position deviation and the reference shape of the corresponding simulated tetrahedral mesh model are found;
  • the elastic object deformation motion model is established according to the material attribute coefficient when the position deviation is minimum and the reference shape of the corresponding simulated tetrahedral mesh model.
  • a simulated tetrahedral mesh model for simulating the deformation motion of the elastic object is established according to the static point cloud, including:
  • a simulated tetrahedral mesh model for simulating the deformation motion of the elastic object is established.
  • Each vertex of the static surface mesh model is linearly interpolated with the center of gravity coordinates of each tetrahedral space of the simulated tetrahedral mesh model.
  • the driving simulation tetrahedral mesh model tracks the dynamic point cloud sequence, and the tracking deformation motion sequence of the simulated tetrahedral mesh model is obtained, including:
  • a virtual external force is applied to each node of the simulated tetrahedral mesh model to drive each node of the simulated tetrahedral mesh model to track the corresponding dynamic point cloud sequence, and the deformed simulated tetrahedral mesh model is obtained.
  • the difference is less than the preset value, including: the distance between each vertex of the static surface mesh model and the dynamic point cloud sequence corresponding to the deformation of the elastic object is less than a preset distance, or the attraction between them is less than Preset attraction.
  • a reference shape of the simulated tetrahedral mesh model corresponding to the current material property coefficient of the elastic object is obtained, including:
  • the initial value estimate of the material property coefficient of the elastic object includes:
  • the initial value of the material property coefficient of the elastic object is determined according to the matching degree of the vibration frequency corresponding to the first main mode and the vibration frequency of the actual collected data.
  • the following target equation is used to solve the shape corresponding to the minimum residual force residual as the reference shape of the simulated tetrahedral mesh model:
  • R is the rotation matrix
  • K is the stiffness matrix
  • x s is the static shape of the elastic object
  • X is the reference shape of the simulated tetrahedral network model
  • M is the mass of the elastic object
  • g is the gravitational acceleration
  • the Jacobian matrix is:
  • f is the virtual elastic force applied to each node of the simulated tetrahedral mesh model
  • X ij is the position of the i-th node in the j-direction of the reference shape of the simulated tetrahedral mesh model
  • x s is the elastic object
  • the static shape X is the reference shape of the simulated tetrahedral mesh model
  • V is the volume of each tetrahedral element of the simulated tetrahedral mesh model
  • E is the constant matrix associated with the material properties of the elastic object
  • R is the simulation A rigid rotation matrix for each tetrahedral element of a tetrahedral mesh model.
  • the material property coefficients of the elastic object are estimated, including:
  • the material attribute coefficients of other nodes of the simulated tetrahedral mesh model are obtained.
  • the embodiment of the present invention relies on the data-driven method to track and analyze the point cloud sequence of the elastic object deformation motion, specifically First, a static point cloud and a dynamic point cloud sequence of the elastic object are collected; secondly, a simulated tetrahedral mesh model for simulating the deformation motion of the elastic object is established according to the static point cloud; The simulated tetrahedral mesh model is tracked to track the dynamic point cloud sequence, and the tracking deformation motion sequence of the simulated tetrahedral mesh model is obtained; then, the material property coefficient of the elastic object and the corresponding simulated tetrahedral mesh are iterative
  • the reference shape of the model the following operations are performed in each iteration cycle: obtaining the reference shape of the simulated tetrahedral mesh model corresponding to the current material property coefficient of the elastic object; according to the current elastic object material property coefficient and the corresponding simulated tetrahedral mesh model
  • the reference shape drives the simulated tetrahedral mesh model to simulate the deformation motion of the elastic object from the same initial deformation, and obtains a simulated deformation motion sequence of the simulated tetrahedral mesh model; and calculates the simulated deformation motion sequence and tracking Position deviation of the deformation motion sequence; updating the material attribute coefficient of the elastic object along a direction in which the position deviation is reduced; until the material attribute coefficient of the minimum position deviation and the reference shape of the corresponding simulated tetrahedral mesh model are found; finally, according to Material attribute coefficient and corresponding when the position deviation is the smallest Reference tetrahedral shape simulation model, elastic deformation of the object to establish a motion model.
  • FIG. 1 is a schematic flow chart of a modeling method for deformation motion of an elastic body in an embodiment of the present invention
  • FIG. 2 is a schematic diagram of a point cloud sequence of an elastic object collected in an embodiment of the present invention
  • FIG. 3 is a schematic diagram of a static surface mesh model of an elastic object established in an embodiment of the present invention.
  • FIG. 4 is a schematic diagram of a simulated tetrahedral mesh model in an embodiment of the present invention.
  • Figure 5 is a schematic view showing the position of the plant model control point placed in the axial direction in the embodiment of the present invention
  • Fig. 6 is a schematic view showing the modeling of the deformation motion of the elastic body in one embodiment of the present invention.
  • FIG. 1 is a schematic flow chart of a modeling method for deformation motion of an elastic object according to an embodiment of the present invention. As shown in FIG. 1, the method includes the following steps:
  • Step 101 Collect a static point cloud of the elastic object and a dynamic point cloud sequence during the deformation motion
  • Step 102 Establish a simulated tetrahedral mesh model for simulating the deformation motion of the elastic object according to the static point cloud;
  • Step 103 Driving the simulated tetrahedral mesh model to track the dynamic point cloud sequence, and obtaining a tracking deformation motion sequence of the simulated tetrahedral mesh model;
  • Step 104 Iteratively estimating the material attribute coefficient of the elastic object and the reference shape of the corresponding simulated tetrahedral mesh model; performing the following operations in each iteration period: obtaining the simulated tetrahedral mesh model corresponding to the current material attribute coefficient of the elastic object Reference shape; according to the current elastic object material property coefficient and the corresponding reference shape of the simulated tetrahedral mesh model, the simulated tetrahedral mesh model is driven to simulate the elastic object deformation motion from the same initial deformation, and the simulated tetrahedral mesh model is obtained.
  • the simulated deformation motion sequence calculating the positional deviation of the simulated deformation motion sequence and the tracking deformation motion sequence; updating the material property coefficient of the elastic object along the direction in which the position deviation is reduced; until the material attribute coefficient and the corresponding value when the position deviation is found to be the smallest Simulating the reference shape of the tetrahedral mesh model;
  • Step 105 Establish a deformation motion model of the elastic object according to the material attribute coefficient when the position deviation is minimum and the reference shape of the corresponding simulated tetrahedral mesh model.
  • the modeling method of the deformation motion of the elastic object firstly collects the static point cloud and the dynamic point cloud sequence of the elastic object; secondly, according to the static point cloud, the simulation tetrahedron network for simulating the deformation motion of the elastic object is established. Then, the driving simulation tetrahedral mesh model is used to track the dynamic point cloud sequence, and the tracking deformation motion sequence of the simulated tetrahedral mesh model is obtained. Then, the material property coefficient of the elastic object and the corresponding simulated tetrahedral mesh model are iteratively estimated.
  • each iteration cycle performs the following operations: obtaining a reference shape of the simulated tetrahedral mesh model corresponding to the current material property coefficient of the elastic object; according to the current elastic object material property coefficient and the corresponding simulated tetrahedral mesh model
  • the reference shape drives the simulated tetrahedral mesh model to simulate the deformation motion of the elastic object from the same initial deformation, and obtains the simulated deformation motion sequence of the simulated tetrahedral mesh model; calculates the positional deviation of the simulated deformation motion sequence and the tracking deformation motion sequence; More along the direction that reduces the positional deviation
  • the reference shape establishes a deformation motion model of the elastic object.
  • Step 101 we use three Kinect combinations to collect the static point cloud sequence of the elastic object and the dynamic point cloud sequence during the deformation motion, and FIG. 2 is a schematic diagram of the point cloud sequence of the collected elastic object.
  • Step 101 is a process of collecting a static point cloud and a dynamic point cloud sequence of an elastic object.
  • a simulated tetrahedral mesh model for simulating the deformation motion of the elastic object is established, which may include:
  • a static surface mesh model of the elastic object is established; since the static surface mesh model is very accurate, it can also be called a static fine surface mesh model, and FIG. 3 is a static surface mesh model;
  • FIG. 4 is a schematic diagram of the simulated tetrahedral mesh model
  • Each vertex of the static surface mesh model is linearly interpolated with the center of gravity coordinates of each tetrahedral space of the simulated tetrahedral mesh model.
  • the static surface mesh model in Figure 3 has 15,368 vertices, corresponding to the simulated tetrahedral mesh model, with 9,594 nodes.
  • This step 102 is a process of establishing a simulated tetrahedral mesh model for simulating the deformation motion of an elastic object.
  • the mesh of the fine surface mesh model is then passed to the volume data generation tool TETGEN, which derives a relatively coarse tetrahedral mesh for physical simulation, which is used as a template to track the point cloud sequence (as shown in Figure 4).
  • TETGEN volume data generation tool
  • an embedding strategy is adopted here, and the vertices on the fine mesh can be represented by the coordinates of the center of gravity of the tetrahedral space, and the two are linear interpolation relationships.
  • driving the simulated tetrahedral mesh model to track the dynamic point cloud sequence, and obtaining the tracking deformation motion sequence of the simulated tetrahedral mesh model may include the following steps:
  • a virtual external force is applied to each node of the simulated tetrahedral mesh model to drive each node of the simulated tetrahedral mesh model to track the corresponding dynamic point cloud sequence, and the deformed simulated tetrahedral mesh model is obtained.
  • the difference is less than the preset value, which may include: the distance between each vertex of the static surface mesh model and the dynamic point cloud sequence corresponding to the deformation of the elastic object is less than a preset distance, or an attraction between them Less than the preset attraction.
  • the above step 103 is a process of driving the simulation tetrahedral mesh model deformation tracking elastic body deformation motion.
  • the step 103 is a physics-based probability tracking method, and the motion tracking needs to process the point cloud data with noise, and also considers the problems of occlusion, rapid motion, and large deformation. Therefore, we transform the deformation motion tracking into a maximum posterior probability (MAP) problem, and use the expected maximum method (EM) for iterative solution: (E step) to find the optimal correspondence according to the current point cloud and node position; (M step The mobile node position is such that the above correspondence is the maximum likelihood estimation.
  • MAP maximum posterior probability
  • EM expected maximum method
  • Our task is to find the shape matching s based on a given point cloud c.
  • the mesh deformation tracking matching point cloud can be expressed as a maximum a posteriori probability problem:
  • the EM algorithm is used for solving.
  • p refers to the probability distribution
  • the second term of the above formula reflects the potential energy of the deformed object model, so it can be optimized by physical simulation. We add virtual power to each node:
  • ⁇ is a proportional coefficient of the virtual force
  • the external force drives the mesh deformation to match the point cloud shape.
  • equation of motion is:
  • a co-rotating linear finite element model is used to simulate the deformed object
  • M is the mass matrix
  • ⁇ and ⁇ are the two coefficients corresponding to Rayleigh damping
  • R is the rotation.
  • the matrix is obtained by polar decomposition of the deformation gradient
  • K is the stiffness matrix
  • x is the deformed shape of the simulated tetrahedral mesh model
  • X is the reference shape of the tetrahedral mesh model
  • fext is the external force.
  • the physical simulation part is packaged in a third-party library VEGA FEM. In order to speed up the calculation, a nesting strategy is adopted here.
  • Each iterative calculation maps the external force to a tracking tetrahedral model with a small number of nodes to simulate the tracking motion, and then interpolates the node displacement back to the vertex of the static surface mesh model.
  • the motion tracking process is to apply an external force simulation until the EM iteration converges.
  • the traditional technical scheme generally uses a non-rigid registration algorithm to match the template grid and the point cloud sequence frame by frame. This method firstly does not have a fast calculation speed; secondly, it has high requirements for point cloud data. It is impossible to deal with large deformation motions and situations with more noise interference items; the matching results again will result in a mesh topology because there is no fusion physical constraint.
  • step 103 the technical solution of the present application uses a probability-based tracking algorithm to find the probability correspondence between the mesh vertex and the point cloud, and drives the mesh deformation motion driven by the physics engine, and the matching result fusion physical Constraints (satisfying the constraints of the above equations), the mesh topology will not be out of shape, with faster running speed, better tracking effect, and stronger robustness.
  • the above step 103 can make the point cloud data of the elastic object deformation motion model have small noise, and can also solve the problems of occlusion, rapid motion and large deformation.
  • a reference shape of the simulated tetrahedral mesh model corresponding to the current material attribute coefficient of the elastic object is obtained, including:
  • Shape X Where E is the Young's modulus, ⁇ is the Poisson's coefficient, and ⁇ and ⁇ are the two coefficients corresponding to Rayleigh damping.
  • the solution problem can be expressed as a spatiotemporal optimization problem.
  • the following objective equation F measures the positional deviation of the simulation and tracking sequences:
  • the following target equation is used to solve the shape corresponding to the minimum residual force residual as the reference shape of the simulated tetrahedral mesh model:
  • R is the rotation matrix
  • K is the stiffness matrix
  • x s is the static shape of the elastic object
  • X is the reference shape of the simulated tetrahedral network model
  • M is the mass of the elastic object
  • g is the gravitational acceleration
  • the reference shape and the static shape of the object model should be distinguished.
  • the reference shape should be unaffected by the gravity factor. Otherwise, the most obvious distortion phenomenon is when the object starts from the static shape.
  • a significant shape changes Take the deformation movement of the plant. For example, the leaves of the plant are affected by gravity. There is a slight movement that pulls down first. Then the simulation movement is inevitably unreal. Because the inventor considered this problem and solved the above objective equation.
  • the resultant residual force is applied to the simulated tetrahedral network model to obtain the reference shape of the simulated tetrahedral network model. Applying this method, a simulation model of the deformation motion of the elastic body with strong realism can be generated.
  • the Jacobian matrix is:
  • f is the virtual elastic force applied to each node of the simulated tetrahedral mesh model
  • X ij is the position of the i-th node in the j-direction of the reference tetrahedral mesh model reference shape
  • x s is the elastic object Static shape
  • X is the reference shape of the simulated tetrahedral mesh model
  • V is the volume of each tetrahedral element of the simulated tetrahedral mesh model
  • E is the constant matrix related to the material properties of the elastic object
  • R is the simulated four sides.
  • V is the volume of the tetrahedral element and E is a 6 ⁇ 6 constant matrix, which is related to the elastic properties of the material.
  • F is a matrix
  • the 6 ⁇ 12 matrix B depends only on X, and the interior is composed of B m , where B ij is used to represent the elements of the i-th row and the j-th column in B m :
  • k q represents the qth column of K, which is the partial differential of the qth column of X, and is actually the component of the i-th term of node j, which is represented by X ij here.
  • the material property coefficients of the elastic object are estimated, including:
  • the material attribute coefficients of other tetrahedral elements of the simulated tetrahedral mesh model are obtained.
  • t is the number of the frame
  • k is the node number.
  • F(p) we use a gradient-free simplex optimization method.
  • control points we place the control points in the axis direction and the interpolation weights are represented by the normalized axis distance.
  • the control points are assumed to be assigned to different locations according to the assumption, and the interpolation weights are controlled by the dual harmonic function.
  • the initial value estimate of the material property coefficient of the elastic object includes:
  • the initial value of the material property coefficient of the elastic object is determined according to the matching degree of the vibration frequency corresponding to the first main mode and the vibration frequency of the actual collected data.
  • main mode refers to the vibration mode corresponding to the minimum feature value.
  • the objective function F(p) since the objective function F(p) usually contains a plurality of local minimum values. Therefore, it is important to provide a suitable input parameter for the final successful optimization solution, that is, to give an optimal value when estimating the material attribute coefficient for the first time, which is beneficial to subsequent rapid calculation and improve the efficiency of modeling.
  • a novel strategy that uses modal analysis and coordinate descent methods to obtain appropriate initial material property coefficient values.
  • Young's modulus E affects vibration frequency during eigenvalue decomposition Intuitively understand that the softer the material, the smaller the vibration frequency of the object.
  • the coordinate descent method is used to sequentially update each material parameter. During each iteration, we search linearly based on the sensitivity of variable changes to changes in the value of the objective function.
  • FIG. 6 is a schematic diagram showing the modeling implementation of the deformation motion of the elastic body according to an embodiment of the present invention.
  • the static shape and the dynamic point cloud motion sequence of the elastic object are first collected; then the core of the system is The alternate iterative optimization strategy performs the deformation motion tracking and parameter estimation part in turn, and the results of each iteration can significantly improve the other part of the effect.
  • the system uses a physics-based probability tracking algorithm to drive the mesh deformation to register each frame point cloud, and output the position of each frame of the tetrahedral mesh.
  • the parameter estimation section simultaneously optimizes the estimated material property coefficient, damping coefficient, and object model reference shape.
  • the reference shape here refers to the original shape of the object model without any external force including gravity, and therefore the object reference shape and the static shape affected by gravity should be different.
  • the partitioning strategy is adopted. After the initial estimated physical parameters are given, the current model reference shape is solved according to the static balance equation, and then the forward simulation is performed using the set of materials and the reference shape data to obtain the same initial deformation condition. The sequence of motion; and the difference between the shape of the motion sequence and the tracking result is calculated as a criterion for evaluating the parameters of the group. The best value is found after multiple iterations; finally, the system generates a simulation model that enables realistic interaction.

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Mathematical Physics (AREA)
  • Data Mining & Analysis (AREA)
  • Computational Mathematics (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Pure & Applied Mathematics (AREA)
  • General Engineering & Computer Science (AREA)
  • Algebra (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • Computing Systems (AREA)
  • Operations Research (AREA)
  • Processing Or Creating Images (AREA)

Abstract

一种弹性物体变形运动的建模方法，该方法包括：采集弹性物体的静态点云和变形运动过程中的动态点云序列（101）；建立用于仿真弹性物体变形运动的仿真四面体网格模型（102）；驱动仿真四面体网格模型跟踪动态点云序列，得到跟踪变形运动序列（103）；迭代估计弹性物体的材质属性系数和对应的参考形状；每个迭代周期均执行以下操作：获得当前材质属性系数对应的参考形状；根据该系数和参考形状，驱使仿真四面体网格模型从相同的初始形变下变形运动，得到仿真变形运动序列；计算仿真变形运动序列与跟踪变形运动序列的位置偏差；沿着位置偏差减小的方向更新材质属性系数；直到找到位置偏差最小时材质属性系数和对应的仿真四面体网格模型的参考形状（104）；根据位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状，建立弹性物体变形运动模型（105）。通过上述建模方法，可以建立逼真的弹性物体变形运动模型。

Description

弹性物体变形运动的建模方法 技术领域
本发明涉及仿真建模技术领域，特别涉及一种弹性物体变形运动的建模方法。
背景技术
目前，近几十年来，计算机图形学技术得到了长足的发展，其研究成果在影视游戏、虚拟仿真、设计制造等领域都得到了广泛的应用。除了处理精细几何模型和渲染逼真效果，采用恰当物体数学物理模型，生成与控制真实感的运动是一个待继续深入研究的问题。
具有物理真实感的运动生成与控制是影视制作、动画领域不可或缺的重要组成部分。但是传统方法中数学模型的过度简化和参数不准确严重损害了仿真结果的精确性，阻碍了这一技术在实际工业领域的广泛应用。
传统的建模方法主要是：借助力传感设备在物体的不同位置施加强度已知的外力，然后通过位置跟踪设备，获得不同外力条件下稳定时物体的形状变化，最后利用机器学习、概率统计等理论方法建立应力与应变的关系曲线，从而获得被测量对象的数学物理模型进行仿真。然而，应用现有弹性物体运动建模方法所建运动模型不够逼真。
发明内容
本发明实施例提供了一种弹性物体变形运动的建模方法，用以建立逼真的弹性物体变形运动模型，该方法包括：
采集弹性物体的静态点云和变形运动过程中的动态点云序列；
根据静态点云，建立用于仿真弹性物体变形运动的仿真四面体网格模型；
驱动仿真四面体网格模型跟踪动态点云序列，得到仿真四面体网格模型的跟踪变形运动序列；
迭代估计弹性物体的材质属性系数和对应的仿真四面体网格模型的参考形状；每个迭代周期均执行以下操作：获得弹性物体的当前材质属性系数对应的仿真四面体网格模型的参考形状；根据当前弹性物体材质属性系数和对应的仿真四面体网格模型的参考形状，驱使仿真四面体网格模型从相同的初始形变下，仿真弹性物体变形运动，得到仿真四面体网格模型的仿真变形运动序列；计算仿真变形运动序列与跟踪变形运动序列的位 置偏差；沿着使得位置偏差减小的方向更新弹性物体的材质属性系数；直到找到位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状；
根据位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状，建立弹性物体变形运动模型。
在一个实施例中，根据静态点云，建立用于仿真弹性物体变形运动的仿真四面体网格模型，包括：
根据静态点云，建立弹性物体的静态表面网格模型；
根据静态表面网格模型，建立用于仿真弹性物体变形运动的仿真四面体网格模型；
静态表面网格模型的每个顶点与仿真四面体网格模型的每个四面体空间重心坐标是线性插值关系。
在一个实施例中，驱动仿真四面体网格模型跟踪动态点云序列，得到仿真四面体网格模型的跟踪变形运动序列，包括：
找到弹性物体变形后所有动态点云序列与仿真四面体网格模型的所有节点的最大概率对应关系；
根据最大概率对应关系，向仿真四面体网格模型的每个节点施加虚拟外力，驱使仿真四面体网格模型的每个节点跟踪对应的动态点云序列，得到变形后的仿真四面体网格模型的每个节点的位置；
根据线性插值关系和变形后的仿真四面体网格模型的每个节点的位置，找到变形后的静态表面网格模型的每个顶点位置；
计算变形后的静态表面网格模型的每个顶点位置与弹性物体变形后对应的动态点云序列位置之间的差异，直到差异小于预设值时，得到当前点云数据对应的仿真四面体网格模型变形后的跟踪变形运动序列。
在一个实施例中，差异小于预设值包括：静态表面网格模型的每个顶点与弹性物体变形后对应的动态点云序列之间的距离小于预设距离，或它们之间的吸引力小于预设吸引力。
在一个实施例中，获得弹性物体的当前材质属性系数对应的仿真四面体网格模型的参考形状，包括：
验证当前材质属性系数和对应的仿真四面体网格模型的参考形状是否满足物理静力平衡方程。
在一个实施例中，弹性物体的材质属性系数的初值估计，包括：
从弹性物体变形运动时的多个变形模态中选择第一个主要的模态；
确定第一个主要的模态对应的振动频率；
根据第一个主要的模态对应的振动频率与实际采集数据振动频率的匹配程度，确定弹性物体的材质属性系数的初值。
在一个实施例中，求解弹性物体的仿真四面体网格模型的参考形状时，采用如下目标方程求解合力残差最小时对应的形状作为仿真四面体网格模型的参考形状：
Figure PCTCN2015088126-appb-000001
其中，R是旋转矩阵，K是刚度矩阵，xs是弹性物体的静态形状，X是仿真四面体网络模型的参考形状，M为弹性物体的质量，g为重力加速度；
驱使仿真四面体网格模型变形时，计算施加给每个节点的虚拟弹性力相对于参考形状的雅克比矩阵，雅克比矩阵为：
Figure PCTCN2015088126-appb-000002
其中，f为施加给仿真四面体网格模型的每个节点的虚拟弹性力，Xij为仿真四面体网格模型的参考形状中第i个节点在j方向上的位置，xs是弹性物体的静态形状，X是仿真四面体网格模型的参考形状，V是仿真四面体网格模型的每个四面体元素的体积，E是与弹性物体的材质属性相关的常量矩阵，R是代表仿真四面体网格模型的每个四面体元素的刚性旋转矩阵。
在一个实施例中，估计弹性物体的材质属性系数，包括：
从仿真四面体网格模型的不同位置选出多个节点作为控制点；
根据弹性物体的材质分布属性，估计出不同位置的控制点的不同材质属性系数；
根据线性插值算法和不同位置的控制点的不同材质属性系数，求出仿真四面体网格模型的其他节点的材质属性系数。
与传统方法中，借助力传感设备在物体的不同位置施加强度已知的外力，然后通过位置跟踪设备获得不同外力条件下稳定时物体的形状变化，最后利用机器学习、概率统计等理论方法建立应力与应变的关系曲线，从而获得被测量对象的数学物理模型，进行建模的方法相比较，本发明实施例依赖数据驱动的方法对弹性物体变形运动的点云序列进行跟踪分析建模，具体地，首先，采集弹性物体的静态点云和动态点云序列；其次，根据所述静态点云，建立用于仿真弹性物体变形运动的仿真四面体网格模型；接着，驱 动所述仿真四面体网格模型跟踪所述动态点云序列，得到所述仿真四面体网格模型的跟踪变形运动序列；然后，迭代估计弹性物体的材质属性系数和对应的仿真四面体网格模型的参考形状；每个迭代周期均执行以下操作：获得弹性物体的当前材质属性系数对应的仿真四面体网格模型的参考形状；根据当前弹性物体材质属性系数和对应的仿真四面体网格模型的参考形状，驱使所述仿真四面体网格模型从相同的初始形变下，仿真弹性物体变形运动，得到所述仿真四面体网格模型的仿真变形运动序列；计算所述仿真变形运动序列与跟踪变形运动序列的位置偏差；沿着使得位置偏差减小的方向更新弹性物体的材质属性系数；直到找到位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状；最后，根据位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状，建立弹性物体变形运动模型。通过上述建模方法，可以建立逼真的弹性物体变形运动模型。
附图说明
为了更清楚地说明本发明实施例中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。在附图中：
图1是本发明实施例中弹性物体变形运动的建模方法的流程示意图；
图2是本发明实施例中采集到的弹性物体的点云序列的示意图；
图3是本发明实施例中建立的弹性物体静态表面网格模型的示意图；
图4是本发明实施例中仿真四面体网格模型的示意图；
图5是本发明实施例中对植物模型控制点安放在轴方向位置的示意图；
图6是本发明一个实施例中弹性物体变形运动的建模实施时的示意图。
具体实施方式
为使本发明实施例的目的、技术方案和优点更加清楚明白，下面结合附图对本发明实施例做进一步详细说明。在此，本发明的示意性实施例及其说明用于解释本发明，但并不作为对本发明的限定。
图1是本发明实施例中弹性物体变形运动的建模方法的流程示意图，如图1所示，该方法包括如下步骤：
步骤101：采集弹性物体的静态点云和变形运动过程中的动态点云序列；
步骤102：根据静态点云，建立用于仿真弹性物体变形运动的仿真四面体网格模型；
步骤103：驱动仿真四面体网格模型跟踪动态点云序列，得到仿真四面体网格模型的跟踪变形运动序列；
步骤104：迭代估计弹性物体的材质属性系数和对应的仿真四面体网格模型的参考形状；每个迭代周期均执行以下操作：获得弹性物体的当前材质属性系数对应的仿真四面体网格模型的参考形状；根据当前弹性物体材质属性系数和对应的仿真四面体网格模型的参考形状，驱使仿真四面体网格模型从相同的初始形变下，仿真弹性物体变形运动，得到仿真四面体网格模型的仿真变形运动序列；计算仿真变形运动序列与跟踪变形运动序列的位置偏差；沿着使得位置偏差减小的方向更新弹性物体的材质属性系数；直到找到位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状；
步骤105：根据位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状，建立弹性物体变形运动模型。
本发明实施例提供的弹性物体变形运动的建模方法，首先，采集弹性物体的静态点云和动态点云序列；其次，根据静态点云，建立用于仿真弹性物体变形运动的仿真四面体网格模型；接着，驱动仿真四面体网格模型跟踪动态点云序列，得到仿真四面体网格模型的跟踪变形运动序列；然后，迭代估计弹性物体的材质属性系数和对应的仿真四面体网格模型的参考形状；每个迭代周期均执行以下操作：获得弹性物体的当前材质属性系数对应的仿真四面体网格模型的参考形状；根据当前弹性物体材质属性系数和对应的仿真四面体网格模型的参考形状，驱使仿真四面体网格模型从相同的初始形变下，仿真弹性物体变形运动，得到仿真四面体网格模型的仿真变形运动序列；计算仿真变形运动序列与跟踪变形运动序列的位置偏差；沿着使得位置偏差减小的方向更新弹性物体的材质属性系数；直到找到位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状；最后，根据位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状，建立弹性物体变形运动模型。通过上述建模方法，可以建立逼真的弹性物体变形运动模型。
下面就对本发明实施例中提到的各步骤进行详细描述：
在上述步骤101中，我们使用三个Kinect组合，来采集弹性物体的静态点云序列和变形运动过程中的动态点云序列，图2即为采集到的弹性物体的点云序列的示意图，该步骤101即为采集弹性物体的静态点云和动态点云序列的过程。
在上述步骤102中，根据静态点云，建立用于仿真弹性物体变形运动的仿真四面体网格模型，可以包括：
根据静态点云，建立弹性物体的静态表面网格模型；由于该静态表面网格模型十分精确，又可称为静态精细表面网格模型，图3所示即为静态表面网格模型；
根据静态表面网格模型，建立用于仿真弹性物体变形运动的仿真四面体网格模型；图4即为仿真四面体网格模型的示意图；
静态表面网格模型的每个顶点与仿真四面体网格模型的每个四面体空间重心坐标是线性插值关系。
图3中静态表面网格模型有15368个顶点，对应到仿真四面体网格模型上，有9594个节点。
该步骤102即为建立用于仿真弹性物体变形运动的仿真四面体网格模型的过程。
具体实施时，首先，创建物体的精细表面网格模型(如图3所示)。然后将精细表面网格模型的网格传递给体数据生成工具TETGEN，导出用来物理仿真的相对粗糙的四面体网格，它将作为模板用以追踪点云序列(如图4所示)。为了从仿真的四面体元素中得到精细网格的位置信息，这里采用一种嵌入策略，精细网格上的顶点均可用所在四面体空间重心坐标表示，两者是线性插值关系。
在上述步骤103中，驱动仿真四面体网格模型跟踪动态点云序列，得到仿真四面体网格模型的跟踪变形运动序列，可以包括如下步骤：
找到弹性物体变形后所有动态点云序列与仿真四面体网格模型的所有节点的最大概率对应关系；
根据最大概率对应关系，向仿真四面体网格模型的每个节点施加虚拟外力，驱使仿真四面体网格模型的每个节点跟踪对应的动态点云序列，得到变形后的仿真四面体网格模型的每个节点的位置；
根据线性插值关系和变形后的仿真四面体网格模型的每个节点的位置，找到变形后的静态表面网格模型的每个顶点位置；
计算变形后的静态表面网格模型的每个顶点位置与弹性物体变形后对应的动态点云序列位置之间的差异，直到差异小于预设值时，得到当前点云数据对应的仿真四面体网格模型变形后的跟踪变形运动序列。
在一个实施例中，差异小于预设值可以包括：静态表面网格模型的每个顶点与弹性物体变形后对应的动态点云序列之间的距离小于预设距离，或它们之间的吸引力小于预设吸引力。
上述步骤103即为驱动仿真四面体网格模型变形跟踪弹性物体变形运动的过程。
具体实施时，该步骤103为基于物理的概率跟踪方法，运动跟踪需要处理带噪声的点云数据，同时还要考虑遮挡、快速运动和大幅度变形的问题。因此我们将变形运动跟踪转化成一个最大后验概率(MAP)问题，利用期望最大方法(EM)进行迭代求解：(E步)根据当前点云和节点位置寻找最优的对应关系；(M步)移动节点位置使得上述对应关系为最大似然估计。
下面详细说明上述跟踪的过程。
假定c＝c1:N,1≤n≤N表示采集到的点云中点的坐标，s＝s1:k,1≤k≤K表示网格节点的位置。我们的任务就是基于给定的点云c求出形状匹配的s。点云中的点和网格上点的关系并不知晓，而是作为一个隐变量Zkn，表明观测点Cn可能来自于节点sk，假设Cn符合正态分布，具体为：cn～N(sk,∑k)，协方差矩阵为∑k＝σ2I，其中，σ为方差，I为单位矩阵。网格变形跟踪匹配点云可以表示成一个最大后验概率的问题：
Figure PCTCN2015088126-appb-000003
此时，运用EM算法进行求解。在E部分，基于对隐变量p(Zkn|s,c)(p指代的是概率分布)的期望，找到总对数联合概率最低边界log p(s,c)；在M部分，通过调整四面体网格的顶点位置来最大化上一步的最低边界：
Figure PCTCN2015088126-appb-000004
上式的第二项反映了变形物体模型的势能，因此可以采取物理仿真方式进行优化求解。我们通过向每个节点添加虚拟力：
Figure PCTCN2015088126-appb-000005
其中，η为虚拟力的比例系数；
作为外力驱动网格变形去匹配点云形状。这里运动方程为：
Figure PCTCN2015088126-appb-000006
这里采用共旋转的线性有限元模型来仿真变形物体，M是质量矩阵，D＝αM+βK是瑞利阻尼(Rayleigh damping)矩阵，α和β是瑞利阻尼对应的两个系数，R是旋转矩阵，通过对变形梯度做极分解得到，K是刚度矩阵，x是仿真四面体网格模型变形后的形状，X仿真四面体网格模型的参考形状，fext是外力合力。物理仿真部分采用第三方库VEGA FEM封装而成。为了加速计算，这里采用嵌套策略，每次迭代计算将外力映射给节点数量较少的跟踪四面体模型来仿真跟踪运动，然后把节点位移再插值回静态表面网格模型顶点上。运动跟踪过程即为施加外力仿真直至EM迭代收敛。
对于上述提到的运动跟踪部分，传统的技术方案一般采用非刚性配准算法来逐帧匹配模板网格和点云序列，这种方法首先运算速度不够快；其次对点云数据有较高要求，无法处理大变形运动和存在较多噪声干扰项的情形；再次匹配结果因为没有融合物理约束，网格拓扑会走样。而通过上述关于步骤103的记载，本申请技术方案采用基于概率的跟踪算法，找到网格顶点与点云的概率对应关系，并在物理引擎的驱动下驱使网格变形运动，同时匹配结果融合物理约束(满足上述的各个公式方程的约束)，网格拓扑不会走样，有更快的运行速度，更好的跟踪效果，更强的鲁棒性。上述步骤103可以使得建立弹性物体变形运动模型的点云数据噪声小，同时还可以解决遮挡、快速运动和大幅度变形的问题。
在上述步骤104中，获得弹性物体的当前材质属性系数对应的仿真四面体网格模型的参考形状，包括：
验证当前材质属性系数和对应的仿真四面体网格模型的参考形状是否满足物理静力平衡方程。
具体实施时，上述步骤104包括参数估计和参考形状的优化部分，在该部分中，我们通过数据驱动的方法估计未知的静态材质属性系数p＝(E,ν,α,β)和相应的参考形状X。其中E是杨氏模量，ν是泊松系数，α和β是瑞利阻尼对应的两个系数。求解问题可以表示成时空优化问题，下面的目标方程F测量了仿真和跟踪序列的位置偏差：
Figure PCTCN2015088126-appb-000007
其中，
Figure PCTCN2015088126-appb-000008
是跟踪的输出结果，xt是仿真结果在第t帧的位置。这个时空优化问题高维度，非线性且目标函数是非凸函数，传统方法不能有效地解决。因此我们提出了一种 新颖的分治策略交替地迭代优化X和p，即通过该方法可以寻找出最优的材质属性系数和参考形状(即位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状)。
在算法的每次迭代中，我们首先执行参考形状的优化估计，保证当前的材质属性系数p和其对应的参考形状X满足物理静态平衡的约束，这样该步骤融合了物理约束。接着我们采用标准的向下搜索方法沿着减少轨迹偏差的方向去更新材质属性系数p。我们循环执行该策略直到足够收敛。
在一个实施例中，求解弹性物体的仿真四面体网格模型的参考形状时，采用如下目标方程求解合力残差最小时对应的形状作为仿真四面体网格模型的参考形状：
Figure PCTCN2015088126-appb-000009
其中，R是旋转矩阵，K是刚度矩阵，xs是弹性物体的静态形状，X是仿真四面体网络模型的参考形状，M为弹性物体的质量，g为重力加速度；
具体实施时，为了重建可信的仿真运动，物体模型的参考形状和静态形状应该区分开来，参考形状应该是不受重力因素影响的，否则最明显的失真现象就是当物体从静止形状开始仿真时，会发生明显的形状改变。拿植物的变形运动举例来说，植物叶子都会受重力影响，有一个先向下耷拉一下的轻微动作，那么仿真运动必然不真实，由于发明人考虑到了这个问题，并且采用了上述目标方程，求解出的合力残差施加给仿真四面体网络模型，得到仿真四面体网络模型的参考形状。应用该方法，可以生成真实感强的弹性物体变形运动的仿真模型。
具体的推导过程为：由于物体本身重力是唯一导致这两种形状差别的原因，我们去除加速度和速度项，简化运动方程来优化求解模型的参考形状，用xs表示物体静态形状，那么简化后的力平衡方程为：
RK(RTxs-X)＝Mg；
于是为了得到最小的合力残差，有如下优化目标：
Figure PCTCN2015088126-appb-000010
我们应用当前合力残差作为虚拟外力施加给仿真模型，让其通过仿真变形来更新模型参考形状。仿真持续运行直到合力残差足够小并且仿真达到稳定状态。这种方式相比于传统优化方法的好处在于其强健和快速，并且本身变形满足物理约束。
驱使仿真四面体网格模型变形时，计算施加给每个节点的虚拟弹性力相对于参考形状的雅克比矩阵，雅克比矩阵为：
Figure PCTCN2015088126-appb-000011
其中，f为施加给仿真四面体网格模型的每个节点的虚拟弹性力，Xij为仿真四面体网格模型参考形状中第i个节点在j方向上的位置，xs是弹性物体的静态形状，X是仿真四面体网格模型的参考形状，V是仿真四面体网格模型的每个四面体元素的体积，E是与弹性物体的材质属性相关的常量矩阵，R是代表仿真四面体网格模型的每个四面体元素的刚性旋转矩阵。
即，为实现快速稳定的仿真以解决上述优化问题，我们采用隐式求解方法。这是需要计算弹性力相对于参考形状的雅克比矩阵
Figure PCTCN2015088126-appb-000012
下面给出详细的推导过程：拟定
Figure PCTCN2015088126-appb-000013
是四面体网格在未变形时节点的位置，
Figure PCTCN2015088126-appb-000014
是参考形状的节点位置。共旋转的线性有限元模型的弹性力大小为：
f＝VRBTE(RTx-X)；
其中，V是四面体元素的体积，E是一个6×6的常量矩阵，跟材质的弹性属性相关。矩阵R＝diag(R,R,R,R)是一个分块对角矩阵，R是代表元素的刚性旋转，通过对模型的变形梯度做极分解得到，即F＝RS，S是对称矩阵。对于一个四面体元素，F是矩阵
Figure PCTCN2015088126-appb-000015
左上方3×3子块矩阵，其中，
Figure PCTCN2015088126-appb-000016
Bm＝Vm -1都为4×4的矩阵。
6×12的矩阵B仅依赖于X，内部由Bm构成，这里用Bij表示Bm中第i行第j列的元素：
Figure PCTCN2015088126-appb-000017
我们要求的雅克比矩阵为：
Figure PCTCN2015088126-appb-000018
kq代表K的第q列，即为X第q列的偏微分，实际上也是就是节点j第i项组成部分，这里用Xij表示。应用链式法则，我们有：
Figure PCTCN2015088126-appb-000019
这里单引号(')代表偏微分
Figure PCTCN2015088126-appb-000020
是
Figure PCTCN2015088126-appb-000021
第q项的标准基。下面是具体每一小项的求导结果：
Figure PCTCN2015088126-appb-000022
关于R'的计算更加复杂，还是应用链式法则，可以得到：
Figure PCTCN2015088126-appb-000023
最后，关于上式右边第二项的求导，我们有：
Figure PCTCN2015088126-appb-000024
Figure PCTCN2015088126-appb-000025
我们采用滴水观音模型的参考形状。应用该参考形状进行重力的变形仿真可以得到非常准确的静态形状，与采集到的静态形状基本一致。对于参数估计部分，现有的技术一般认为物体模型的参考形状已知或将参考形状等同于静态形状，这样分析的结果在简化后必然不精确。本技术方案通过上述技术方案，能同时优化估计出物体模型的物理参数和对应的参考形状。
在一个实施例中，估计弹性物体的材质属性系数，包括：
从仿真四面体网格模型的不同位置选出多个节点作为控制点；
根据弹性物体的材质分布属性，估计出不同位置的控制点的不同材质属性系数；
根据线性插值算法和不同位置的控制点的不同材质属性系数，求出仿真四面体网格模型的其他四面体元素的材质属性系数。
具体实施时，我们采用上述方法，针对弹性物体的不同位置的不同材质，给不同位置的材质分别估算材质属性系数，那么就可以生成非均一材质的数学物理模型，使得所建仿真模型更加的真实。
具体地，我们需要弹性物体的真实运动轨迹去估计其弹性材质参数和阻尼系数。给模型赋予合适的材质分布，并给予同样的初始位置条件，就可能通过仿真重建一样的运动。考虑到参考形状可以通过上述方法分开计算得到，之前的目标方程可以重写为：
Figure PCTCN2015088126-appb-000026
t是帧的编号，k是节点编号。为了最小化F(p)，我们采用无梯度的单纯形优化方法求解。另外，由于单一材质的模型不能很好的真实还原物体变形运动，我们引入控制点的概念来解决此问题，通过指定控制点不同的材质属性，模型其他节点通过线性插值也得到不均匀分布的材质属性。如图5所示，对植物模型，我们将控制点安放在轴方向，插值权重用归一化的轴距离表示。对于生成的恐龙模型，控制点根据设想认为指定在不同位置，插值权重采用双调和函数控制。
在一个实施例中，弹性物体的材质属性系数的初值估计，包括：
从弹性物体变形运动时的多个变形模态中选择第一个主要的模态；
确定第一个主要的模态对应的振动频率；
根据第一个主要的模态对应的振动频率与实际采集数据振动频率的匹配程度，确定弹性物体的材质属性系数的初值。
上述“主要的模态”是指最小特征值对应的振动模态。
具体实施时，由于目标函数F(p)通常包含多个局部最小值。因此提供一个合适的输入参数对最终成功优化求解至关重要，即在首次估计材质属性系数时，给出一个最佳值，那么有利于后续快速的计算，提高建模的效率。具体地，我们提出了一种新颖的策略，利用模态分析和坐标下降法来得到合适的初始材质属性系数值。
在模态分析中，小的形变置换
Figure PCTCN2015088126-appb-000027
被线性化成u＝Φz,其中Φ＝[Φ1,Φ2,...,Φk]的每一列代表了一个变形模态。其可以通过广义特征值分解求得Kφi＝λiMφi。一般将较小特征值对应的特征向量作为基建立降维的模态坐标系
Figure PCTCN2015088126-appb-000028
对应每个模态的自然频率。
在特征值分解过程中，杨氏模量E会影响振动频率
Figure PCTCN2015088126-appb-000029
直观理解，材质越软的物体振动频率越小。
我们将真实位移投影到相应模态下从而得到真实的振动频率。如果估计的杨氏模量越接近基准真实值，那么频率差距应该越小。这里采用坐标下降法法去顺序更新每一个材质参数。在每次迭代过程中，我们基于变量变化对目标函数值变化的敏感程度线性搜索。
我们另外还用两个真实模型(杯垫和衣架)来验证我们的发明测量精确度，我们测得的杨氏模量分别为7.0e6和5.6e6。结果和真实实物分别做了静态受力和动态变形的对比。实验结果验证了发明的可行性，测量结果精确度较高。
图6即为本发明一个实施例中弹性物体变形运动的建模实施的示意图，如图6所示，首先采集构建弹性物体的静态形状和动态点云运动序列；接着是本***的核心，一种交替迭代的优化策略轮流执行变形运动跟踪和参数估计部分，每次迭代运行结果能显著改进另一部分效果。将静态精细模型作为跟踪的模板网格，***采用基于物理的概率跟踪算法去驱使网格变形配准每一帧点云，并输出四面体网格每一帧节点位置。下一步，参数估计部分同时优化估计材质属性系数，阻尼系数和物体模型参考形状。这里参考形状指的是物体模型不受任何外力包括重力作用时的原始形状，因此物体参考形状和受重力影响的静止形状应有所差别。在优化部分采取分治策略，给定初始估计物理参数后，根据静力平衡方程求解出当前模型参考形状，然后使用这组材质和参考形状数据进行正向仿真，得到在相同初始变形条件下的运动序列；并计算运动序列形状与跟踪结果之间的差异，作为评价该组参数的标准。多次迭代后找到最佳值；最终，***生成一套能进行真实感交互的仿真模型。
本发明实施例提供的弹性物体变形运动的建模方法可以达到的有益技术效果为：
1)、实验采集平台容易搭建，设备平价，无需额外成本。
2)、运动重建满足物理约束，能处理大幅度变形，对噪声和数据确实有较强的容错性和鲁棒性。
3)、参数估计快速高效，利用分治策略同时计算物理参数和参考形状使得参数估计和后续仿真更加准确逼真。
以上仅为本发明的优选实施例而已，并不用于限制本发明，对于本领域的技术人员来说，本发明实施例可以有各种更改和变化。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (8)

  1. 一种弹性物体变形运动的建模方法，其特征在于，包括：
    采集弹性物体的静态点云和变形运动过程中的动态点云序列；
    根据所述静态点云，建立用于仿真弹性物体变形运动的仿真四面体网格模型；
    驱动所述仿真四面体网格模型跟踪所述动态点云序列，得到所述仿真四面体网格模型的跟踪变形运动序列；
    迭代估计弹性物体的材质属性系数和对应的仿真四面体网格模型的参考形状；每个迭代周期均执行以下操作：获得弹性物体的当前材质属性系数对应的仿真四面体网格模型的参考形状；根据当前弹性物体材质属性系数和对应的仿真四面体网格模型的参考形状，驱使所述仿真四面体网格模型从相同的初始形变下，仿真弹性物体变形运动，得到所述仿真四面体网格模型的仿真变形运动序列；计算所述仿真变形运动序列与跟踪变形运动序列的位置偏差；沿着使得位置偏差减小的方向更新弹性物体的材质属性系数；直到找到位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状；
    根据位置偏差最小时的材质属性系数和对应的仿真四面体网格模型的参考形状，建立弹性物体变形运动模型。
  2. 如权利要求1所述的弹性物体变形运动的建模方法，其特征在于，根据所述静态点云，建立用于仿真弹性物体变形运动的仿真四面体网格模型，包括：
    根据所述静态点云，建立弹性物体的静态表面网格模型；
    根据所述静态表面网格模型，建立用于仿真弹性物体变形运动的仿真四面体网格模型；
    所述静态表面网格模型的每个顶点与所述仿真四面体网格模型的每个四面体空间重心坐标是线性插值关系。
  3. 如权利要求2所述的弹性物体变形运动的建模方法，其特征在于，驱动所述仿真四面体网格模型跟踪所述动态点云序列，得到所述仿真四面体网格模型的跟踪变形运动序列，包括：
    找到弹性物体变形后所有动态点云序列与所述仿真四面体网格模型的所有节点的最大概率对应关系；
    根据所述最大概率对应关系，向所述仿真四面体网格模型的每个节点施加虚拟外力，驱使所述仿真四面体网格模型的每个节点跟踪对应的动态点云序列，得到变形后的所述仿真四面体网格模型的每个节点的位置；
    根据所述线性插值关系和变形后的所述仿真四面体网格模型的每个节点的位置，找到变形后的静态表面网格模型的每个顶点位置；
    计算变形后的静态表面网格模型的每个顶点位置与弹性物体变形后对应的动态点云序列位置之间的差异，直到差异小于预设值时，得到当前点云数据对应的仿真四面体网格模型变形后的跟踪变形运动序列。
  4. 如权利要求3所述的弹性物体变形运动的建模方法，其特征在于，所述差异小于预设值包括：静态表面网格模型的每个顶点与弹性物体变形后对应的动态点云序列之间的距离小于预设距离，或它们之间的吸引力小于预设吸引力。
  5. 如权利要求1所述的弹性物体变形运动的建模方法，其特征在于，获得弹性物体的当前材质属性系数对应的仿真四面体网格模型的参考形状，包括：
    验证当前材质属性系数和对应的仿真四面体网格模型的参考形状是否满足物理静力平衡方程。
  6. 如权利要求1所述的弹性物体变形运动的建模方法，其特征在于，弹性物体的材质属性系数的初值估计，包括：
    从弹性物体变形运动时的多个变形模态中选择第一个主要的模态；
    确定所述第一个主要的模态对应的振动频率；
    根据第一个主要的模态对应的振动频率与实际采集数据振动频率的匹配程度，确定弹性物体的材质属性系数的初值。
  7. 如权利要求1所述的弹性物体变形运动的建模方法，其特征在于，求解弹性物体的仿真四面体网格模型的参考形状时，采用如下目标方程求解合力残差最小时对应的形状作为仿真四面体网格模型变形的参考形状：
    Figure PCTCN2015088126-appb-100001
    其中，R是旋转矩阵，K是刚度矩阵，xs是弹性物体的静态形状，X是仿真四面体网络模型的参考形状，M为弹性物体的质量，g为重力加速度；
    驱使仿真四面体网格模型变形时，计算施加给每个节点的虚拟弹性力相对于仿真参考形状的雅克比矩阵，雅克比矩阵为：
    Figure PCTCN2015088126-appb-100002
    其中，f为施加给仿真四面体网格模型的每个节点的虚拟弹性力，Xij为仿真四面体网格模型参考形状中第i个节点在j方向上的位置，xs是弹性物体的静态形状，X是仿真四面体网格模型的参考形状，V是仿真四面体网格模型的每个四面体元素的体积，E是与弹性物体的材质属性相关的常量矩阵，R是代表仿真四面体网格模型的每个四面体元素的刚性旋转矩阵。
  8. 如权利要求1所述的弹性物体变形运动的建模方法，其特征在于，估计弹性物体的材质属性系数，包括：
    从仿真四面体网格模型的不同位置选出多个节点作为控制点；
    根据弹性物体的材质分布属性，估计出不同位置的控制点的不同材质属性系数；
    根据线性插值算法和不同位置的控制点的不同材质属性系数，求出仿真四面体网格模型的其他四面体元素的材质属性系数。
PCT/CN2015/088126 2015-08-26 2015-08-26 弹性物体变形运动的建模方法 WO2017031718A1 (zh)

Priority Applications (2)

Application Number Priority Date Filing Date Title
PCT/CN2015/088126 WO2017031718A1 (zh) 2015-08-26 2015-08-26 弹性物体变形运动的建模方法
US15/392,688 US10394979B2 (en) 2015-08-26 2016-12-28 Method and device for elastic object deformation modeling

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
PCT/CN2015/088126 WO2017031718A1 (zh) 2015-08-26 2015-08-26 弹性物体变形运动的建模方法

Related Child Applications (1)

Application Number Title Priority Date Filing Date
US15/392,688 Continuation US10394979B2 (en) 2015-08-26 2016-12-28 Method and device for elastic object deformation modeling

Publications (1)

Publication Number Publication Date
WO2017031718A1 true WO2017031718A1 (zh) 2017-03-02

Family

ID=58099480

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/CN2015/088126 WO2017031718A1 (zh) 2015-08-26 2015-08-26 弹性物体变形运动的建模方法

Country Status (2)

Country Link
US (1) US10394979B2 (zh)
WO (1) WO2017031718A1 (zh)

Cited By (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109670271A (zh) * 2019-01-15 2019-04-23 上海交通大学 基于敏感点多点增量成形的大型薄壁加筋板的校形方法
CN111539126A (zh) * 2020-05-14 2020-08-14 山东慧泰智能科技有限公司 一种化工储罐蒸汽云***事故模拟分析方法和***
CN114347029A (zh) * 2022-01-11 2022-04-15 大连理工大学 一种用于气动软体机器人快速模拟的模型降阶方法
CN114741745A (zh) * 2022-06-10 2022-07-12 浙江理工大学 一种基于纱线动力学的纺织品成型模拟预测方法

Families Citing this family (14)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
EP3136293A1 (en) * 2015-08-28 2017-03-01 Thomson Licensing Method and device for processing an image of pixels, corresponding computer program product and computer readable medium
CN111159906B (zh) * 2019-12-31 2023-06-16 上海理工大学 基于adams二次开发的汽车悬架动态模型搭建方法
US11597092B1 (en) 2020-03-26 2023-03-07 Amazon Technologies, Ine. End-of-arm tool with a load cell
CN111460743B (zh) * 2020-04-07 2022-08-09 西安交通大学 一种数据驱动的仿生推进水动力性能参数采集方法
CN111597720B (zh) * 2020-05-19 2022-06-10 上海工程技术大学 一种考虑空间大变形的介入手术导丝建模方法及其应用
CN112036029A (zh) * 2020-05-22 2020-12-04 江西星盈科技有限公司 电池组模组膨胀变形的仿真方法
US11534924B1 (en) * 2020-07-21 2022-12-27 Amazon Technologies, Inc. Systems and methods for generating models for automated handling of vehicles
US11534915B1 (en) 2020-08-05 2022-12-27 Amazon Technologies, Inc. Determining vehicle integrity based on observed behavior during predetermined manipulations
CN112233248B (zh) * 2020-10-19 2023-11-07 广东省计量科学研究院（华南国家计量测试中心） 基于三维点云的表面平整度的检测方法、***及介质
CN112529010B (zh) * 2020-12-04 2023-03-24 浙江大学计算机创新技术研究院 一种基于在线局部特征提取的点云识别方法
CN112528522A (zh) * 2020-12-29 2021-03-19 福州智永信息科技有限公司 一种轻量级柔体效果实现方法及终端
CN113910211B (zh) * 2021-11-06 2023-06-16 中国船舶集团有限公司第七一六研究所 考虑自重的协作机器***性变形建模方法及终端
CN114182488B (zh) * 2021-12-16 2023-11-21 海信冰箱有限公司 滚筒洗衣机仿真方法、装置、计算机可读介质及洗衣机
CN115659643B (zh) * 2022-10-26 2023-09-26 西南交通大学 一种铁路轨道高聚物弹性元件动静刚度比的评价方法

Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101216950A (zh) * 2008-01-16 2008-07-09 浙江大学 一种基于线性微分算子的弹性形变模拟方法
CN101634618A (zh) * 2009-08-27 2010-01-27 浙江大学 一种三维弹性模量成像方法
CN103942377A (zh) * 2014-04-11 2014-07-23 浙江大学 一种面向弹性物体制造的逆向形状设计方法
CN105069826A (zh) * 2015-08-26 2015-11-18 中国科学院深圳先进技术研究院 弹性物体变形运动的建模方法

Family Cites Families (42)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5524636A (en) * 1992-12-21 1996-06-11 Artann Corporation Dba Artann Laboratories Method and apparatus for elasticity imaging
US5988862A (en) * 1996-04-24 1999-11-23 Cyra Technologies, Inc. Integrated system for quickly and accurately imaging and modeling three dimensional objects
US6456289B1 (en) * 1999-04-23 2002-09-24 Georgia Tech Research Corporation Animation system and method for a animating object fracture
FR2820534B1 (fr) * 2001-02-08 2003-06-27 Inst Nat Rech Inf Automat Dispositif de simulation de deformations de materiaux, notamment de tissus corporels mous
JP4635392B2 (ja) * 2001-08-09 2011-02-23 コニカミノルタホールディングス株式会社 ３次元物体の表面形状モデリング装置、および、プログラム
US7373284B2 (en) * 2004-05-11 2008-05-13 Kimberly-Clark Worldwide, Inc. Method of evaluating the performance of a product using a virtual environment
US7346382B2 (en) * 2004-07-07 2008-03-18 The Cleveland Clinic Foundation Brain stimulation models, systems, devices, and methods
US7542034B2 (en) * 2004-09-23 2009-06-02 Conversion Works, Inc. System and method for processing video images
US20080259073A1 (en) * 2004-09-23 2008-10-23 Conversion Works, Inc. System and method for processing video images
US20080246836A1 (en) * 2004-09-23 2008-10-09 Conversion Works, Inc. System and method for processing video images for camera recreation
US7493243B2 (en) * 2004-12-27 2009-02-17 Seoul National University Industry Foundation Method and system of real-time graphical simulation of large rotational deformation and manipulation using modal warping
US7650266B2 (en) * 2005-05-09 2010-01-19 Nvidia Corporation Method of simulating deformable object using geometrically motivated model
US8004517B1 (en) * 2005-06-24 2011-08-23 Geomagic, Inc. Methods, apparatus and computer program products that model three-dimensional surface structures
JP2007034724A (ja) 2005-07-27 2007-02-08 Glory Ltd 画像処理装置、画像処理方法および画像処理プログラム
EP2081494B1 (en) * 2006-11-16 2018-07-11 Vanderbilt University System and method of compensating for organ deformation
EP2618102A2 (en) * 2006-11-21 2013-07-24 Mantisvision Ltd. 3d geometric modeling and 3d video content creation
EP2006803A1 (en) * 2007-06-19 2008-12-24 Agfa HealthCare NV Method of segmenting anatomic entities in 3D digital medical images
CN100583160C (zh) 2007-10-19 2010-01-20 北京航空航天大学 一种基于细节编码及重构的物理变形方法
US8384714B2 (en) * 2008-05-13 2013-02-26 The Board Of Trustees Of The Leland Stanford Junior University Systems, methods and devices for motion capture using video imaging
WO2010004466A1 (en) * 2008-07-10 2010-01-14 C-True Ltd. Three dimensional mesh modeling
US9245382B2 (en) * 2008-10-04 2016-01-26 Microsoft Technology Licensing, Llc User-guided surface reconstruction
CN101425188A (zh) 2008-12-03 2009-05-06 中国人民解放军国防科学技术大学 带刚体核的广义弹簧振子形变仿真方法
DE102009007935B4 (de) * 2009-02-06 2011-06-30 Lewandowski, Angela, 30159 Vorrichtung und Verfahren für die künstliche Verlängerung von Nägeln
EP2317348B1 (en) * 2009-10-30 2014-05-21 Services Pétroliers Schlumberger Method for building a depositional space corresponding to a geological domain
US8565909B2 (en) * 2010-02-24 2013-10-22 Disney Enterprises, Inc. Fabrication of materials with desired characteristics from base materials having determined characteristics
US9135738B2 (en) * 2011-05-02 2015-09-15 Disney Enterprises, Inc. Efficient elasticity for character skinning
US10019440B2 (en) * 2011-05-27 2018-07-10 Adobe Systems Incorporated Methods and apparatus for three-dimensional (3D) sketching
US9235928B2 (en) * 2012-01-24 2016-01-12 University Of Southern California 3D body modeling, from a single or multiple 3D cameras, in the presence of motion
EP2674913B1 (en) * 2012-06-14 2014-07-23 Softkinetic Software Three-dimensional object modelling fitting & tracking.
US9477798B1 (en) * 2012-06-20 2016-10-25 Synopsys, Inc. Moving mesh system and method for finite element/finite volume simulations
WO2014025791A1 (en) * 2012-08-06 2014-02-13 The University Of North Carolina At Chapel Hill Simulation-based estimation of elasticity parameters and use of same for non-invasive cancer detection and cancer staging
EP2880628A1 (en) * 2012-08-06 2015-06-10 Vanderbilt University Enhanced method for correcting data for deformations during image guided procedures
US10642945B2 (en) * 2012-10-24 2020-05-05 New York University Structural weak spot analysis
US20160000518A1 (en) * 2013-02-11 2016-01-07 Neomedz Sàrl Tracking apparatus for tracking an object with respect to a body
CN104933755B (zh) * 2014-03-18 2017-11-28 华为技术有限公司 一种静态物体重建方法和***
WO2015187483A1 (en) * 2014-06-05 2015-12-10 Geocosm, LLC Predicting sediment and sedimentary rock properties
CN104063894B (zh) * 2014-06-13 2017-02-22 中国科学院深圳先进技术研究院 点云三维模型重建方法及***
CN104050720B (zh) * 2014-06-24 2017-02-15 中国科学院深圳先进技术研究院 点云曲面重建方法及***
CA2923917A1 (en) * 2015-03-17 2016-09-17 Roel Vertegaal Flexible display for a mobile computing device
US9842411B2 (en) * 2015-04-03 2017-12-12 Disney Enterprises, Inc. Geometric multigrid on incomplete linear octrees for simulating deformable animated characters
US9767572B2 (en) * 2015-05-01 2017-09-19 Raytheon Company Systems and methods for 3D point cloud processing
WO2017113260A1 (zh) * 2015-12-30 2017-07-06 中国科学院深圳先进技术研究院 一种三维点云模型重建方法及装置

Patent Citations (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN101216950A (zh) * 2008-01-16 2008-07-09 浙江大学 一种基于线性微分算子的弹性形变模拟方法
CN101634618A (zh) * 2009-08-27 2010-01-27 浙江大学 一种三维弹性模量成像方法
CN103942377A (zh) * 2014-04-11 2014-07-23 浙江大学 一种面向弹性物体制造的逆向形状设计方法
CN105069826A (zh) * 2015-08-26 2015-11-18 中国科学院深圳先进技术研究院 弹性物体变形运动的建模方法

Cited By (7)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN109670271A (zh) * 2019-01-15 2019-04-23 上海交通大学 基于敏感点多点增量成形的大型薄壁加筋板的校形方法
CN109670271B (zh) * 2019-01-15 2023-03-07 上海交通大学 基于敏感点多点增量成形的大型薄壁加筋板的校形方法
CN111539126A (zh) * 2020-05-14 2020-08-14 山东慧泰智能科技有限公司 一种化工储罐蒸汽云***事故模拟分析方法和***
CN114347029A (zh) * 2022-01-11 2022-04-15 大连理工大学 一种用于气动软体机器人快速模拟的模型降阶方法
CN114347029B (zh) * 2022-01-11 2023-09-12 大连理工大学 一种用于气动软体机器人快速模拟的模型降阶方法
CN114741745A (zh) * 2022-06-10 2022-07-12 浙江理工大学 一种基于纱线动力学的纺织品成型模拟预测方法
CN114741745B (zh) * 2022-06-10 2022-09-16 浙江理工大学 一种基于纱线动力学的纺织品成型模拟预测方法

Also Published As

Publication number Publication date
US20170109465A1 (en) 2017-04-20
US10394979B2 (en) 2019-08-27

Similar Documents

Publication Publication Date Title
WO2017031718A1 (zh) 弹性物体变形运动的建模方法
CN105069826B (zh) 弹性物体变形运动的建模方法
Liang et al. Differentiable cloth simulation for inverse problems
CN105843073B (zh) 一种基于气动力不确定降阶的机翼结构气动弹性稳定性分析方法
US7616204B2 (en) Method of simulating dynamic objects using position based dynamics
CN104268943A (zh) 一种基于欧拉-拉格朗日耦合方法的流体仿真方法
CN110289104B (zh) 软组织按压和形变恢复的模拟方法
Sun et al. Topology optimization based on level set for a flexible multibody system modeled via ANCF
US20230061175A1 (en) Real-Time Simulation of Elastic Body
US20230286144A1 (en) Real-time vibration-suppression control for robotic systems
JP2023544175A (ja) メッシュ表現およびグラフニューラルネットワークを使用した物理的環境のシミュレーション
CN113591324A (zh) 一种软体机器人仿真方法、装置、电子设备及存储介质
CN104463934B (zh) 一种“质点‑弹簧”***驱动的点集模型动画自动生成方法
CN103839280B (zh) 一种基于视觉信息的人体姿态跟踪方法
JP2022184829A (ja) ３ｄ形状最適化のためのディープパラメータ化
CN112800690B (zh) 一种基于群智能优化算法下的水下折展机构参数优化方法
CN112308952B (zh) 模仿给定视频中人物动作的3d角色动作生成***和方法
CN105807093A (zh) 一种基于粒子图像测速技术的加速度测量方法及装置
CN105095555A (zh) 一种速度场的无散平滑处理方法及装置
JP4612989B2 (ja) メカトロニクスシステムのシミュレーション方法
CN116416407A (zh) 衣物仿真中碰撞和自碰撞处理方法、***、设备及介质
Huang et al. A survey on fast simulation of elastic objects
CN112862957B (zh) 一种基于约束投影的gpu并行试衣仿真方法
CN101655832B (zh) 一种基于标量场梯度的物理变形方法
Cheng Algorithm of CAD surface generation for complex pipe model in Industry 4.0 background

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 15901989

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

32PN Ep: public notification in the ep bulletin as address of the adressee cannot be established

Free format text: NOTING OF LOSS OF RIGHTS PURSUANT TO RULE 112(1) EPC (EPO FORM 1205A DATED 22.06.2018)

122 Ep: pct application non-entry in european phase

Ref document number: 15901989

Country of ref document: EP

Kind code of ref document: A1