Reputation: 739
I have a data set with 4 columns: x,y,z, and value, let's say:
x y z value
0 0 0 0
0 1 0 0
0 2 0 0
1 0 0 0
1 1 0 1
1 2 0 1
2 0 0 0
2 1 0 0
2 2 0 0
I would like to calculate the center of mass CM = (x_m,y_m,z_m)
of all values. In the present example, I would like to see (1,1.5,0)
as output.
I thought this must be a trivial problem, but I can't find a solution to it in the internet. scipy.ndimage.measurements.center_of_mass
seems to be the right thing, but unfortunately, the function always returns two values (instead of 3). In addition, I can't find any documentation on how to set up an ndimage
from an array: Would I use a numpy array N of shape (9,4)
? Would then N[:,0] be the x-coordinate?
Any help is highly appreciated.
Upvotes: 11
Views: 50138
Reputation: 1694
Why did ndimage.measurements.center_of_mass
not give the expected result?
The key is in how the input data masses
was represented by an array of 4-tuples (x, y, z, value)
# x y z value
[[0, 0, 0, 0],
[0, 1, 0, 0],
[0, 2, 0, 0],
[1, 0, 0, 0],
[1, 1, 0, 1],
[1, 2, 0, 1],
[2, 0, 0, 0],
[2, 1, 0, 0],
[2, 2, 0, 0]]
The array masses
here represents the 3-D position and weights of each mass.
Note however that this python array structure is only a 2-D array. It's shape is (9, 4).
The input you need to pass to ndimage to get the expected result is a 3-D array containing zeros everywhere and the weight of each mass at the appropriate coordinates within the array, like this:
from scipy import ndimage
import numpy
masses = numpy.zeros((3, 3, 1))
# x y z value
masses[1, 1, 0] = 1
masses[1, 2, 0] = 1
CM = ndimage.measurements.center_of_mass(masses)
# x y z
# (1.0, 1.5, 0.0)
Which is exactly the expected output.
Note the limitation of this solution (and the ndimage library) is it requires non-negative integer coordinates. Also will not be efficient for large and/or sparse volumes because each "pixel" of the ndimage needs to be instantiated in memory.
Upvotes: 1
Reputation: 434
Another option is to use the scipy center of mass:
from scipy import ndimage
import numpy
masses = numpy.array([[0, 0, 0, 0],
[0, 1, 0, 0],
[0, 2, 0, 0],
[1, 0, 0, 0],
[1, 1, 0, 1],
[1, 2, 0, 1],
[2, 0, 0, 0],
[2, 1, 0, 0],
[2, 2, 0, 0]])
ndimage.measurements.center_of_mass(masses)
Upvotes: 5
Reputation: 25813
How about:
# x y z value
table = np.array([[ 5. , 1.3, 8.3, 9. ],
[ 6. , 6.7, 1.6, 5.9],
[ 9.1, 0.2, 6.2, 3.7],
[ 2.2, 2. , 6.7, 4.6],
[ 3.4, 5.6, 8.4, 7.3],
[ 4.8, 5.9, 5.7, 5.8],
[ 3.7, 1.1, 8.2, 2.2],
[ 0.3, 0.7, 7.3, 4.6],
[ 8.1, 1.9, 7. , 5.3],
[ 9.1, 8.2, 3.3, 5.3]])
def com(xyz, mass):
mass = mass.reshape((-1, 1))
return (xyz * mass).mean(0)
print(com(table[:, :3], table[:, 3]))
Upvotes: 2
Reputation: 2926
The simplest way I can think of is this: just find an average of the coordinates of mass components weighted by each component's contribution.
import numpy
masses = numpy.array([[0, 0, 0, 0],
[0, 1, 0, 0],
[0, 2, 0, 0],
[1, 0, 0, 0],
[1, 1, 0, 1],
[1, 2, 0, 1],
[2, 0, 0, 0],
[2, 1, 0, 0],
[2, 2, 0, 0]])
nonZeroMasses = masses[numpy.nonzero(masses[:,3])] # Not really necessary, can just use masses because 0 mass used as weight will work just fine.
CM = numpy.average(nonZeroMasses[:,:3], axis=0, weights=nonZeroMasses[:,3])
Upvotes: 16