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 <[email protected]> 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 _______________________________________________ debian-med-commit mailing list [email protected] http://lists.alioth.debian.org/cgi-bin/mailman/listinfo/debian-med-commit
