1

Online Estimation of Point-based Volumetric
Stiffness Model Using Joint Torque Sensors
Shaoxiong Yao1 and Kris Hauser1

Abstract—This paper introduces a method for learning theforce response of novel, heterogeneous deformable objects directlyfrom robot sensor data without prior knowledge. We present a
new volumetric stiffness field representation and a point-based
collision dynamics model that lends itself to rapid estimation
using force or torque measurements. By estimating the stiffnessof each colliding point using an extended Kalman filter, ourmethod can update a stiffness field over a grid of 8 million cells
at 17Hz. Experiments show that for the task of predicting theforce response from a robot with joint torque sensors repeatedly
touching parts of an artificial plant, our method is more accurate
than other learning methods by a wide margin.

I. I NTRODUCTION
Robots deployed in domestic [1] and agricultural [2] environments must frequently interact with deformable objects likeclothes and plants. Although it is relatively well known how
to plan paths and manipulation trajectories if a deformableobject model is assumed [3, 4], it is far less clear how toestimate a model of the deformable object in the first place.
Various deformable object models have been used in robotics
applications [5]. Material parameter identification methods areavailable to match observed deformation or force/torque data
to FEM [6] or spring-mass models [7]. However, they arecomputationally expensive and assume a given mesh topology.
Tactile models based on machine learning methods predict theforce response at varying touch locations, and avoid strong
assumptions about the mesh structure [8, 9]. However, theseprior works assume a roughly planar object, discrete touches,and coherence in the object’s force response across the object’ssurface.
Inspired by the problem of plant manipulation, our methodaddresses deformation modeling for highly heterogeneous
objects and continuous movement through a volume. We
propose a volumetric stiffness field (VSF) model where theobject’s potential resistance to movement is represented bypoints anchored to cells on a voxel grid. The stiffness ofeach Hookean spring is treated as an unknown distribution
and estimated by our algorithm. As the robot arm sweeps
across voxels, a fast point-based simulator estimates how thepoints are displaced by the arm, and an extended-Kalman filter
estimator updates the stiffness distribution as a function offorce or torque measurements. Similar probabilistic methods
have been used to estimate the locations of rigid contacts
from torque measurements [10, 11], and similar volumetricmethods have been used to map the shape of rigid objects
1 : S. Yao and K. Hauser are with the Department of Computer Science, University of Illinois at Urbana-Champaign, IL, USA. {syao16,kkhauser}@illinois.edu

Robot sensors

Jointangles 𝑞𝑞 𝑡𝑡

Jointtorques 𝜏𝜏 𝑡𝑡

Recursiveestimator

Contact simulator

Simulated point 𝑝𝑝𝑖𝑖
displacements 𝑢𝑢𝑖𝑖

Joint torques
prediction 𝜏𝜏̂ 𝑡𝑡
Estimation
update Δ𝐾𝐾𝑖𝑖

Volumetric stiffness 𝐾𝐾𝑖𝑖
HardK(N/m)

SoftFig. 1. Our volumetric stiffness field (VSF) estimation system consists of:(1) a contact simulator that simulates the displacement of points pushed by
the robot arm, (2) a VSF containing estimated stiffness distributions, (3) arecursive state estimator that updates the VSF beliefs.

from contact [12, 13]. Compared to these methods, the VSF
model has two important features:• Independence between points allows efficient realtime
estimation;• Flexibility to represent objects with heterogeneous stiffness, e.g. plants with hard branches and soft leaves.
We test our method in a scenario in which a robot arm
touches an artificial plant. Our method updates the belief ofa volumetric stiffness field with 8 million voxels at 17 Hz.
The estimated stiffness field is used to predict joint torques
when the same location is touched multiple times, and despitethe simplicity of the point-based dynamics model, our method
is approximately 5 times more accurate than learning-based
approaches. Moreover, we develop a stiffness extrapolation
method using color correlation, and the estimated stiffness
fields qualitatively identify stiffer regions like branches fromsofter regions like leaves. These results show our method has
applications to tactile modeling of general deformable objects.
II. M ETHODThe VSF models a deformable object as a set of points{pi }i=1..N from a region W expressed in a world frame. Whenthe robot arm sweeps through W , a subset of points will bepushed and displaced by the arm. We denote the indices ofcontact points as Wc ⊂ {1, .., N } and the displacement ofpoint pi at time t as uti ∈ R3 . The model makes three key
assumptions on the displacements and force responses:

