Merge pull request #206 from sangfrois/split_utility

Add the possibility to split multi-run physiological recordings

.gitignore

@@ -126,3 +126,4 @@ dmypy.json
 
 docs/generated/
 .vscode/
+phys2bids/tests/data/*

docs/howto.rst

@@ -279,6 +279,74 @@ By looking at this figure, we can work out that we need a smaller threshold in o
     Tip: Time 0 is the time of first trigger
     ------------------------------------------------
 
+Splitting your input file into multiple run output files
+^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
+
+If your file contains more than one (f)MRI acquisition (or runs), you can provide multiple values to ``-ntp`` in order to get multiple ``.tsv.gz`` outputs. If the TR of the entire session is consistent (i.e. all your acquisitions have the same TR), then you can specify one value after ``-tr``.
+
+By specifying the number of timepoints in each acquisition, ``phys2bids`` will recursively cut the input file by detecting the first trigger of the entire session and the ones after the number of timepoints you specified.
+
+.. code-block:: shell
+
+    phys2bids -in two_scans_samefreq_all.txt -chtrig 2 -ntp 536 398 -tr 1.2 -thr 2
+
+Now, instead of counting the trigger timepoints once, ``physbids`` will check the trigger channel recursively with all the values listed in ``-ntp``. The logger will inform you about the number of timepoints left at each iteration.
+
+.. code-block:: shell
+
+    INFO:phys2bids.physio_obj:Counting trigger points
+    INFO:phys2bids.physio_obj:The trigger is in channel 2
+    INFO:phys2bids.physio_obj:The number of timepoints found with the manual threshold of 2.0000 is 934
+    INFO:phys2bids.physio_obj:Checking number of timepoints
+    INFO:phys2bids.physio_obj:Found just the right amount of timepoints!
+    WARNING:phys2bids.slice4phys:
+    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    phys2bids will split the input file according to the given -tr and -ntp arguments
+    ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    INFO:phys2bids.physio_obj:Counting trigger points
+    INFO:phys2bids.physio_obj:The trigger is in channel 2
+    INFO:phys2bids.physio_obj:The number of timepoints found with the manual threshold of 2.0000 is 934
+    INFO:phys2bids.physio_obj:Checking number of timepoints
+    WARNING:phys2bids.physio_obj:Found 398 timepoints more than expected!
+    Assuming extra timepoints are at the end (try again with a more liberal thr)
+    INFO:phys2bids.slice4phys:
+    --------------------------------------------------------------
+    Slicing between 0.0 seconds and 961.381 seconds
+    --------------------------------------------------------------
+    INFO:phys2bids.physio_obj:Counting trigger points
+    INFO:phys2bids.physio_obj:The trigger is in channel 2
+    INFO:phys2bids.physio_obj:The number of timepoints found with the manual threshold of 2.0000 is 400
+    INFO:phys2bids.physio_obj:Checking number of timepoints
+    WARNING:phys2bids.physio_obj:Found 2 timepoints more than expected!
+    Assuming extra timepoints are at the end (try again with a more liberal thr)
+    INFO:phys2bids.slice4phys:
+    --------------------------------------------------------------
+    Slicing between 952.381 seconds and 1817.96 seconds
+    --------------------------------------------------------------
+    INFO:phys2bids.viz:Plot trigger
+    INFO:phys2bids.viz:Plot trigger
+    INFO:phys2bids.phys2bids:Preparing 2 output files.
+    INFO:phys2bids.phys2bids:Exporting files for run 1 freq 1000.0
+
+The logger also notifies you about the slicing points used (the first always being from the beginning of session, until the specified number of timepoints after the first trigger). The user can also check the resulting slice by looking at the plot of the trigger channel for each run. Each slice is adjusted with a padding after the last trigger. Such padding can be specified while calling ``phys2bids`` with ``-pad``. If nothing is specified, the default value of 9 seconds will be used. This padding is also applied at the beginning (before the first trigger of the run) of the 2nd to last run.
+
+What if I have multiple acquisition types ?
+*******************************************
+
+The user can also benefit from this utility when dealing with multiple **acquisition types** such as different functional scans with different TRs. Like ``-ntp``, ``-tr`` can take multiple values. **Though, if more than one value is specified, they require the same amount of values**. The idea is simple: if you only have one acquisition type, the one ``-tr`` input you gave will be broadcast through all runs, but if you have different acquisition types, you have to list all of them in order.
+
+.. warning::
+    There are currently no ``multi-run tutorial files`` available along with the package (under ``phys2bids/tests/data``). Although, you can visit `phys2bids OSF <https://osf.io/3txqr/files/>`_ storage to access a LabChart physiological recording with multiple fMRI acquisitions. Find it under ``labchart/chicago``.
+
+.. note::
+    **Why would I have more than one fMRI acquisition in the physiological recording?**
+
+    The idea is to reduce human error and have a good padding around your fMRI scan!
+
+    Synchronization between start of both fMRI and physiological acquisitions can be difficult, so it is safer to have as few physiological recordings as possible, with multiple imaging sequences.
+
+    Moreover, if you want to correct for recording/physiological delays, you will love *that bit of recorded information* around your fMRI scan!
+
 Generating outputs in BIDs format
 #################################
 

phys2bids/bids.py

@@ -1,5 +1,5 @@
-import os
 import logging
+import os
 from csv import reader
 from pathlib import Path
 
@@ -25,6 +25,9 @@
                 # ampere: electric current
                 'ampere': 'A', 'amp': 'A', 'amps': 'A',
                 # second: time and hertzs: frequency
+                # siemens: electric conductance (e.g. EDA)
+                'siemens': 'S',
+                # second: time and hertzs
                 '1/hz': 's', '1/hertz': 's', 'hz': 'Hz',
                 '1/s': 'Hz', '1/second': 'Hz', '1/seconds': 'Hz',
                 '1/sec': 'Hz', '1/secs': 'Hz', 'hertz': 'Hz',
@@ -96,7 +99,7 @@ def bidsify_units(orig_unit):
     return orig_unit
 
 
-def use_heuristic(heur_file, sub, ses, filename, outdir, record_label=''):
+def use_heuristic(heur_file, sub, ses, filename, outdir, run='', record_label=''):
     """
     Import and use the heuristic specified by the user to rename the file.
 
@@ -130,7 +133,7 @@ def use_heuristic(heur_file, sub, ses, filename, outdir, record_label=''):
 
     # Initialise a dictionary of bids_keys that has already "recording"
     bids_keys = {'sub': '', 'ses': '', 'task': '', 'acq': '', 'ce': '',
-                 'dir': '', 'rec': '', 'run': '', 'recording': record_label}
+                 'dir': '', 'rec': '', 'run': run, 'recording': record_label}
 
     # Start filling bids_keys dictionary and path with subject and session
     if sub.startswith('sub-'):
@@ -150,7 +153,7 @@ def use_heuristic(heur_file, sub, ses, filename, outdir, record_label=''):
 
     # Load heuristic and use it to fill dictionary
     heur = utils.load_heuristic(heur_file)
-    bids_keys.update(heur.heur(Path(filename).stem))
+    bids_keys.update(heur.heur(Path(filename).stem, run))
 
     # If bids_keys['task'] is still empty, stop the program
     if not bids_keys['task']:

phys2bids/cli/run.py

@@ -1,7 +1,6 @@
 # -*- coding: utf-8 -*-
-"""
-Parser for phys2bids
-"""
+"""Parser for phys2bids."""
+
 
 import argparse
 
@@ -10,7 +9,7 @@
 
 def _get_parser():
     """
-    Parses command line inputs for this function
+    Parse command line inputs for this function.
 
     Returns
     -------
@@ -27,8 +26,8 @@ def _get_parser():
     required.add_argument('-in', '--input-file',
                           dest='filename',
                           type=str,
-                          help='The name of the file containing physiological data, with or '
-                               'without extension.',
+                          help='The name of the file containing physiological '
+                               'data, with or without extension.',
                           required=True)
     optional.add_argument('-info', '--info',
                           dest='info',
@@ -61,11 +60,6 @@ def _get_parser():
                                'in the current folder. Edit the heur_ex.py file in '
                                'heuristics folder.',
                           default=None)
-    # optional.add_argument('-hdir', '--heur-dir',
-    #                       dest='heurdir',
-    #                       type=str,
-    #                       help='Folder containing heuristic file.',
-    #                       default='.')
     optional.add_argument('-sub', '--subject',
                           dest='sub',
                           type=str,
@@ -94,26 +88,39 @@ def _get_parser():
                           default=None)
     optional.add_argument('-ntp', '--numtps',
                           dest='num_timepoints_expected',
+                          nargs='*',
                           type=int,
-                          help='Number of expected timepoints (TRs). '
-                               'Default is 0. Note: the estimation of when the '
-                               'neuroimaging acquisition started cannot take place '
-                               'with this default.',
-                          default=0)
+                          help='Number of expected trigger timepoints (TRs). '
+                               'Default is None. Note: the estimation of beggining of '
+                               'neuroimaging acquisition cannot take place with this default. '
+                               'If you\'re running phys2bids on a multi-run recording, '
+                               'give a list of each expected ntp for each run.',
+                          default=None)
     optional.add_argument('-tr', '--tr',
                           dest='tr',
+                          nargs='*',
                           type=float,
                           help='TR of sequence in seconds. '
-                               'Default is 0 second.',
-                          default=0)
+                               'If you\'re running phys2bids on a multi-run recording, '
+                               'you can give a list of each expected ntp for each run, '
+                               'or just one TR if it is consistent throughout the session.',
+                          default=None)
     optional.add_argument('-thr', '--threshold',
                           dest='thr',
                           type=float,
                           help='Threshold to use for trigger detection. '
-                               'If "ntp" and "TR" are specified, phys2bids automatically computes '
-                               'a threshold to detect the triggers. Use this parameter to set it '
-                               'manually',
-                               default=None)
+                               'If "ntp" and "TR" are specified, phys2bids '
+                               'automatically computes a threshold to detect '
+                               'the triggers. Use this parameter to set it manually. '
+                               'This parameter is necessary for multi-run recordings. ',
+                          default=None)
+    optional.add_argument('-pad', '--padding',
+                          dest='pad',
+                          type=float,
+                          help='Padding in seconds used around a single run '
+                               'when separating multi-run session files. '
+                               'Default is 9 seconds.',
+                          default=9)
     optional.add_argument('-chnames', '--channel-names',
                           dest='ch_name',
                           nargs='*',

phys2bids/heuristics/heur_euskalibur.py

@@ -1,7 +1,7 @@
 import fnmatch
 
 
-def heur(physinfo):
+def heur(physinfo, run=''):
     """
     Set of if .. elif statements to fill BIDS names.
 

phys2bids/heuristics/heur_test_acq.py

@@ -1,7 +1,7 @@
 import fnmatch
 
 
-def heur(physinfo):
+def heur(physinfo, run=''):
     """
     Set of if .. elif statements to fill BIDS names.
 

phys2bids/heuristics/heur_test_multifreq.py

@@ -1,7 +1,7 @@
 import fnmatch
 
 
-def heur(physinfo):
+def heur(physinfo, run=''):
     """
     Set of if .. elif statements to fill BIDS names.
 

phys2bids/heuristics/heur_tutorial.py

@@ -1,7 +1,7 @@
 import fnmatch
 
 
-def heur(physinfo):
+def heur(physinfo, run=''):
     """
     Set of if .. elif statements to fill BIDS names.
 

phys2bids/interfaces/acq.py

@@ -1,9 +1,7 @@
 #!/usr/bin/env python3
 # -*- coding: utf-8 -*-
 
-"""
-phys2bids interface for acqknowledge files.
-"""
+"""phys2bids interface for acqknowledge files."""
 import logging
 import warnings
 

phys2bids/interfaces/txt.py

@@ -1,9 +1,7 @@
 #!/usr/bin/env python3
 # -*- coding: utf-8 -*-
 
-"""
-phys2bids interface for txt files.
-"""
+"""phys2bids interface for txt files."""
 
 import logging
 from collections import Counter