I would like to use Python to get the path coordinates of the shortest (geodesic) path on a surface in 3D between two points. The solution is somehow already shown in the accepted answer in this post, however, I am not able to reproduce the result of Fedor, since I don't understand the format that is expected by meshlib to provide the point locations.
The geodesic path calculation is done with the computeGeodesicPath function which expects MeshTriPoints as input. Unfortunately, I don't understand how this format works and how I can convert to it from cartesian point coordinates. Here is an example for a simple cylindrical 3D surface, but I want it to be able to work also on more complex surfaces.
import numpy as np
import meshlib.mrmeshnumpy as mnp
def cyl2cart(rho,phi,z):
return rho*np.cos(phi),rho*np.sin(phi),z
# cylinder coordinates
N = 101
radius = 5
rho = radius*np.ones(N)
phi = np.linspace(0,2*np.pi,37) # 10 degree steps
z = np.linspace(-10,10,N)
# surface grid coordinates
x = np.outer(rho,np.cos(phi))
y = np.outer(rho,np.sin(phi))
z = np.outer(z,np.ones(len(phi)))
# create mesh
mesh = mnp.meshFromUVPoints(x, y, z)
# define two points on the surface
# point1
rhop = radius
phip = -10/180.*np.pi
zp = -3
xp1,yp1,zp1 = cyl2cart(rhop, phip, zp)
# point2
rhop = radius
phip = 60/180.*np.pi
zp = 8
xp2,yp2,zp2 = cyl2cart(rhop, phip, zp)
# This part is not working because I don't understand how to provide the point coordinates
# p1 = mpy.MeshTriPoint(??)
# p2 = mpy.MeshTriPoint(??)
# mpy.computeGeodesicPath(mesh, p1, p2)
#%% optional plot with plotly
import plotly.graph_objects as go
# extract numpy arrays
verts = mnp.getNumpyVerts(mesh)
faces = mnp.getNumpyFaces(mesh.topology)
# prepare data for plotly
vertsT = np.transpose(verts)
facesT = np.transpose(faces)
# draw
fig = go.Figure(data=[
go.Mesh3d(
x = vertsT[0],
y = vertsT[1],
z = vertsT[2],
i = facesT[0],
j = facesT[1],
k = facesT[2]
)
])
fig.add_trace(go.Scatter3d(x=[xp1,xp2],y=[yp1,yp2],z=[zp1,zp2],mode='markers',marker={'size':8,'symbol':'circle'}))
fig.update_layout(scene=dict(aspectmode='data'))
fig.show(renderer='browser')
mpy.MeshTriPoint
is a point that is attached to mesh object. To get mpy.MeshTriPoint
from cartesian point, project world point to mesh and find corresponding MeshTriPoint:
from meshlib import mrmeshpy as mm
mesh = mm.makeSphere( mm.SphereParams() ) # simple mesh for example
startPoint = mm.Vector3f(-1,-1,-1) # cartesian start coordinate
stopPoint = mm.Vector3f(1,1,1) # cartesian stop coordinate
# find corresponding MeshTriPoints
startMtp = mm.findProjection( startPoint, mesh).mtp
stopMtp = mm.findProjection( stopPoint, mesh).mtp
path = mm.computeGeodesicPath(mesh,startMtp,stopMtp,mm.GeodesicPathApprox.DijkstraBiDir)
# print path
for p in path:
pc = mesh.edgePoint(p)
print(pc.x,pc.y,pc.z)