FEM tutorial: Issues to fix

[ftadel]

1. SimNIBS: Remeshing the SimNIBS model with iso2mesh doesn't work. Here is the error message:

    FEM> 1. scalp: C:\Work\Protocols\TutorialFem\anat\Subject01\tess_scalp.mat
    FEM> 2. skull: C:\Work\Protocols\TutorialFem\anat\Subject01\tess_skull.mat
    FEM> 3. white: C:\Work\Protocols\TutorialFem\anat\Subject01\tess_white.mat
    FEM> 4.  gray: C:\Work\Protocols\TutorialFem\anat\Subject01\tess_gray.mat
    FEM> 5.   csf: C:\Work\Protocols\TutorialFem\anat\Subject01\tess_csf.mat
     
    generating tetrahedral mesh from closed surfaces ...
    creating volumetric mesh from a surface mesh ...

    ***************************************************************************
    ** Error: Line 117:  surf2mesh (line 117)
    ** Tetgen command failed
    **
    ** Call stack:
    ** >surf2mesh.m at 117
    ** >process_fem_mesh.m>Compute at 442

2. Pre-processing: The PSD screen capture does not seem to be from this file (not the same time definition)
3. Averaging: The averaging screen capture shows a weighted average: this is confusing => Uncheck the "weighted" box
4. Noise covariance:
   • The name of the epochs in the screen capture is wrong ("2" and not "STI101")
   • The is a mention to a DC offset removal that was not applied when importing the epochs
   • The time window doesn't match between the text and the screen capture.
