wip:add functionality for variable rois

tests fail for pointing mode sample data (but not for the miniscan).
the shapes and datatypes etc all match up, but the signature values are
different. This is expected because it seems, at least for the pointing
mode, the values in the TDMS file are written as
(roi1_cycle0,...,roi10_cycle0, roi1_cycle1,...,roi10_cycle1,
roi1_cycle2,...) and we are reading the first cycles_per_trial values
for roi0. But I don't understand why the miniscan test passes?

However, this code does create a reasonable-looking variable roi file for real life
data (more tests needed in the future for this though)

The failure is related to the order of things inside the TDMS file, and
the code executed in the two inner loops (channel and roi) of
_write_roi_data. But I am struggling to figure out how one would do this
well.

Maybe we need to handle legacy and variable roi code completely
separately, although it would be more elegant if we didn't have to do
this.

src/silverlabnwb/nwb_file.py

@@ -702,6 +702,7 @@ def read_functional_data(self, folder_path):
         gains = self.imaging_info.gains
         # Iterate over ROIs, which are nested inside each imaging plane section
         all_rois = {}
+        all_roi_dimensions = {}
         seg_iface = self.nwb_file.processing['Acquired_ROIs'].get("ImageSegmentation")
         for plane_name, plane in seg_iface.plane_segmentations.items():
             self.log('  Defining ROIs for plane {}', plane_name)
@@ -711,14 +712,18 @@ def read_functional_data(self, folder_path):
             # https://github.com/NeurodataWithoutBorders/pynwb/issues/673
             for roi_num, roi_ind in self.roi_mapping[plane_name].items():
                 roi_name = 'ROI_{:03d}'.format(roi_num)
+                roi_dimensions = plane[roi_ind, 'dimensions']
+                if 'pixels_per_miniscan' in self.roi_mapping.keys() and 'num_lines' in self.roi_mapping.keys():
+                    all_roi_dimensions[roi_num] = [int(plane[roi_ind, 'pixels_per_miniscan']), int(plane[roi_ind, 'num_lines'])]
+                else:
+                    all_roi_dimensions[roi_num] = plane[roi_ind, 'dimensions']
                 if roi_num not in all_rois.keys():
                     all_rois[roi_num] = {}
                 for ch, channel in {'A': 'Red', 'B': 'Green'}.items():
                     # Set zero data for now; we'll read the real data later
                     # TODO: The TDMS uses 64 bit floats; we may not really need that precision!
                     # The exported data seems to be rounded to unsigned ints. Issue #15.
-                    roi_dimensions = plane[roi_ind, 'dimensions']
-                    data_shape = np.concatenate((roi_dimensions, [num_times]))[::-1]
+                    data_shape = np.concatenate((all_roi_dimensions[roi_num], [num_times]))[::-1]
                     data = np.zeros(data_shape, dtype=np.float64)
                     # Create the timeseries object and fill in standard metadata
                     ts_name = 'ROI_{:03d}_{}'.format(roi_num, channel)
@@ -728,7 +733,7 @@ def read_functional_data(self, folder_path):
                     data_attrs['format'] = 'raw'
                     pixel_size_in_m = self.imaging_info.field_of_view / 1e6 / self.imaging_info.frame_size
                     data_attrs['field_of_view'] = roi_dimensions * pixel_size_in_m
-                    data_attrs['imaging_plane'] = plane.imaging_plane
+                    data_attrs['imaging_plane'] = self.roi_reader.get_roi_imaging_plane(roi_num, plane_name, self)
                     data_attrs['pmt_gain'] = gains[channel]
                     data_attrs['scan_line_rate'] = 1 / self.cycle_time
                     # TODO The below are not supported, so will require an extension
@@ -761,13 +766,7 @@ def read_functional_data(self, folder_path):
                     #     series_ref_in_epoch.make_group('timeseries', ts)
         # We need to write the zero-valued timeseries before editing them!
         self._write()
-        # The shape to put the TDMS data in for more convenient indexing
-        # TODO Are the roi_dimensions always the same across ROIs? (it seems that
-        # this was the implication from the previous version of the code, as it
-        # always used the last value of roi_dimensions - but that may be a bug?)
-        ch_data_shape = np.concatenate((roi_dimensions,
-                                        [len(all_rois), cycles_per_trial]))[::-1]
-        self._write_roi_data(all_rois, len(trials), cycles_per_trial, ch_data_shape, folder_path)
+        self._write_roi_data(all_rois, len(trials), cycles_per_trial, all_roi_dimensions, folder_path)
 
     def add_custom_silverlab_data(self, include_opto=True):
         metadata_class = get_class('SilverLabMetaData', 'silverlab_extended_schema')
@@ -795,7 +794,7 @@ def add_custom_silverlab_data(self, include_opto=True):
         self._write()
 
     def _write_roi_data(self, all_rois, num_trials, cycles_per_trial,
-                        ch_data_shape, folder_path):
+                        all_roi_dimensions_pixels, folder_path):
         """Edit the NWB file directly to add the real ROI data."""
         with h5py.File(self.nwb_path, 'a') as out_file:
             # Iterate over trials, reading data from the TDMS file for each
