fishyness in fwd['sol'] vs. fwd['info'] ordering of channels

test_fwd_ordering_gist.py

from os import path as op
import mne
from mne.datasets import testing
import numpy as np
import matplotlib.pyplot as plt

sample_path = op.join(testing.data_path(download=False), 'MEG', 'sample')
raw_fname = op.join(sample_path, 'sample_audvis_trunc_raw.fif')
trans_fname = op.join(sample_path, 'sample_audvis_trunc-trans.fif')

sample_bem_path = op.join(testing.data_path(download=False), 'subjects',
                          'sample', 'bem')
bem_fname = op.join(sample_bem_path, 'sample-320-320-320-bem-sol.fif')
src_fname = op.join(sample_bem_path, 'sample-oct-2-src.fif')

raw = mne.io.Raw(raw_fname, preload=True, proj=True)

# reorder info-struct, don't worry about the data as it won't be used!
# new order: EOG, EEG, MEG, STI
new_order = [375] + range(315, 375) + range(0,306) + range(306, 315)

# Are these the only two that need reordering?
raw.info['ch_names'] = [raw.info['ch_names'][n] for n in new_order]
raw.info['chs'] = [raw.info['chs'][n] for n in new_order]

fwd = mne.make_forward_solution(raw.info, trans_fname, src_fname,
                                bem_fname, eeg=True, mindist=5.0)
fwd = mne.convert_forward_solution(fwd, surf_ori=True)

# Following slightly modified from mne/minimimum_norm/inverse.py#_prepare_forward
fwd_ch_names = [c['ch_name'] for c in fwd['info']['chs']]
ch_names = [c['ch_name'] for c in raw.info['chs']
        if ((c['ch_name'] not in raw.info['bads'])
             and (c['ch_name'] in fwd_ch_names))]

# But the row_names of the solution are ALWAYS: MEG first, then EEG
fwd_row_names = fwd['sol']['row_names']

# So the first 306 rows of the gain matrix are MEG, then 60 EEG (minus bads)
# Make a square identity matrix to illustrate the point
gain = np.eye(fwd['sol']['data'].shape[0])

fig, axs = plt.subplots(nrows = 1, ncols = 3)
axs[0].imshow(gain, origin='lower', cmap='Greens', interpolation='nearest')
axs[0].set_title('Original ordering')

# Reorder gain using the ordering found above from fwd['info']['chs']
fwd_idx = [fwd_ch_names.index(name) for name in ch_names]
gain_fwd_ch_names = gain[fwd_idx]
axs[1].imshow(gain_fwd_ch_names, origin='lower', cmap='Greens', interpolation='nearest')
axs[1].set_title('''Ordering from fwd['info']['chs']''')

# Reorder gain using the ordering found above from fwd['sol']['row_names']
fwd_idx = [fwd_row_names.index(name) for name in ch_names]
gain_fwd_row_names = gain[fwd_idx]
axs[2].imshow(gain_fwd_row_names, origin='lower', cmap='Greens', interpolation='nearest')
axs[2].set_title('''Ordering from fwd['sol']['row_names']''')


plt.show()