5. Inverse model: Redo screen captures with correct average file name.
6. Dipole scanning:
   • Redo screen capture to show only the LINKS (we can't select the shared kernels in Process1, as the text said before I edited it)
   • The dipole scanning process is executed over 40-60ms, but the results are illustrated at 21ms...
   • On my end, peak is at 22ms, not 21ms
   • If you scan the dipoles at multiple time points, you need to explain in the tutorial how to navigate in time and/or overlay multiple time points. Otherwise, simply scan at one specific time sample (eg. 21ms), this is simpler to display and post-process.
   • I find the second figure confusing. It should be split in two: in this section, show only the merge+renaming of the dipoles. Then make another figure to illustrate the dipoles display, and move it to the following section (P20/N20)
7. Comparison: DTI vs isotropic:
   • Add text to explain that you change the coloring of the dipoles "by group"
   • Explain how to change center the slices of the MRI 3D on a dipole (click on the dipole, the click on the "View/3D" button in the toolbar of the Dipole info panel)
   • "When we compare between the EEG and the MEG, there is a difference of 20mm between the localization between the MEG and the EEG": It looks less than 2cm in your screen captures... It would be worth illustrating how to compute this correctly: click on each dipole to see their locations (add screen capture of the dipole info panel), and compute the norm of the difference of the two vectors (sqrt((x1-x2)^2+(y1-y2)^2+(z1-z2)^2)) from the Matlab command line to get the exact distance (screen capture of the Matlab command window).
   • You need to explain the granularity of these results: you are using the cortex as the source space, the output of the dipole scanning can only be one of the vertices of the cortex surface. The two methods lead to the selection a selection of difference vertices, but it doesn't mean that there is "an error of XX millimeters" between the two methods.
   • "as well as a difference in the orientation in the coronal view": No, this difference in orientation is clearly visible in the three views. You can also compute the angle between the two from the orientation vector read in the dipoles info panel.
8. Comparison: MEG vs EEG:
   • "If we compare the difference between the EEG and MEG dipoles": Do you merge again the dipoles in a different way for this? If so, you need to tell the user to repeat the operation (possibly, this should be added in the previous section "Dipole scanning", as this is where you did the other merge). => Or you can keep the section but skip the figure , as it doesn't bring any additional interesting information (left and right figures look too similar) - I'm not sure what you can see here that you can't see in the in the previous figure.
9. Comparison: FEM vs BEM:
   • I get different plots than yours for EEG FEM DTI vs. EEG BEM:
     
     
     
10. DUNEuro advanced options: Merged into the DUNEuro tutorial: https://neuroimage.usc.edu/brainstorm/Tutorials/Duneuro#DUNEuro_options:_Advanced_.5BTODO.5D
    => Must be edited, most options are not documented.
11. Removed text (too drafty): "Dipoles outside of the grid for the MEG: remove this dipole or replace it with its spatial neighbors (not the neighbor in the list) ==> not easy to handle write small tuto? "
12. Review associated FEM tutorials (Francois):
    • https://neuroimage.usc.edu/brainstorm/Tutorials/Duneuro => DONE
      • Documentation of all the advanced options missing:
    • https://neuroimage.usc.edu/brainstorm/Tutorials/FemTensors => DONE
    • https://neuroimage.usc.edu/brainstorm/Tutorials/FemMesh => DONE
      • Iso2mesh: Removed link "From our test", which was pointing at https://neuroimage.usc.edu/brainstorm/reviewMEEGForward?highlight=(duneuro), which does not exist.
    • Juan, Takfarinas: Please have a look to make sure that it still looks like what you wanted.

[ftadel]

@juangpc @tmedani

I also posted a draft for the processing script:
https://neuroimage.usc.edu/brainstorm/Tutorials/FemMedianNerve#Scripting
https://github.com/brainstorm-tools/brainstorm3/blob/master/toolbox/script/tutorial_fem.m

Please use try to use it to reproduce your computations:
it will allow to automate the rebuild of everything when we change a parameter in the script or when we update Brainstorm or Duneuro.

For addressing the issues: please post your comments as you work on them, and I will edit the first comment of this thread to keep track of what is left in the todo list.

Thanks!

[ftadel]

SimNIBS: Remeshing the SimNIBS model with iso2mesh doesn't work. Here is the error message:

Run through the simnibs part myself and everything seems to be working.

The SimNIBS process works perfectly.
What I couldn't get to work is the extraction of the surfaces + remeshing with iso2mesh:
https://neuroimage.usc.edu/brainstorm/Tutorials/FemMedianNerve#Remesh_with_Iso2mesh

If this works for you, we need to understand the differences between our two setups.

• Is you Brainstorm up-to-date? (if this is the case, we should be using the same version of iso2mesh)
• What version of SimNIBS are you running?
• Could you try with the default Brainstorm fiuducials, as illustrated in the tutorial, to avoid any possible bias due to different coordinate systems?

With this information my results locate the fiducials at these locations:
Should we include these captures in the tutorial?
Should we include or propose these coordinates in the tutorial?

In the tutorial I used as a starting point, it was requested to "Refine the registration with head points". This makes accurate localization of the fiducials useless. => "Because we will use the digitized head shape of the subject to refine the MRI/MEG registration, we don't need the position of the fiducials to be very accurate."

If you prefer using accurate fiducial positions, we need to disable the automatic registration, and add several pages of instructions to the tutorial to position correctly these fiducials (including the screen capture and positions you shared, description of how you defined these points, etc.)

Last option: we include lots of information about the registration (accurate fiducials) just for teaching purposes, AND discard it all by running the automatic registration. However, the main topic of this tutorial is not coregistration, therefore I'd rather avoid adding sub-branches of discussions that are not related with the topic of interest - this would end-up being only "noise" for people interested primarily in FEM modeling. Coregistration is discussed in earlier tutorials, that the user already read...

[tmedani]

Hi Francois,
I started reviewing all these points one by one.
Some figures can be different indeed since I have tried many fif files and I did not update the images yet.

Can I start editing the tuto or are you still working on it?

[juangpc]

Regarding iso2mesh remeshing:
Reviewing the code, it fails here: system() returns a floating point value which is quite strange actually. Not sure if iso2mesh people know about this.

Anyhow, I've tried different combination of surfaces and the remeshing works ok pretty much with any combination of them other than using them all. So I think it is worth trying to modify slightly simnibs values. I'm testing right now. I'll keep posting.

[ftadel]

Some figures can be different indeed since I have tried many fif files and I did not update the images yet.
Can I start editing the tuto or are you still working on it?

I am mostly done with all the tutorials: Duneuro, Median nerve, Tensors.
You can go ahead and edit them.

I'll work today on the FemMesh tutorial. I hope I'll be done before you start your workday :)

[ftadel]

@fangq
We're having issues remeshing with iso2mesh the meshes generated by SimNIBS: https://neuroimage.usc.edu/brainstorm/Tutorials/FemMedianNerve#Remesh_with_Iso2mesh
Error reported in the first post, and documented along the thread.
In case you have suggestions...

