Created
November 16, 2017 03:33
-
-
Save kingjr/b437f363d9e6e9c6f12d1ab3c5dd1870 to your computer and use it in GitHub Desktop.
failed attempt: fields removes sources
This file contains bidirectional Unicode text that may be interpreted or compiled differently than what appears below. To review, open the file in an editor that reveals hidden Unicode characters.
Learn more about bidirectional Unicode characters
from mayavi import mlab | |
import numpy as np | |
import mne | |
from mne.datasets import sample | |
from mne.minimum_norm import make_inverse_operator, apply_inverse | |
mne.set_log_level(False) | |
data_path = sample.data_path() | |
subjects_dir = data_path + '/subjects' | |
raw_fname = data_path + '/MEG/sample/sample_audvis_filt-0-40_raw.fif' | |
raw = mne.io.read_raw_fif(raw_fname) # already has an average reference | |
events = mne.find_events(raw, stim_channel='STI 014') | |
event_id = dict(visual_left=3, visual_right=4) | |
raw.info['bads'] = ['MEG 2443', 'EEG 053'] | |
picks = mne.pick_types(raw.info, meg=True, eeg=False, eog=True, exclude='bads') | |
epochs = mne.Epochs(raw, events, event_id, tmin=-.2, tmax=.5, proj=True, | |
picks=picks, baseline=(None, 0), | |
reject=dict(grad=4000e-13, mag=4e-12, eog=150e-6)) | |
evoked = epochs.average() | |
trans_fname = data_path + ('/MEG/sample/sample_audvis_raw-trans.fif') | |
maps = mne.make_field_map(evoked, trans_fname, subject='sample', | |
subjects_dir=subjects_dir, n_jobs=1) | |
# source | |
noise_cov = mne.compute_covariance(epochs, tmax=0., method='empirical') | |
del epochs | |
fname_fwd = data_path + '/MEG/sample/sample_audvis-meg-oct-6-fwd.fif' | |
fwd = mne.read_forward_solution(fname_fwd) | |
fwd = mne.convert_forward_solution(fwd, surf_ori=True) | |
# Restrict forward solution as necessary for MEG | |
fwd = mne.pick_types_forward(fwd, meg=True, eeg=False) | |
# make an MEG inverse operator | |
inverse_operator = make_inverse_operator(evoked.info, fwd, noise_cov, | |
loose=0.2, depth=0.8) | |
stc = apply_inverse(evoked, inverse_operator, lambda2=1. / 3. ** 2, | |
method="dSPM", pick_ori=None) | |
# Plot source estimates | |
brain = stc.plot(surface='pial', hemi='both', subjects_dir=subjects_dir, | |
clim=dict(kind='value', lims=[2, 4, 15]), background='w', | |
initial_time=.100, time_unit='s', colormap='Oranges') | |
fig = mlab.gcf() | |
children = fig.children | |
# plot sensors | |
fig = mne.viz.plot_alignment(evoked.info, trans_fname, fig=fig, | |
coord_frame='mri', surfaces=[], | |
subject='sample', subjects_dir=subjects_dir, | |
meg='sensors') | |
# convert sources estimates into meters (if done before sensors, doesn't work?) | |
for child in children: | |
loc = np.array(child.data.points.data) | |
child.data.points.data = loc / 1e3 | |
mlab.view(19, 82, .0006) | |
# attempt to plot fields | |
surf = maps[0]['surf'] | |
map_data = maps[0]['data'] | |
time_idx = np.where(evoked.times > .170)[0][0] | |
pick = mne.pick_types(evoked.info, meg=True, eeg=False, ref_meg=False) | |
data = np.dot(map_data, evoked.data[pick, time_idx]) | |
mesh = mne.viz._3d._create_mesh_surf(surf, fig, scalars=data) | |
cont = mlab.pipeline.contour_surface(mesh) |
Sign up for free
to join this conversation on GitHub.
Already have an account?
Sign in to comment