Hello,
I would like to port a matlab library which provides functions for rigid body mechanics such as operations on homogeneous matrices (in SE(3)), twists (in se(3)) and so.
In matlab, the library worked on 3d matrices: n homoneous matrices were stacked along the 3d dimension. This speeded up computations on multiple matrices, at the price that native matrix opertators such H1 * H2 did not work (it is not defined for 3d matrices).
In this spirit, in numpy a set of rotation matrices could be built in the following way:
def rotx(theta): """ SE(3) matrices corresponding to a rotation around x-axis. Theta is a 1-d array """ costheta = np.cos(theta) sintheta = np.sin(theta) H = np.zeros((theta.size,4,4)) H[:,0,0] = 1 H[:,3,3] = 1 H[:,1,1] = costheta H[:,2,2] = costheta H[:,2,1] = sintheta H[:,1,2] = sintheta return H
I'm now seeking advices regarding an implementation with numpy (it's my first time with numpy).
- Is there any difference between a 3d array and a 1-d array of 2-d arrays (seems not)
- would you use a 3d-array or a list of 2d arrays or anything else to store these matrices ?
- Is there a way to work easily on multiple matrices of the same shape ? (in single-instruction/multiple-data spirit) for instance if A.shape==(3,2,4) and B.shape=(3,4,1), what is the best way to compute C for i = (0,1,2) C[i,:,:] = np.dot(A[i,:,:],B[i,:,:])
- is it a good idea/not to hard to subclass ndarray with a HomogenousMatrix class ? (one could then redefine inv() for instance)
- Is my implementation of rotx efficient ?
Thank you for any help, and please feel free to send me back reading anything I missed.
regards
On Fri, Nov 21, 2008 at 2:57 AM, Sébastien Barthélemy barthelemy@crans.orgwrote:
Hello,
I would like to port a matlab library which provides functions for rigid body mechanics such as operations on homogeneous matrices (in SE(3)), twists (in se(3)) and so.
In matlab, the library worked on 3d matrices: n homoneous matrices were stacked along the 3d dimension. This speeded up computations on multiple matrices, at the price that native matrix opertators such H1
- H2 did not work (it is not defined for 3d matrices).
In this spirit, in numpy a set of rotation matrices could be built in the following way:
def rotx(theta): """ SE(3) matrices corresponding to a rotation around x-axis. Theta is a 1-d array """ costheta = np.cos(theta) sintheta = np.sin(theta) H = np.zeros((theta.size,4,4)) H[:,0,0] = 1 H[:,3,3] = 1 H[:,1,1] = costheta H[:,2,2] = costheta H[:,2,1] = sintheta H[:,1,2] = sintheta return H
I'm now seeking advices regarding an implementation with numpy (it's my first time with numpy).
- Is there any difference between a 3d array and a 1-d array of 2-d
arrays (seems not)
- would you use a 3d-array or a list of 2d arrays or anything else to
store these matrices ?
- Is there a way to work easily on multiple matrices of the same shape
? (in single-instruction/multiple-data spirit) for instance if A.shape==(3,2,4) and B.shape=(3,4,1), what is the best way to compute C for i = (0,1,2) C[i,:,:] = np.dot(A[i,:,:],B[i,:,:])
There is work going on along these lines called generalized ufuncs, but it isn't there yet. The easiest current workaround is to use lists of matrices or arrays. Or you can get clever and do the matrix multiplications yourself using the current numpy operators, i.e.
C = (A[:,:,:,newaxis]*B[:,newaxis,:,:]).sum(axis=-2)
But this uses more memory than the list of matrices approach.
Chuck
2008/11/21 Sébastien Barthélemy barthelemy@crans.org:
In this spirit, in numpy a set of rotation matrices could be built in the following way:
def rotx(theta): """ SE(3) matrices corresponding to a rotation around x-axis. Theta is a 1-d array """ costheta = np.cos(theta) sintheta = np.sin(theta) H = np.zeros((theta.size,4,4)) H[:,0,0] = 1 H[:,3,3] = 1 H[:,1,1] = costheta H[:,2,2] = costheta H[:,2,1] = sintheta H[:,1,2] = sintheta return H
Btw, you can just create an array of these elements directly:
np.array([[costheta, -sintheta, 0], [sintheta, costheta , 0], [0 , 0 , 1]])
Regards Stéfan
def rotx(theta): """ SE(3) matrices corresponding to a rotation around x-axis. Theta is a 1-d array """ costheta = np.cos(theta) sintheta = np.sin(theta) H = np.zeros((theta.size,4,4)) H[:,0,0] = 1 H[:,3,3] = 1 H[:,1,1] = costheta H[:,2,2] = costheta H[:,2,1] = sintheta H[:,1,2] = sintheta return H
Btw, you can just create an array of these elements directly:
np.array([[costheta, -sintheta, 0], [sintheta, costheta , 0], [0 , 0 , 1]])
Are you sure ? Here it reports ValueError: setting an array element with a sequence. probably because theta, sintheta and costheta are 1-d arrays of n>1 elements.
However, I'll drop them and use list of matrices as Charles suggested.
Thanks to both if you.
2008/11/24 Sébastien Barthélemy barthelemy@crans.org:
Are you sure ? Here it reports ValueError: setting an array element with a sequence. probably because theta, sintheta and costheta are 1-d arrays of n>1 elements.
Sorry, I missed that detail.
Cheers Stéfan