FYI We have just released a tutorial based on Brain2mesh:
https://neuroimage.usc.edu/brainstorm/Tutorials/FemTensors#FEM_mesh
(please wait for a few days before reviewing the linked tutorial named "FEM mesh generation", what is online is an old draft that I'll update soon)

[juangpc]

Hi @fangq,
The error we're experiencing is puzzling. Let me give a bit of info.

We have a pair of coregistered T1 and T2 images which we mesh with SimNIBS . We want to remesh them with Iso2mesh. For that we extract the surfaces of whyte, gray, csf, skull and scalp. And we try to call Iso2mesh to remesh them following the tutorial ftadel is mentioning.

I see a problem in a system() call but it might not be related to our actual issue. I'll explain.

When calling Iso2mesh, this call to system() returns status as a floating point value (-1.0737e+09), instead of an integer which is weird, and the cmdout is empty. Since the execution of tetgen call is defectuous, appart from the error check, readtetgen also fails a few lines after and thus surf2mesh fails.

There are quite a few questions in Matlabcentral regarding system() not working properly showing a behavior very similar to what I see in surf2mesh (like here). (or see... 1 2 3 4 ...). Apparently, system() modifies the path variable of the terminal it opens which can lead to collision of different binary files being executed. However, in surf2mesh the system call to tetgen includes the executable's full path, so the environment variable shouldn't affect. I've also checked with dependency walker and there doesn't seem to be a dll wrongly linked (some people have reported that libiomp5md.dll can cause problems, and in deed I found it in Anaconda, in SimNIBS and in my own path, however I couldn't solve any issues). So maybe this is not the cause of our problems.

We call mergemesh here with the surfaces and it runs great (the bemMerge cell can be downloaded here). A few lines after that, when calling surf2mesh here it fails. However if we repeat the process but selecting different combinations of two, three or four of the 5 surfaces it works. So, like I mentioned at the beginning, I don't know if these two problems are related but the behavior is weird. Plenty of free memory on my system.

Hope it helps.

[tmedani]

Hi @ftadel @juangpc ,

I'm investigating this issue regarding the re-mesh with iso2mesh,
it was working a few weeks ago, I'm checking what is changed since then.

It seems to be related to version of the tetgen,
the last argument in the sur2mesh function, 'tetgen1.5' ,
but I'm not sure 100%, something else is happening.

I will give you update asap.

Thanks

[fangq]

@ftadel @juangpc and @tmedani, the status checking was added recently - I haven't encountered strange status code like this before, but happy to take a look.

can you set a break-point at the system call, evaluate the string parameter feed to system, which is basically a shell command, copy paste that command in a terminal and run it, what do you see?

if it still looks strange, please send me all input variables and full command to that surf2mesh call and I will take a look.

[juangpc]

can you set a break-point at the system call, evaluate the string parameter feed to system, which is basically a shell command, copy paste that command in a terminal and run it, what do you see?

If you do that, the command executes with normal output (not empty) and exit status (not floating point).

 "C:\Users\<<username>>\.brainstorm\plugins\iso2mesh\iso2mesh\bin\tetgen1.5.exe" -A -q1.414a1e-08  "C:\Users\<<username>>\.brainstorm\tmp\post_vmesh.poly"

This is the output:

Opening C:\Users\j\.brainstorm\tmp\post_vmesh.poly.
Delaunizing vertices...
Delaunay seconds:  3.526
Creating surface mesh ...
Surface mesh seconds:  0.669
Constrained Delaunay...

However the readtetgen call that attempts to read the result of tetgen still fails.

[fangq]

the above log looks like incomplete - likely tetgen crashed in the middle, can you type echo $? and see if it has a non-zero return value?

my suspicion is that this is related to 32bit/64bit executable issue, which was fixed last year: fangq/iso2mesh#11

if you run iso2mesh on a 64bit windows machine, the tetgen executable (mcpath('tetgen1.5')) s2m runs should be tetgen1.5_x86-64.exe, instead of tetgen1.5.exe.

the part that handles the exe search fallback is here
https://github.com/fangq/iso2mesh/blob/master/mcpath.m#L46-L52

PS: this issue was also reported previously here: fangq/iso2mesh#43

[fangq]

on a side note, I'd like to understand the rationales for a user to choose the combination of SimNIBS+iso2mesh remeshing - what is the remeshing step trying to achieve?

on a related question - from the tutorial SimNIBS meshing can take 4-5 hours and the mesh is very dense, why not use brain2mesh? any issue for brain2mesh's outputs? I am asking because I want to see what aspect that brain2mesh can be further improved (at least matching the output quality of SimNIBS), while giving full control on mesh density.

[juangpc]

the above log looks like incomplete - likely tetgen crashed in the middle, can you type echo $? and see if it has a non-zero return value?

I get no return value. Not sure if it is my terminal though.

if you run iso2mesh on a 64bit windows machine, the tetgen executable (mcpath('tetgen1.5')) s2m runs should be tetgen1.5_x86-64.exe, instead of tetgen1.5.exe.

Could it be that we're downloading an outdated version? tetgen1.5_x86-64.exe is not present in the release v1.9.2. And the mcpath.m Brainstorm users are running doesn't search for x64 executables.

Brainstorm downloads from Iso2mesh's releases the following one: https://github.com/fangq/iso2mesh/releases/tag/v1.9.2
Which was generated on Feb 13 2020. The part that handles 64bit exe was committed a week later fangq/iso2mesh@ddc9a16

[juangpc]

If I manually clone Iso2mesh's updated repo and install it for Brainstorm to use it. I'm able to get a different output. system() seems to be now working fine. So it probably was that x32, x64 confusion. Thank you for that.

The error I get now is:

I can remesh no problem all surfaces except the scalp. So it seems there might be some problem there.

[juangpc]

any issue for brain2mesh's outputs?

Isn't it necessary to have matlab's Image processing toolbox installed?

[tmedani]

Hi @fangq and thanks for your feedback and advice.

PS: this issue was also reported previously here: fangq/iso2mesh#43

indeed we have discussed this issue some time ago but it appears in some cases that I'm also investigating.

on a side note, I'd like to understand the rationales for a user to choose the combination of SimNIBS+iso2mesh remeshing - what is the remeshing step trying to achieve?

This is related to two main reasons,

• The first is that for the FEM solver of the EEG/MEG leadfield with the Duneuro fails in some cases, the recommendation of the DUNEuro team is to avoid the hole in the mesh. They are working to improve this problem, but it may take some time. The SimNibs outputs has in some case holes (air cavities) on the mesh. The reconstruction the FEM mesh from the outer surfaces of the tissue was the first solution. - An other solution is to mesh the air cavities and count it as a new material. SimNibs team implemented this solution but we did not included it yet in Brainstorm. This require a lot of change that we will add soon.

• The second reason to reconstruct the mesh from surfaces is for the sEEG/ECOG application. In these modelities the area of the interest is within the brain tissues (sources and electrodes), therefore it's not important to keep the outer layers in the simulation, we build the model only by using the inner tissues (wm, gm and csf). This allow us to increase the mesh resolution within the volume of the brain and around the intracranial electrodes (finer resolution around the electrodes).

[ftadel]

Brainstorm downloads from Iso2mesh's releases the following one: https://github.com/fangq/iso2mesh/releases/tag/v1.9.2
Which was generated on Feb 13 2020. The part that handles 64bit exe was committed a week later fangq/iso2mesh@ddc9a16

@fangq Could you recompile a newer version, as you did last year?
We'd only need the "-allinone.zip" package.
Thanks!

[fangq]

Re @ftadel: @fangq Could you recompile a newer version, as you did last year?
We'd only need the "-allinone.zip" package.

done, see
https://github.com/fangq/iso2mesh/releases/tag/v1.9.6

Re @juangpc: I can remesh no problem all surfaces except the scalp

from the error message, looks like the input surface contains self-intersections, there must be something in the earlier steps that has caused this. please check

Re @juangpc: Isn't it necessary to have matlab's Image processing toolbox installed?

indeed - the only thing I need are imdilate/imfill, see fangq/brain2mesh#10, still haven't located a portable replacement

The first is that for the FEM solver of the EEG/MEG leadfield with the Duneuro fails in some cases,

did you mean brain2mesh created surfaces have holes? that would be unexpected, I would love to get a reproducer and test

it's not important to keep the outer layers in the simulation, we build the model only by using the inner tissues (wm, gm and csf)

in brain2mesh, only wm/gm are required, all other segments (csf, bone, scalp) are optional, see

fangq/brain2mesh@3e88bd8
https://github.com/fangq/brain2mesh/blob/master/brain2mesh.m#L32

again, I am aware of the output quality differences between SimNIBS and SPM seg derived brain2mesh meshes, just want to if there is anything easy to fix

[ftadel]

@fangq: done, see: https://github.com/fangq/iso2mesh/releases/tag/v1.9.6

Thanks!
I updated the link in Brainstorm (this will update iso2mesh automatically after Brainstorm is updated):
96f1494

@juangpc @tmedani It solves the remeshing issue with the tutorial! (at least on my computer)
https://neuroimage.usc.edu/brainstorm/Tutorials/FemMedianNerve

@fangq: did you mean brain2mesh created surfaces have holes? that would be unexpected, I would love to get a reproducer and test

No, brain2mesh fills entirely the head volume.
The "problem" is related only with SimNIBS, which doesn't mesh some of the air cavities, causing DUNEuro to crash in some cases (I haven't seen this myself, in the current tutorials, it works).

SimNIBS (air cavities in the anterior part of the skull):

Surfaces extracted and remeshed with Iso2mesh (no more holes):

@fangq: again, I am aware of the output quality differences between SimNIBS and SPM seg derived brain2mesh meshes, just want to if there is anything easy to fix

I can't judge the quality of the tetrahedral mesh.

From a more general point of view in Brainstorm source estimation: SimNIBS has the advantage of running a full CAT12 pipeline, which generates good quality triangular meshes of the grey matter.
However, we can run an updated version of CAT12 as a Brainstorm plugin, with more tailored options for generating anatomical atlases. Therefore SimNIBS has no real advantage, as we have an option for a 100% Matlab/Brainstorm-based pipeline, which runs faster and with better-quality similar pial triangular meshes.

@juangpc @tmedani Maybe you can elaborate on the choice of recommending SimNIBS over Brain2mesh+CAT12?

[tmedani]

Hi all,
Sorry for being slow in this discussion, I was on vacation, and then had a lot of emails and accumulated deadlines.

Thanks, @fangq, and @ftadel for solving the meshing issues

@juangpc @tmedani It solves the remeshing issue with the tutorial! (at least on my computer)
https://neuroimage.usc.edu/brainstorm/Tutorials/FemMedianNerve

I will check it on my computer as well asap.

The "problem" is related only with SimNIBS, which doesn't mesh some of the air cavities, causing DUNEuro to crash in some cases (I haven't seen this myself, in the current tutorials, it works).

