Finally, I managed to plot the mesh using the native Makie functions. The key point was to convert Meshes.SimpleMesh
to the coordinates and faces recognizable by Makie.mesh
function. Though this conversion is present in the source code of Meshes
, it is buried inside a sophisticated (for me) code. I found it in its original form here.
The final solution is the following:
import GLMakie as mak
import Interpolations
import MarchingCubes
import Meshes
# Prepare data:
x = collect(range(-10, 10, 101))
y = collect(range(-10, 10, 101))
z = collect(range(-10, 10, 101))
G = [zi <= 0 || sqrt(xi^2+yi^2+zi^2) <= 5 ? 1.0 : 0.0 for xi=x, yi=y, zi=z]
F = [cos(sqrt(xi^2+zi^2)) for xi=x, yi=y, zi=z]
# Extract the isosurface of G for a given isovalue:
isovalue = 1
mc = MarchingCubes.MC(G, Int; x, y, z)
MarchingCubes.march(mc, isovalue)
# Prepare the mesh using Meshes:
mesh_points = @. Meshes.Point3(mc.vertices)
mesh_triangles = @. Meshes.connect(Tuple(mc.triangles), Meshes.Triangle)
mesh = Meshes.SimpleMesh(mesh_points, mesh_triangles)
# Extract coordinates and faces recognizable by GLMakie.mesh:
mesh_vertices = Meshes.vertices(mesh)
xyz = @. collect(Meshes.coordinates(mesh_vertices))
xyz = transpose(reduce(hcat, xyz))
mesh_elements = Meshes.elements(Meshes.topology(mesh))
faces = @. collect(Meshes.indices(mesh_elements))
faces = transpose(reduce(hcat, faces))
# Convert F to colors by interpolating it onto the coordinates of the mesh vertices:
itp = Interpolations.linear_interpolation((x, y, z), F)
Fcolors = [itp(xyz[i,1], xyz[i,2], xyz[i,3]) for i=1:size(xyz,1)]
# Visualize:
mak.set_theme!(lightposition=mak.Vec3f(1,1,1), ambient=mak.RGBf(1,1,1))
fig = mak.Figure()
ax = mak.Axis3(fig[1,1]; aspect=:data)
mak.mesh!(ax, xyz, faces; color=Fcolors)
mak.display(fig)
# mak.save("out.png", fig)
Thank you everyone for your help!