Custom MEG sensor configurations for simulation

💬 Support & Discussions
1

Hello,

Is there a way to create custom MEG montages? The functionality of MNE is amazing and I would like to use some of the preprocessing tools (SSS, etc.) for use with simulation of different sensor geometries.

Some more details: I would like to experiment with different sensor geometries and noise profiles using my own simulation test bed. The hope is to create custom Raw objects and apply MNE preprocessing tools to this data.

I found this post discussing custom OPM sensor montages, but this problem seems to be solved by MNE adding support for this persons OPM system.

I also found mne.use_coil_def which seems to be able to use a custom coil file. Is there any guidance on how to format this file for use with my simulation?

Any advice or help would be greatly appreciated.

Thanks in advance,
Jordan

2

Just want to make sure you saw this example where they use a custom coil_def.dat:

https://mne.tools/stable/auto_examples/time_frequency/source_power_spectrum_opm.html#sphx-glr-auto-examples-time-frequency-source-power-spectrum-opm-py

You can for example look at the custom coil def they supply if you download the dataset. Some details on the format are here:

https://mne.tools/stable/overview/implementation.html?highlight=coil_def#implemented-coil-geometries

Here is a code snippet I’ve used for simulating OPM sensors:

Example code 
data = pymatreader.read_mat('python.mat')
data['R'] = data['R'].T
n_chs = len(data['R'])
info = mne.create_info(n_chs, 1000., 'mag')
info['dev_head_t'] = mne.transforms.Transform('meg', 'head', np.eye(4))
# FIFF.FIFFV_COIL_QUSPIN_ZFOPM_MAG   = 8001
# FIFF.FIFFV_COIL_QUSPIN_ZFOPM_MAG2  = 8002
# FIFF.FIFFV_COIL_FIELDLINE_OPM_MAG_GEN1    = 8101
# FIFF.FIFFV_COIL_KERNEL_OPM_MAG_GEN1    = 8201
for ii, ch in enumerate(info['chs']):
    ch['coil_type'] = 9999
    ch['loc'][:] = np.concatenate(
        [data[key][ii] for key in ('R', 'ex', 'ey', 'ez')])
trans = mne.transforms.Transform('head', 'mri', np.eye(4))

...

d = data['d']
dd = np.sqrt(3 / 5) * d
w = np.zeros(9)
w[0] = 4 * d ** 2 * (16 / 81)
w[1:5] = 4 * d ** 2 * (25 / 324)
w[5:] = 4 * d ** 2 * (10 / 81)

fname_def = 'coil_def_cube.dat'
with open(fname_def, 'w') as fid:
    fid.write('# Custom\n')
    fid.write(f'1 9999 2 9 3e-03 0.000e+00 "9-Point Cubature"\n')
    for ww, (xw, yw) in zip(w, signs):
        x = xw * dd
        y = yw * dd
        z = 0.
        fid.write(
            f'  {2 * ww:0.4f} {x:0.3f} {y:0.3f} {z:0.3f} 0.000 0.000 1.000\n')

...

        with mne.forward.use_coil_def(fname_def):
            mne.viz.plot_alignment(
                info, trans,
                surfaces=surfaces,
                show_axes=True, bem=bem, meg='sensors', src=src,
                fig=fig)

...

    with mne.forward.use_coil_def(fname_def):
        fwd = mne.make_forward_solution(info, trans, src, bem)

This does some wrong stuff like equating the head, MEG, and MRI coordinate frames, but it should at least get you started, and plot_alignment should allow you to verify things look reasonable, etc.

3

Hi Eric,

Thank you for responding to my message. I had not seen that example that you sent over. It was helpful. I have been trying to hack together to get things working on my end, but haven’t had much success. To summarize what I’ve done:

  1. I was able to make a custom info structure with my simulated data.
  2. Made a custom coil_def file.

Most of the MNE functionality works except for things that require information about sensor geometry. If I try to plot_alignment or use SSS I get a similar error. I’m thinking that it has to do with the coil_def file being specified incorrectly. I’ve dropped a coil def below. Does this look generally correct? I have 16 OPMs arranged on a flat array.

# Custom
1 9999 1 9 3e-03 0.000e+00 "flat array"
 0.0625 -0.088 0.086 0.000 0.000 0.000 1.000
 0.0625 -0.029 0.086 0.000 0.000 0.000 1.000
 0.0625 0.029 0.086 0.000 0.000 0.000 1.000
 0.0625 0.085 0.086 0.000 0.000 0.000 1.000
 0.0625 -0.085 0.029 0.000 0.000 0.000 1.000
 0.0625 -0.029 0.029 0.000 0.000 0.000 1.000
 0.0625 0.029 0.029 0.000 0.000 0.000 1.000
 0.0625 0.088 0.029 0.000 0.000 0.000 1.000
 0.0625 -0.088 -0.029 0.000 0.000 0.000 1.000
 0.0625 -0.029 -0.029 0.000 0.000 0.000 1.000
 0.0625 0.029 -0.029 0.000 0.000 0.000 1.000
 0.0625 0.088 -0.029 0.000 0.000 0.000 1.000
 0.0625 -0.088 -0.086 0.000 0.000 0.000 1.000
 0.0625 -0.029 -0.086 0.000 0.000 0.000 1.000
 0.0625 0.029 -0.086 0.000 0.000 0.000 1.000
 0.0625 0.088 -0.086 0.000 0.000 0.000 1.000

