streamline

Author: sappelhoff

README.md

@@ -28,16 +28,15 @@ functions to project the 3D locations to 2D space and plot them.
 
 ## How to work with it
 
-- `git clone` the repository (or download as .zip and unpack)
+- `git clone` the repository (or download as `.zip` and unpack)
 - `cd eeg_positions`
 - Using your python environment of choice, install the package and its
   dependencies locally using `pip install -e .`
 - Run the tests using `pytest` (you might have to `pip install pytest` first)
-- Calculate and plot electrodes by calling `python eeg_positions.py` in the
-  `eeg_positions/eeg_positions` directory
+- Calculate and plot electrodes by calling `python eeg_positions/calc_positions.py`
 - Check out `contour_labels.py` for the order how electrodes are computed
-- ... and see `utils.py` for the `find_point_at_fraction` function that is
-  the core of the computations.
+- ... and see `utils.py` for the `find_point_at_fraction` function that is the
+  core of the computations.
 
 ## References
 
@@ -46,6 +45,19 @@ functions to project the 3D locations to 2D space and plot them.
   112(4), 713-719. doi:
   [10.1016/S1388-2457(00)00527-7](https://www.biosemi.com/publications/pdf/Oostenveld2001b.pdf)
 
+## Acknowledgements
+
+Explicit thanks to:
+
+- Robert Oostenveld for writing his [blog post](http://robertoostenveld.nl/electrode/)
+  on electrodes
+- Ed Williams for the helpful correspondence and discussions about
+  "intermediate points on a great circle" (see his
+  [aviation formulary](http://www.edwilliams.org/avform.htm#Intermediate))
+- "N. Bach" and "Nominal Animal" who helped me to figure out the math for
+  the `find_point_at_fraction` function (see this
+  [math.stackexchange.com post](https://math.stackexchange.com/questions/2800845/find-intermediate-points-on-small-circle-of-a-sphere/2805204#2805204))
+
 ---
 
 # Examples
@@ -89,13 +101,13 @@ montage.plot()
 
 ## Interactively viewing 3D coordinates
 
-reproduce by running `python eeg_positions.py`
+reproduce by running `python eeg_positions/calc_positions.py`
 
 ![img: coordinate system](./images/3d_view.png)
 
 ## Projections to 2D
 
-reproduce by running `python eeg_positions.py`
+reproduce by running `python eeg_positions/calc_positions.py`
 
 ### 10-20 system
 ![img: coordinate system](./images/1020.png)
@@ -115,8 +127,10 @@ reproduce by running `python eeg_positions.py`
 ### 3D Axes and Cartesian Coordinate System
 
 - Imagine the x-axis pointing roughly towards the viewer with increasing values
-- The y-axis is orthogonal to the x-axis, pointing to the right of the viewer with increasing values
-- The z-axis is orthogonal to the xy-plane and pointing vertically up with increasing values
+- The y-axis is orthogonal to the x-axis, pointing to the right of the viewer
+  with increasing values
+- The z-axis is orthogonal to the xy-plane and pointing vertically up with
+  increasing values
 
 ![img: coordinate system](./images/coords_cartesian.png)
 
@@ -132,15 +146,16 @@ sphere:
 
 ## Cartesian Coordinates
 
-- The left preauricular point = (-1, 0, 0) ... coincides with T9
-- The right preauricular point = (1, 0, 0) ... coincides with T10
-- The nasion = (0, 1, 0) ... coincides with Nz
-- The inion = (0, -1, 0) ... coincides with Iz
-- The vertex = (0, 0, 1) ... coincides with Cz
+- The left preauricular point = `(-1, 0, 0)` ... coincides with T9
+- The right preauricular point = `(1, 0, 0)` ... coincides with T10
+- The nasion = `(0, 1, 0)` ... coincides with Nz
+- The inion = `(0, -1, 0)` ... coincides with Iz
+- The vertex = `(0, 0, 1)` ... coincides with Cz
 
-Hence, the equator of the sphere goes through T9, T10, and Nz.
+Hence, the equator of the sphere goes through T9, T10, Nz and Iz.
 
-**Note that these are ASSUMPTIONS**. It would be equally valid to assume the equator going through T7, T8, and Fpz.
+**Note that these are ASSUMPTIONS**. It would be equally valid to assume the
+equator going through T7, T8, Fpz, and Oz.
 
 # EEG Electrode Position Data
 

eeg_positions/__init__.py

@@ -1 +1,3 @@
 """Initialize eeg_positions."""
+
+__version__ = '1.0.0'

eeg_positions/eeg_positions.py renamed to eeg_positions/calc_positions.py

@@ -8,9 +8,11 @@
 
 import pandas as pd
 
-from utils import (get_xyz, find_point_at_fraction, plot_spherical_head,
-                   plot_2d_head, stereographic_projection)
-from contour_labels import all_contours, system1020, system1010, system1005
+from eeg_positions.utils import (get_xyz, find_point_at_fraction,
+                                 plot_spherical_head, plot_2d_head,
+                                 stereographic_projection)
+from eeg_positions.contour_labels import (ALL_CONTOURS, SYSTEM1020, SYSTEM1010,
+                                          SYSTEM1005)
 
 
 if __name__ == '__main__':
@@ -32,7 +34,7 @@
 
     # Calculate all positions
     # -----------------------
-    for contour in all_contours:
+    for contour in ALL_CONTOURS:
 
         if len(contour) == 21:
             midpoint_idx = 10
@@ -49,15 +51,15 @@
 
         # Calculate all other points at fractions of distance
         # see `contour_labels.py` and `test_contour_labels.py`
-        other_points = {}
+        other_ps = {}
         for i, label in enumerate(contour):
-            other_points[label] = find_point_at_fraction(p1,
-                                                         p2,
-                                                         p3,
-                                                         f=i/(len(contour)-1))
+            other_ps[label] = find_point_at_fraction(p1,
+                                                     p2,
+                                                     p3,
+                                                     frac=i/(len(contour)-1))
 
         # Append to data frame
-        tmp = pd.DataFrame.from_dict(other_points, orient='index')
+        tmp = pd.DataFrame.from_dict(other_ps, orient='index')
         tmp.columns = ['x', 'y', 'z']
         tmp['label'] = tmp.index
         df = df.append(tmp, ignore_index=True, sort=True)
@@ -71,19 +73,19 @@
     fname_template = os.path.join(fpath, '..', 'data', 'standard_{}.tsv')
 
     # First in 3D, then in 2D for each system
-    for system, fmt in zip([system1020, system1010, system1005],
+    for system, fmt in zip([SYSTEM1020, SYSTEM1010, SYSTEM1005],
                            ['1020', '1010', '1005']):
 
         idx = df.label.isin(system)
         system_df = df.loc[idx, :]
-        system_df.sort_values(by='label', inplace=True)
+        system_df = system_df.sort_values(by='label')
         system_df.to_csv(fname_template.format(fmt), sep='\t', na_rep='n/a',
                          index=False, float_format='%.4f')
 
         # Now in 2D using stereographic projection
-        xs, ys = stereographic_projection(system_df.values[:, 1],
-                                          system_df.values[:, 2],
-                                          system_df.values[:, 3])
+        xs, ys = stereographic_projection(system_df.to_numpy()[:, 1],
+                                          system_df.to_numpy()[:, 2],
+                                          system_df.to_numpy()[:, 3])
         system_df = system_df.loc[:, ['label', 'x', 'y']]
         system_df['x'] = xs
         system_df['x'] = ys
@@ -108,11 +110,11 @@
         # 2D
         fig2, ax2 = plot_2d_head()
 
-        xs, ys = stereographic_projection(df.x, df.y, df.z)
+        xs, ys = stereographic_projection(df['x'], df['y'], df['z'])
 
         ax2.scatter(xs, ys, marker='.', color='r')
 
-        for lab, x, y in zip(list(df.label), xs, ys):
+        for lab, x, y in zip(list(df['label']), xs, ys):
             ax2.annotate(lab, xy=(x, y), fontsize=5)
 
         ax2.set_title('standard_{}'.format(system))

eeg_positions/contour_labels.py

@@ -114,7 +114,7 @@
                'POO8']
 
 # List of lists. Note the specific ordering.
-all_contours = [sagittal_Nz_Cz_Iz,
+ALL_CONTOURS = [sagittal_Nz_Cz_Iz,
                 horizontal_Nz_T10_Iz, horizontal_Nz_T9_Iz,
                 coronal_T9_Cz_T10,
                 horizontal_NFpz_T10h_OIz, horizontal_NFpz_T9h_OIz,
@@ -125,17 +125,20 @@
                 contour_PO, contour_POO]
 
 # Defining which electrodes belong into which standard system
-system1020 = ['Fp1', 'Fpz', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8',
+SYSTEM1020 = ['Fp1', 'Fpz', 'Fp2', 'F7', 'F3', 'Fz', 'F4', 'F8',
               'T7', 'C3', 'Cz', 'C4', 'T8', 'P7', 'P3', 'Pz', 'P4', 'P8',
               'O1', 'Oz', 'O2']
 
-system1010 = system1020 + ['Nz', 'AF7', 'AFz', 'AF8', 'F9', 'F5', 'F1', 'F2',
+SYSTEM1010 = SYSTEM1020 + ['Nz', 'AF7', 'AFz', 'AF8', 'F9', 'F5', 'F1', 'F2',
                            'F6', 'F10', 'FT9', 'FT7', 'FC5', 'FC3', 'FC1',
                            'FCz', 'FC2', 'FC4', 'FC6', 'FT8', 'FT10', 'T9',
                            'C5', 'C1', 'C2', 'C6', 'T10', 'TP7', 'CP5', 'CP3',
                            'CP1', 'CPz', 'CP2', 'CP4', 'CP6', 'TP8', 'P9',
                            'P5', 'P1', 'P2', 'P6', 'P10', 'PO9', 'PO7', 'POz',
                            'PO8', 'PO10', 'I1', 'Iz', 'I2']
 
-system1005 = list(set([label for contour in all_contours
-                       for label in contour]))
+SYSTEM1005 = list()
+for contour in ALL_CONTOURS:
+    for label in contour:
+        if label not in SYSTEM1005:
+            SYSTEM1005.append(label)

eeg_positions/tests/test_contour_labels.py

@@ -1,21 +1,17 @@
 """Test whether the contour labels are complete."""
 
-from eeg_positions.contour_labels import (all_contours, system1020, system1010,
-                                          system1005)
+from eeg_positions.contour_labels import (ALL_CONTOURS, SYSTEM1020, SYSTEM1010,
+                                          SYSTEM1005)
 
 
 def test_contour_lengths():
     """Check that we have 17 or 21 electrode per contour."""
-    for contour in all_contours:
+    for contour in ALL_CONTOURS:
         assert len(contour) in [17, 21]
 
 
 def test_system_labels():
     """Check if systems have the correct number of labels."""
-    assert len(system1005) == 345
-    assert len(system1020) == 21
-    assert len(system1010) == 71
-
-
-if __name__ == '__main__':
-    print(all_contours)
+    assert len(SYSTEM1005) == 345
+    assert len(SYSTEM1020) == 21
+    assert len(SYSTEM1010) == 71

eeg_positions/tests/test_utilities.py

@@ -52,16 +52,16 @@ def test_find_point_at_fraction():
     p1 = (1., 0., 0.)
     p2 = (0., 0., 1.)
     p3 = (-1., 0., 0.)
-    p = find_point_at_fraction(p1, p2, p3, f=0.)
-    assert p == p1
-    p = find_point_at_fraction(p1, p2, p3, f=1.)
-    assert p == p3
-    p = find_point_at_fraction(p1, p2, p3, f=0.5)
-    assert p == p2
+    point = find_point_at_fraction(p1, p2, p3, frac=0.)
+    assert point == p1
+    point = find_point_at_fraction(p1, p2, p3, frac=1.)
+    assert point == p3
+    point = find_point_at_fraction(p1, p2, p3, frac=0.5)
+    assert point == p2
 
     # Assert error when equal points
     with pytest.raises(ValueError):
-        find_point_at_fraction(p1, p1, p3, 0.5)
+        find_point_at_fraction(p1, p1, p3, frac=0.5)
 
 
 def test_stereographic_projection():
@@ -73,7 +73,7 @@ def test_stereographic_projection():
     df = pd.DataFrame(data)
 
     # We know where Cz and Nz should be in 2D:
-    xs, ys = stereographic_projection(df.x, df.y, df.z)
+    xs, ys = stereographic_projection(df['x'], df['y'], df['z'])
     np.testing.assert_allclose(xs, np.array((0, 0)))
     np.testing.assert_allclose(ys, np.array((0, 1)))
 

eeg_positions/utils.py

@@ -7,28 +7,28 @@
 from mpl_toolkits.mplot3d import Axes3D  # noqa: F401
 
 
-def find_point_at_fraction(p1, p2, p3, f):
+def find_point_at_fraction(p1, p2, p3, frac):
     """Find a point on an arc spanned by three points.
 
     Given three points `p1`, `p2` and `p3` on a sphere with origin (0, 0, 0),
-    find the coordinates of a point `p` at a fraction `f` of the overall
+    find the coordinates of a `point` at a fraction `frac` of the overall
     distance on an arc spanning from `p1` over `p2` to `p3`. Given this
-    assumption, for fractions of zero, `p` will equal `p1`; for fractions of
-    one, `p` will equal `p3`; and for fractions of one half, `p` will equal
-    `p2` [1]_.
+    assumption, for fractions of zero, `point` will equal `p1`; for fractions
+    of one, `point` will equal `p3`; and for fractions of one half, `point`
+    will equal `p2` [1]_.
 
     Parameters
     ----------
     p1, p2, p3 : tuple
         Each tuple containing x, y, z cartesian coordinates.
-    f : float
+    frac : float
         Fraction of distance from `p1` to `p3` over p2` at which
         to find coordinates of `p`.
 
     Returns
     -------
-    p : tuple
-     Containing x, y, z cartesian coordinates.
+    point : tuple
+        Containing x, y, z cartesian coordinates.
 
     Notes
     -----
@@ -45,13 +45,16 @@ def find_point_at_fraction(p1, p2, p3, f):
 
     Examples
     --------
-    >>> find_point_at_fraction((1., 0., 0.), (0., 0., 1.), (-1., 0., 0.), f=0.)
+    >>> p1 = (1., 0., 0.)
+    >>> p2 = (0., 0., 1.)
+    >>> p3 = (-1., 0., 0.)
+    >>> find_point_at_fraction(p1, p2, p3, frac=0.)
     (1., 0., 0.)
-    >>> find_point_at_fraction((1., 0., 0.), (0., 0., 1.), (-1., 0., 0.), f=.5)
+    >>> find_point_at_fraction(p1, p2, p3, frac=.5)
     (0., 0., 1.)
-    >>> find_point_at_fraction((1., 0., 0.), (0., 0., 1.), (-1., 0., 0.), f=1.)
+    >>> find_point_at_fraction(p1, p2, p3, frac=1.)
     (-1., 0., 0.)
-    >>> find_point_at_fraction((1., 0., 0.), (0., 0., 1.), (-1., 0., 0.), f=.3)
+    >>> find_point_at_fraction(p1, p2, p3, frac=.3)
     (0.5878, 0.0, 0.809)
 
     """
@@ -126,14 +129,14 @@ def find_point_at_fraction(p1, p2, p3, f):
         theta = theta + 2*np.pi
 
     # Now calculate coordinates at fraction
-    x = xc + xu * np.cos(f*theta) + xv*np.sin(f*theta)
-    y = yc + yu * np.cos(f*theta) + yv*np.sin(f*theta)
-    z = zc + zu * np.cos(f*theta) + zv*np.sin(f*theta)
+    x = xc + xu * np.cos(frac*theta) + xv*np.sin(frac*theta)
+    y = yc + yu * np.cos(frac*theta) + yv*np.sin(frac*theta)
+    z = zc + zu * np.cos(frac*theta) + zv*np.sin(frac*theta)
 
     # Round to 4 decimals and collect the points in tuple
-    p = np.asarray((x, y, z))
-    p = tuple(p.round(decimals=4))
-    return p
+    point = np.asarray((x, y, z))
+    point = tuple(point.round(decimals=4))
+    return point
 
 
 # Convenient helper function to access xyz coordinates from df