simtk.org

Dear Opensim community,

I am trying to extract a body (in this example the tibia_r) angular velocity using the SimBodyEngine on Python 2.7 without going through the Analyze Tool.
- The model is the Gait2354
- The states file consists of joint coordinates extracted from IK
- I realise the state variable derivatives
- the COM_global is correctly computed
- However the value obtained for the w_local is all zeroes

Has anyone faced this issue or am I missing a crucial step ?

Thank you very much for your help!

Best Regards

This is the Python code that I am running:

##############################################
model = opensim.Model('Gait2354.osim')
bodyset = model.getBodySet()
bodypart = bodyset.get('tibia_r')
states = 'IK_joint_coordinates.mot'

### COM in local
COM_local = opensim.Vec3(0)
COM_local = bodypart.getMassCenter()
ground_body = model.getGround()

MotStorage = opensim.Storage(states)
n_coord = model.getNumCoordinates()

COM1 = []
COM2 = []
COM3 = []
W1 = []
W2 = []
W3 = []

if MotStorage.isInDegrees():
r = radians
else:
r = lambda x: x

## Realize the state
for n_frame in range (0 , MotStorage.getSize()):

## initialize state
state = model.initSystem()

## loop through coordinate set
for i in range(n_coord):
name = model.getCoordinateSet().get(i).getName()
index = MotStorage.getStateIndex(name)
vector = MotStorage.getStateVector(n_frame).getData().get(index)
model.updCoordinateSet().get(i).setValue(state,r(vector))

## generate vector for COM in global
COM_global = opensim.Vec3(0)
w_local = opensim.Vec3(0)

## transform COM to global
model.getSimbodyEngine().transformPosition(state, bodypart, COM_local, ground_body, COM_global)

## realize to the acceleration stage
model.computeStateVariableDerivatives(state)

## get the angular velocity
model.getSimbodyEngine().getAngularVelocityBodyLocal(state, bodypart, w_local);

COM1.append([COM_global.get(0)])
COM2.append([COM_global.get(1)])
COM3.append([COM_global.get(2)])
W1.append([w_local.get(0)])
W2.append([w_local.get(1)])
W3.append([w_local.get(2)])
COM_glob = np.hstack((COM1,COM2,COM3))
W_loc = np.hstack((W1,W2,W3))

Maybe the initializeState call is causing the problem. You should initialize the state only once before the first for loop.

Dear Dimitar,

Thank you for taking the time to reply to my questions. Unfortunately even if the state initialisation takes place before the first loop, the results for the angular velocity do not change and remain zero. When trying to compute afterwards:
model.getSimbodyEngine().getAcceleration(state, bodypart, COM_global, acc)

I do manage to extract some values for the acceleration but the angular velocity always remains zero. The angular acceleration also delivers zero values.

It looks like you are setting the positions of the coordinates but leaving the velocities set to zero. I would expect to see something like https://simtk.org/api_docs/opensim/api_ ... dd7a8febd3

Also, in general, you should enforce the position-level constraints only once all coordinates have been set. Please see the comments for the Coordinate::setValue() method: https://simtk.org/api_docs/opensim/api_ ... f2e88d7c42.

Dear Thomas,

Thank you very much for taking the time to reply! Indeed, this was missing. After setting the speed value, I am able to extract the angular velocity. Important thing is to differentiate the coordinates. Here are the few lines of code I added if someone is interested in the same operations:

speed = np.diff(vector,0)
model.updCoordinateSet().get(i).setValue(state,r(vector))
model.updCoordinateSet().get(i).setSpeedValue(state,r(speed))

Hello,

Thanks for this post. I am trying to get whole body CoM (center of mass) during a dynamic motion using Maltab Scripts.


I have done the following in Matlab and then I did Body Kinematics Analysis to get CoM position. However, the outcomes for CoM are different and Maltab Script takes a long time to calculate CoM. Here is what I used in Maltab. I wonder if you could help me understand where I am making the mistake. Thank you.

import org.opensim.modeling.*
model_folder = '\OpenSim\4.0\Code\Matlab';
model = Model(fullfile(model_folder, 'subject01_simbody.osim'));
state = model.initSystem();

motion_storage = Storage("W1_ik.mot");MotStorage = motion_storage;
motion_storage.getSize();
n_coord = model.getNumCoordinates();

% ## Realize the state
for n_frame = 0:MotStorage.getSize()-1

% ## initialize state
% state = model.initSystem();

for i=0:n_coord-1

name = model.getCoordinateSet().get(i).getName();
index = MotStorage.getStateIndex(name);
vector = MotStorage.getStateVector(n_frame).getData().get(index);
model.updCoordinateSet().get(i).setValue(state,vector);
model.computeStateVariableDerivatives(state);

end

p = model.calcMassCenterPosition(state);
% COMLoc(i+1,:,n_frame+1) = str2num(char(ArrayDouble.getValuesFromVec3(p)));
COMLoc(n_frame+1, = str2num(char(ArrayDouble.getValuesFromVec3(p)));

end

You can calculate the total CoM directly from the model:

https://simtk.org/api_docs/opensim/api_ ... 8a90aaf4cd

Thanks Dimitar! I think I am using the same method you kindly shared. My problem were the following:

1. Why the results are different between this script (that I shared before here) and BodyKinematics analysis in OpenSim GUI?
2. Also, scripting takes a few second to calculate CoM while BodyKinematics analysis is much faster.

I appreciate any ideas.

Cheers,
Hossein

I think the problem is in the setValue, because vector variable is probably in degrees, while Simbody/OpenSim operates in radiance.

The execution time is different because the analyze tool is compiled and optimized for speed, whereas the Matlab script file is possibly not or in some case uses JIT compilation to optimize portions of code.