2(i) Hookean force response. Each point has a virtual springwith stiffness Ki ∈ R+ connecting it to its rest positionpi . The force response Fit at time t is:
Fit = −Ki · uti(1)(ii) Quasi-static, stick-slip motion. We assume the robot
arm moves slowly enough such that points move quasistatically, and points in contact can either stick to the armor slip along the surface of the arm;(iii) Independent motion. We assume there is no interaction
between points, and the belief over of the stiffness
coefficient of each point is independent (P (Ki , Kj ) =P (Ki )P (Kj ) for i ̸= j);
This last assumption is the most restrictive compared to massspring and FEM models, but as we shall see it allows us toperform estimation of the stiffness field extremely quickly andis a reasonable approximation for some objects. Also note that
we make no distinction between empty and occupied space;with a properly estimated VSF, a point pi in empty spaceshould have stiffness Ki = 0.
At each time step, we sequentially step the contact simulator, query the VSF, and update the VSF using a recursiveestimator. Firstly, the contact simulator will update the contact
set Wct and displacements uti given the new arm configuration
q t . Then, using uti and Kit , the volumetric stiffness field willpredict force response Fit for i ∈ Wc . The joint torques τ twill be a low-dimensional linear observation of Fit as well as
Kit . Finally, the recursive estimator will update the belief ofKit .
Contact simulator: The contact simulator receives therobot arm’s current configuration q t at every time step. It
also maintains an internal state of contact points Wct−1
and displacements {ut−1}i∈Wct−1 in the previous time step.
iBecause of the independence assumption (iii), the simulator
will update the displacement of each point in parallel. Thesimulator simulates three effects: (1) contact detection by
distance threshold δ, (2) stick motion with friction coefficientµ, (3) slip motion in tangential direction with fixed step size λin each iteration. Here δ, µ and λ are the hyper-parameters of
the simulator. The simulator performs constant-time updates toeach point and the time complexity is bounded by O(|Wct |).
We assume the robot arm moves quasi-statically (ii), so
the external torques τ can be obtained by subtracting thegravitational torques from total torques. The relation between
joint torques and force responses is given by:Xτt =J(q t , pi + uti )⊤ Fit ,(2)i∈Wctwhere J(q , pi + uti ) is the contact Jacobian matrix that maps
joint velocity to the translational velocity at pi + uti .

(a) Configuration 1(b) Configuration 2

Fig. 2. Two configurations of the robot arm touching the same location withthe same end-effector pose.

Because of the independence assumption (iii), the belief ofeach κi is independent, which means P (κi , κj ) = P (κi ) ·P (κj ) if i ̸= j. It means we assume a diagonal covariancematrix with no stiffness estimation correlation between different points.
To perform a Bayesian update of the belief, we need atransition model of κi and an observation model of τ t . We
assume a static transition model with growing uncertainty:
κt+1= κti + ϵκ ,i(4)where ϵκ ∼ N (0, σκ2 ) and σκ2 is a pre-defined parameter.
We also assume there is independent Gaussian noise ϵτ ∼N (0, στ2 I) in torque measurements at every time step. Combining equations 1 and 2, we haveX
τt = −J(q t , pi + uti )⊤ uti |κti | + ϵτ
(5)i∈WctBecause the absolute value function makes the observation
function non-linear, we use an extended-Kalman filter (EKF)
to estimate stiffness coefficients. The observation matrix Ht
can be derived by linearizing 5 and each column of Ht equalsto ∂τ t /∂κti ,∂τ t= −J(q t , pi + uti )⊤ uti · sgn(κti )∂κti(6)

In practice, we find the EKF has history forgetting problem
where torque prediction to previous configurations is inconsistent to previous measurements. We maintain a buffer of
EKF observation matrices H t and observed τ t . At each step,after we update the filter with current observation, we replay afixed number of randomly sampled past measurements and donormal Kalman update for each sample. We find this technique
can mitigate the problem that filter forgets past observations,which also improves torque prediction accuracy on a testsequence.

t

