Obtain CM position and velocity during a simulation

[Brian Umberger]

Hi All,

Using the Matlab interface to the API, it is possible to obtain the current velocity (or position) of the center of mass of any segment in a model during a simulation using code such as this:

massCenter = Vec3(0.0,0.0,0.0);
velocity = Vec3(0.0,0.0,0.0);
bodySet = osimModel.getBodySet();
bodySet.get('segment_name').getMassCenter(massCenter);
simbodyEngine = osimModel.getSimbodyEngine();
simbodyEngine.getVelocity(osimState,osimModel.getBodySet().get('segment_name'),massCenter,velocity);

My question is how best to obtain the position and velocity of the whole body center of mass during a simulation (something besides querying all of the body segments and then calculating a mass-weighted average). I found a nice C++ example that accomplished this using getMultibodySystem, getMatterSubsystem, and calcSystemMassCenterLocationInGround, but I could not get anything like that to work in Matlab. Any suggestions or insights are greatly appreciated!

Brian

[jimmy d]

Hi Brian,
You could use the class 'bodykinematics' and its subclass 'setRecordCenterOfMass'. The Doxygen documentation can be found here; https://simtk.org/api_docs/opensim/api_ ... atics.html

Let us know if that helps,
-james

[Ajay Seth]

The easiest way to get the center of mass info for the whole model, is to ask the model.

SimTK::Vec3 Model::calcMassCenterPosition(const SimTK::State &s) const
SimTK::Vec3 Model::calcMassCenterVelocity(const SimTK::State &s) const
SimTK::Vec3 Model::calcMassCenterAcceleration(const SimTK::State &s) const

Note, when calling for mass center acceleration of the model you must have computed the model acceleration already or you will get an exception. This can be done from the model by calling:

void Model::computeStateVariableDerivatives(const SimTK::State &s) const

Please let us know if you run into any issues. Cheers.

[Brian Umberger]

Ajay and James,

Thanks so much for the responses. I will give this a try and let you know if I encounter any problems.

Best,
Brian

[Brian Umberger]

Last month, I asked about obtaining whole-body center of mass kinematics. This is in the context of running a forward simulation to evaluate the objective function for an optimization problem. I have one case where I only need to know the CM kinematics at the initial time and the final time. Ajay’s suggestion to use calcMassCenterPosition/Velocity worked perfectly in that case (Thanks, Ajay!)

I also have a case where, immediately after a simulation runs, I want to access the full time series of the CM kinematics to find certain extreme values (e.g., minimum position, maximum velocity). James’ suggestion to add a BodyKinematics analysis to the model seems to be the way to go for this case (Thanks, James!) I am able to get the results of the BodyKinematics analysis in an output file using bodyKinematics.printResults, but I am unclear how to get the results into a Matlab variable so that I can examine the CM position or velocity without having to load the file generated by printResults. In other words, I don’t want to write out the file, I’d rather access the values directly in my Matlab code.

I saw that the BodyKinematics class has functions called getPositionStorage and getVelocityStorage. I also noticed that the Storage class has functions such as getDataColumn, but I could not figure out from the doxygen documentation how to use these functions to get the CM kinematics into one or more Matlab variables (assuming they are indeed the right functions to use).

Any suggestions are greatly appreciated!

Thanks,
Brian

[shayan moradkhani]

Last month, I asked about obtaining whole-body center of mass kinematics. This is in the context of running a forward simulation to evaluate the objective function for an optimization problem. I have one case where I only need to know the CM kinematics at the initial time and the final time. Ajay’s suggestion to use calcMassCenterPosition/Velocity worked perfectly in that case (Thanks, Ajay!)

I also have a case where, immediately after a simulation runs, I want to access the full time series of the CM kinematics to find certain extreme values (e.g., minimum position, maximum velocity). James’ suggestion to add a BodyKinematics analysis to the model seems to be the way to go for this case (Thanks, James!) I am able to get the results of the BodyKinematics analysis in an output file using bodyKinematics.printResults, but I am unclear how to get the results into a Matlab variable so that I can examine the CM position or velocity without having to load the file generated by printResults. In other words, I don’t want to write out the file, I’d rather access the values directly in my Matlab code.

I saw that the BodyKinematics class has functions called getPositionStorage and getVelocityStorage. I also noticed that the Storage class has functions such as getDataColumn, but I could not figure out from the doxygen documentation how to use these functions to get the CM kinematics into one or more Matlab variables (assuming they are indeed the right functions to use).

Any suggestions are greatly appreciated!

Thanks,
Brian

hi,
in case you could figure out how to handle your posted issue, i would appreciate it very much helping me as well.
or if anyone has any idea about the same thing quoted here, that would be helpful to share.

