failed attempt: fields removes sources

plot_source_helmet.py

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
fig = mlab.figure(bgcolor=(0., 0., 0.), fgcolor=(1., 1., 1.), size=(800, 800))

brain = stc.plot(surface='white', hemi='both', subjects_dir=subjects_dir,
                 clim=dict(kind='value', lims=[2, 4, 15]), background='w',
                 initial_time=.100, time_unit='s', colormap='rocket',
                 foreground='k', figure=fig)

# convert sources estimates into meters
for child in fig.children[:2]:  # next two are shared, skip
    loc = np.array(child.data.points.data)
    child.data.points.data = loc / 1e3
    child.update()

# plot sensors
mne.viz.plot_alignment(evoked.info, trans_fname, fig=fig,
                       coord_frame='mri',  surfaces=['head-dense'],
                       subject='sample', subjects_dir=subjects_dir,
                       meg='sensors')
fig.children[-2].children[0].children[0].actor.property.opacity = 0.33

# 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)

mlab.view(19, 82, 0.6)  # m