Recursive estimator: Our goal is to have a probabilistic
estimation of Ki ∈ R+ . Because the stiffness is non-negative,
we cannot use a Gaussian belief which has support for negativevalues. We reparametrize the stiffness as observation of ahidden state Ki = |κi |, and κi ∈ R is unconstrained. Thebelief of each κi is Gaussian:κi ∈ N (µκi , σκ2 i )(3)III. E XPERIMENTS
We perform a real-world experiment using a Kinova Gen3
robot arm touching an artificial plant. We consider W as abox centered at the arm base with dimensions 2.4m × 2.4m ×2.4m, which fully contains reachable region by the arm. We
uniformly discretize this box into a grid with 2003 points andthe distance between neighbor points is 1.2cm.

3TABLE I
T ORQUE PREDICTION STATISTICS
MLP

0.50 ± 0.39
1.41 ± 0.39
0.15 ± 0.12
0.43 ± 0.21

2.45 ± 2.89
8.50 ± 2.89
0.76 ± 0.45
1.59 ± 1.18

GP6.04 ± 3.10 5.16 ± 4.81
7.73 ± 3.60 11.98 ± 3.65
1.48 ± 1.84 1.30 ± 1.66
6.32 ± 3.64 5.35 ± 5.03

1412B. Stiffness extrapolation
Although our technique builds a stiffness map in a relatively
uninformed and localized manner, we consider how a small
number of touches can be used to generalize across theunexplored volume. Certainly, touching the object all over toobtain a complete stiffness map would be time-consuming. We
examine a possible approach that uses correlations in visualfeatures to perform this generalization. Here we take onlya single depth image of the artificial plant before touches
happen. The left part of Fig. 4 shows raw stiffness map
obtained from 8 touches to similar locations on the plant.
The configuration of the arm is shown as transparent in thefigure. The raw stiffness map has spurious stiffness estimation||A. Stiffness estimation for torque prediction
We first evaluate whether an estimated VSF can accurately
predict joint torques when the same location is touched but
a different configuration. Fig. 2 shows the robot arm touches
the plant with the same end-effector pose but different joint
configurations. We touch the plant from 11 different configurations and collect 400 joint angles and torques measurementsin each touch. To evaluate our method, we take a k-foldcross validation approach. We use each trajectory as the test
sequence and estimate a volumetric stiffness field from theremaining 10 trajectories. Note that we do not update thestiffness estimation while evaluating on the test sequence.
We compare our method to three baseline methods. A
two-layer multilayer perceptron(MLP) with 128 hidden units
trained with Adam[14] optimizer. Gaussian Process(GP) regression predicts joint torques at each joint from joint angles
using an RBF kernel. k-nearest neighbors (k-NN) regression
using N = 50 nearest neighbors and an RBF kernel weights.
Evaluation metrics: We denote predicted torques as τ̂ andthe ground truth torques as τ . We measured both the absolute
torque prediction error ∥τ̂ − τ ∥2 and the relative torqueprediction error ∥τ̂ − τ ∥2 /∥τ ∥2 . The relative torque prediction
error is recorded when ∥τ ∥2 ≥ 1 Nm, which indicates the arm
has touched the plant. We compute the average and maximum
error in each sequence and report the mean and standard
deviation of error over the 11 sequences in Table I.
Results: Our method has consistently smaller prediction
error compared to the baseline methods as shown in Table I.
By explicitly modeling a stiffness field, we can accurately
predict joint torques when the robot arm touches the samelocation. Baseline methods that only use similarity in joint
angles are insufficient to predict joint torques. We also plot
the changes of torque prediction error in one test trajectoryin Fig. 3. We may see the GP method has incorrect non-zero
torques prediction when the arm is not in contact.

VSFMLPGPk-NN

16k-NN||2 (Nm)

Avg ∥τ̂ − τ ∥2 (N/m)
Max ∥τ̂ − τ ∥2 (N/m)
Avg ∥τ̂ − τ ∥2 /∥τ ∥2Max ∥τ̂ − τ ∥2 /∥τ ∥2VSF

10864
20050

100 150 200 250 300 350 400
Time step

Fig. 3. Torque prediction error by different methods in one test trajectory.

Zoom in view
touched region

Spurious
stiffness
estimation

Touched region

Fig. 4. Visualizing qualitative stiffness field estimation. The color bar at
the top indicates the predicted stiffness. Left figure shows raw estimation ofpoints with stiffness greater than 2 N/m. There is spurious stiffness estimation
because of estimation ambiguity from limited number of touches. Right figure
shows stiffness prediction to every point in the plant using color correlation.