Here is the code I am using:

# create an MNE structure based on our data
info = mne.create_info(ch_names = experiment.magdata.columns.to_list(),
                           ch_types = ["mag"] * 16, 
                           sfreq = experiment.sampling_freq)
# coordinate transformation between head and array - equating for now. Not correct
info['dev_head_t'] = mne.transforms.Transform('meg', 'head', np.eye(4))
trans = mne.transforms.Transform('head', 'mri', np.eye(4))

# add sensor geom to ch field
for ii, ch in enumerate(info['chs']):
    ikey = experiment.magdata.columns[ii]
    ch['coil_type'] = 9999
    tmp_loc = np.hstack((
        np.array(experiment.sensor_dict[ikey]['location']),
        np.array(experiment.sensor_dict[ikey]['rot_mat']).flatten(),
    ))
    # print(tmp_loc)
    ch['loc'][:] = tmp_loc

# make custom coil file
fname_def = 'coil_def_flat_array.dat'
with open(fname_def, 'w') as fid:
    # write header
    fid.write('# Custom\n')
    fid.write(f'1 9999 1 9 3e-03 0.000e+00 "flat array"\n')
    # iterate through our sensor_dict to build 
    for ikey in experiment.sensor_dict.keys():
        # stich locations stuff togethter
        tmp_loc = np.hstack(
            np.array(experiment.sensor_dict[ikey]['location']),
            )
        # print(tmp_loc)
        fid.write(
            f' {1/16:0.4f} {tmp_loc[0]:0.3f} {tmp_loc[1]:0.3f} {tmp_loc[2]:0.3f} 0.000 0.000 1.000\n'
        ) 

# plot sensor montage (generally looks right)
_=mne.viz.plot_sensors(info)

# cast into mne raw 
data_rawformat = mne.io.RawArray(experiment.magdata.to_numpy().T, info=info)

When I run plot_alignment or mne.preprocessing.maxwell_filter() I get an assertion error

# use other brains
data_path = mne.datasets.opm.data_path()
subject = 'OPM_sample'
subjects_dir = op.join(data_path, 'subjects')
trans = mne.transforms.Transform('head', 'mri', np.eye(4))
# use custom coil def to check alignment
with mne.forward.use_coil_def(fname_def):
    mne.viz.plot_alignment(
        info, trans, dig = False,
        subject=subject,
        subjects_dir=subjects_dir,
        surfaces=[], meg = ['sensors'],coord_frame="meg"
    )

Error below:

---------------------------------------------------------------------------
AssertionError                            Traceback (most recent call last)
/Users/jordan.crivelli.decker/Carnival-Data/notebooks/custom_mne_arrangement.ipynb Cell 10' in <cell line: 3>()
      2 surfaces=dict(brain=0.4, outer_skull=0.6, head=None)
      3 with mne.forward.use_coil_def(fname_def):
----> 4     mne.viz.plot_alignment(
      5         info, trans, dig = False,
      6         subject=subject,
      7         subjects_dir=subjects_dir,
      8         surfaces=[], meg = ['sensors'],coord_frame="meg"
      9     )

File <decorator-gen-171>:12, in plot_alignment(info, trans, subject, subjects_dir, surfaces, coord_frame, meg, eeg, fwd, dig, ecog, src, mri_fiducials, bem, seeg, fnirs, show_axes, dbs, fig, interaction, verbose)

File /opt/anaconda3/envs/mne/lib/python3.9/site-packages/mne/viz/_3d.py:752, in plot_alignment(info, trans, subject, subjects_dir, surfaces, coord_frame, meg, eeg, fwd, dig, ecog, src, mri_fiducials, bem, seeg, fnirs, show_axes, dbs, fig, interaction, verbose)
    749 # plot sensors (NB snapshot_brain_montage relies on the last thing being
    750 # plotted being the sensors, so we need to do this after the surfaces)
    751 if picks.size > 0:
--> 752     _plot_sensors(renderer, info, to_cf_t, picks, meg, eeg, fnirs,
    753                   warn_meg, head_surf, 'm')
    755 if src is not None:
    756     atlas_ids, colors = read_freesurfer_lut()

File /opt/anaconda3/envs/mne/lib/python3.9/site-packages/mne/viz/_3d.py:1151, in _plot_sensors(renderer, info, to_cf_t, picks, meg, eeg, fnirs, warn_meg, head_surf, units, sensor_opacity, orient_glyphs, scale_by_distance, project_points, surf, check_inside, nearest)
   1149 """Render sensors in a 3D scene."""
...
     91 desc = line[desc_start:desc_end]
     92 vals = np.fromstring(line[:desc_start].strip(),
     93                      dtype=float, sep=' ')

AssertionError:

Thank you very much for your help so far!

Best,
Jordan