@@ -811,12 +810,17 @@ def _write_roi_data(self, all_rois, num_trials, cycles_per_trial,
                     # TODO: Consider precision: the round() here is to match the exported data...
                     ch_data = np.round(tdms_file.channel_data('Functional Imaging Data',
                                                               'Channel {} Data'.format(ch)))
-                    ch_data = ch_data.reshape(ch_data_shape)
                     # Copy each ROI's data into the NWB
+                    offset_from_previous_rois = 0
                     for roi_num, data_paths in all_rois.items():
+                        roi_shape = all_roi_dimensions_pixels[roi_num]
+                        roi_ch_data = ch_data[offset_from_previous_rois:
+                                              offset_from_previous_rois + cycles_per_trial*roi_shape[0]*roi_shape[1]]
+                        roi_ch_data = roi_ch_data.reshape(np.concatenate((roi_shape, [cycles_per_trial]))[::-1])
                         channel_path = out_file[data_paths[channel]]
-                        channel_path[time_segment, ...] = ch_data[:, roi_num - 1, ...]
-        # Update our reference to the NWB file, since it's now out of sync
+                        channel_path[time_segment, ...] = roi_ch_data
+                        offset_from_previous_rois = offset_from_previous_rois + cycles_per_trial*roi_shape[0]*roi_shape[1]
+                        # Update our reference to the NWB file, since it's now out of sync
         # We need to keep a reference to the IO object, as the file contents are
         # not read until needed
         self.nwb_io = NWBHDF5IO(self.nwb_path, 'r')
@@ -985,7 +989,7 @@ def add_rois(self, roi_path):
         organised by ROI number and channel name, so we can iterate there. Issue #16.
         """
         self.log('Loading ROI locations from {}', roi_path)
-        roi_data = self.roi_reader.read_roi_table(roi_path)
+        self.roi_data = self.roi_reader.read_roi_table(roi_path)
         module = self.nwb_file.create_processing_module(
             'Acquired_ROIs',
             'ROI locations and acquired fluorescence readings made directly by the AOL microscope.')
@@ -996,10 +1000,10 @@ def add_rois(self, roi_path):
         self.custom_silverlab_dict['imaging_mode'] = self.mode.name
         if self.mode is Modes.pointing:
             # Sanity check that each ROI is a single pixel
-            assert np.all(roi_data.num_pixels == 1)
+            assert np.all(self.roi_data.num_pixels == 1)
         # Figure out which plane each ROI is in
-        assert (roi_data['z_start'] == roi_data['z_stop']).all()  # Planes are flat in Z
-        grouped = roi_data.groupby('z_start', sort=False)
+        assert (self.roi_data['z_start'] == self.roi_data['z_stop']).all()  # Planes are flat in Z
+        grouped = self.roi_data.groupby('z_start', sort=False)
         # Iterate over planes and define ROIs
         self.roi_mapping = {}  # mapping from ROI ID to row index (used to look up ROIs)
         for plane_z, roi_group in grouped:
@@ -1030,9 +1034,10 @@ def add_rois(self, roi_path):
                 if self.mode is Modes.pointing:
                     assert row.num_pixels == 1, 'Unexpectedly large ROI in pointing mode'
                     num_x_pixels = num_y_pixels = 1
-                assert row.num_pixels == num_x_pixels * num_y_pixels, (
-                    'ROI is not rectangular: {} != {} * {}'.format(
-                        row.num_pixels, num_x_pixels, num_y_pixels))
+                # FIXME: the assertion below needs to take the resolution into account.
+                # assert row.num_pixels == num_x_pixels * num_y_pixels, (
+                #     'ROI is not rectangular: {} != {} * {}'.format(
+                #         row.num_pixels, num_x_pixels, num_y_pixels))
                 # Record the ROI dimensions for ease of lookup when adding functional data
                 dimensions = np.array([num_x_pixels, num_y_pixels], dtype=np.int32)
                 for i in range(row.num_pixels):

src/silverlabnwb/rois.py

@@ -52,6 +52,11 @@ def get_reader(cls, header):
     def __init__(self):
         pass
 
+    @classmethod
+    def get_roi_imaging_plane(cls, roi_number, plane_name, nwb_file):
+        """Get the imaging plane belonging to a specific ROI"""
+        return nwb_file.nwb_file.processing['Acquired_ROIs'].get("ImageSegmentation")[plane_name].imaging_plane
+
     @property
     def columns(self):
         column_descriptions = {
@@ -115,5 +120,26 @@ def __init__(self):
 
 class RoiReaderv300Variable(RoiReaderv300):
     """A reader for LabView version 3.0.0, supporting variable shape ROIs."""
-    # FIXME This should do something different when we ask to get the plane for a ROI.
-    pass
+
+    def get_roi_imaging_plane(self, roi_number, plane_name, nwb_file):
+        """
+        Gets the imaging plane for a variable size ROI.
+
+        Overrides base method for variable size ROIs to create (if necessary, it may have been created previously
+        for a different optical channel) and return an appropriately spatially calibrated imaging plane. This is
+        needed because variable size ROIs need to have their own imaging plane, unlike fixed sized ROIs,
+        which can share an imaging plane with other ROIs. """
+        roi_row_idx = [index for index, value in enumerate(nwb_file.roi_data.roi_index) if value == roi_number]
+        assert len(roi_row_idx) == 1
+        roi_row_idx = roi_row_idx[0]
+        resolution = nwb_file.roi_data.resolution[roi_row_idx]
+        z_plane = nwb_file.nwb_file.processing['Acquired_ROIs'].get("ImageSegmentation")[plane_name].imaging_plane
+        new_plane_name = z_plane.name + '_ROI_' + str(roi_number)
+        if new_plane_name not in nwb_file.nwb_file.imaging_planes.keys():
+            nwb_file.add_imaging_plane(
+                name=new_plane_name,
+                description='Imaging plane for variable size ROI nr. ' + str(roi_number),
+                origin_coords=z_plane.origin_coords,
+                grid_spacing=z_plane.grid_spacing*resolution
+            )
+        return nwb_file.nwb_file.imaging_planes[new_plane_name]