rois.to_svg() fix

cortex/align.py

@@ -200,3 +200,76 @@ def autotweak(subject, xfmname):
         print('Saved transform as (%s, %s)'%(subject, xfmname+'_auto'))
     finally:
         shutil.rmtree(cache)
+
+
+def automatic_bbregister(subject, xfmname, reference, noclean=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
+    transform called `xfmname` will be automatically stored in the database.
+
+    It's good practice to open up this transform afterward in the manual aligner and check how it worked.
+    Do that using the following (with the same `subject` and `xfmname` used here, no need for `reference`):
+    > align.manual(subject, xfmname)
+
+    Parameters
+    ----------
+    subject : str
+        Subject identifier.
+    xfmname : str
+        String identifying the transform to be created.
+    reference : str
+        Path to a nibabel-readable image that will be used as the reference for this transform.
+        Usually, this is a single (3D) functional data volume.
+    noclean : bool, optional
+        If True intermediate files will not be removed from /tmp (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 shlex
+    import shutil
+    import tempfile
+    import subprocess as sp
+
+    from .database import db
+    from .xfm import Transform
+    from .options import config
+
+    fsl_prefix = config.get("basic", "fsl_prefix")
+    schfile = os.path.join(os.path.split(os.path.abspath(__file__))[0], "bbr.sch")
+
+    retval = None
+    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('Running bbregister')
+        cmd = 'bbregister --s {subject} --mov {absref} --bold --init-fsl --reg {cache}/register.dat --fslmat {cache}/out.mat'
+        cmd = cmd.format(cache=cache, absref=absreference, subject=subject)
+
+        if sp.call(cmd, shell=True) != 0:
+            raise IOError('Error calling bbregister')
+
+        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')
+
+    finally:
+        if not noclean:
+            shutil.rmtree(cache)
+        else:
+            retval = cache
+
+    return retval

cortex/polyutils.py

@@ -43,6 +43,55 @@ def __init__(self, pts, polys):
         self._rlfac_solvers = dict()
         self._nLC_solvers = dict()
 
+    @property
+    def poly_graph(self):
+        """NetworkX undirected graph representing polygons of this Surface.
+        """
+        import networkx as nx
+        from collections import defaultdict
+        edges = defaultdict(list)
+        for ii,(a,b,c) in enumerate(self.polys):
+            edges[frozenset([a,b])].append(ii)
+            edges[frozenset([a,c])].append(ii)
+            edges[frozenset([b,c])].append(ii)
+
+        #nedges = len(edges)
+        #ii,jj = np.vstack(edges.values()).T
+        #polymat = sparse.coo_matrix((np.ones((nedges,)), (ii, jj)), shape=[len(self.polys)]*2)
+        polygraph = nx.Graph()
+        polygraph.add_edges_from(((p[0], p[1], dict(verts=k)) for k,p in edges.iteritems()))
+        return polygraph
+
+    def get_boundary(self, vertices, remove_danglers=False):
+        """Return interior and exterior boundary sets for `vertices`.
+        If `remove_danglers` is True vertices in the internal boundary with
+        only one neighbor in the internal boundary will be moved to the external
+        boundary.
+        """
+        if not len(vertices):
+            return [], []
+
+        import networkx as nx
+
+        # Use networkx to get external boundary
+        external_boundary = set(nx.node_boundary(self.graph, vertices))
+
+        # Find adjacent vertices to get inner boundary
+        internal_boundary = set.union(*[set(self.graph[v].keys())
+                                        for v in external_boundary]).intersection(set(vertices))
+
+        if remove_danglers:
+            ingraph = self.graph.subgraph(internal_boundary)
+            danglers = [n for n,d in ingraph.degree().items() if d==1]
+            while danglers:
+                internal_boundary -= set(danglers)
+                external_boundary |= set(danglers)
+
+                ingraph = self.graph.subgraph(internal_boundary)
+                danglers = [n for n,d in ingraph.degree().items() if d<2]
+
+        return list(internal_boundary), list(external_boundary)
+
     @property
     @_memo
     def ppts(self):

cortex/rois.py

@@ -7,11 +7,12 @@
 import networkx as nx
 
 from db import surfs
+
+from dataset import db
 from svgroi import get_roipack, _make_layer, _find_layer, parser
 from lxml import etree
-from dataset import VertexData
-from polyutils import Surface, boundary_edges
-from utils import get_curvature, add_roi
+from polyutils import Surface
+
 import quickflat
 
 class ROIpack(object):
@@ -29,9 +30,12 @@ def load_roifile(self):
             print("ROI file %s doesn't exist.." % self.roifile)
             return
 
-        # Create basic VertexData to avoid expensive initialization..
-        empty = VertexData(None, self.subject)
-
+        #Used to create basic VertexData to avoid expensive initialization..
+        #empty = VertexData(None, self.subject)
+
+        #hack just gets curvature VertexData to prevent error from old way
+        empty = db.get_surfinfo(self.subject)
+
         # Load ROIs from file
         if self.roifile.endswith("npz"):
             roidata = np.load(self.roifile)
@@ -79,7 +83,7 @@ def to_svg(self, open_inkscape=False, filename=None):
         # Add default layers
         # Add curvature
         from matplotlib import cm
-        curv = VertexData(np.hstack(get_curvature(self.subject)), self.subject)
+        curv = db.get_surfinfo(self.subject)
         fp = cStringIO.StringIO()
         curvim = quickflat.make_png(fp, curv, height=1024, with_rois=False, with_labels=False,
                                     with_colorbar=False, cmap=cm.gray,recache=True)