# Visualize 3D polygon of brain regions given points representing exterior border
**URL:** https://discourse.julialang.org/t/visualize-3d-polygon-of-brain-regions-given-points-representing-exterior-border/101247
**Category:** Visualization
**Tags:** makie
**Created:** [July 6, 2023, 9:15am UTC](https://discourse.julialang.org/t/visualize-3d-polygon-of-brain-regions-given-points-representing-exterior-border/101247 "2023-07-06T09:15:15Z")
**Posts on this page:** 4
**Page:** 1
### Author:  [@vandenman](https://discourse.julialang.org/u/vandenman)
#### Post date: [July 6, 2023, 9:15am UTC](https://discourse.julialang.org/t/visualize-3d-polygon-of-brain-regions-given-points-representing-exterior-border/101247/1 "2023-07-06T09:15:15Z")
I’m trying to create a 3d plot of some brain regions. However, the data I have is not a mesh but a series of 3d points that indicate the _exterior border_ of the regions. How can I plot this with GLMakie?
Below is the code I’m using. The two data files, `ROI_MNI_V4_List.mat` and `ROI_MNI_V4_Border.mat` can be obtained from [EBRAINS - Knowledge Graph Search](https://search.kg.ebrains.eu/instances/Dataset/f8758eda-483e-45fe-8a88-a1fc806dde18) and then going to `Get data`.
```julia
using GLMakie
import Colors, GeometryBasics
import MAT # for reading matlab data files
# [5ae59095] Colors v0.12.10
# [e9467ef8] GLMakie v0.8.6
# [5c1252a2] GeometryBasics v0.4.7
# [23992714] MAT v0.10.5
# wherever you downloaded the aal atlas
root = joinpath(pwd(), "aal_for_SPM12")
roi_list = MAT.matread(joinpath(root, "ROI_MNI_V4_List.mat"))
roi_border = MAT.matread(joinpath(root, "ROI_MNI_V4_Border.mat"))
distinct_ids = unique(vec(Int.(roi_list["ROI"]["ID"])))
ids = eachindex(distinct_ids)
ids_indices = [findall(==(distinct_ids[id]), vec(roi_border["BORDER_V"])) for id in ids]
# a vector with the boundary regions for each id
border_points = [
GeometryBasics.Point3f.(eachcol(view(roi_border["BORDER_XYZ"], :, indices)))
for indices in ids_indices
]
mesh(border_points[1]) # visualize only the first region for now
```
The code above sort of plots something okay, but the mesh is not a volume but a disjoint set of triangles.
I tried a bunch of things, but I also noticed that e.g., `GeometryBasics.mesh` fails completely for this type of data and returns an empty Mesh:
```julia
julia> GeometryBasics.mesh(border_points[1])
Mesh{3, Float32, Triangle}
```
I figured that the border points essentially represent a 3d polygon, but that doesn’t help much either:
```julia
julia> polygon = GeometryBasics.Polygon(border_points[1])
julia> GeometryBasics.mesh(polygon)
ERROR: MethodError: no method matching earcut_triangulate(::Vector{Vector{Point{3, Float32}}})
Closest candidates are:
earcut_triangulate(::Vector{Vector{Point2{Float64}}})
@ GeometryBasics ~/.julia/packages/GeometryBasics/Du43H/src/triangulation.jl:181
earcut_triangulate(::Vector{Vector{Point{2, Float32}}})
@ GeometryBasics ~/.julia/packages/GeometryBasics/Du43H/src/triangulation.jl:189
earcut_triangulate(::Vector{Vector{Point2{Int64}}})
@ GeometryBasics ~/.julia/packages/GeometryBasics/Du43H/src/triangulation.jl:197
...
Stacktrace:
[1] faces(polygon::GeometryBasics.Polygon{3, Float32, Point{3, Float32}, GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}, Vector{GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}}})
@ GeometryBasics ~/.julia/packages/GeometryBasics/Du43H/src/triangulation.jl:224
[2] decompose(#unused#::Type{GeometryBasics.NgonFace{3, GeometryBasics.OffsetInteger{-1, UInt32}}}, primitive::GeometryBasics.Polygon{3, Float32, Point{3, Float32}, GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}, Vector{GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}}})
@ GeometryBasics ~/.julia/packages/GeometryBasics/Du43H/src/interfaces.jl:103
[3] decompose_triangulate_fallback(primitive::GeometryBasics.Polygon{3, Float32, Point{3, Float32}, GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}, Vector{GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}}}; pointtype::Type, facetype::Type)
@ GeometryBasics ~/.julia/packages/GeometryBasics/Du43H/src/meshes.jl:91
[4] decompose_triangulate_fallback
@ ~/.julia/packages/GeometryBasics/Du43H/src/meshes.jl:89 [inlined]
[5] mesh(polygon::GeometryBasics.Polygon{3, Float32, Point{3, Float32}, GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}, Vector{GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}}}; pointtype::Type, facetype::Type, normaltype::Nothing)
@ GeometryBasics ~/.julia/packages/GeometryBasics/Du43H/src/meshes.jl:161
[6] mesh(polygon::GeometryBasics.Polygon{3, Float32, Point{3, Float32}, GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}, Vector{GeometryBasics.LineString{3, Float32, Point{3, Float32}, Base.ReinterpretArray{GeometryBasics.Line{3, Float32}, 1, Tuple{Point{3, Float32}, Point{3, Float32}}, GeometryBasics.TupleView{Tuple{Point{3, Float32}, Point{3, Float32}}, 2, 1, Vector{Point{3, Float32}}}, false}}}})
@ GeometryBasics ~/.julia/packages/GeometryBasics/Du43H/src/meshes.jl:158
[7] top-level scope
```
For what it’s worth, if I visualize all regions then this does produce some rather artsy figures:
but as soon as you look at a different angle it’s clear that everything is messed up 😅
Does anybody have any ideas on how to fix this?
---
### Author:  [@johnh](https://discourse.julialang.org/u/johnh)
#### Post date: [July 6, 2023, 10:10am UTC](https://discourse.julialang.org/t/visualize-3d-polygon-of-brain-regions-given-points-representing-exterior-border/101247/2 "2023-07-06T10:10:24Z")
A quick search returns this
> **[3d](https://plotly.com/julia/3d-surface-plots/)**
>
> Detailed examples of 3D Surface Plots including changing color, size, log axes, and more in Julia.
---
### Author:  [@sdanisch](https://discourse.julialang.org/u/sdanisch)
#### Post date: [July 6, 2023, 10:38am UTC](https://discourse.julialang.org/t/visualize-3d-polygon-of-brain-regions-given-points-representing-exterior-border/101247/3 "2023-07-06T10:38:51Z")
I don’t think surface in plotly will do this, since it’s not a surface which is described by the points.
The simplest I could find is to compute the convex hull (which is probably not optimal, not sure what other algorithms there are to compute the triangles out of an unordered pointcloud).
This works:
```julia
function compute_hull(points3d::AbstractMatrix)
hull = chull(collect(points3d))
tris = [GLTriangleFace(hull.simplices[i, :]...) for i in 1:size(hull.simplices, 1)]
return normal_mesh(Makie.to_vertices(points3d), tris)
end
# a vector with the boundary regions for each id
border_meshes = [
compute_hull(view(roi_border["BORDER_XYZ"], :, indices)')
for indices in ids_indices
]
mesh(border_meshes; color=1:length(border_meshes), colormap=(:viridis, 0.5), transparency=true) # visualize only the first region for now
```
---
### Author:  [@vandenman](https://discourse.julialang.org/u/vandenman)
#### Post date: [July 6, 2023, 11:25am UTC](https://discourse.julialang.org/t/visualize-3d-polygon-of-brain-regions-given-points-representing-exterior-border/101247/4 "2023-07-06T11:25:19Z")
Thanks that works perfectly!
To save anybody trying to replicate this some time, the function `chull` comes from the library `QHull`.