Permalink
Browse files

Merge branch 'master' into add_freesurfer_wm_segmentation

  • Loading branch information...
anwarnunez committed Jan 30, 2019
2 parents aa787e9 + 7073deb commit b5ab4fad5fe2b8d8491af7b7357081dd2a330534
@@ -48,3 +48,8 @@ Citation
If you use pycortex in published work, please cite the [pycortex paper](http://dx.doi.org/10.3389/fninf.2015.00023):

_Gao JS, Huth AG, Lescroart MD and Gallant JL (2015) Pycortex: an interactive surface visualizer for fMRI. Front. Neuroinform. 9:23. doi: 10.3389/fninf.2015.00023_

Getting help
-----------
Please post on [NeuroStars](https://neurostars.org/) with the tag `pycortex` to
ask questions about how to use Pycortex.
@@ -82,7 +82,99 @@ def view_callback(aligner):

return m

def automatic(subject, xfmname, reference, noclean=False, bbrtype="signed", pre_flirt_args=''):
def fs_manual(subject, xfmname, output_name="register.lta", wm_color="blue", pial_color="red", noclean=False):
"""Open Freesurfer FreeView GUI for manually aligning/adjusting a functional
volume to the cortical surface for `subject`. This creates a new transform
called `xfmname`. The name of a nibabel-readable file (e.g. NIfTI) should be
supplied as `reference`. This image will be copied into the database.
IMPORTANT: This function assumes that the resulting .lta file is saved as:
"{default folder chosen by FreeView (should be /tmp/fsalign_xxx)}/{output_name}".
Parameters
----------
subject : str
Subject identifier.
xfmname : str
The name of the transform to be modified.
output_name : str
The name of the .lta file generated after FreeView editing.
wm_color : str | "blue"
Color of the white matter surface. Default is "blue". This can
also be adjusted in the FreeView GUI.
pial_color : str | "red"
Color of the pial surface. Default is "red". This can also be adjusted
the FreeView GUI.
noclean : boolean | False
If True, intermediate files will not be removed from /tmp/fsalign_xxx
(this is useful for debugging things), and the returned value will be
the name of the temp directory. Default False.
Returns
-------
Nothing unless noclean is true.
"""

import subprocess as sp
import tempfile
import shutil
from .xfm import Transform
from .database import db

retval = None

try:
try:
cache = tempfile.mkdtemp(prefix="fsalign_")
sub_xfm = db.get_xfm(subject, xfmname)

# if masks have been cached, quit! user must remove them by hand
from glob import glob
if len(glob(db.get_paths(subject)['masks'].format(xfmname=xfmname, type='*'))):
print('Refusing to overwrite existing transform %s because there are cached masks. Delete the masks manually if you want to modify the transform.' % xfmname)
raise ValueError('Exiting...')
except IOError:
print("Transform does not exist!")

# Load transform-relevant things
reference = sub_xfm.reference.get_filename()
xfm_dir = os.path.dirname(reference)
_ = sub_xfm.to_freesurfer(os.path.join(cache, "register.dat"), subject) # Transform in freesurfer .dat format

# Command for FreeView and run
cmd = ("freeview -v $SUBJECTS_DIR/{sub}/mri/orig.mgz "
"{ref}:reg={reg} "
"-f $SUBJECTS_DIR/{sub}/surf/lh.white:edgecolor={wmc} $SUBJECTS_DIR/{sub}/surf/rh.white:edgecolor={wmc} "
"$SUBJECTS_DIR/{sub}/surf/lh.pial:edgecolor={pialc} $SUBJECTS_DIR/{sub}/surf/rh.pial:edgecolor={pialc}")
cmd = cmd.format(sub=subject, ref=reference, reg=os.path.join(cache, "register.dat"),
wmc=wm_color, pialc=pial_color)

# Run and save transform when user is done editing
if sp.call(cmd, shell=True) != 0:
raise IOError("Problem with FreeView!")
else:
# Convert transform into .dat format
reg_dat = os.path.join(cache, os.path.splitext(output_name)[0] + ".dat")
cmd = "lta_convert --inlta {inlta} --outreg {regdat}"
cmd = cmd.format(inlta=os.path.join(cache, output_name), regdat=reg_dat)
if sp.call(cmd, shell=True) != 0:
raise IOError("Error converting lta into dat!")

# Save transform to pycortex
xfm = Transform.from_freesurfer(reg_dat, reference, subject)
db.save_xfm(subject, xfmname, xfm.xfm, xfmtype='coord', reference=reference)
print("saved xfm")

finally:
if not noclean:
shutil.rmtree(cache)
else:
retval = cache

return retval


def automatic(subject, xfmname, reference, noclean=False, bbrtype="signed", pre_flirt_args='', use_fs_bbr=False):
"""Create an automatic alignment using the FLIRT boundary-based alignment (BBR) from FSL.
If `noclean`, intermediate files will not be removed from /tmp. The `reference` image and resulting
@@ -111,6 +203,11 @@ def automatic(subject, xfmname, reference, noclean=False, bbrtype="signed", pre_
The 'bbrtype' argument that is passed to FLIRT.
pre_flirt_args : str, optional
Additional arguments that are passed to the FLIRT pre-alignment step (not BBR).
use_fs_bbr : bool, optional
If True will use freesurfer bbregister instead of FSL BBR.
save_dat : bool, optional
If True, will save the register.dat file from freesurfer bbregister into
freesurfer's $SUBJECTS_DIR/subject/tmp.
Returns
-------
@@ -132,28 +229,41 @@ def automatic(subject, xfmname, reference, noclean=False, bbrtype="signed", pre_
try:
cache = tempfile.mkdtemp()
absreference = os.path.abspath(reference)
raw = db.get_anat(subject, type='raw').get_filename()
bet = db.get_anat(subject, type='brainmask').get_filename()
wmseg = db.get_anat(subject, type='whitematter').get_filename()
#Compute anatomical-to-epi transform
print('FLIRT pre-alignment')
cmd = '{fslpre}flirt -in {epi} -ref {bet} -dof 6 {pre_flirt_args} -omat {cache}/init.mat'.format(
fslpre=fsl_prefix, cache=cache, epi=absreference, bet=bet, pre_flirt_args=pre_flirt_args)
if sp.call(cmd, shell=True) != 0:
raise IOError('Error calling initial FLIRT')

print('Running BBR')
# Run epi-to-anat transform (this is more stable than anat-to-epi in FSL!)
cmd = '{fslpre}flirt -in {epi} -ref {raw} -dof 6 -cost bbr -wmseg {wmseg} -init {cache}/init.mat -omat {cache}/out.mat -schedule {schfile} -bbrtype {bbrtype}'
cmd = cmd.format(fslpre=fsl_prefix, cache=cache, raw=bet, wmseg=wmseg, epi=absreference, schfile=schfile, bbrtype=bbrtype)
if sp.call(cmd, shell=True) != 0:
raise IOError('Error calling BBR flirt')
if use_fs_bbr:
print('Running freesurfer BBR')
cmd = 'bbregister --s {sub} --mov {absref} --init-fsl --reg {cache}/register.dat --t1'
cmd = cmd.format(sub=subject, absref=absreference, cache=cache)

if sp.call(cmd, shell=True) != 0:
raise IOError('Error calling freesurfer BBR!')

xfm = Transform.from_freesurfer(os.path.join(cache, "register.dat"), absreference, subject)
else:
raw = db.get_anat(subject, type='raw').get_filename()
bet = db.get_anat(subject, type='brainmask').get_filename()
wmseg = db.get_anat(subject, type='whitematter').get_filename()
#Compute anatomical-to-epi transform
print('FLIRT pre-alignment')
cmd = '{fslpre}flirt -in {epi} -ref {bet} -dof 6 {pre_flirt_args} -omat {cache}/init.mat'.format(
fslpre=fsl_prefix, cache=cache, epi=absreference, bet=bet, pre_flirt_args=pre_flirt_args)
if sp.call(cmd, shell=True) != 0:
raise IOError('Error calling initial FLIRT')


print('Running BBR')
# Run epi-to-anat transform (this is more stable than anat-to-epi in FSL!)
cmd = '{fslpre}flirt -in {epi} -ref {raw} -dof 6 -cost bbr -wmseg {wmseg} -init {cache}/init.mat -omat {cache}/out.mat -schedule {schfile} -bbrtype {bbrtype}'
cmd = cmd.format(fslpre=fsl_prefix, cache=cache, raw=bet, wmseg=wmseg, epi=absreference, schfile=schfile, bbrtype=bbrtype)
if sp.call(cmd, shell=True) != 0:
raise IOError('Error calling BBR flirt')

x = np.loadtxt(os.path.join(cache, "out.mat"))
# Pass transform as FROM epi TO anat; transform will be inverted
# back to anat-to-epi, standard direction for pycortex internal
# storage by from_fsl
xfm = Transform.from_fsl(x,absreference,raw)

x = np.loadtxt(os.path.join(cache, "out.mat"))
# Pass transform as FROM epi TO anat; transform will be inverted
# back to anat-to-epi, standard direction for pycortex internal
# storage by from_fsl
xfm = Transform.from_fsl(x,absreference,raw)
# Save as pycortex 'coord' transform
xfm.save(subject,xfmname,'coord')
print('Success')
@@ -18,6 +18,7 @@ def brainmask(outfile, subject):
raw = db.get_anat(subject, type='raw').get_filename()
print('Brain masking anatomical...')
cmd = '{fsl_prefix}bet {raw} {bet} -B -v'.format(fsl_prefix=fsl_prefix, raw=raw, bet=outfile)
print("Calling: %s"%cmd)
assert sp.call(cmd, shell=True) == 0, "Error calling fsl-bet"

def whitematter(outfile, subject, do_voxelize=False):
@@ -28,7 +28,7 @@ def _call_blender(filename, code, blender_path=default_blender):
"""
with tempfile.NamedTemporaryFile() as tf:
print("In new named temp file: %s"%tf.name)
startcode=_base_imports
startcode = _base_imports
endcode = "\nbpy.ops.wm.save_mainfile(filepath='{fname}')".format(fname=filename)
cmd = "{blender_path} -b {fname} -P {tfname}".format(blender_path=blender_path, fname=filename, tfname=tf.name)
if not os.path.exists(filename):
@@ -165,6 +165,15 @@ def fs_cut(fname, subject, hemi, freesurfer_subject_dir=None):
def write_patch(bname, pname, mesh="hemi"):
"""Write out the mesh 'mesh' in the blender file 'bname' into patch file 'pname'
This is a necessary step for flattening the surface in freesurfer
Parameters
----------
bname : str
blender file name that contains the mesh
pname : str
name of patch file to be saved
mesh : str
name of mesh in blender file
"""
p = xdrlib.Packer()
p.pack_string(pname.encode())
@@ -71,7 +71,8 @@ def add_vcolor(hemis, mesh=None, name='color'):
else:
for i, j in enumerate(loopidx):
vcolor.data[i].color = color[j]
else: #older blender version, need to iterate faces instead
else:
# older blender version, need to iterate faces instead
print("older blender found...")
if not isinstance(color[0], (list, tuple)):
for i in range(len(mesh.faces)):
@@ -100,6 +101,18 @@ def add_shapekey(shape, name=None):
return key

def write_patch(filename, pts, edges=None):
"""Writes a patch file that is readable by freesurfer.
Parameters
----------
filename : name for patch to write. Should be of the form
<subject>.flatten.3d
pts : array-like
points in the mesh
edges : array-like
edges in the mesh
"""
if edges is None:
edges = set()

@@ -111,8 +124,16 @@ def write_patch(filename, pts, edges=None):
else:
fp.write(struct.pack('>i3f', i+1, *pt))

def save_patch(fname, mesh='hemi'):
"""Saves patch to file that can be read by pycortex"""
def _get_pts_edges(mesh):
"""Function called within blender to get non-cut vertices & edges
Operates on a mesh object within an open instance of blender.
Parameters
----------
mesh : str
name of mesh to cut
"""
if isinstance(mesh, str):
mesh = D.meshes[mesh]

@@ -147,10 +168,12 @@ def save_patch(fname, mesh='hemi'):
for i, vert in enumerate(mesh.vertices):
if vert.select:
smore.add(i)

bpy.ops.object.mode_set(mode='EDIT')
bpy.ops.mesh.select_all(action='DESELECT')
bpy.ops.object.mode_set(mode='OBJECT')
# Leave cuts (+ area around them) selected.
# Uncomment the next lines to revert to previous behavior
# (deselecting everything)
# bpy.ops.object.mode_set(mode='EDIT')
# bpy.ops.mesh.select_all(action='DESELECT')
# bpy.ops.object.mode_set(mode='OBJECT')

fverts = set()
if hasattr(mesh, "polygons"):
@@ -166,6 +189,11 @@ def save_patch(fname, mesh='hemi'):
edges = mwall_edge | (smore - seam)
verts = fverts - seam
pts = [(v, D.shape_keys['Key'].key_blocks['inflated'].data[v].co) for v in verts]
return verts, pts, edges

def save_patch(fname, mesh='hemi'):
"""Saves patch to file that can be read by freesurfer"""
verts, pts, edges = _get_pts_edges(mesh)
write_patch(fname, pts, edges)

def read_xdr(filename):
@@ -174,6 +202,7 @@ def read_xdr(filename):
pts = u.unpack_array(p.unpack_double)
polys = u.unpack_array(p.unpack_uint)
return pts, polys

def write_xdr(filename, pts, polys):
with open(filename, "wb") as fp:
p = xdrlib.Packer()
@@ -5,16 +5,25 @@ fsl_prefix = fsl5.0-

[dependency_paths]
# The following specify paths to the binary executable files
# for blender and inkscape (which are external dependencies of
# pycortex). "inkscape" and "blender" alone should work on Linux,
# since the binaries should already be on the system path.
# for external dependencies of pycortex. The default values
# (just "inkscape" and "blender") should work on Linux,
# since the binaries should already be on the system path if
# blender and inkscape are installed.
# For MacOS, unless you have done some manual configuration on
# your own, these will need to be something like:
# blender = /Applications/Blender/blender.app/Contents/MacOS/blender
# inkscape = /Applications/Inkscape.app/Contents/MacOS/Inkscape
# Windows is not currently supported.
# Windows is not currently supported.
inkscape = inkscape
blender = blender
# SLIM and meshlab are optional dependencies for pycortex, that
# provide a fast-and-dirty alternative to flattening brains with
# freesurfer (SLIM), and an alternative way to visualize
# intermediate stages of cortical mesh & flatmap generatioN (meshlab)
# Download SLIM code here: https://github.com/MichaelRabinovich/Scalable-Locally-Injective-Mappings
# Download meshal here: http://www.meshlab.net/
slim = None
meshlab = None

[mayavi_aligner]
line_width = 1
Oops, something went wrong.

0 comments on commit b5ab4fa

Please sign in to comment.