in doxygen i found "SimTK::Vec3 OpenSim::Model::calcMassCenterVelocity ( const SimTK::State & s) const"
but i cant get it correct in matlab; i typed:"v=foot1.getMultibodySystem().calcMassCenterVelocity(s);" and i got "Undefined function 'calcMassCenterVelocity' for input arguments of type'org.opensim.modeling.SWIGTYPE_p_SimTK__MultibodySystem'

or i found the following:
velocity=Vec3();
masscenter=Vec3();
simbodyEngine = foot1.getSimbodyEngine();
SimbodyEngine.getVelocity(s, body, masscenter, velocity);

but when i run it:"No method 'getVelocity' with matching signature found for class 'org.opensim.modeling.SimbodyEngine'." for this i checked the doxygen and arguments seemed to be correct:"irtual void getVelocity (const SimTK::State &s, const OpenSim::Body &aBody, const SimTK::Vec3 &aPoint, SimTK::Vec3 &rVel) const "

for my optimization, i need to get the velocity at each time-interval to do the rest of the studies. i have already used pointkinematics to store the kinematics data of CoM, but is there any other way to get the CoM velocity through the simulation and not after the simulation in matlab?
regards
Shayan

[jimmy d]

Code: Select all

import org.opensim.modeling.*

model = Model('gait2392_simbody.osim');
s = model.initSystem()

% returns Vec3()
p = model.calcMassCenterPosition(s);
v = model.calcMassCenterVelocity(s);
a = model.calcMassCenterAcceleration(s);

You should go through the Dynamic Walking example.

[Monica Armengol]

Hi James,

Thanks for your answer. I am a beginner using APIs in opensim. How do you update the states "s" to be able to know how position, velocity and acceleration change over time?. I imagine you need to insert your code into a loop and link it to IK/ID, but I don't know how to do that.

Thanks,

Monica

[Hossein Mokhtarzadeh]

Hi all,

I am sorry this may be long.

I am not sure where I am making the mistake in getting whole body CoM position during a simulation (by loading Inverse Kinematics results).

I am trying two methods. My problem is that I cannot get the same results as I do with BodyKinematics analyses.

1. Using "StatesTrajectory.createFromStatesStorage"

import org.opensim.modeling.*
model_folder = '.\OpenSim\4.0\Code\Matlab';
model = Model(fullfile(model_folder, 'subject01_simbody.osim'));
cd('.\IKResults');
model.initSystem();
filename = "W101_ik.mot";
sto = Storage(filename); % sto = STOFileAdapter.read(filename);
sto.isInDegrees() % is it in Degree
sto.setInDegrees(0) % so make it Radians I guess
states = StatesTrajectory.createFromStatesStorage(model, sto, true);
p = model.calcMassCenterPosition(states);

Using this first method, I get an error in using calcMassCenterPosition method. I guess I need to initialise states at some point but not sure how! Am I right?


2. By loading states (joint angles at each frames) and then using model.calcMassCenterPosition(state)

import org.opensim.modeling.*
model_folder = '.\4.0\Code\Matlab';
model = Model(fullfile(model_folder, 'subject01_simbody.osim'));
state = model.initSystem();
cd('.\IKResults');
motion_storage = Storage("W301_ik.mot");
motion_storage.getSize();
n_coord = model.getNumCoordinates();

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

for i=0:n_coord-1

name = model.getCoordinateSet().get(i).getName();
index = motion_storage.getStateIndex(name);
vector = motion_storage.getStateVector(n_frame).getData().get(index);
model.updCoordinateSet().get(i).setValue(state,degtorad(vector));
% speed = diff(vector);
% model.updCoordinateSet().get(i).setSpeedValue(state,speed)
end

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

This method gives me quite small CoM. I am assuming the outcome of calcMassCenterPosition is in global coordinate. but this seems to be in local coordinate system.

I would greatly appreciate any comments.

[Umer Huzaifa]

Hossein,

I may have figured out using your Method 1 (thanks for sharing it!). Here is how I modified it:

Code: Select all

sto=Storage('leg_states.sto')
states = StatesTrajectory.createFromStorage(model, sto, true)   % my model was in the workspace so not defining that separately
Straj = StatesTrajectory()        % Probably dont need it but I am keeping it
states = Straj.createFromStatesStorage(model, sto, true)
model.calcMassCenterPosition(states.get(9))      % get() function is the key here. states itself is a collection of the states so cant apply
model.calcMassCenterPosition(states.get(12))    % directly to the model
model.calcMassCenterPosition(states.get(199))
model.calcMassCenterPosition(states.get(20))
model.calcMassCenterPosition(states.get(0))

I myself am not quite done with it yet because the states are coming out to be NANs. But I am pretty sure the rest of it is now in place.

Hope this helps anyone else still struggling with it.