[med-svn] [python-mne] 09/353: ENH : adding write_label function + function that converts an stc to a label
Yaroslav Halchenko
debian at onerussian.com
Fri Nov 27 17:24:24 UTC 2015
This is an automated email from the git hooks/post-receive script.
yoh pushed a commit to tag 0.4
in repository python-mne.
commit 44593b8d06da6e53381a237174e2464d7fb2032b
Author: Alexandre Gramfort <alexandre.gramfort at inria.fr>
Date: Fri Nov 11 13:00:25 2011 -0500
ENH : adding write_label function + function that converts an stc to a label
---
mne/__init__.py | 3 +-
mne/label.py | 103 +++++++++++++++++++++++++++++++++++++++++++++++-
mne/source_estimate.py | 2 +
mne/tests/test_label.py | 41 ++++++++++++++++++-
4 files changed, 145 insertions(+), 4 deletions(-)
diff --git a/mne/__init__.py b/mne/__init__.py
index 45ee7fb..0136e9a 100644
--- a/mne/__init__.py
+++ b/mne/__init__.py
@@ -11,7 +11,8 @@ from .source_estimate import read_stc, write_stc, SourceEstimate, morph_data, \
from .surface import read_bem_surfaces, read_surface
from .source_space import read_source_spaces
from .epochs import Epochs
-from .label import label_time_courses, read_label, label_sign_flip
+from .label import label_time_courses, read_label, label_sign_flip, \
+ write_label, stc_to_label
from .misc import parse_config, read_reject_parameters
from .transforms import transform_coordinates
from .proj import read_proj
diff --git a/mne/label.py b/mne/label.py
index 13320fe..2675d65 100644
--- a/mne/label.py
+++ b/mne/label.py
@@ -1,7 +1,9 @@
+import os
import numpy as np
from scipy import linalg
-from .source_estimate import read_stc
+from .source_estimate import read_stc, mesh_edges
+from .surface import read_surface
def read_label(filename):
@@ -46,6 +48,36 @@ def read_label(filename):
return label
+def write_label(filename, label):
+ """Write a FreeSurfer label
+
+ Parameters
+ ----------
+ filename : string
+ Path to label file to produce.
+ label : dict
+ The label structure.
+ """
+ if not filename.endswith('lh.label') or not filename.endswith('rh.label'):
+ filename += '-' + label['hemi'] + '.label'
+
+ print 'Saving label to : %s' % filename
+
+ fid = open(filename, 'wb')
+ n_vertices = len(label['vertices'])
+ data = np.zeros((n_vertices, 5), dtype=np.float)
+ data[:, 0] = label['vertices']
+ data[:, 1:4] = 1e3 * label['pos']
+ data[:, 4] = label['values']
+ fid.write("#%s\n" % label['comment'])
+ fid.write("%d\n" % n_vertices)
+ for d in data:
+ fid.write("%d %f %f %f %f\n" % tuple(d))
+
+ print '[done]'
+ return label
+
+
def label_time_courses(labelfile, stcfile):
"""Extract the time courses corresponding to a label file from an stc file
@@ -119,3 +151,72 @@ def label_sign_flip(label, src):
# Comparing to the direction of the first right singular vector
flip = np.sign(np.dot(ori, Vh[:, 0]))
return flip
+
+
+def stc_to_label(stc, src, smooth=5):
+ """Compute a label from the non-zero sources in an stc object.
+
+ Parameters
+ ----------
+ stc : SourceEstimate
+ The source estimates
+ src : list of dict or string
+ The source space over which are defined the source estimates.
+ If it's a string it should the subject name (e.g. fsaverage).
+
+ Returns
+ -------
+ labels : list of dict
+ The generated labels. One per hemisphere containing sources.
+ """
+ from scipy import sparse
+
+ if not stc.is_surface():
+ raise ValueError('SourceEstimate should be surface source '
+ 'estimates')
+
+ if isinstance(src, str):
+ if 'SUBJECTS_DIR' in os.environ:
+ subjects_dir = os.environ['SUBJECTS_DIR']
+ else:
+ raise ValueError('SUBJECTS_DIR environment variable not set')
+ surf_path_from = os.path.join(subjects_dir, src, 'surf')
+ rr_lh, tris_lh = read_surface(os.path.join(surf_path_from,
+ 'lh.white'))
+ rr_rh, tris_rh = read_surface(os.path.join(surf_path_from,
+ 'rh.white'))
+ rr = [rr_lh, rr_rh]
+ tris = [tris_lh, tris_rh]
+ else:
+ if len(src) != 2:
+ raise ValueError('source space should contain the 2 hemispheres')
+ tris = [src[0]['tris'], src[1]['tris']]
+ rr = [1e3 * src[0]['rr'], 1e3 * src[1]['rr']]
+
+ labels = []
+ cnt = 0
+ for hemi, this_vertno, this_tris, this_rr in \
+ zip(['lh', 'rh'], stc.vertno, tris, rr):
+ if len(this_vertno) == 0:
+ continue
+ e = mesh_edges(this_tris)
+ e.data[e.data == 2] = 1
+ n_vertices = e.shape[0]
+ this_data = stc.data[cnt:cnt + len(this_vertno)]
+ cnt += len(this_vertno)
+ e = e + sparse.eye(n_vertices, n_vertices)
+ idx_use = this_vertno[np.any(this_data, axis=1)]
+ for k in range(smooth):
+ e_use = e[:, idx_use]
+ data1 = e_use * np.ones(len(idx_use))
+ idx_use = np.where(data1)[0]
+
+ label = dict()
+ label['comment'] = 'Label from stc'
+ label['vertices'] = idx_use
+ label['pos'] = this_rr[idx_use]
+ label['values'] = np.ones(len(idx_use))
+ label['hemi'] = hemi
+ labels.append(label)
+
+ return labels
diff --git a/mne/source_estimate.py b/mne/source_estimate.py
index 5ce0f1f..1f6a315 100644
--- a/mne/source_estimate.py
+++ b/mne/source_estimate.py
@@ -128,6 +128,8 @@ class SourceEstimate(object):
np.arange(self.data.shape[1]))
self.vertno = [vl['vertices']]
else: # surface source spaces
+ if fname.endswith('-lh.stc') or fname.endswith('-rh.stc'):
+ fname = fname[:-7]
lh = read_stc(fname + '-lh.stc')
rh = read_stc(fname + '-rh.stc')
self.data = np.r_[lh['data'], rh['data']]
diff --git a/mne/tests/test_label.py b/mne/tests/test_label.py
index 846c999..bd4af93 100644
--- a/mne/tests/test_label.py
+++ b/mne/tests/test_label.py
@@ -1,21 +1,58 @@
+import os
import os.path as op
+from numpy.testing import assert_array_almost_equal
from nose.tools import assert_true
from ..datasets import sample
-from .. import label_time_courses
+from .. import label_time_courses, read_label, write_label, stc_to_label, \
+ SourceEstimate, read_source_spaces
+
examples_folder = op.join(op.dirname(__file__), '..', '..', 'examples')
data_path = sample.data_path(examples_folder)
stc_fname = op.join(data_path, 'MEG', 'sample', 'sample_audvis-meg-lh.stc')
label = 'Aud-lh'
label_fname = op.join(data_path, 'MEG', 'sample', 'labels', '%s.label' % label)
+src_fname = op.join(data_path, 'MEG', 'sample',
+ 'sample_audvis-eeg-oct-6p-fwd.fif')
def test_label_io_and_time_course_estimates():
- """Test IO for STC files
+ """Test IO for label + stc files
"""
values, times, vertices = label_time_courses(label_fname, stc_fname)
assert_true(len(times) == values.shape[1])
assert_true(len(vertices) == values.shape[0])
+
+
+def test_label_io():
+ """Test IO of label files
+ """
+ label = read_label(label_fname)
+ write_label('foo', label)
+ label2 = read_label('foo-lh.label')
+
+ for key in label.keys():
+ if key in ['comment', 'hemi']:
+ assert_true(label[key] == label2[key])
+ else:
+ assert_array_almost_equal(label[key], label2[key], 5)
+
+
+def test_stc_to_label():
+ """Test stc_to_label
+ """
+ src = read_source_spaces(src_fname)
+ stc = SourceEstimate(stc_fname)
+ os.environ['SUBJECTS_DIR'] = op.join(data_path, 'subjects')
+ labels1 = stc_to_label(stc, src='sample', smooth=3)
+ labels2 = stc_to_label(stc, src=src, smooth=3)
+ assert_true(len(labels1) == len(labels2))
+ for l1, l2 in zip(labels1, labels2):
+ for key in l1.keys():
+ if key in ['comment', 'hemi']:
+ assert_true(l1[key] == l1[key])
+ else:
+ assert_array_almost_equal(l1[key], l2[key], 4)
--
Alioth's /usr/local/bin/git-commit-notice on /srv/git.debian.org/git/debian-med/python-mne.git
More information about the debian-med-commit
mailing list