[ODE] Work done by torque in joints

[Krister Wolff]

Hi,
In my application, I set Torque and Desired Velocity values to joints in 
order to control the motions of the corresponding bodies. E.g.
dJointSetHingeParam(jointID, dParamVel, v);
dJointSetHingeParam(jointID, dParamFMax, 10.0);
where 'v' is some suitable value computed by my controller.
For some reason, I need to compute the Work done by Torque in each joint 
every time step (generally, mechanical work is the path integral of the 
force => work in circular plane motion can be obtained as 'fabs(Torque) 
* delta_phi', where phi is the joint angle).
Currently, I do this by multiplying the change of a joint angle over one 
time step with the dParamFMax value. However, this gives only a ruff 
estimate of the work done by the force (torque) over one timestep. 
Simply because that value, dParamFMax, is not applied in every time 
step, it is an upper bound! Most likely, many times the applied force 
could be much smaller than the dParamFMax value...
(however, one can assume that if delta_phi is 'large', the applied force 
was probably 'large', and vice versa...)
Right?
So, how to get the correct value of the applied force in every time step 
out of ode? I've been looking into the dJointFeedback structure, but I'm 
not really shure if this should give me what I want. I find its 
documentation a bit 'short', and further it doesn't seem to distinguish 
between axis1 and axis2 of the Universal joint (I'm using both Unversal 
and Hinge), which is necessary for me. I need to know the explicit 
Torque along both axis1 and axis2, as well as delta_phi of axis1 and 
axis2(that one is of course solved already), in order to correctly 
compute the work for an Universal joint. The mentioned torque values are 
present (?) in the simulation, but how to explicitly read them out?
Has anyone an idea how to do this?
Best,
Krister Wolff