brainglobe / brainglobe-heatmap

Rendering anatomical heatmaps with brainrender and matplotlib
https://brainglobe.info
MIT License
34 stars 10 forks source link

Slicer crashes at specific planes #11

Closed LeonardoLupori closed 11 months ago

LeonardoLupori commented 2 years ago

Hi!
I'm creating various 2D planes out of a render by slicing a scene and I noticed a bizarre behavior. Slicing the scene at very specific planes makes the software stall. It does not raise any error but it never completes running. After trying in a notebook it makes the kernel die.

after a bit of trial and error, I narrowed it down to an example code that on my machine reproduces this behavior consistently. Specifically, it does not like the plane 11147 but it does like 11155, a few micrometers behind. This is just an example, but it happens randomly with very few planes, but it is still annoying when slicing a set of planes obtained with np.linspace since you never know whether the planes you want to slice will work or not.

I solved it by just changing the plane a bit, but I figured maybe it's worth investigating since it does not even throw an error.

import bgheatmaps as bgh
from bgheatmaps.slicer import Slicer

heatmap = bgh.heatmap(
    {'VISp':10,'fiber tracts':5},      # Dummy data
    position=5000,
    orientation="frontal",
    thickness=100,
    title="bugReport",
    cmap='Blues',
    vmin=0,
    vmax=10,
    format="2D",
)
heatmap.scene.close()

blessedPlane = 11155
cursedPlane = 11147

print(f"Slicing blessed plane...")
s = Slicer(blessedPlane, "frontal", 100, heatmap.scene.root)
projected, _ = s.get_structures_slice_coords(
    heatmap.regions_meshes, heatmap.scene.root
)
print("Done")

print(f"Slicing cursed plane...")
s = Slicer(cursedPlane, "frontal", 100, heatmap.scene.root)
projected, _ = s.get_structures_slice_coords(
    heatmap.regions_meshes, heatmap.scene.root
)
print("Done")

I'm on:

FedeClaudi commented 2 years ago

Hey, thanks for reporting this. I'm guessing it might complain when the plane is "tangent" to one of the brain regions. Is there anything special about the planes that you can see (e.g. using Planner)? What happens if you change which brainregions you're using? Do you know where in the code the error comes up?

LeonardoLupori commented 2 years ago

In the previous code, the structure 'fiber tracts' seems to be the one responsible for this behavior.

I also sliced my collection of areas (about 250 areas) and, if I specifically remove fiber tracts from the structures set I don't get this behavior anymore. This of course is not a comprehensive test, since I didn't test all the possible planes, but it seems that this is the structure that is causing the most problems.

I could not notice anything in particular about these planes, but by debugging the script I noticed that it stalls not at the instantiation of the Slicer object, but when the function get_structures_slice_coords() is called, so it must be something inside that.

FedeClaudi commented 2 years ago

Thanks for the update and the extra info. It sounds like there might be something funky in the mesh for fiber tracts - but hopefully it's not one that people would be using much!

carlocastoldi commented 1 year ago

I think these issues are relevant: https://github.com/marcomusy/vedo/issues/696 https://discourse.vtk.org/t/vtkintersectionpolydatafilter-crashing/9428

carlocastoldi commented 1 year ago

I managed to make it work reliably by not using bg-heatmaps.Slicer and backporting a binding function to vtx from vedo 2023.4.6.

Do you want me to open a PR so that we can continue discussing it over there?

adamltyson commented 11 months ago

Hi @LeonardoLupori, I think this should be fixed now following #21. The slightly updated example below doesn't cause any errors for me. Would you mind testing main sometime? A new release will be out soon.

import brainglobe_heatmap as bgh
from brainglobe_heatmap.slicer import Slicer

heatmap = bgh.heatmap(
    {'VISp':10,'fiber tracts':5},      # Dummy data
    position=5000,
    orientation="frontal",
    thickness=100,
    title="bugReport",
    cmap='Blues',
    vmin=0,
    vmax=10,
    format="2D",
)

blessedPlane = 11155
cursedPlane = 11147

print(f"Slicing blessed plane...")
s = Slicer(blessedPlane, "frontal", 100, heatmap.scene.root)
projected, _ = s.get_structures_slice_coords(
    heatmap.regions_meshes, heatmap.scene.root
)
print("Done")

print(f"Slicing cursed plane...")
s = Slicer(cursedPlane, "frontal", 100, heatmap.scene.root)
projected, _ = s.get_structures_slice_coords(
    heatmap.regions_meshes, heatmap.scene.root
)
print("Done")