in the empty space, because of limited number of touches fromsimilar configurations.
For touched points in the point cloud we look up thestiffness values of nearby touched voxels, creating a datasetthat maps RGB value to estimated stiffness. For all untouchedpoints, we use k-NN regression to predict the stiffness valuefrom the color. We use k = 50 neighbors and an RBF kernel.
From Fig. 4 right, we can see even with this relatively naïveapproach, the stiffness predictions align with our intuition thatbranches should be stiffer than leaves. In our future work, weare interested in exploring how other features can be used tobetter generalize stiffness estimation across the whole object.
IV. C ONCLUSION
In this paper, we introduce VSF, a new deformable object
stiffness model which can be estimated online using an EKF.
Our real world experiments show an estimated VSF canbe used to predict joint torques when the same location istouched. We also show the VSF estimated from a small
number of touches can be extrapolated to untouched regions
using color correlation.

4R EFERENCES
[1]
[2]
[3]
[4]

[5]
[6]

[7]
[8]

[9]
[10]
[11]
[12]
[13]
[14]

C.-A. Smarr et al., “Domestic robots for older adults: Attitudes,preferences, and potential,” Int. J. Social Rob., vol. 6, pp. 229–247,
2014.
E. Kayacan, Z. Zhang, and G. Chowdhary, “Embedded high precisioncontrol and corn stand counting algorithms for an ultra-compact 3d
printed field robot,” in Robotics: Science and Systems, 2018.
A. Jain, M. D. Killpack, A. Edsinger, and C. C. Kemp, “Reachingin clutter with whole-arm tactile sensing,” Int. J. Robotics Research,
vol. 32, no. 4, pp. 458–482, 2013.
B. Frank, C. Stachniss, R. Schmedding, M. Teschner, and W. Burgard,“Learning object deformation models for robot motion planning,”
Robotics and Autonomous Systems, vol. 62, no. 8, pp. 1153–1174,
2014.
V. E. Arriola-Rios, P. Guler, F. Ficuciello, D. Kragic, B. Siciliano, andJ. L. Wyatt, “Modeling of deformable objects for robotic manipulation:
A tutorial and review,” Frontiers in Robotics and AI, vol. 7, 2020.
B. Frank, R. Schmedding, C. Stachniss, M. Teschner, and W. Burgard,“Learning the elasticity parameters of deformable objects with amanipulation robot,” in IEEE/RSJ Int. Conf. Intel. Robots and Systems(IROS), 2010, pp. 1877–1883.
V. E. Arriola-Rios and J. L. Wyatt, “A multimodal model of object
deformation under robotic pushing,” IEEE Transactions on Cognitiveand Developmental Systems, vol. 9, no. 2, pp. 153–169, 2017.
S. Caccamo, P. Güler, H. Kjellström, and D. Kragic, “Active perception
and modeling of deformable surfaces using gaussian processes andposition-based dynamics,” IEEE-RAS Int. Conf. Humanoid Robots,
pp. 530–537, 2016.
Y. Zhu, K. Lu, and K. Hauser, “Semi-empirical simulation of learned
force response models for heterogeneous elastic objects,” in IEEE Int.
Conf. Robotics and Automation (ICRA), 2020, pp. 1646–1652.
L. Manuelli and R. Tedrake, “Localizing external contact using proprioceptive sensors: The contact particle filter,” in IEEE/RSJ Int. Conf.
Intel. Robots and Systems (IROS), 2016, pp. 5062–5069.
T. Pang and J. Umenberger, “Identifying external contacts from joint
torque measurements on serial robotic arms and its limitations,” IEEE
Int. Conf. Robotics and Automation, pp. 6476–6482, 2021.
B. Saund, S. Choudhury, S. S. Srinivasa, and D. Berenson, “Theblindfolded robot: A bayesian approach to planning with contact
feedback,” in ISRR, 2019.
J. Bimbo, A. S. Morgan, and A. M. Dollar, “Force-based simultaneous
mapping and object reconstruction for robotic manipulation,” IEEE
Robotics and Automation Letters, vol. 7, pp. 4749–4756, 2022.
D. P. Kingma and J. Ba, “Adam: A method for stochastic optimization,” in 3rd International Conference on Learning Representations,
ICLR 2015, San Diego, CA, USA, May 7-9, 2015, Conference Track
Proceedings, Y. Bengio and Y. LeCun, Eds., 2015. [Online]. Available:
http://arxiv.org/abs/1412.6980.