marcomusy / vedo

A python module for scientific analysis of 3D data based on VTK and Numpy
https://vedo.embl.es
MIT License
2.03k stars 264 forks source link

Vectorized creation of ellipsoids #943

Open ttsesm opened 11 months ago

ttsesm commented 11 months ago

Hi @marcomusy,

is it possible to create multiple ellipsoids in a vectorized form at once without using a for loop?

I tried to use following code snippet, but it fails:

# Create a plotter
plotter = vd.Plotter()

# Number of ellipsoids to create
num_ellipsoids = 3

a_values = np.random.uniform(1.0, 2.0, size=(num_ellipsoids, 3, 3)).reshape(-1,3,3)
positions = np.random.uniform(-3, 3, size=(num_ellipsoids, 3))

# ellipsoids = vd.Ellipsoid(positions, a_values[0,:,:].reshape(1,-1), a_values[1,:,:].reshape(1,-1), a_values[2,:,:].reshape(1,-1))
#
# # Add all ellipsoids to the plot
# plotter += ellipsoids

# Create and add ellipsoids to the plot
for i in range(num_ellipsoids):
    ellipsoid = vd.Ellipsoid(positions[i], a_values[i, 0, :], a_values[i, 1, :], a_values[i, 2, :])
    plotter.add(ellipsoid)

# Show the plot
plotter.show()

Using the for loop above works fine, but if I try to pass the axes values all together I am getting the following error:

 ellipsoids = vd.Ellipsoid(positions, a_values[0,:,:].reshape(1,-1), a_values[1,:,:].reshape(1,-1), a_values[2,:,:].reshape(1,-1))
  File "/home/Development/......../lib/python3.9/site-packages/vedo/shapes.py", line 2888, in __init__
    angle = np.arcsin(np.dot(axis1, axis2))
ValueError: shapes (1,9) and (1,9) not aligned: 9 (dim 1) != 1 (dim 0)

Probably you could replace the np.dot with its vectorized form?

marcomusy commented 11 months ago

Hi @ttsesm sorry for the late reply! No, there is no way to do that. If you have very many ellipsoids (or any other mesh object) you can add them to a python list and them merge to treat them as a single (though disconnected) Mesh by using merge():

ellis = merge(my_list_of_ellipsoids)
ttsesm commented 11 months ago

I see, thanks. Well the idea was to avoid the creation of the ellipsoids with the for loop in first place which takes forever if I go to a number of a couple thousand ellipsoids that I have now.

The merge() might help to plot more ellipsoids since I've noticed that once I was plotting any number above 5k ellipsoids the plotting window was becoming kind of non-responsive.

marcomusy commented 11 months ago

This takes ~2.7s on my pc:

from vedo import *

centers = np.random.rand(5000,3)*50
#create ellipsoids

ellis = []
for center in centers:
    ell = Ellipsoid(center, res=12)
    ell.scale(np.random.rand(3))
    ellis.append(ell)
ellis = merge(ellis).lighting("off")
# show(ellis, axes=1).close()

Screenshot from 2023-10-10 15-57-15

you can lower a bit the resolution to 12 i guess?

ttsesm commented 11 months ago

Sure in this number of ellipsoids my machine works quite fine as well, but for my object I need to render < 100k ellipsoids :-p

marcomusy commented 11 months ago

I suggest you only render your ellipsoids in a smaller region of interest by clicking on it. Having 100k objects on the screen sounds a bit weird to me..

ttsesm commented 11 months ago

I suggest you only render your ellipsoids in a smaller region of interest by clicking on it. Having 100k objects on the screen sounds a bit weird to me..

Hahahaha... I know indeed it is. I am actually trying to plot the ellipsoids from a 3D Gassian Splatting model which are then used to render an image, since each ellipsoid contains also color information. The more the ellipsoids you have the higher quality of image you will get after rasterization.

So in principle, if you plot your ellipsoids this should look like a point cloud but instead of points you have ellipsoids in different shape, orientation and size. If you use only the center point of each ellipsoid then you can have a "normal" point cloud.

ttsesm commented 10 months ago

Hi Marco,

I will ask it here, so that not opening a new issue. Is it possible to create an ellipsoid where its surface split by a grid with cells of equal (or roughly equal) area?

I guess if you have something like in image bellow where you split the ellipsoid in n number of rings and then for each ring create a number of facets that should have more or less the same area. Δλ should be a bit more tricky to address though :thinking:

image

marcomusy commented 10 months ago

Yes - this is a starting point, but you still need compute the functional relation between the step and the circle heights that make the areas equal (but that should be simple geometry)

from vedo import *

circles = []
a3 = 1
for h in np.arange(0,a3, 0.1):
    z = cos(h/a3) ### need to find the right 
    r = h         ###  relationship between z and h
    circle = Circle(r=r, res=20).z(z)
    circle.lw(2).wireframe().color('k')
    circles.append(circle)

# create the mesh by merging the ribbon strips
rbs = []
for i in range(len(circles) - 1):
    rb = Ribbon(circles[i], circles[i+1], closed=True, res=[20,1])
    rb.lw(1).c("blue9")
    rbs.append(rb)
# mesh = merge(rbs).clean()

show(rbs, axes=1, viewup="z")

also keep in mind that if you want quad cells the 4 points must stay on the same plane (which is indeed the case in the above example)

ttsesm commented 10 months ago

Thanks Marco,

I will have a look on your snippet and I'll update you.