Repository: bids-standard/pybv

Issue:
potential indexing issue
Could it be that we have an indexing issue in our code?

BrainVision uses 1-based indexing for writing the events, whereas in Python and MNE-Python, we use 0-based indexing.

Now imagine a user passing `events`: https://github.com/bids-standard/pybv/blob/cf8aa9d6904149d68ad9d9500315db5e785664a2/pybv/io.py#L46-L50

These events will most likely be 0-indexed.

However, later down in the code -when we write a VMRK file, we translate these indices one to one, disregarding that BrainVision used 1-based indexing. See: https://github.com/bids-standard/pybv/blob/cf8aa9d6904149d68ad9d9500315db5e785664a2/pybv/io.py#L163-L167

Later on in the testing, we actually see this issue: https://github.com/bids-standard/pybv/blob/cf8aa9d6904149d68ad9d9500315db5e785664a2/pybv/tests/test_bv_writer.py#L95-L96

note the comment. I may be wrong ... but it seems to me we should handle the indexing better

cc @tstenner @choldgraf 

Discussion:
codecov[bot]: # [Codecov](https://codecov.io/gh/bids-standard/pybv/pull/29?src=pr&el=h1) Report
> Merging [#29](https://codecov.io/gh/bids-standard/pybv/pull/29?src=pr&el=desc) into [master](https://codecov.io/gh/bids-standard/pybv/commit/ae5ab65ab2819a3cb1e3643e631b1bba5ba198fe?src=pr&el=desc) will **increase** coverage by `0.4%`.
> The diff coverage is `100%`.

[![Impacted file tree graph](https://codecov.io/gh/bids-standard/pybv/pull/29/graphs/tree.svg?width=650&token=L7VvFoE1SK&height=150&src=pr)](https://codecov.io/gh/bids-standard/pybv/pull/29?src=pr&el=tree)

```diff
@@           Coverage Diff            @@
##           master     #29     +/-   ##
========================================
+ Coverage   94.89%   95.3%   +0.4%     
========================================
  Files           3       3             
  Lines         196     213     +17     
========================================
+ Hits          186     203     +17     
  Misses         10      10
```


| [Impacted Files](https://codecov.io/gh/bids-standard/pybv/pull/29?src=pr&el=tree) | Coverage Δ | |
|---|---|---|
| [pybv/io.py](https://codecov.io/gh/bids-standard/pybv/pull/29/diff?src=pr&el=tree#diff-cHlidi9pby5weQ==) | `92.7% <100%> (+0.39%)` | :arrow_up: |
| [pybv/tests/test\_bv\_writer.py](https://codecov.io/gh/bids-standard/pybv/pull/29/diff?src=pr&el=tree#diff-cHlidi90ZXN0cy90ZXN0X2J2X3dyaXRlci5weQ==) | `100% <100%> (ø)` | :arrow_up: |

------

[Continue to review full report at Codecov](https://codecov.io/gh/bids-standard/pybv/pull/29?src=pr&el=continue).
> **Legend** - [Click here to learn more](https://docs.codecov.io/docs/codecov-delta)
> `Δ = absolute <relative> (impact)`, `ø = not affected`, `? = missing data`
> Powered by [Codecov](https://codecov.io/gh/bids-standard/pybv/pull/29?src=pr&el=footer). Last update [ae5ab65...b06d01d](https://codecov.io/gh/bids-standard/pybv/pull/29?src=pr&el=lastupdated). Read the [comment docs](https://docs.codecov.io/docs/pull-request-comments).

sappelhoff: Thanks for the review @choldgraf All tests pass, I am merging this now and preparing a new release :-)

Patch:
diff --git a/docs/changelog.rst b/docs/changelog.rst
index 4e671ea..48ace17 100644
--- a/docs/changelog.rst
+++ b/docs/changelog.rst
@@ -17,10 +17,12 @@ current
 
 Changelog
 ~~~~~~~~~
+Add measurement date parameter to public API, by `Stefan Appelhoff`_ (`#29 <https://github.com/bids-standard/pybv/pull/29>`_)
 Add binary format parameter to public API by `Tristan Stenner`_ (`#22 <https://github.com/bids-standard/pybv/pull/22>`_)
 
 Bug
 ~~~
+fix bug with events indexing. VMRK events are now correctly written with 1-based indexing, by `Stefan Appelhoff`_ (`#29 <https://github.com/bids-standard/pybv/pull/29>`_)
 fix bug with events that only have integer codes of length less than 3, by `Stefan Appelhoff`_ (`#26 <https://github.com/bids-standard/pybv/pull/26>`_)
 
 0.0.2
diff --git a/pybv/io.py b/pybv/io.py
index ba6615f..ff9a9ae 100644
--- a/pybv/io.py
+++ b/pybv/io.py
@@ -11,9 +11,11 @@
 import codecs
 import os
 import os.path as op
+import datetime
+
 import numpy as np
 
-from . import __version__
+from pybv import __version__
 
 # ascii as future formats
 supported_formats = {
@@ -26,7 +28,7 @@ supported_orients = {'multiplexed'}
 
 def write_brainvision(data, sfreq, ch_names, fname_base, folder_out,
                       events=None, resolution=1e-7, scale_data=True,
-                      fmt='binary_float32'):
+                      fmt='binary_float32', meas_date=None):
     """Write raw data to BrainVision format.
 
     Parameters
@@ -45,9 +47,9 @@ def write_brainvision(data, sfreq, ch_names, fname_base, folder_out,
         The folder where output files will be saved.
     events : ndarray, shape (n_events, 2)
         Events to write in the marker file. This array has two columns.
-        The first column is the index of each event (corresponding to the
-        "time" dimension of the data array). The second column is a number
-        associated with the "type" of event.
+        The first column is the zero-based index of each event (corresponding
+        to the "time" dimension of the data array). The second column is a
+        number associated with the "type" of event.
     resolution : float | ndarray
         The resolution **in volts** in which you'd like the data to be stored.
         By default, this will be 1e-7, or .1 microvolts. Since data is stored
@@ -64,6 +66,11 @@ def write_brainvision(data, sfreq, ch_names, fname_base, folder_out,
     fmt : str
         Binary format the data should be written as. Valid choices are
         'binary_float32' (default) and 'binary_int16'.
+    meas_date : datetime.datetime | str | None
+        The measurement date of the data specified as a datetime.datetime
+        object. Alternatively, can be a string in the format:
+        "YYYYMMDDhhmmssuuuuuu". "u" stands for microseconds. If None, defaults
+        to '00000000000000000000'.
     """
     # Create output file names/paths
     if not op.isdir(folder_out):
@@ -106,8 +113,22 @@ def write_brainvision(data, sfreq, ch_names, fname_base, folder_out,
 
     _chk_fmt(fmt)
 
+    # measurement date
+    if not isinstance(meas_date, (str, datetime.datetime, type(None))):
+        raise ValueError('`meas_date` must be of type str, datetime.datetime, '
+                         'or None but is of type '
+                         '"{}"'.format(type(meas_date)))
+    elif meas_date is None:
+        meas_date = '00000000000000000000'
+    elif isinstance(meas_date, datetime.datetime):
+        meas_date = meas_date.strftime('%Y%m%d%H%M%S%f')
+    elif not (meas_date.isdigit() and len(meas_date) == 20):
+        raise ValueError('Got a str for `meas_date`, but it was not formatted '
+                         'as expected. Please supply a str in the format: '
+                         '"YYYYMMDDhhmmssuuuuuu".')
+
     # Write output files
-    _write_vmrk_file(vmrk_fname, eeg_fname, events)
+    _write_vmrk_file(vmrk_fname, eeg_fname, events, meas_date)
     _write_vhdr_file(vhdr_fname, vmrk_fname, eeg_fname, data, sfreq,
                      ch_names, orientation='multiplexed', format=fmt,
                      resolution=resolution)
@@ -136,7 +157,7 @@ def _chk_multiplexed(orientation):
     return orientation == 'multiplexed'
 
 
-def _write_vmrk_file(vmrk_fname, eeg_fname, events):
+def _write_vmrk_file(vmrk_fname, eeg_fname, events, meas_date):
     """Write BrainvVision marker file."""
     with codecs.open(vmrk_fname, 'w', encoding='utf-8') as fout:
         print(r'Brain Vision Data Exchange Marker File, Version 1.0', file=fout)  # noqa: E501
@@ -149,9 +170,10 @@ def _write_vmrk_file(vmrk_fname, eeg_fname, events):
         print(r'[Marker Infos]', file=fout)
         print(r'; Each entry: Mk<Marker number>=<Type>,<Description>,<Position in data points>,', file=fout)  # noqa: E501
         print(r';             <Size in data points>, <Channel number (0 = marker is related to all channels)>', file=fout)  # noqa: E501
+        print(r';             <Date (YYYYMMDDhhmmssuuuuuu)>', file=fout)
         print(r'; Fields are delimited by commas, some fields might be omitted (empty).', file=fout)  # noqa: E501
         print(r'; Commas in type or description text are coded as "\1".', file=fout)  # noqa: E501
-        print(r'Mk1=New Segment,,1,1,0,0', file=fout)
+        print(r'Mk1=New Segment,,1,1,0,{}'.format(meas_date), file=fout)
 
         if events is None or len(events) == 0:
             return
@@ -165,11 +187,12 @@ def _write_vmrk_file(vmrk_fname, eeg_fname, events):
         twidth = twidth if twidth > 3 else 3
         tformat = 'S{:>' + str(twidth) + '}'
 
-        for ii, irow in enumerate(range(len(events)), start=2):
-            i_ix = events[irow, 0]
+        for marker_number, irow in enumerate(range(len(events)), start=2):
+            i_ix = events[irow, 0] + 1  # BrainVision uses 1-based indexing
             i_val = events[irow, 1]
             print(r'Mk{}=Stimulus,{},{},1,0'
-                  .format(ii, tformat.format(i_val), i_ix), file=fout)
+                  .format(marker_number, tformat.format(i_val), i_ix),
+                  file=fout)
 
 
 def _optimize_channel_unit(resolution):


Test patch:
diff --git a/pybv/tests/test_bv_writer.py b/pybv/tests/test_bv_writer.py
index b5ff9e9..222029b 100644
--- a/pybv/tests/test_bv_writer.py
+++ b/pybv/tests/test_bv_writer.py
@@ -11,6 +11,8 @@ import os
 import os.path as op
 from shutil import rmtree
 from tempfile import mkdtemp
+from time import gmtime
+from datetime import datetime
 
 import pytest
 
@@ -75,13 +77,31 @@ def test_bv_bad_format():
     rmtree(tmpdir)
 
 
-def test_bv_writer_oi_cycle():
+@pytest.mark.parametrize("meas_date,match",
+                         [(1, '`meas_date` must be of type str, datetime'),
+                          ('', 'Got a str for `meas_date`, but it was'),
+                          ('1973', 'Got a str for `meas_date`, but it was')])
+def test_bad_meas_date(meas_date, match):
+    """Test that bad measurement dates raise errors."""
+    tmpdir = _mktmpdir()
+    with pytest.raises(ValueError, match=match):
+        write_brainvision(data, sfreq, ch_names, fname, tmpdir,
+                          meas_date=meas_date)
+
+    rmtree(tmpdir)
+
+
+@pytest.mark.parametrize("meas_date",
+                         [('20000101120000000000'),
+                          (datetime(2000, 1, 1, 12, 0, 0, 0))])
+def test_bv_writer_oi_cycle(meas_date):
     """Test that a write-read cycle produces identical data."""
     tmpdir = _mktmpdir()
 
     # Write, then read the data to BV format
     write_brainvision(data, sfreq, ch_names, fname, tmpdir, events=events,
-                      resolution=np.power(10., -np.arange(10)))
+                      resolution=np.power(10., -np.arange(10)),
+                      meas_date=meas_date)
     vhdr_fname = op.join(tmpdir, fname + '.vhdr')
     raw_written = mne.io.read_raw_brainvision(vhdr_fname, preload=True)
     # Delete the first annotation because it's just marking a new segment
@@ -92,8 +112,8 @@ def test_bv_writer_oi_cycle():
     # sfreq
     assert sfreq == raw_written.info['sfreq']
 
-    # Event timing should be within one index of originals
-    assert_allclose(events[:, 0], events_written[:, 0], 1)
+    # Event timing should be exactly the same
+    assert_array_equal(events[:, 0], events_written[:, 0])
     assert_array_equal(events[:, 1], events_written[:, 2])
     # Should be 2 unique event types
     assert len(event_id) == 2
@@ -104,6 +124,11 @@ def test_bv_writer_oi_cycle():
     # channels
     assert ch_names == raw_written.ch_names
 
+    # measurement date, we do not test microsecs
+    unix_seconds = raw_written.info['meas_date'][0]
+    time_str = ('{:04}{:02}{:02}{:02}{:02}{:02}'.format(*gmtime(unix_seconds)))
+    assert time_str == '20000101120000'  # 1st of Jan, 2000 at 12:00 and 0 secs
+
     rmtree(tmpdir)