Deprecate read_montage and class Montage

doc/changes/latest.inc

@@ -21,13 +21,19 @@ Changelog
 
 - Add :func:`mne.channels.compute_dev_head_t` to compute Device-to-Head transformation from a montage by `Joan Massich`_ and `Alex Gramfort`_.
 
+- Add :func:`mne.channels.read_dig_fif` to read digitization coordinates from ``.fif`` files by `Joan Massich`_ and `Alex Gramfort`_.
+
+- Add :func:`mne.channels.read_dig_egi` to read digitization coordinates from EGI ``.xml`` files by `Joan Massich`_ and `Alex Gramfort`_.
+
 - Add :func:`mne.channels.read_dig_polhemus_isotrak` and :func:`mne.channels.read_polhemus_fastscan` to read Polhemus data by `Joan Massich`_
 
-- Add support for making epochs with duplicated events, by allowing three policies: "error" (default), "drop", or "merge" in :class:`mne.Epochs` by `Stefan Appelhoff`_
+- Add :func:`mne.channels.read_dig_captrack` to read BrainVision CapTrak (BVCT) digitization coordinate files by `Stefan Appelhoff`_ and `Joan Massich`_
 
 - Add :func:`mne.channels.make_dig_montage` to create :class:`mne.channels.DigMontage` objects out of np.arrays by `Joan Massich`_
 
-- Add support for reading in BrainVision CapTrak (BVCT) digitization coordinate files in :func:`mne.channels.read_dig_montage` by `Stefan Appelhoff`_
+- Add :func:`mne.channels.read_dig_eeglab` to read :class:`mne.channels.DigMontage` from EEGLAB ``.loc``, ``.locs``, and ``.eloc`` files by `Joan Massich`_ and `Alex Gramfort`_.
+
+- Add support for making epochs with duplicated events, by allowing three policies: "error" (default), "drop", or "merge" in :class:`mne.Epochs` by `Stefan Appelhoff`_
 
 - Allow :meth:`mne.Annotations.crop` to support negative ``tmin`` and ``tmax`` by `Joan Massich`_
 
@@ -178,6 +184,8 @@ Bug
 API
 ~~~
 
+- Deprecate ``mne.channels.Montage`` class, ``mne.channels.read_montage`` and ``mne.channels.read_dig_montage`` function by `Joan Massich`_.
+
 - Deprecate passing ``Montage``, ``str`` as montage parameter in :meth:`mne.io.Raw.set_montage` by `Joan Massich`_.
 
 - Deprecate ``set_dig`` parameter in :meth:`mne.io.Raw.set_montage` and ``update_ch_names`` in ``mne.io.RawEEGLAB.set_montage`` when using :class:`mne.channels.DigMontage` as  by `Joan Massich`_.
@@ -194,10 +202,6 @@ API
 
 - New boolean parameter ``show_scrollbars`` for :meth:`mne.io.Raw.plot`, :meth:`mne.Epochs.plot`, and :meth:`mne.preprocessing.ICA.plot_sources` (and associated functions) that allows hiding the scrollbars and buttons for a "zen mode" data browsing experience. When the plot window has focus, zen mode can be toggled by pressing :kbd:`z`, by `Daniel McCloy`_.
 
-- Deprecate passing ``np.arrays`` as ``hsp``, ``hpi`` or ``elp`` parameters to ``mne.channels.read_dig_montage`` by `Joan Massich`_.
-
-- Deprecate ``fif`` and ``egi`` parameters in ``mne.channels.read_dig_montage`` by `Joan Massich`_.
-
 - Deprecate ``mne.evoked.grand_average`` in favor of :func:`mne.grand_average` (which works on both :class:`~mne.Evoked` and :class:`~mne.time_frequency.AverageTFR`) by `Daniel McCloy`_
 
 - Deprecate ``exclude`` parameter in :func:`mne.viz.plot_ica_sources` and :meth:`mne.preprocessing.ICA.plot_sources`, instead always use the ``exclude`` attribute of the ICA object by `Daniel McCloy`_.

doc/python_reference.rst

@@ -323,6 +323,7 @@ Projections:
    make_eeg_layout
    make_grid_layout
    make_standard_montage
+   read_standard_montage
    find_ch_connectivity
    read_ch_connectivity
    equalize_channels

mne/channels/__init__.py

@@ -9,7 +9,9 @@
                       get_builtin_montages, make_dig_montage,
                       read_dig_egi, read_dig_captrack, read_dig_fif,
                       read_dig_polhemus_isotrak, read_polhemus_fastscan,
-                      compute_dev_head_t, make_standard_montage)
+                      compute_dev_head_t, make_standard_montage,
+                      read_standard_montage, read_dig_hpts
+                      )
 from .channels import (equalize_channels, rename_channels, fix_mag_coil_types,
                        read_ch_connectivity, _get_ch_type,
                        find_ch_connectivity, make_1020_channel_selections)
