I have two numpy arrays - one containing coordinates of 3D points and the second with triangles composed from those points. I need to calculate top surface 2D area and volume of this triangles (tetrahedrons?) with base height of the lowest point.
Minimal working example of what I got. It work OK, but is slow.
import numpy as np
pts = np.array([[744, 547, 695], [784, 511, 653], [779, 546, 746], [784, 489, 645], [834, 423, 614], [619, 541, 598]])
trs = np.array([[1, 2, 3], [2, 3, 4], [3, 4, 5], [4, 5, 0], [5, 0, 1]])
h_min = np.min(pts[:, 2])
print(f"H min : {h_min} m", )
a = 0
v = 0
for tr in trs:
p1, p2, p3 = tr
x1, y1, z1 = pts[p1]
x2, y2, z2 = pts[p2]
x3, y3, z3 = pts[p3]
area2d = abs(0.5 * (((x2 - x1) * (y3 - y1)) - ((x3 - x1) * (y2 - y1))))
h_mean = (z1 + z2 + z3) / 3 - h_min
v_tr = area2d * h_mean
a += area2d
v += v_tr
print(f"AREA : {a} m2")
print(f"VOLUME: {v} m3")
The problem is that really those arrays contains millions of points and triangles and this is calculating way too long. I have found method called numpy vectorization, but have no idea how to make it work in my case. Can anyone explain to me if it is even possible? Fast volume calculation is a must. Area - would be great. Thank you!
You could use vectorized operations. Way efficient compared to the method you have
[x1, x2, x3], [y1,y2,y3], [z1,z2,z3] = pts[trs].T
area2d = abs(0.5 * (((x2 - x1) * (y3 - y1)) - ((x3 - x1) * (y2 - y1))))
a = area2d.sum()
v = area2d @ ( (z1 + z2 + z3) / 3 - h_min)
print(f"AREA : {a} m2")
print(f"VOLUME: {v} m3")