@ftadel do you mean that Duneuro works even on the model with holes?
for both modalities EEG/MEG?
this is interesting to investigate because it was not working on my case for MEG.
I have already reported this issue to the Duneuro dev team.

@juangpc @tmedani Maybe you can elaborate on the choice of recommending SimNIBS over Brain2mesh+CAT12?

We have already integrated brain2mesh within brainstorm and tested in some case, however since brain2mesh need a segmented volume as an input, we need first to segment the MRI, and from our testings, we have observed bad segmentation outputs in some cases which lead to a bad mesh.

However since simnibs has the full process with integrated methods that correct the segmentation, so in most cases, it has a better mesh, this was the main reason why we have recommended SimNibs for now.

So, if we have a robust method for segmenting the whole, then of course brain2mesh is the straightforward method to recommend since it's Matlab based, and it can be much faster.

[ftadel]

do you mean that Duneuro works even on the model with holes? for both modalities EEG/MEG?
this is interesting to investigate because it was not working on my case for MEG.

Actually, I realize that I haven't let the MEG FEM run on my old laptop until the end, so I don't know...
EEG works, MEG maybe not.
If this is only the MEG that crashes, could you please add this information (after testing it) in the tutorial, in the section about the iso2mesh remeshing?
Thanks

[ftadel]

@tmedani
Is there anything left to do from this issue?

[tmedani]

Hi @ftadel,

Thanks for checking on that.
yes, there are some remaining points, and then I need to update the tutorial.

I jumped to something else, but I will be back to that and complete it.