@@ -26,6 +28,7 @@
     'read_ch_connectivity', 'read_dig_captrack', 'read_dig_egi',
     'read_dig_fif', 'read_dig_montage', 'read_dig_polhemus_isotrak',
     'read_layout', 'read_montage', 'read_polhemus_fastscan',
+    'read_standard_montage', 'read_dig_hpts',
 
     # Helpers
     'rename_channels', 'make_1020_channel_selections',

mne/channels/_dig_montage_utils.py

@@ -122,7 +122,10 @@ def _read_dig_montage_egi(
         elif kind == 1:
             dig_ch_pos['EEG %03d' %
                        (len(dig_ch_pos.keys()) + 1)] = coordinates
-            # XXX: we should do something with this (ref and eeg get mixed)
+            # XXX: The EGI reader needs to be fixed with this code here.
+            # As a reference channel it should be called EEG000 or
+            # REF to follow the conventions. I should be:
+            # dig_ch_pos['REF'] = coordinates
 
         # Fiducials
         elif kind == 2:

mne/channels/_standard_montage_utils.py

@@ -7,6 +7,7 @@
 import numpy as np
 
 from functools import partial
+import xml.etree.ElementTree as ElementTree
 
 from .montage import make_dig_montage
 from ..transforms import _sph_to_cart
@@ -22,73 +23,53 @@
 
 def _egi_256(head_size):
     fname = op.join(MONTAGE_PATH, 'EGI_256.csd')
-    # Label, Theta, Phi, Radius, X, Y, Z, off sphere surface
-    options = dict(comments='//',
-                   dtype=(_str, 'f4', 'f4', 'f4', 'f4', 'f4', 'f4', 'f4'))
-    ch_names, _, _, _, xs, ys, zs, _ = _safe_np_loadtxt(fname, **options)
-    pos = np.stack([xs, ys, zs], axis=-1)
-
-    # Fix pos to match Montage code
-    pos *= head_size / np.median(np.linalg.norm(pos, axis=1))
+    montage = _read_csd(fname, head_size)
+    ch_pos = montage._get_ch_pos()
 
     # For this cap, the Nasion is the frontmost electrode,
     # LPA/RPA we approximate by putting 75% of the way (toward the front)
     # between the two electrodes that are halfway down the ear holes
-    nasion = pos[ch_names.index('E31')]
-    lpa = (0.75 * pos[ch_names.index('E67')] +
-           0.25 * pos[ch_names.index('E94')])
-    rpa = (0.75 * pos[ch_names.index('E219')] +
-           0.25 * pos[ch_names.index('E190')])
+    nasion = ch_pos['E31']
+    lpa = 0.75 * ch_pos['E67'] + 0.25 * ch_pos['E94']
+    rpa = 0.75 * ch_pos['E219'] + 0.25 * ch_pos['E190']
 
-    return make_dig_montage(
-        ch_pos=OrderedDict(zip(ch_names, pos)),
+    fids_montage = make_dig_montage(
         coord_frame='unknown', nasion=nasion, lpa=lpa, rpa=rpa,
     )
 
+    montage += fids_montage  # add fiducials to montage
+
+    return montage
+
 
 def _easycap(basename, head_size):
     fname = op.join(MONTAGE_PATH, basename)
-    options = dict(skip_header=1, dtype=(_str, 'i4', 'i4'))
-    ch_names, theta, phi = _safe_np_loadtxt(fname, **options)
+    # ignore existing fiducials to adjust to mne head coord frame
+    fid_names = None
+    montage = _read_theta_phi_in_degrees(fname, head_size, fid_names)
 
-    radii = np.full(len(phi), head_size)
-    pos = _sph_to_cart(np.stack(
-        [radii, np.deg2rad(phi), np.deg2rad(theta)],
-        axis=-1,
-    ))
-    nasion = np.concatenate([[0],  pos[ch_names.index('Fpz'), 1:]])
-    nasion *= head_size / np.linalg.norm(nasion)
-    lpa = np.mean([pos[ch_names.index('FT9')],
-                   pos[ch_names.index('TP9')]], axis=0)
+    ch_pos = montage._get_ch_pos()
+
+    nasion = np.concatenate([[0], ch_pos['Fpz'][1:]])
+    lpa = np.mean([ch_pos['FT9'],
+                   ch_pos['TP9']], axis=0)
     lpa *= head_size / np.linalg.norm(lpa)  # on sphere
-    rpa = np.mean([pos[ch_names.index('FT10')],
-                   pos[ch_names.index('TP10')]], axis=0)
+    rpa = np.mean([ch_pos['FT10'],
+                   ch_pos['TP10']], axis=0)
     rpa *= head_size / np.linalg.norm(rpa)
 
-    return make_dig_montage(
-        ch_pos=OrderedDict(zip(ch_names, pos)),
+    fids_montage = make_dig_montage(
         coord_frame='unknown', nasion=nasion, lpa=lpa, rpa=rpa,
     )
 
+    montage += fids_montage  # add fiducials to montage
 
-def _hydrocel(basename, head_size):
-    fid_names = ('FidNz', 'FidT9', 'FidT10')
-    fname = op.join(MONTAGE_PATH, basename)
-    options = dict(dtype=(_str, 'f4', 'f4', 'f4'))
-    ch_names, xs, ys, zs = _safe_np_loadtxt(fname, **options)
+    return montage
 
-    pos = np.stack([xs, ys, zs], axis=-1)
-    ch_pos = OrderedDict(zip(ch_names, pos))
-    nasion, lpa, rpa = [ch_pos.pop(n) for n in fid_names]
-    scale = head_size / np.median(np.linalg.norm(pos, axis=-1))
-    for value in ch_pos.values():
-        value *= scale
-    nasion *= scale
-    lpa *= scale
-    rpa *= scale
 
-    return make_dig_montage(ch_pos=ch_pos, coord_frame='unknown',
-                            nasion=nasion, rpa=rpa, lpa=lpa)
+def _hydrocel(basename, head_size):
+    fname = op.join(MONTAGE_PATH, basename)
+    return _read_sfp(fname, head_size)
 
 
 def _str_names(ch_names):
@@ -103,22 +84,9 @@ def _safe_np_loadtxt(fname, **kwargs):
 
 
 def _biosemi(basename, head_size):
-    fid_names = ('Nz', 'LPA', 'RPA')
     fname = op.join(MONTAGE_PATH, basename)
-    options = dict(skip_header=1, dtype=(_str, 'i4', 'i4'))
-    ch_names, theta, phi = _safe_np_loadtxt(fname, **options)
-
-    radii = np.full(len(phi), head_size)
-    pos = _sph_to_cart(np.stack(
-        [radii, np.deg2rad(phi), np.deg2rad(theta)],
-        axis=-1,
-    ))
-
-    ch_pos = OrderedDict(zip(ch_names, pos))
-    nasion, lpa, rpa = [ch_pos.pop(n) for n in fid_names]
-
-    return make_dig_montage(ch_pos=ch_pos, coord_frame='unknown',
-                            nasion=nasion, lpa=lpa, rpa=rpa)
+    fid_names = ('Nz', 'LPA', 'RPA')
+    return _read_theta_phi_in_degrees(fname, head_size, fid_names)
 
 
 def _mgh_or_standard(basename, head_size):
@@ -193,3 +161,161 @@ def _mgh_or_standard(basename, head_size):
     'standard_primed': partial(_mgh_or_standard,
                                basename='standard_primed.elc'),
 }
+
+
+def _read_sfp(fname, head_size):
+    """Read .sfp BESA/EGI files."""
+    # fname has been already checked
+    fid_names = ('FidNz', 'FidT9', 'FidT10')
+    options = dict(dtype=(_str, 'f4', 'f4', 'f4'))
+    ch_names, xs, ys, zs = _safe_np_loadtxt(fname, **options)
+
+    pos = np.stack([xs, ys, zs], axis=-1)
+    ch_pos = OrderedDict(zip(ch_names, pos))
+    # no one grants that fid names are there.
+    nasion, lpa, rpa = [ch_pos.pop(n, None) for n in fid_names]
+
+    if head_size is not None:
+        scale = head_size / np.median(np.linalg.norm(pos, axis=-1))
+        for value in ch_pos.values():
+            value *= scale
+        nasion = nasion * scale if nasion is not None else None
+        lpa = lpa * scale if lpa is not None else None
+        rpa = rpa * scale if rpa is not None else None
+
+    return make_dig_montage(ch_pos=ch_pos, coord_frame='unknown',
+                            nasion=nasion, rpa=rpa, lpa=lpa)
+
+
+def _read_csd(fname, head_size):
+    # Label, Theta, Phi, Radius, X, Y, Z, off sphere surface
+    options = dict(comments='//',
+                   dtype=(_str, 'f4', 'f4', 'f4', 'f4', 'f4', 'f4', 'f4'))
+    ch_names, _, _, _, xs, ys, zs, _ = _safe_np_loadtxt(fname, **options)
+    pos = np.stack([xs, ys, zs], axis=-1)
+
+    if head_size is not None:
+        pos *= head_size / np.median(np.linalg.norm(pos, axis=1))
+
+    return make_dig_montage(
+        ch_pos=OrderedDict(zip(ch_names, pos)),
+    )
+
+
+def _read_elc(fname, head_size):
+    """Read .elc files.
+
+    Parameters
+    ----------
+    fname : str
+        File extension is expected to be '.elc'.
+    head_size : float | None
+        The size of the head in [m]. If none, returns the values read from the
+        file with no modification.
+
+    Returns
+    -------
+    montage : instance of DigMontage
+        The montage in [m].
+    """
+    fid_names = ('Nz', 'LPA', 'RPA')
+
+    ch_names_, pos = [], []
+    with open(fname) as fid:
+        # _read_elc does require to detect the units. (see _mgh_or_standard)
+        for line in fid:
+            if 'UnitPosition' in line:
+                units = line.split()[1]
+                scale = dict(m=1., mm=1e-3)[units]
+                break
+        else:
+            raise RuntimeError('Could not detect units in file %s' % fname)
+        for line in fid:
+            if 'Positions\n' in line:
+                break
+        pos = []
+        for line in fid:
+            if 'Labels\n' in line:
+                break
+            pos.append(list(map(float, line.split())))
+        for line in fid:
+            if not line or not set(line) - {' '}:
+                break
+            ch_names_.append(line.strip(' ').strip('\n'))
+
+    pos = np.array(pos) * scale
+    if head_size is not None:
+        pos *= head_size / np.median(np.linalg.norm(pos, axis=1))
+
+    ch_pos = OrderedDict(zip(ch_names_, pos))
+    nasion, lpa, rpa = [ch_pos.pop(n, None) for n in fid_names]
+
+    return make_dig_montage(ch_pos=ch_pos, coord_frame='unknown',
+                            nasion=nasion, lpa=lpa, rpa=rpa)
+
+
+def _read_theta_phi_in_degrees(fname, head_size, fid_names):
+    options = dict(skip_header=1, dtype=(_str, 'i4', 'i4'))
+    ch_names, theta, phi = _safe_np_loadtxt(fname, **options)
+
+    radii = np.full(len(phi), head_size)
+    pos = _sph_to_cart(np.stack(
+        [radii, np.deg2rad(phi), np.deg2rad(theta)],
+        axis=-1,
+    ))
+
+    nasion, lpa, rpa = None, None, None
+    if fid_names is not None:
+        ch_pos = OrderedDict(zip(ch_names, pos))
+        nasion, lpa, rpa = [ch_pos.pop(n, None) for n in fid_names]
+
+    return make_dig_montage(ch_pos=ch_pos, coord_frame='unknown',
+                            nasion=nasion, lpa=lpa, rpa=rpa)
+
+
+def _read_elp_besa(fname, head_size):
+    # This .elp is not the same as polhemus elp. see _read_isotrak_elp_points
+    dtype = np.dtype('S8, S8, f8, f8, f8')
+    try:
+        data = np.loadtxt(fname, dtype=dtype, skip_header=1)
+    except TypeError:
+        data = np.loadtxt(fname, dtype=dtype, skiprows=1)
+
+    ch_names = data['f1'].astype(str).tolist()
+    az = data['f2']
+    horiz = data['f3']
+    radius = np.abs(az / 180.)
+    az = np.deg2rad(np.array([h if a >= 0. else 180 + h
+                              for h, a in zip(horiz, az)]))
+    pol = radius * np.pi
+    rad = data['f4'] / 100
+    pos = _sph_to_cart(np.array([rad, az, pol]).T)
+
+    if head_size is not None:
+        pos *= head_size / np.median(np.linalg.norm(pos, axis=1))
+
+    return make_dig_montage(ch_pos=OrderedDict(zip(ch_names, pos)))
+
+
+def _read_brainvision(fname, head_size, unit):
+    # 'BrainVision Electrodes File' format
+    # Based on BrainVision Analyzer coordinate system: Defined between
+    # standard electrode positions: X-axis from T7 to T8, Y-axis from Oz to
+    # Fpz, Z-axis orthogonal from XY-plane through Cz, fit to a sphere if
+    # idealized (when radius=1), specified in millimeters
+    if unit not in ['auto', 'mm']:
+        raise ValueError('`unit` must be "auto" or "mm" for .bvef files.')
+    root = ElementTree.parse(fname).getroot()
+    ch_names = [s.text for s in root.findall("./Electrode/Name")]
+    theta = [float(s.text) for s in root.findall("./Electrode/Theta")]
+    pol = np.deg2rad(np.array(theta))
+    phi = [float(s.text) for s in root.findall("./Electrode/Phi")]
+    az = np.deg2rad(np.array(phi))
+    rad = [float(s.text) for s in root.findall("./Electrode/Radius")]
+    rad = np.array(rad)  # specified in mm
+    pos = _sph_to_cart(np.array([rad, az, pol]).T)
+
+    if head_size is not None:
+        pos *= head_size / np.median(np.linalg.norm(pos, axis=1))
+
+    return make_dig_montage(ch_pos=OrderedDict(zip(ch_names, pos)))