Skip to content

suite2p.io package#

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

BinaryFile #

Source code in suite2p/io/binary.py 
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
class BinaryFile:

    def __init__(self, Ly, Lx, filename, n_frames=None,
                 dtype="int16", write=False):
        """
        Open or create a Suite2p binary file backed by a memory-mapped numpy array.

        Parameters
        ----------
        Ly : int
            Height of each frame in pixels.
        Lx : int
            Width of each frame in pixels.
        filename : str
            Path to the binary file to read from or write to.
        n_frames : int, optional
            Number of frames. Required when creating a new file for writing.
            Inferred from file size when reading.
        dtype : str, optional (default "int16")
            Data type of each pixel value.
        write : bool, optional (default False)
            If True, open the file for reading and writing. If False, open read-only.
        """
        self.Ly = Ly
        self.Lx = Lx
        self.filename = filename
        self.dtype = dtype
        self.write = write

        if write and n_frames is None and not os.path.exists(self.filename):
            raise ValueError(
                "need to provide number of frames n_frames when writing file")
        elif not write:
            n_frames = self.n_frames
        shape = (n_frames, self.Ly, self.Lx)
        if write:
            mode = "r+" if os.path.exists(self.filename) else "w+"
        else:
            mode = "r"
        self.file = np.memmap(self.filename, mode=mode, dtype=self.dtype, shape=shape)
        self._index = 0
        self._can_read = True

    @staticmethod
    def convert_numpy_file_to_suite2p_binary(from_filename, to_filename):
        """
        Convert a numpy file (.npy, .npz) to a Suite2p binary file.

        Parameters
        ----------
        from_filename : str
            Path to the source numpy file to convert.
        to_filename : str
            Path to the binary file that will be created.
        """
        np.load(from_filename).tofile(to_filename)

    @property
    def nbytesread(self):
        """Number of bytes per frame (fixed for a given file)."""
        return np.int64(2 * self.Ly * self.Lx)

    @property
    def nbytes(self):
        """Total number of bytes in the file."""
        return os.path.getsize(self.filename)

    @property
    def n_frames(self):
        """Total number of frames in the file."""
        return int(self.nbytes // self.nbytesread)

    @property
    def shape(self):
        """
        Return the dimensions of the data in the file.

        Returns
        -------
        n_frames : int
            Number of frames.
        Ly : int
            Height of each frame in pixels.
        Lx : int
            Width of each frame in pixels.
        """
        return self.n_frames, self.Ly, self.Lx

    @property
    def size(self):
        """
        Return the total number of pixels across all frames.

        Returns
        -------
        size : int
            Product of n_frames * Ly * Lx.
        """
        return np.prod(np.array(self.shape).astype(np.int64))

    def close(self):
        """
        Closes the file.
        """
        self.file._mmap.close()

    def __enter__(self):
        return self

    def __exit__(self, exc_type, exc_val, exc_tb):
        self.close()

    def __setitem__(self, *items):
        indices, data = items
        if data.dtype != "int16":
            self.file[indices] = np.minimum(data, 2**15 - 2).astype("int16")
        else:
            self.file[indices] = data

    def __getitem__(self, *items):
        indices, *crop = items
        return self.file[indices]

    def sampled_mean(self):
        """
        Compute the mean image from up to 1000 evenly spaced frames.

        Returns
        -------
        mean_img : numpy.ndarray
            Mean image of shape (Ly, Lx), averaged over the sampled frames.
        """
        n_frames = self.n_frames
        nsamps = min(n_frames, 1000)
        inds = np.linspace(0, n_frames, 1 + nsamps).astype(np.int64)[:-1]
        frames = self.file[inds].astype(np.float32)
        return frames.mean(axis=0)

    @property
    def data(self):
        """
        Return all frames in the file.

        Returns
        -------
        frames : numpy.ndarray
            All frame data as an array of shape (n_frames, Ly, Lx).
        """
        return self.file[:]

    def write_tiff(self, fname, range_dict={}):
        """
        Write binary file contents to a TIFF file, optionally cropped by frame/y/x ranges.

        Parameters
        ----------
        fname : str
            Output TIFF file path.
        range_dict : dict, optional
            Dictionary with optional keys "frame_range", "y_range", "x_range", each
            a tuple of (start, stop) indices. Defaults to the full extent of each
            dimension.
        """
        n_frames, Ly, Lx = self.shape
        frame_range, y_range, x_range = (0,n_frames), (0, Ly), (0, Lx)
        with TiffWriter(fname, bigtiff=True) as f:
            # Iterate through current data and write each frame to a tiff
            # All ranges should be Tuples(int,int)
            if 'frame_range' in range_dict:
                frame_range = range_dict['frame_range']
            if 'x_range' in range_dict:
                x_range = range_dict['x_range']
            if 'y_range' in range_dict:
                y_range = range_dict['y_range']
            logger.info('Frame Range: {}, y_range: {}, x_range{}'.format(frame_range, y_range, x_range))
            for i in range(frame_range[0], frame_range[1]):
                curr_frame = np.floor(self.file[i, y_range[0]:y_range[1], x_range[0]:x_range[1]]).astype(np.int16)
                f.write(curr_frame, contiguous=True)
        logger.info('Tiff has been saved to {}'.format(fname))

data property #

data

Return all frames in the file.

Returns:

Name Type Description
frames ndarray

All frame data as an array of shape (n_frames, Ly, Lx).

n_frames property #

n_frames

Total number of frames in the file.

nbytes property #

nbytes

Total number of bytes in the file.

nbytesread property #

nbytesread

Number of bytes per frame (fixed for a given file).

shape property #

shape

Return the dimensions of the data in the file.

Returns:

Name Type Description
n_frames int

Number of frames.
Ly int

Height of each frame in pixels.
Lx int

Width of each frame in pixels.

size property #

size

Return the total number of pixels across all frames.

Returns:

Name Type Description
size int

Product of n_frames * Ly * Lx.

__init__ #

__init__(Ly, Lx, filename, n_frames=None, dtype='int16', write=False)

Open or create a Suite2p binary file backed by a memory-mapped numpy array.

Parameters:

Name Type Description Default
Ly int

Height of each frame in pixels.

 required
Lx int

Width of each frame in pixels.

 required
filename str

Path to the binary file to read from or write to.

 required
n_frames int

Number of frames. Required when creating a new file for writing. Inferred from file size when reading.

 None
dtype str, optional (default "int16")

Data type of each pixel value.

 'int16'
write bool, optional (default False)

If True, open the file for reading and writing. If False, open read-only.

 False
Source code in suite2p/io/binary.py 
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
def __init__(self, Ly, Lx, filename, n_frames=None,
             dtype="int16", write=False):
    """
    Open or create a Suite2p binary file backed by a memory-mapped numpy array.

    Parameters
    ----------
    Ly : int
        Height of each frame in pixels.
    Lx : int
        Width of each frame in pixels.
    filename : str
        Path to the binary file to read from or write to.
    n_frames : int, optional
        Number of frames. Required when creating a new file for writing.
        Inferred from file size when reading.
    dtype : str, optional (default "int16")
        Data type of each pixel value.
    write : bool, optional (default False)
        If True, open the file for reading and writing. If False, open read-only.
    """
    self.Ly = Ly
    self.Lx = Lx
    self.filename = filename
    self.dtype = dtype
    self.write = write

    if write and n_frames is None and not os.path.exists(self.filename):
        raise ValueError(
            "need to provide number of frames n_frames when writing file")
    elif not write:
        n_frames = self.n_frames
    shape = (n_frames, self.Ly, self.Lx)
    if write:
        mode = "r+" if os.path.exists(self.filename) else "w+"
    else:
        mode = "r"
    self.file = np.memmap(self.filename, mode=mode, dtype=self.dtype, shape=shape)
    self._index = 0
    self._can_read = True

close #

close()

Closes the file.

Source code in suite2p/io/binary.py 
114
115
116
117
118
def close(self):
    """
    Closes the file.
    """
    self.file._mmap.close()

convert_numpy_file_to_suite2p_binary staticmethod #

convert_numpy_file_to_suite2p_binary(from_filename, to_filename)

Convert a numpy file (.npy, .npz) to a Suite2p binary file.

Parameters:

Name Type Description Default
from_filename str

Path to the source numpy file to convert.

 required
to_filename str

Path to the binary file that will be created.

 required
Source code in suite2p/io/binary.py 
57
58
59
60
61
62
63
64
65
66
67
68
69
@staticmethod
def convert_numpy_file_to_suite2p_binary(from_filename, to_filename):
    """
    Convert a numpy file (.npy, .npz) to a Suite2p binary file.

    Parameters
    ----------
    from_filename : str
        Path to the source numpy file to convert.
    to_filename : str
        Path to the binary file that will be created.
    """
    np.load(from_filename).tofile(to_filename)

sampled_mean #

sampled_mean()

Compute the mean image from up to 1000 evenly spaced frames.

Returns:

Name Type Description
mean_img ndarray

Mean image of shape (Ly, Lx), averaged over the sampled frames.
Source code in suite2p/io/binary.py 
137
138
139
140
141
142
143
144
145
146
147
148
149
150
def sampled_mean(self):
    """
    Compute the mean image from up to 1000 evenly spaced frames.

    Returns
    -------
    mean_img : numpy.ndarray
        Mean image of shape (Ly, Lx), averaged over the sampled frames.
    """
    n_frames = self.n_frames
    nsamps = min(n_frames, 1000)
    inds = np.linspace(0, n_frames, 1 + nsamps).astype(np.int64)[:-1]
    frames = self.file[inds].astype(np.float32)
    return frames.mean(axis=0)

write_tiff #

write_tiff(fname, range_dict={})

Write binary file contents to a TIFF file, optionally cropped by frame/y/x ranges.

Parameters:

Name Type Description Default
fname str

Output TIFF file path.

 required
range_dict dict

Dictionary with optional keys "frame_range", "y_range", "x_range", each a tuple of (start, stop) indices. Defaults to the full extent of each dimension.

 {}
Source code in suite2p/io/binary.py 
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
def write_tiff(self, fname, range_dict={}):
    """
    Write binary file contents to a TIFF file, optionally cropped by frame/y/x ranges.

    Parameters
    ----------
    fname : str
        Output TIFF file path.
    range_dict : dict, optional
        Dictionary with optional keys "frame_range", "y_range", "x_range", each
        a tuple of (start, stop) indices. Defaults to the full extent of each
        dimension.
    """
    n_frames, Ly, Lx = self.shape
    frame_range, y_range, x_range = (0,n_frames), (0, Ly), (0, Lx)
    with TiffWriter(fname, bigtiff=True) as f:
        # Iterate through current data and write each frame to a tiff
        # All ranges should be Tuples(int,int)
        if 'frame_range' in range_dict:
            frame_range = range_dict['frame_range']
        if 'x_range' in range_dict:
            x_range = range_dict['x_range']
        if 'y_range' in range_dict:
            y_range = range_dict['y_range']
        logger.info('Frame Range: {}, y_range: {}, x_range{}'.format(frame_range, y_range, x_range))
        for i in range(frame_range[0], frame_range[1]):
            curr_frame = np.floor(self.file[i, y_range[0]:y_range[1], x_range[0]:x_range[1]]).astype(np.int16)
            f.write(curr_frame, contiguous=True)
    logger.info('Tiff has been saved to {}'.format(fname))

BinaryFileCombined #

Source code in suite2p/io/binary.py 
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
class BinaryFileCombined:

    def __init__(self, LY, LX, Ly, Lx, dy, dx, read_filenames):
        """
        Open multiple Suite2p binary files for combined reading across ROIs/planes.

        Stitches multiple binary files into a single full-frame view by placing each
        ROI at its (dy, dx) offset within a canvas of size (LY, LX).

        Parameters
        ----------
        LY : int
            Height of the full combined frame in pixels.
        LX : int
            Width of the full combined frame in pixels.
        Ly : numpy.ndarray
            Array of per-ROI frame heights.
        Lx : numpy.ndarray
            Array of per-ROI frame widths.
        dy : numpy.ndarray
            Array of y-offsets for placing each ROI in the full frame.
        dx : numpy.ndarray
            Array of x-offsets for placing each ROI in the full frame.
        read_filenames : list of str
            Paths to the binary files to read, one per ROI.
        """
        self.LY = LY
        self.LX = LX
        self.Ly = Ly
        self.Lx = Lx
        self.dy = dy
        self.dx = dx
        self.read_filenames = read_filenames

        self.read_files = [
            BinaryFile(ly, lx, read_filename)
            for (ly, lx, read_filename) in zip(self.Ly, self.Lx, self.read_filenames)
        ]
        n_frames = np.zeros(len(self.read_files))
        for i, rf in enumerate(self.read_files):
            n_frames[i] = rf.n_frames
        assert (n_frames == n_frames[0]).sum() == len(self.read_files)
        self._index = 0
        self._can_read = True

    def __enter__(self):
        return self

    def __exit__(self, exc_type, exc_val, exc_tb):
        self.close()

    def close(self):
        """Close all underlying binary files."""
        for n in range(len(self.read_files)):
            self.read_files[n].close()

    @property
    def nbytes(self):
        """Total number of bytes per ROI file, as a numpy array."""
        nbytes = np.zeros(len(self.read_files), np.int64)
        for i, read_file in enumerate(self.read_files):
            nbytes[i] = read_file.nbytes
        return nbytes

    @property
    def n_frames(self):
        """Total number of frames (from the first ROI file)."""
        return self.read_files[0].n_frames

    def __getitem__(self, *items):
        indices, *crop = items
        data0 = self.read_files[0][indices]
        data_all = np.zeros((data0.shape[0], self.LY, self.LX), "int16")
        for n, read_file in enumerate(self.read_files):
            if n > 0:
                data0 = self.read_file[indices]
            data_all[:, self.dy[n]:self.dy[n] + self.Ly[n],
                     self.dx[n]:self.dx[n] + self.Lx[n]] = data0

        return data_all

n_frames property #

n_frames

Total number of frames (from the first ROI file).

nbytes property #

nbytes

Total number of bytes per ROI file, as a numpy array.

__init__ #

__init__(LY, LX, Ly, Lx, dy, dx, read_filenames)

Open multiple Suite2p binary files for combined reading across ROIs/planes.

Stitches multiple binary files into a single full-frame view by placing each ROI at its (dy, dx) offset within a canvas of size (LY, LX).

Parameters:

Name Type Description Default
LY int

Height of the full combined frame in pixels.

 required
LX int

Width of the full combined frame in pixels.

 required
Ly ndarray

Array of per-ROI frame heights.

 required
Lx ndarray

Array of per-ROI frame widths.

 required
dy ndarray

Array of y-offsets for placing each ROI in the full frame.

 required
dx ndarray

Array of x-offsets for placing each ROI in the full frame.

 required
read_filenames list of str

Paths to the binary files to read, one per ROI.

 required
Source code in suite2p/io/binary.py 
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
def __init__(self, LY, LX, Ly, Lx, dy, dx, read_filenames):
    """
    Open multiple Suite2p binary files for combined reading across ROIs/planes.

    Stitches multiple binary files into a single full-frame view by placing each
    ROI at its (dy, dx) offset within a canvas of size (LY, LX).

    Parameters
    ----------
    LY : int
        Height of the full combined frame in pixels.
    LX : int
        Width of the full combined frame in pixels.
    Ly : numpy.ndarray
        Array of per-ROI frame heights.
    Lx : numpy.ndarray
        Array of per-ROI frame widths.
    dy : numpy.ndarray
        Array of y-offsets for placing each ROI in the full frame.
    dx : numpy.ndarray
        Array of x-offsets for placing each ROI in the full frame.
    read_filenames : list of str
        Paths to the binary files to read, one per ROI.
    """
    self.LY = LY
    self.LX = LX
    self.Ly = Ly
    self.Lx = Lx
    self.dy = dy
    self.dx = dx
    self.read_filenames = read_filenames

    self.read_files = [
        BinaryFile(ly, lx, read_filename)
        for (ly, lx, read_filename) in zip(self.Ly, self.Lx, self.read_filenames)
    ]
    n_frames = np.zeros(len(self.read_files))
    for i, rf in enumerate(self.read_files):
        n_frames[i] = rf.n_frames
    assert (n_frames == n_frames[0]).sum() == len(self.read_files)
    self._index = 0
    self._can_read = True

close #

close()

Close all underlying binary files.

Source code in suite2p/io/binary.py 
261
262
263
264
def close(self):
    """Close all underlying binary files."""
    for n in range(len(self.read_files)):
        self.read_files[n].close()

temporary_pointer #

temporary_pointer(file)

Context manager that saves and restores a file's pointer position.

Parameters:

Name Type Description Default
file file object

An open file with tell() and seek() methods.

 required
Source code in suite2p/io/binary.py 
195
196
197
198
199
200
201
202
203
204
205
206
207
@contextmanager
def temporary_pointer(file):
    """
    Context manager that saves and restores a file's pointer position.

    Parameters
    ----------
    file : file object
        An open file with tell() and seek() methods.
    """
    orig_pointer = file.tell()
    yield file
    file.seek(orig_pointer)

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

dcimg_to_binary #

dcimg_to_binary(dbs, settings, reg_file, reg_file_chan2)

finds dcimg files and writes them to binaries

Parameters:

Name Type Description Default
settings

"nplanes", "data_path", "save_path", "save_folder", "fast_disk", "nchannels", "keep_movie_raw", "look_one_level_down"

 required

Returns:

Type Description
settings : dictionary of first plane

settings["reg_file"] or settings["raw_file"] is created binary assigns keys "Ly", "Lx", "tiffreader", "first_tiffs", "nframes", "meanImg", "meanImg_chan2"
Source code in suite2p/io/dcam.py 
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
def dcimg_to_binary(dbs, settings, reg_file, reg_file_chan2):
    """finds dcimg files and writes them to binaries

    Parameters
    ----------
    settings: dictionary
        "nplanes", "data_path", "save_path", "save_folder", "fast_disk",
        "nchannels", "keep_movie_raw", "look_one_level_down"

    Returns
    -------
        settings : dictionary of first plane
            settings["reg_file"] or settings["raw_file"] is created binary
            assigns keys "Ly", "Lx", "tiffreader", "first_tiffs",
            "nframes", "meanImg", "meanImg_chan2"
    """

    t0 = time.time()
    # # copy settings to list where each element is settings for each plane
    # settings1 = utils.init_settings(settings)

    # # open all binary files for writing
    # # look for dcimg in all requested folders
    # settings1, fs, reg_file, reg_file_chan2 = utils.find_files_open_binaries(settings1, False)
    # settings = settings1[0]
    fs = dbs[0]["file_list"]

    # loop over all dcimg files
    iall = 0
    ik = 0

    for file_name in fs:
        # open dcimg
        dcimg_file = dcimg.DCIMGFile(file_name)

        nplanes = 1
        nchannels = 1
        nframes = dcimg_file.shape[0]

        iblocks = np.arange(0, nframes, dbs[0]["batch_size"])
        if iblocks[-1] < nframes:
            iblocks = np.append(iblocks, nframes)

        if nchannels > 1:
            nfunc = dbs[0]["functional_chan"] - 1
        else:
            nfunc = 0

        # loop over all frames
        for ichunk, onset in enumerate(iblocks[:-1]):
            offset = iblocks[ichunk + 1]
            im_p = dcimg_file[onset:offset, :, :]
            im2mean = im_p.mean(axis=0).astype(np.float32) / len(iblocks)
            for ichan in range(nchannels):
                nframes = im_p.shape[0]
                im2write = im_p[:]
                Ly, Lx = im_p.shape[1], im_p.shape[2]
                for j in range(0, nplanes):
                    if iall == 0:
                        dbs[j]["meanImg"] = np.zeros((Ly, Lx), np.float32)
                        if nchannels > 1:
                            dbs[j]["meanImg_chan2"] = np.zeros((Ly, Lx), np.float32)
                        dbs[j]["nframes"] = 0
                        dbs[j]["Ly"], dbs[j]["Lx"] = Ly, Lx
                    if ichan == nfunc:
                        dbs[j]["meanImg"] += np.squeeze(im2mean)
                        reg_file[j].write(
                            bytearray(im2write[:].astype("uint16")))
                    else:
                        dbs[j]["meanImg_chan2"] += np.squeeze(im2mean)
                        reg_file_chan2[j].write(
                            bytearray(im2write[:].astype("uint16")))

                    dbs[j]["nframes"] += im2write.shape[0]


            ik += nframes
            iall += nframes

        dcimg_file.close()

        # write settings files
    # write settings files
    for db in dbs:
        db["meanImg"] /= db["nframes"]
        if nchannels > 1:
            db["meanImg_chan2"] /= db["nframes"]
        np.save(db["db_path"], db)
        np.save(db["settings_path"], settings)

    # close all binary files and write settings files
    for j in range(0, nplanes):
        reg_file[j].close()
        if nchannels > 1:
            reg_file_chan2[j].close()
    return dbs

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

h5py_to_binary #

h5py_to_binary(dbs, settings, reg_file, reg_file_chan2)

Read HDF5 files and write interleaved plane/channel data to binary files.

Iterates over all HDF5 files listed in dbs[0]["file_list"], de-interleaves planes and channels, and writes each plane's frames to the corresponding binary file. Supports 3D, 4D, and 5D HDF5 datasets (higher-dimensional data is flattened to frames x Ly x Lx before de-interleaving).

Parameters:

Name Type Description Default
dbs list of dict

Database dictionaries for each plane. Must contain keys "file_list", "nplanes", "nchannels", "batch_size", "h5py_key", and "functional_chan". Updated in-place with "Ly", "Lx", "nframes", "nframes_per_folder", "meanImg", and "meanImg_chan2".

 required
settings dict

Suite2p settings dictionary, saved alongside each plane's database.

 required
reg_file list of file objects

Opened binary files for writing each plane's functional channel data.

 required
reg_file_chan2 list of file objects

Opened binary files for writing each plane's second channel data (used only when nchannels > 1).

 required

Returns:

Name Type Description
dbs list of dict

Updated database dictionaries with image dimensions, frame counts, and mean images populated.
Source code in suite2p/io/h5.py 
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
def h5py_to_binary(dbs, settings, reg_file, reg_file_chan2):
    """
    Read HDF5 files and write interleaved plane/channel data to binary files.

    Iterates over all HDF5 files listed in `dbs[0]["file_list"]`, de-interleaves
    planes and channels, and writes each plane's frames to the corresponding binary
    file. Supports 3D, 4D, and 5D HDF5 datasets (higher-dimensional data is
    flattened to frames x Ly x Lx before de-interleaving).

    Parameters
    ----------
    dbs : list of dict
        Database dictionaries for each plane. Must contain keys "file_list",
        "nplanes", "nchannels", "batch_size", "h5py_key", and "functional_chan".
        Updated in-place with "Ly", "Lx", "nframes", "nframes_per_folder",
        "meanImg", and "meanImg_chan2".
    settings : dict
        Suite2p settings dictionary, saved alongside each plane's database.
    reg_file : list of file objects
        Opened binary files for writing each plane's functional channel data.
    reg_file_chan2 : list of file objects
        Opened binary files for writing each plane's second channel data
        (used only when nchannels > 1).

    Returns
    -------
    dbs : list of dict
        Updated database dictionaries with image dimensions, frame counts, and
        mean images populated.
    """
    if not HAS_H5PY:
        raise ImportError("h5py is required for this file type, please 'pip install h5py'")

    nplanes = dbs[0]["nplanes"]
    nchannels = dbs[0]["nchannels"]
    h5list = dbs[0]["file_list"]

    keys = dbs[0]["h5py_key"]
    if isinstance(keys, str):
        keys = [keys]
    iall = 0
    for j in range(nplanes):
        dbs[j]["nframes_per_folder"] = np.zeros(len(h5list), np.int32)

    for ih5, h5 in enumerate(h5list):
        with h5py.File(h5, "r") as f:
            # if h5py data is 5D or 4D instead of 3D, assume that
            # data = (nchan x) (nframes x) nplanes x pixels x pixels
            # 5D/4D data is flattened to process the same way as interleaved data
            for key in keys:
                hdims = f[key].ndim
                # keep track of the plane identity of the first frame (channel identity is assumed always 0)
                ncp = nplanes * nchannels
                nbatch = ncp * math.ceil(dbs[0]["batch_size"] / ncp)
                nframes_all = f[key].shape[
                    0] if hdims == 3 else f[key].shape[0] * f[key].shape[1]
                nbatch = min(nbatch, nframes_all)
                nfunc = dbs[0]["functional_chan"] - 1 if nchannels > 1 else 0
                # loop over all tiffs
                ik = 0
                while 1:
                    if hdims == 3:
                        irange = np.arange(ik, min(ik + nbatch, nframes_all), 1)
                        if irange.size == 0:
                            break
                        im = f[key][irange, :, :]
                    else:
                        irange = np.arange(
                            ik / ncp, min(ik / ncp + nbatch / ncp, nframes_all / ncp),
                            1).astype(int)
                        if irange.size == 0:
                            break
                        im = f[key][irange, ...]
                        if im.ndim == 5 and im.shape[0] == nchannels:
                            im = im.transpose((1, 0, 2, 3, 4))
                        # flatten to frames x pixels x pixels
                        im = np.reshape(im, (-1, im.shape[-2], im.shape[-1]))
                    nframes = im.shape[0]
                    if type(im[0, 0, 0]) == np.uint16:
                        im = im / 2
                    for j in range(0, nplanes):
                        if iall == 0:
                            dbs[j]["meanImg"] = np.zeros((im.shape[1], im.shape[2]),
                                                          np.float32)
                            if nchannels > 1:
                                dbs[j]["meanImg_chan2"] = np.zeros(
                                    (im.shape[1], im.shape[2]), np.float32)
                            dbs[j]["nframes"] = 0
                        i0 = nchannels * ((j) % nplanes)
                        im2write = im[np.arange(int(i0) +
                                                nfunc, nframes, ncp), :, :].astype(
                                                    np.int16)
                        reg_file[j].write(bytearray(im2write))
                        dbs[j]["meanImg"] += im2write.astype(np.float32).sum(axis=0)
                        if nchannels > 1:
                            im2write = im[np.arange(int(i0) + 1 -
                                                    nfunc, nframes, ncp), :, :].astype(
                                                        np.int16)
                            reg_file_chan2[j].write(bytearray(im2write))
                            dbs[j]["meanImg_chan2"] += im2write.astype(
                                np.float32).sum(axis=0)
                        dbs[j]["nframes"] += im2write.shape[0]
                        dbs[j]["nframes_per_folder"][ih5] += im2write.shape[0]
                    ik += nframes
                    iall += nframes

    # update dbs with image dimensions and mean images
    do_registration = settings["run"]["do_registration"]
    for db in dbs:
        db["Ly"] = im2write.shape[1]
        db["Lx"] = im2write.shape[2]
        if not do_registration:
            db["yrange"] = np.array([0, db["Ly"]])
            db["xrange"] = np.array([0, db["Lx"]])
        db["meanImg"] /= db["nframes"]
        if nchannels > 1:
            db["meanImg_chan2"] /= db["nframes"]
        # Save db and settings to each plane folder
        np.save(db["db_path"], db)
        np.save(db["settings_path"], settings)

    return dbs

VideoReader #

Uses cv2 to read video files

Source code in suite2p/io/movie.py 
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
class VideoReader:
    """ Uses cv2 to read video files """
    def __init__(self, filenames: list):
        """ Uses cv2 to open video files and obtain their details for reading

        Parameters
        ------------
        filenames : int
            list of video files
        """
        cumframes = [0]
        containers = []
        Ly = []
        Lx = []
        for f in filenames:  # for each video in the list
            cap = cv2.VideoCapture(f)
            containers.append(cap)
            Lx.append(int(cap.get(cv2.CAP_PROP_FRAME_WIDTH)))
            Ly.append(int(cap.get(cv2.CAP_PROP_FRAME_HEIGHT)))
            cumframes.append(cumframes[-1] + int(cap.get(cv2.CAP_PROP_FRAME_COUNT)))
        cumframes = np.array(cumframes).astype(int)
        Ly = np.array(Ly)
        Lx = np.array(Lx)
        if (Ly==Ly[0]).sum() < len(Ly) or (Lx==Lx[0]).sum() < len(Lx):
            raise ValueError("videos are not all the same size in y and x")
        else:
            Ly, Lx = Ly[0], Lx[0]

        self.filenames = filenames
        self.cumframes = cumframes 
        self.n_frames = cumframes[-1]
        self.Ly = Ly
        self.Lx = Lx
        self.containers = containers
        self.fs = containers[0].get(cv2.CAP_PROP_FPS)

    def close(self) -> None:
        """
        Closes the video files
        """
        for i in range(len(self.containers)):  # for each video in the list
            cap = self.containers[i]
            cap.release()

    def __enter__(self):
        return self

    def __exit__(self, exc_type, exc_val, exc_tb):
        self.close()

    @property
    def shape(self) -> Tuple[int, int, int]:
        """
        The dimensions of the data in the file

        Returns
        -------
        n_frames: int
            The number of frames
        Ly: int
            The height of each frame
        Lx: int
            The width of each frame
        """
        return self.n_frames, self.Ly, self.Lx

    def get_frames(self, cframes):
        """
        read frames "cframes" from videos

        Parameters
        ------------
        cframes : np.array
            start and stop of frames to read, or consecutive list of frames to read
        """
        cframes = np.maximum(0, np.minimum(self.n_frames - 1, cframes))
        cframes = np.arange(cframes[0], cframes[-1] + 1).astype(int)
        # find which video the frames exist in (ivids is length of cframes)
        ivids = (cframes[np.newaxis, :] >= self.cumframes[1:, np.newaxis]).sum(axis=0)
        nk = 0
        im = np.zeros((len(cframes), self.Ly, self.Lx), "uint8")
        for n in np.unique(ivids):  # for each video in cumframes
            cfr = cframes[ivids == n]
            start = cfr[0] - self.cumframes[n]
            end = cfr[-1] - self.cumframes[n] + 1
            nt0 = end - start
            capture = self.containers[n]
            if int(capture.get(cv2.CAP_PROP_POS_FRAMES)) != start:
                capture.set(cv2.CAP_PROP_POS_FRAMES, start)
            fc = 0
            ret = True
            while fc < nt0 and ret:
                ret, frame = capture.read()
                if ret:
                    im[nk + fc] = cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY)
                else:
                    logger.info("img load failed, replacing with prev..")
                    im[nk + fc] = im[nk + fc - 1]
                fc += 1
            nk += nt0
        return im

shape property #

shape

The dimensions of the data in the file

Returns:

Name Type Description
n_frames int

The number of frames
Ly int

The height of each frame
Lx int

The width of each frame

__init__ #

__init__(filenames)

Uses cv2 to open video files and obtain their details for reading

Parameters:

Name Type Description Default
filenames int

list of video files

 required
Source code in suite2p/io/movie.py 
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
def __init__(self, filenames: list):
    """ Uses cv2 to open video files and obtain their details for reading

    Parameters
    ------------
    filenames : int
        list of video files
    """
    cumframes = [0]
    containers = []
    Ly = []
    Lx = []
    for f in filenames:  # for each video in the list
        cap = cv2.VideoCapture(f)
        containers.append(cap)
        Lx.append(int(cap.get(cv2.CAP_PROP_FRAME_WIDTH)))
        Ly.append(int(cap.get(cv2.CAP_PROP_FRAME_HEIGHT)))
        cumframes.append(cumframes[-1] + int(cap.get(cv2.CAP_PROP_FRAME_COUNT)))
    cumframes = np.array(cumframes).astype(int)
    Ly = np.array(Ly)
    Lx = np.array(Lx)
    if (Ly==Ly[0]).sum() < len(Ly) or (Lx==Lx[0]).sum() < len(Lx):
        raise ValueError("videos are not all the same size in y and x")
    else:
        Ly, Lx = Ly[0], Lx[0]

    self.filenames = filenames
    self.cumframes = cumframes 
    self.n_frames = cumframes[-1]
    self.Ly = Ly
    self.Lx = Lx
    self.containers = containers
    self.fs = containers[0].get(cv2.CAP_PROP_FPS)

close #

close()

Closes the video files

Source code in suite2p/io/movie.py 
51
52
53
54
55
56
57
def close(self) -> None:
    """
    Closes the video files
    """
    for i in range(len(self.containers)):  # for each video in the list
        cap = self.containers[i]
        cap.release()

get_frames #

get_frames(cframes)

read frames "cframes" from videos

Parameters:

Name Type Description Default
cframes array

start and stop of frames to read, or consecutive list of frames to read

 required
Source code in suite2p/io/movie.py 
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
def get_frames(self, cframes):
    """
    read frames "cframes" from videos

    Parameters
    ------------
    cframes : np.array
        start and stop of frames to read, or consecutive list of frames to read
    """
    cframes = np.maximum(0, np.minimum(self.n_frames - 1, cframes))
    cframes = np.arange(cframes[0], cframes[-1] + 1).astype(int)
    # find which video the frames exist in (ivids is length of cframes)
    ivids = (cframes[np.newaxis, :] >= self.cumframes[1:, np.newaxis]).sum(axis=0)
    nk = 0
    im = np.zeros((len(cframes), self.Ly, self.Lx), "uint8")
    for n in np.unique(ivids):  # for each video in cumframes
        cfr = cframes[ivids == n]
        start = cfr[0] - self.cumframes[n]
        end = cfr[-1] - self.cumframes[n] + 1
        nt0 = end - start
        capture = self.containers[n]
        if int(capture.get(cv2.CAP_PROP_POS_FRAMES)) != start:
            capture.set(cv2.CAP_PROP_POS_FRAMES, start)
        fc = 0
        ret = True
        while fc < nt0 and ret:
            ret, frame = capture.read()
            if ret:
                im[nk + fc] = cv2.cvtColor(frame, cv2.COLOR_BGR2GRAY)
            else:
                logger.info("img load failed, replacing with prev..")
                im[nk + fc] = im[nk + fc - 1]
            fc += 1
        nk += nt0
    return im

movie_to_binary #

movie_to_binary(dbs, settings, reg_file, reg_file_chan2)

finds movie files and writes them to binaries

Parameters:

Name Type Description Default
dbs list of dict

Database dictionaries for each plane. Must contain keys "file_list", "nplanes", "nchannels", "batch_size", "h5py_key", and "functional_chan". Updated in-place with "Ly", "Lx", "nframes", "nframes_per_folder", "meanImg", and "meanImg_chan2".

 required
settings dict

Suite2p settings dictionary, saved alongside each plane's database.

 required
reg_file list of file objects

Opened binary files for writing each plane's functional channel data.

 required
reg_file_chan2 list of file objects

Opened binary files for writing each plane's second channel data (used only when nchannels > 1).

 required

Returns:

Name Type Description
dbs list of dict

Updated database dictionaries with image dimensions, frame counts, and mean images populated.
Source code in suite2p/io/movie.py 
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
def movie_to_binary(dbs, settings, reg_file, reg_file_chan2):
    """  finds movie files and writes them to binaries

    Parameters
    ----------
    dbs : list of dict
        Database dictionaries for each plane. Must contain keys "file_list",
        "nplanes", "nchannels", "batch_size", "h5py_key", and "functional_chan".
        Updated in-place with "Ly", "Lx", "nframes", "nframes_per_folder",
        "meanImg", and "meanImg_chan2".
    settings : dict
        Suite2p settings dictionary, saved alongside each plane's database.
    reg_file : list of file objects
        Opened binary files for writing each plane's functional channel data.
    reg_file_chan2 : list of file objects
        Opened binary files for writing each plane's second channel data
        (used only when nchannels > 1).

    Returns
    -------
    dbs : list of dict
        Updated database dictionaries with image dimensions, frame counts, and
        mean images populated.
    """
    if not HAS_CV2:
        raise ImportError("cv2 is required for this file type, please 'pip install opencv-python-headless'")


    nplanes = dbs[0]["nplanes"]
    nchannels = dbs[0]["nchannels"]
    flist = dbs[0]["file_list"]

    # todo: implement for multiple movies
    # for j in range(nplanes):
    #    dbs[j]["nframes_per_folder"] = np.zeros(len(flist), np.int32)

    ncp = nplanes * nchannels
    nbatch = ncp * int(np.ceil(dbs[0]["batch_size"] / ncp))
    logger.info(flist)
    t0 = time.time()
    with VideoReader(filenames=flist) as vr:
        if settings["fs"]<=0:
            settings["fs"] = vr.fs

        nframes_all = vr.cumframes[-1]
        nbatch = min(nbatch, nframes_all)
        nfunc = settings["functional_chan"] - 1 if nchannels > 1 else 0
        # loop over all video frames
        ik = 0
        while 1:
            irange = np.arange(ik, min(ik + nbatch, nframes_all), 1)
            if irange.size == 0:
                break
            im = vr.get_frames(irange).astype("int16")
            nframes = im.shape[0]
            for j in range(0, nplanes):
                if ik == 0:
                    dbs[j]["meanImg"] = np.zeros((im.shape[1], im.shape[2]),
                                                    np.float32)
                    if nchannels > 1:
                        dbs[j]["meanImg_chan2"] = np.zeros(
                            (im.shape[1], im.shape[2]), np.float32)
                    dbs[j]["nframes"] = 0
                i0 = nchannels * ((j) % nplanes)
                im2write = im[np.arange(int(i0) +
                                        nfunc, nframes, ncp), :, :].astype(
                                            np.int16)
                reg_file[j].write(bytearray(im2write))
                dbs[j]["meanImg"] += im2write.astype(np.float32).sum(axis=0)
                if nchannels > 1:
                    im2write = im[np.arange(int(i0) + 1 -
                                            nfunc, nframes, ncp), :, :].astype(
                                                np.int16)
                    reg_file_chan2[j].write(bytearray(im2write))
                    dbs[j]["meanImg_chan2"] += im2write.astype(
                        np.float32).sum(axis=0)
                dbs[j]["nframes"] += im2write.shape[0]

            ik += nframes
            if ik % (nbatch * 4) == 0:
                logger.info("%d frames of binary, time %0.2f sec." %
                      (ik, time.time() - t0))

    # write settings files
    # update dbs with image dimensions and mean images
    do_registration = settings["run"]["do_registration"]
    for db in dbs:
        db["Ly"] = im2write.shape[1]
        db["Lx"] = im2write.shape[2]
        if not do_registration:
            db["yrange"] = np.array([0, db["Ly"]])
            db["xrange"] = np.array([0, db["Lx"]])
        db["meanImg"] /= db["nframes"]
        if nchannels > 1:
            db["meanImg_chan2"] /= db["nframes"]
        # Save db and settings to each plane folder
        np.save(db["db_path"], db)
        np.save(db["settings_path"], settings)

    return dbs

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

nd2_to_binary #

nd2_to_binary(settings)

finds nd2 files and writes them to binaries

Parameters:

Name Type Description Default
settings

"nplanes", "data_path", "save_path", "save_folder", "fast_disk", "nchannels", "keep_movie_raw", "look_one_level_down"

 required

Returns:

Type Description
settings : dictionary of first plane

settings["reg_file"] or settings["raw_file"] is created binary assigns keys "Ly", "Lx", "tiffreader", "first_tiffs", "nframes", "meanImg", "meanImg_chan2"
Source code in suite2p/io/nd2.py 
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
def nd2_to_binary(settings):
    """finds nd2 files and writes them to binaries

    Parameters
    ----------
    settings: dictionary
        "nplanes", "data_path", "save_path", "save_folder", "fast_disk",
        "nchannels", "keep_movie_raw", "look_one_level_down"

    Returns
    -------
        settings : dictionary of first plane
            settings["reg_file"] or settings["raw_file"] is created binary
            assigns keys "Ly", "Lx", "tiffreader", "first_tiffs",
            "nframes", "meanImg", "meanImg_chan2"
    """

    t0 = time.time()
    # copy settings to list where each element is settings for each plane
    settings1 = utils.init_settings(settings)

    # open all binary files for writing
    # look for nd2s in all requested folders
    settings1, fs, reg_file, reg_file_chan2 = utils.find_files_open_binaries(settings1, False)
    settings = settings1[0]

    # loop over all nd2 files
    iall = 0
    ik = 0
    for file_name in fs:
        # open nd2
        nd2_file = nd2.ND2File(file_name)
        nd2_dims = {k: i for i, k in enumerate(nd2_file.sizes)}

        valid_dimensions = "TZCYX"
        assert set(nd2_dims) <= set(
            valid_dimensions
        ), f"Unknown dimensions {set(nd2_dims)-set(valid_dimensions)} in file {file_name}."

        # Sort the dimensions in the order of TZCYX, skipping the missing ones.
        im = nd2_file.asarray().transpose(
            [nd2_dims[x] for x in valid_dimensions if x in nd2_dims])

        # Expand array to include the missing dimensions.
        for i, dim in enumerate("TZC"):
            if dim not in nd2_dims:
                im = np.expand_dims(im, i)

        nplanes = nd2_file.sizes["Z"] if "Z" in nd2_file.sizes else 1
        nchannels = nd2_file.sizes["C"] if "C" in nd2_file.sizes else 1
        nframes = nd2_file.sizes["T"] if "T" in nd2_file.sizes else 1

        iblocks = np.arange(0, nframes, settings1[0]["batch_size"])
        if iblocks[-1] < nframes:
            iblocks = np.append(iblocks, nframes)

        if nchannels > 1:
            nfunc = settings1[0]["functional_chan"] - 1
        else:
            nfunc = 0

        assert im.max() < 32768 and im.min() >= -32768, "image data is out of range"
        im = im.astype(np.int16)

        # loop over all frames
        for ichunk, onset in enumerate(iblocks[:-1]):
            offset = iblocks[ichunk + 1]
            im_p = np.array(im[onset:offset, :, :, :, :])
            im2mean = im_p.mean(axis=0).astype(np.float32) / len(iblocks)
            for ichan in range(nchannels):
                nframes = im_p.shape[0]
                im2write = im_p[:, :, ichan, :, :]
                for j in range(0, nplanes):
                    if iall == 0:
                        settings1[j]["meanImg"] = np.zeros((im_p.shape[3], im_p.shape[4]),
                                                      np.float32)
                        if nchannels > 1:
                            settings1[j]["meanImg_chan2"] = np.zeros(
                                (im_p.shape[3], im_p.shape[4]), np.float32)
                        settings1[j]["nframes"] = 0
                    if ichan == nfunc:
                        settings1[j]["meanImg"] += np.squeeze(im2mean[j, ichan, :, :])
                        reg_file[j].write(
                            bytearray(im2write[:, j, :, :].astype("int16")))
                    else:
                        settings1[j]["meanImg_chan2"] += np.squeeze(im2mean[j, ichan, :, :])
                        reg_file_chan2[j].write(
                            bytearray(im2write[:, j, :, :].astype("int16")))

                    settings1[j]["nframes"] += im2write.shape[0]
            ik += nframes
            iall += nframes

        nd2_file.close()

    # write settings files
    do_registration = settings1[0]["do_registration"]
    for settings in settings1:
        settings["Ly"] = im.shape[3]
        settings["Lx"] = im.shape[4]
        if not do_registration:
            settings["yrange"] = np.array([0, settings["Ly"]])
            settings["xrange"] = np.array([0, settings["Lx"]])
        np.save(settings["settings_path"], settings)
    # close all binary files and write settings files
    for j in range(0, nplanes):
        reg_file[j].close()
        if nchannels > 1:
            reg_file_chan2[j].close()
    return settings1[0]

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

nwb_to_binary #

nwb_to_binary(settings)

convert nwb file to binary (experimental)

converts single plane single channel nwb file to binary for suite2p processing

Parameters:

Name Type Description Default
settings

requires "nwb_file" key optional keys "nwb_driver", "nwb_series" uses "nplanes", "save_path", "save_folder", "fast_disk", "nchannels", "keep_movie_raw", "look_one_level_down"

 required

Returns:

Type Description
settings : dictionary of first plane

settings["reg_file"] or settings["raw_file"] is created binary assigns keys "Ly", "Lx", "tiffreader", "first_tiffs", "frames_per_folder", "nframes", "meanImg", "meanImg_chan2"
Source code in suite2p/io/nwb.py 
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
def nwb_to_binary(settings):
    """convert nwb file to binary (experimental)

    converts single plane single channel nwb file to binary for suite2p processing

    Parameters
    ----------
    settings: dictionary
        requires "nwb_file" key
        optional keys "nwb_driver", "nwb_series"
        uses "nplanes", "save_path", "save_folder", "fast_disk",
        "nchannels", "keep_movie_raw", "look_one_level_down"

    Returns
    -------
        settings : dictionary of first plane
            settings["reg_file"] or settings["raw_file"] is created binary
            assigns keys "Ly", "Lx", "tiffreader", "first_tiffs",
            "frames_per_folder", "nframes", "meanImg", "meanImg_chan2"

    """

    # force 1 plane 1 chan for now
    settings["nplanes"] = 1
    settings["nchannels"] = 1

    # initialize settings with reg_file and raw_file paths, etc
    settings = utils.init_settings(settings)[0]

    t0 = time.time()
    nplanes = settings["nplanes"]
    nchannels = settings["nchannels"]

    batch_size = settings["batch_size"]
    batch_size = int(nplanes * nchannels * np.ceil(batch_size / (nplanes * nchannels)))

    # open reg_file (and when available reg_file_chan2)
    if "keep_movie_raw" in settings and settings["keep_movie_raw"]:
        reg_file = open(settings["raw_file"], "wb")
        if nchannels > 1:
            reg_file_chan2 = open(settings["raw_file_chan2"], "wb")
    else:
        reg_file = open(settings["reg_file"], "wb")
        if nchannels > 1:
            reg_file_chan2 = open(settings["reg_file_chan2"], "wb")

    nwb_driver = None
    if settings.get("nwb_driver") and isinstance(nwb_driver, str):
        nwb_driver = settings["nwb_driver"]

    with NWBHDF5IO(settings["nwb_file"], "r", driver=nwb_driver) as fio:
        nwbfile = fio.read()

        # get TwoPhotonSeries
        if not settings.get("nwb_series"):
            TwoPhotonSeries_names = []
            for v in nwbfile.acquisition.values():
                if isinstance(v, TwoPhotonSeries):
                    TwoPhotonSeries_names.append(v.name)
            if len(TwoPhotonSeries_names) == 0:
                raise ValueError("no TwoPhotonSeries in NWB file")
            elif len(TwoPhotonSeries_names) > 1:
                raise Warning(
                    "more than one TwoPhotonSeries in NWB file, choosing first one")
            settings["nwb_series"] = TwoPhotonSeries_names[0]

        series = nwbfile.acquisition[settings["nwb_series"]]
        series_shape = nwbfile.acquisition[settings["nwb_series"]].data.shape
        settings["nframes"] = series_shape[0]
        settings["frames_per_file"] = np.array([settings["nframes"]])
        settings["frames_per_folder"] = np.array([settings["nframes"]])
        settings["meanImg"] = np.zeros((series_shape[1], series_shape[2]), np.float32)
        for ik in np.arange(0, settings["nframes"], batch_size):
            ikend = min(ik + batch_size, settings["nframes"])
            im = series.data[ik:ikend]

            # check if uint16
            if im.dtype.type == np.uint16:
                im = (im // 2).astype(np.int16)
            elif im.dtype.type == np.int32:
                im = (im // 2).astype(np.int16)
            elif im.dtype.type != np.int16:
                im = im.astype(np.int16)

            reg_file.write(bytearray(im))
            settings["meanImg"] += im.astype(np.float32).sum(axis=0)

            if ikend % (batch_size * 4) == 0:
                logger.info("%d frames of binary, time %0.2f sec." %
                      (ikend, time.time() - t0))
        gc.collect()

    # write settings files
    settings["Ly"], settings["Lx"] = settings["meanImg"].shape
    settings["yrange"] = np.array([0, settings["Ly"]])
    settings["xrange"] = np.array([0, settings["Lx"]])
    settings["meanImg"] /= settings["nframes"]
    if nchannels > 1:
        settings["meanImg_chan2"] /= settings["nframes"]
    # close all binary files and write settings files
    np.save(settings["settings_path"], settings)
    reg_file.close()
    if nchannels > 1:
        reg_file_chan2.close()

    return settings

read_nwb #

read_nwb(fpath)

read NWB file for use in the GUI

Source code in suite2p/io/nwb.py 
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
def read_nwb(fpath):
    """read NWB file for use in the GUI"""
    with NWBHDF5IO(fpath, "r") as fio:
        nwbfile = fio.read()

        # ROIs
        try:
            rois = nwbfile.processing["ophys"]["ImageSegmentation"][
                "PlaneSegmentation"]["pixel_mask"]
            multiplane = False
        except Exception:
            rois = nwbfile.processing["ophys"]["ImageSegmentation"][
                "PlaneSegmentation"]["voxel_mask"]
            multiplane = True
        stat = []
        for n in range(len(rois)):
            if isinstance(rois[0], np.ndarray):
                stat.append({
                    "ypix":
                        np.array(
                            [rois[n][i][0].astype("int") for i in range(len(rois[n]))]),
                    "xpix":
                        np.array(
                            [rois[n][i][1].astype("int") for i in range(len(rois[n]))]),
                    "lam":
                        np.array([rois[n][i][-1] for i in range(len(rois[n]))]),
                })
            else:
                stat.append({
                    "ypix": rois[n]["x"].astype("int"),
                    "xpix": rois[n]["y"].astype("int"),
                    "lam": rois[n]["weight"],
                })
            if multiplane:
                stat[-1]["iplane"] = int(rois[n][0][-2])
        settings = default_settings()

        if multiplane:
            nplanes = (np.max(np.array([stat[n]["iplane"] for n in range(len(stat))])) +
                       1)
        else:
            nplanes = 1
        stat = np.array(stat)

        # settings with backgrounds
        settings1 = []
        for iplane in range(nplanes):
            settings = default_settings()
            bg_strs = ["meanImg", "Vcorr", "max_proj", "meanImg_chan2"]
            settings["nchannels"] = 1
            for bstr in bg_strs:
                if (bstr in nwbfile.processing["ophys"]["Backgrounds_%d" %
                                                        iplane].images):
                    settings[bstr] = np.array(nwbfile.processing["ophys"]["Backgrounds_%d" %
                                                                     iplane][bstr].data)
                    if bstr == "meanImg_chan2":
                        settings["nchannels"] = 2
            settings["Ly"], settings["Lx"] = settings[bg_strs[0]].shape
            settings["yrange"] = [0, settings["Ly"]]
            settings["xrange"] = [0, settings["Lx"]]
            settings["tau"] = 1.0
            settings["fs"] = nwbfile.acquisition["TwoPhotonSeries"].rate
            settings1.append(settings.copy())

        # fluorescence
        ophys = nwbfile.processing["ophys"]

        def get_fluo(name: str) -> np.ndarray:
            """Extract Fluorescence data."""
            roi_response_series = ophys[name].roi_response_series
            if name in roi_response_series.keys():
                fluo = ophys[name][name].data[:]
            elif "plane0" in roi_response_series.keys():
                for key, value in roi_response_series.items():
                    if key == "plane0":
                        fluo = value.data[:]
                    else:
                        fluo = np.concatenate((fluo, value.data[:]), axis=1)
                fluo = np.transpose(fluo)
            else:
                raise AttributeError(f"Can't find {name} container in {fpath}")
            return fluo

        F = get_fluo("Fluorescence")
        Fneu = get_fluo("Neuropil")
        spks = get_fluo("Deconvolved")

        # cell probabilities
        iscell = [
            ophys["ImageSegmentation"]["PlaneSegmentation"]["iscell"][n]
            for n in range(len(stat))
        ]
        iscell = np.array(iscell)
        probcell = iscell[:, 1]
        iscell_bool = iscell[:, 0].astype("bool")
        # Create redcell as 2-column array for consistency with iscell format
        redcell = np.zeros_like(iscell)
        probredcell = redcell[:, 1]

        if multiplane:
            settings = settings1[0].copy()
            Lx = settings["Lx"]
            Ly = settings["Ly"]
            nX = np.ceil(np.sqrt(settings["Ly"] * settings["Lx"] * len(settings1)) / settings["Lx"])
            nX = int(nX)
            for j in range(len(settings1)):
                settings1[j]["dx"] = (j % nX) * Lx
                settings1[j]["dy"] = int(j / nX) * Ly

            LY = int(np.amax(np.array([settings["Ly"] + settings["dy"] for settings in settings1])))
            LX = int(np.amax(np.array([settings["Lx"] + settings["dx"] for settings in settings1])))
            meanImg = np.zeros((LY, LX))
            max_proj = np.zeros((LY, LX))
            if settings["nchannels"] > 1:
                meanImg_chan2 = np.zeros((LY, LX))

            Vcorr = np.zeros((LY, LX))
            for k, settings in enumerate(settings1):
                xrange = np.arange(settings["dx"], settings["dx"] + settings["Lx"])
                yrange = np.arange(settings["dy"], settings["dy"] + settings["Ly"])
                meanImg[np.ix_(yrange, xrange)] = settings["meanImg"]
                Vcorr[np.ix_(yrange, xrange)] = settings["Vcorr"]
                max_proj[np.ix_(yrange, xrange)] = settings["max_proj"]
                if settings["nchannels"] > 1:
                    if "meanImg_chan2" in settings:
                        meanImg_chan2[np.ix_(yrange, xrange)] = settings["meanImg_chan2"]
                for j in np.nonzero(
                        np.array([stat[n]["iplane"] == k for n in range(len(stat))
                                 ]))[0]:
                    stat[j]["xpix"] += settings["dx"]
                    stat[j]["ypix"] += settings["dy"]
            settings["Vcorr"] = Vcorr
            settings["max_proj"] = max_proj
            settings["meanImg"] = meanImg
            if "meanImg_chan2" in settings:
                settings["meanImg_chan2"] = meanImg_chan2
            settings["Ly"], settings["Lx"] = LY, LX
            settings["yrange"] = [0, LY]
            settings["xrange"] = [0, LX]

        # Compute roi_stats after multiplane coordinates have been adjusted
        stat = roi_stats(stat, settings["Ly"], settings["Lx"], diameter=settings["diameter"])

        # Compute skew after roi_stats (which may filter ROIs)
        dF = F - settings["extraction"]["neuropil_coefficient"] * Fneu
        for n in range(len(stat)):
            stat[n]["skew"] = scipy.stats.skew(dF[n])

    return stat, settings, F, Fneu, spks, iscell, probcell, redcell, probredcell

save_nwb #

save_nwb(save_folder)

convert folder with plane folders to NWB format

Source code in suite2p/io/nwb.py 
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
def save_nwb(save_folder):
    """convert folder with plane folders to NWB format"""

    plane_folders = natsorted([
        Path(f.path)
        for f in os.scandir(save_folder)
        if f.is_dir() and f.name[:5] == "plane"
    ])
    settings1 = [
        np.load(f.joinpath("settings.npy"), allow_pickle=True).item() for f in plane_folders
    ]
    dbs = [
        _load_npy_cross_platform(f.joinpath("db.npy")).item() for f in plane_folders
    ]

    # Load reg_outputs and detect_outputs for background images
    for i, f in enumerate(plane_folders):
        # Merge reg_outputs (contains meanImg, yrange, xrange, etc.)
        reg_path = f.joinpath("reg_outputs.npy")
        if reg_path.exists():
            reg_outputs = np.load(reg_path, allow_pickle=True).item()
            settings1[i] = {**settings1[i], **reg_outputs}

        # Merge detect_outputs (contains Vcorr, max_proj, etc.)
        detect_path = f.joinpath("detect_outputs.npy")
        if detect_path.exists():
            detect_outputs = np.load(detect_path, allow_pickle=True).item()
            settings1[i] = {**settings1[i], **detect_outputs}

        # Add db keys that might be needed (Ly, Lx)
        for key in ["Ly", "Lx"]:
            if key in dbs[i] and key not in settings1[i]:
                settings1[i][key] = dbs[i][key]

    # Get nchannels from the main db or from plane dbs
    nchannels = dbs[0].get("nchannels", 1)

    if NWB and (settings1[0].get("mesoscan") is None):
        if len(settings1) > 1:
            multiplane = True
        else:
            multiplane = False

        settings = settings1[0]
        if "date_proc" in settings:
            session_start_time = settings["date_proc"]
            if not session_start_time.tzinfo:
                session_start_time = session_start_time.astimezone()
        else:
            session_start_time = datetime.datetime.now().astimezone()

        # INITIALIZE NWB FILE
        nwbfile = NWBFile(
            session_description="suite2p_proc",
            identifier=str(dbs[0]["data_path"][0]),
            session_start_time=session_start_time,
        )
        logger.info(nwbfile)

        device = nwbfile.create_device(
            name="Microscope",
            description="My two-photon microscope",
            manufacturer="The best microscope manufacturer",
        )
        optical_channel = OpticalChannel(
            name="OpticalChannel",
            description="an optical channel",
            emission_lambda=500.0,
        )

        imaging_plane = nwbfile.create_imaging_plane(
            name="ImagingPlane",
            optical_channel=optical_channel,
            imaging_rate=settings["fs"],
            description="standard",
            device=device,
            excitation_lambda=600.0,
            indicator="GCaMP",
            location="V1",
            grid_spacing=([2.0, 2.0, 30.0] if multiplane else [2.0, 2.0]),
            grid_spacing_unit="microns",
        )
        # link to external data
        external_data = settings["filelist"] if "filelist" in settings else [""]
        image_series = TwoPhotonSeries(
            name="TwoPhotonSeries",
            dimension=[dbs[0]["Ly"], dbs[0]["Lx"]],
            external_file=external_data,
            imaging_plane=imaging_plane,
            starting_frame=[0 for i in range(len(external_data))],
            format="external",
            starting_time=0.0,
            rate=settings["fs"] * dbs[0]["nplanes"],
        )
        nwbfile.add_acquisition(image_series)

        # processing
        img_seg = ImageSegmentation()
        ps = img_seg.create_plane_segmentation(
            name="PlaneSegmentation",
            description="suite2p output",
            imaging_plane=imaging_plane,
            reference_images=image_series,
        )
        ophys_module = nwbfile.create_processing_module(
            name="ophys", description="optical physiology processed data")
        ophys_module.add(img_seg)

        file_strs = ["F.npy", "Fneu.npy", "spks.npy"]
        file_strs_chan2 = ["F_chan2.npy", "Fneu_chan2.npy"]
        traces, traces_chan2 = [], []
        ncells = np.zeros(len(settings1), dtype=np.int_)
        Nfr = np.array([db["nframes"] for db in dbs]).max()
        for iplane, (settings, db) in enumerate(zip(settings1, dbs)):
            if iplane == 0:
                iscell = np.load(os.path.join(plane_folders[iplane], "iscell.npy"))
                for fstr in file_strs:
                    traces.append(np.load(os.path.join(plane_folders[iplane], fstr)))
                if nchannels > 1:
                    for fstr in file_strs_chan2:
                        traces_chan2.append(
                            np.load(plane_folders[iplane].joinpath(fstr)))
                PlaneCellsIdx = iplane * np.ones(len(iscell))
            else:
                iscell = np.append(
                    iscell,
                    np.load(os.path.join(plane_folders[iplane], "iscell.npy")),
                    axis=0,
                )
                for i, fstr in enumerate(file_strs):
                    trace = np.load(os.path.join(plane_folders[iplane], fstr))
                    if trace.shape[1] < Nfr:
                        fcat = np.zeros((trace.shape[0], Nfr - trace.shape[1]),
                                        "float32")
                        trace = np.concatenate((trace, fcat), axis=1)
                    traces[i] = np.append(traces[i], trace, axis=0)
                if nchannels > 1:
                    for i, fstr in enumerate(file_strs_chan2):
                        traces_chan2[i] = np.append(
                            traces_chan2[i],
                            np.load(plane_folders[iplane].joinpath(fstr)),
                            axis=0,
                        )
                PlaneCellsIdx = np.append(
                    PlaneCellsIdx, iplane * np.ones(len(iscell) - len(PlaneCellsIdx)))

            stat = np.load(os.path.join(plane_folders[iplane], "stat.npy"),
                           allow_pickle=True)
            ncells[iplane] = len(stat)
            for n in range(ncells[iplane]):
                if multiplane:
                    pixel_mask = np.array([
                        stat[n]["ypix"],
                        stat[n]["xpix"],
                        iplane * np.ones(stat[n]["npix"]),
                        stat[n]["lam"],
                    ])
                    ps.add_roi(voxel_mask=pixel_mask.T)
                else:
                    pixel_mask = np.array(
                        [stat[n]["ypix"], stat[n]["xpix"], stat[n]["lam"]])
                    ps.add_roi(pixel_mask=pixel_mask.T)

        ps.add_column("iscell", "two columns - iscell & probcell", iscell)

        rt_region = []
        for iplane, settings in enumerate(settings1):
            if iplane == 0:
                rt_region.append(
                    ps.create_roi_table_region(
                        region=list(np.arange(0, ncells[iplane]),),
                        description=f"ROIs for plane{int(iplane)}",
                    ))
            else:
                rt_region.append(
                    ps.create_roi_table_region(
                        region=list(
                            np.arange(
                                np.sum(ncells[:iplane]),
                                ncells[iplane] + np.sum(ncells[:iplane]),
                            )),
                        description=f"ROIs for plane{int(iplane)}",
                    ))

        # FLUORESCENCE (all are required)
        name_strs = ["Fluorescence", "Neuropil", "Deconvolved"]
        name_strs_chan2 = ["Fluorescence_chan2", "Neuropil_chan2"]

        for i, (fstr, nstr) in enumerate(zip(file_strs, name_strs)):
            for iplane, settings in enumerate(settings1):
                roi_resp_series = RoiResponseSeries(
                    name=f"plane{int(iplane)}",
                    data=np.transpose(traces[i][PlaneCellsIdx == iplane]),
                    rois=rt_region[iplane],
                    unit="lumens",
                    rate=settings["fs"],
                )
                if iplane == 0:
                    fl = Fluorescence(roi_response_series=roi_resp_series, name=nstr)
                else:
                    fl.add_roi_response_series(roi_response_series=roi_resp_series)
            ophys_module.add(fl)

        if nchannels > 1:
            for i, (fstr, nstr) in enumerate(zip(file_strs_chan2, name_strs_chan2)):
                for iplane, settings in enumerate(settings1):
                    roi_resp_series = RoiResponseSeries(
                        name=f"plane{int(iplane)}",
                        data=np.transpose(traces_chan2[i][PlaneCellsIdx == iplane]),
                        rois=rt_region[iplane],
                        unit="lumens",
                        rate=settings["fs"],
                    )

                    if iplane == 0:
                        fl = Fluorescence(roi_response_series=roi_resp_series,
                                          name=nstr)
                    else:
                        fl.add_roi_response_series(roi_response_series=roi_resp_series)

                ophys_module.add(fl)

        # BACKGROUNDS
        # (meanImg, Vcorr and max_proj are REQUIRED)
        bg_strs = ["meanImg", "Vcorr", "max_proj", "meanImg_chan2"]
        for iplane, settings in enumerate(settings1):
            images = Images("Backgrounds_%d" % iplane)
            for bstr in bg_strs:
                if bstr in settings:
                    if bstr == "Vcorr" or bstr == "max_proj":
                        img = np.zeros((settings["Ly"], settings["Lx"]), np.float32)
                        img[
                            settings["yrange"][0]:settings["yrange"][-1],
                            settings["xrange"][0]:settings["xrange"][-1],
                        ] = settings[bstr]
                    else:
                        img = settings[bstr]
                    images.add_image(GrayscaleImage(name=bstr, data=img))

            ophys_module.add(images)

        with NWBHDF5IO(os.path.join(save_folder, "ophys.nwb"), "w") as fio:
            fio.write(nwbfile)
    else:
        logger.info("pip install pynwb OR don't use mesoscope recording")

Copyright © 2023 Yoav Livneh Lab, Authored by Yael Prilutski.

raw_to_binary #

raw_to_binary(ops, use_recorded_defaults=True)

Finds RAW files and writes them to binaries

Parameters:

Name Type Description Default
ops dictionary

"data_path"

 required
use_recorded_defaults bool

Recorded session parameters are used when 'True', otherwise |ops| is expected to contain the following (additional) keys: "nplanes", "nchannels", "fs"

 True

Returns:

Type Description
ops : dictionary of first plane
Source code in suite2p/io/raw.py 
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
def raw_to_binary(ops, use_recorded_defaults=True):

    """ Finds RAW files and writes them to binaries

    Parameters
    ----------
    ops : dictionary
        "data_path"

    use_recorded_defaults : bool
        Recorded session parameters are used when 'True',
        otherwise |ops| is expected to contain the following (additional) keys:
          "nplanes",
          "nchannels",
          "fs"

    Returns
    -------
        ops : dictionary of first plane

    """

    if not HAS_XML:
        raise ImportError("xmltodict is required for RAW file support (pip install xmltodict)")

    # Load raw file configurations
    raw_file_configurations = [_RawFile(path) for path in ops['data_path']]

    # Split ops by captured planes
    ops_paths = _initialize_destination_files(ops, raw_file_configurations, use_recorded_defaults=use_recorded_defaults)

    # Convert all runs in order
    for path in ops['data_path']:
        print(f'Converting raw to binary: `{path}`')
        ops_loaded = [np.load(i, allow_pickle=True)[()] for i in ops_paths]
        _raw2bin(ops_loaded, _RawFile(path))

    # Reload edited ops
    ops_loaded = [np.load(i, allow_pickle=True)[()] for i in ops_paths]

    # Create a mean image with the final number of frames
    _update_mean(ops_loaded)

    # Load & return all ops
    return ops_loaded[0]

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

combined #

combined(save_folder, save=True)

Combine all plane folders in save_folder into a single result file.

Loads per-plane results (stat, F, Fneu, spks, iscell, redcell), shifts ROI coordinates by the tiled offsets, and concatenates them into combined arrays. Multi-plane recordings are arranged to best tile a square. Multi-ROI recordings are arranged by their dx, dy physical localization.

Parameters:

Name Type Description Default
save_folder str

Path to the suite2p output folder containing plane subdirectories (e.g. "plane0", "plane1", ...).

 required
save bool, optional (default True)

If True, save combined results (F.npy, Fneu.npy, spks.npy, stat.npy, db.npy, settings.npy, and optionally Fall.mat) to a "combined" subfolder. If False, only iscell.npy (and redcell.npy) are saved.

 True

Returns:

Name Type Description
stat ndarray

Concatenated ROI statistics across all planes.
db dict

Combined database dictionary with merged mean images and full-frame dimensions.
settings dict

Suite2p settings dictionary.
F ndarray

Combined fluorescence traces of shape (n_cells_total, n_frames).
Fneu ndarray

Combined neuropil traces of shape (n_cells_total, n_frames).
spks ndarray

Combined deconvolved spike traces of shape (n_cells_total, n_frames).
iscell0 ndarray

Binary cell classification labels for each ROI.
iscell1 ndarray

Cell classification probabilities for each ROI.
redcell0 ndarray

Binary red-cell labels for each ROI.
redcell1 ndarray

Red-cell probabilities for each ROI.
hasred bool

Whether red-cell classification data was found.
Source code in suite2p/io/save.py 
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
def combined(save_folder, save=True):
    """
    Combine all plane folders in save_folder into a single result file.

    Loads per-plane results (stat, F, Fneu, spks, iscell, redcell), shifts ROI
    coordinates by the tiled offsets, and concatenates them into combined arrays.
    Multi-plane recordings are arranged to best tile a square. Multi-ROI
    recordings are arranged by their dx, dy physical localization.

    Parameters
    ----------
    save_folder : str
        Path to the suite2p output folder containing plane subdirectories
        (e.g. "plane0", "plane1", ...).
    save : bool, optional (default True)
        If True, save combined results (F.npy, Fneu.npy, spks.npy, stat.npy,
        db.npy, settings.npy, and optionally Fall.mat) to a "combined"
        subfolder. If False, only iscell.npy (and redcell.npy) are saved.

    Returns
    -------
    stat : numpy.ndarray
        Concatenated ROI statistics across all planes.
    db : dict
        Combined database dictionary with merged mean images and full-frame
        dimensions.
    settings : dict
        Suite2p settings dictionary.
    F : numpy.ndarray
        Combined fluorescence traces of shape (n_cells_total, n_frames).
    Fneu : numpy.ndarray
        Combined neuropil traces of shape (n_cells_total, n_frames).
    spks : numpy.ndarray
        Combined deconvolved spike traces of shape (n_cells_total, n_frames).
    iscell0 : numpy.ndarray
        Binary cell classification labels for each ROI.
    iscell1 : numpy.ndarray
        Cell classification probabilities for each ROI.
    redcell0 : numpy.ndarray
        Binary red-cell labels for each ROI.
    redcell1 : numpy.ndarray
        Red-cell probabilities for each ROI.
    hasred : bool
        Whether red-cell classification data was found.
    """
    plane_folders = natsorted([
        f.path for f in os.scandir(save_folder) if f.is_dir() and f.name[:5] == "plane"
    ])
    top_db_path = os.path.join(save_folder, "db.npy")
    if os.path.exists(top_db_path):
        top_db = np.load(top_db_path, allow_pickle=True).item()
        ignore_flyback = set(top_db.get("ignore_flyback") or [])
        if ignore_flyback:
            plane_folders = [f for i, f in enumerate(plane_folders) if i not in ignore_flyback]
    dbs = [
        np.load(os.path.join(f, "db.npy"), allow_pickle=True).item()
        for f in plane_folders
    ]
    settings = np.load(os.path.join(plane_folders[0], "settings.npy"), allow_pickle=True).item()
    if os.path.exists(os.path.join(plane_folders[0], "reg_outputs.npy")):
        dops =  [
        np.load(os.path.join(f, "reg_outputs.npy"), allow_pickle=True).item()
        for f in plane_folders
        ]
        dbs = [{**db, **dop} for db, dop in zip(dbs, dops)]

    if os.path.exists(os.path.join(plane_folders[0], "detect_outputs.npy")):
        dops =  [
        np.load(os.path.join(f, "detect_outputs.npy"), allow_pickle=True).item()
        for f in plane_folders
        ]
        dbs = [{**db, **dop} for db, dop in zip(dbs, dops)]


    dy, dx = compute_dydx(dbs)
    Ly = np.array([db["Ly"] for db in dbs])
    Lx = np.array([db["Lx"] for db in dbs])
    LY = int(np.amax(dy + Ly))
    LX = int(np.amax(dx + Lx))
    meanImg = np.zeros((LY, LX))
    meanImgE = np.zeros((LY, LX))
    logger.info(plane_folders)
    if dbs[0]["nchannels"] > 1:
        meanImg_chan2 = np.zeros((LY, LX))
    if any(["meanImg_chan2_corrected" in db for db in dbs]):
        meanImg_chan2_corrected = np.zeros((LY, LX))
    if any(["max_proj" in db for db in dbs]):
        max_proj = np.zeros((LY, LX))

    Vcorr = np.zeros((LY, LX))
    Nfr = np.amax(np.array([db["nframes"] for db in dbs]))
    ii = 0
    for k, db in enumerate(dbs):
        fpath = plane_folders[k]
        if not os.path.exists(os.path.join(fpath, "stat.npy")):
            continue
        stat0 = np.load(os.path.join(fpath, "stat.npy"), allow_pickle=True)
        xrange = np.arange(dx[k], dx[k] + Lx[k])
        yrange = np.arange(dy[k], dy[k] + Ly[k])
        meanImg[np.ix_(yrange, xrange)] = db["meanImg"]
        if "meanImgE" in db:
            meanImgE[np.ix_(yrange, xrange)] = db["meanImgE"]
        if db["nchannels"] > 1:
            if "meanImg_chan2" in db:
                meanImg_chan2[np.ix_(yrange, xrange)] = db["meanImg_chan2"]
        if "meanImg_chan2_corrected" in db:
            meanImg_chan2_corrected[np.ix_(yrange,
                                           xrange)] = db["meanImg_chan2_corrected"]

        xrange = np.arange(dx[k] + db["xrange"][0], dx[k] + db["xrange"][-1])
        yrange = np.arange(dy[k] + db["yrange"][0], dy[k] + db["yrange"][-1])
        Vcorr[np.ix_(yrange, xrange)] = db["Vcorr"]
        if "max_proj" in db:
            max_proj[np.ix_(yrange, xrange)] = db["max_proj"]
        for j in range(len(stat0)):
            stat0[j]["xpix"] += dx[k]
            stat0[j]["ypix"] += dy[k]
            stat0[j]["med"][0] += dy[k]
            stat0[j]["med"][1] += dx[k]
            stat0[j]["iplane"] = k
        F0 = np.load(os.path.join(fpath, "F.npy"))
        Fneu0 = np.load(os.path.join(fpath, "Fneu.npy"))
        spks0 = np.load(os.path.join(fpath, "spks.npy"))
        iscell0 = np.load(os.path.join(fpath, "iscell.npy"))
        if os.path.isfile(os.path.join(fpath, "redcell.npy")):
            redcell0 = np.load(os.path.join(fpath, "redcell.npy"))
            hasred = True
        else:
            redcell0 = []
            hasred = False
        nn, nt = F0.shape
        if nt < Nfr:
            fcat = np.zeros((nn, Nfr - nt), "float32")
            #logger.info(F0.shape)
            #logger.info(fcat.shape)
            F0 = np.concatenate((F0, fcat), axis=1)
            spks0 = np.concatenate((spks0, fcat), axis=1)
            Fneu0 = np.concatenate((Fneu0, fcat), axis=1)
        if ii == 0:
            F, Fneu, spks, stat, iscell, redcell = F0, Fneu0, spks0, stat0, iscell0, redcell0
        else:
            F = np.concatenate((F, F0))
            Fneu = np.concatenate((Fneu, Fneu0))
            spks = np.concatenate((spks, spks0))
            stat = np.concatenate((stat, stat0))
            iscell = np.concatenate((iscell, iscell0))
            if hasred:
                redcell = np.concatenate((redcell, redcell0))
        ii += 1
        logger.info("appended plane %d to combined view" % k)

    db["meanImg"] = meanImg
    db["meanImgE"] = meanImgE
    if db["nchannels"] > 1:
        db["meanImg_chan2"] = meanImg_chan2
    if "meanImg_chan2_corrected" in db:
        db["meanImg_chan2_corrected"] = meanImg_chan2_corrected
    if "max_proj" in db:
        db["max_proj"] = max_proj
    db["Vcorr"] = Vcorr
    db["Ly"] = LY
    db["Lx"] = LX
    db["xrange"] = [0, db["Lx"]]
    db["yrange"] = [0, db["Ly"]]

    fpath = os.path.join(save_folder, "combined")

    if not os.path.isdir(fpath):
        os.makedirs(fpath)

    db["save_path"] = fpath

    # need to save iscell regardless (required for GUI function)
    np.save(os.path.join(fpath, "iscell.npy"), iscell)
    if hasred:
        np.save(os.path.join(fpath, "redcell.npy"), redcell)
    else:
        redcell = np.zeros_like(iscell)

    if save:
        np.save(os.path.join(fpath, "F.npy"), F)
        np.save(os.path.join(fpath, "Fneu.npy"), Fneu)
        np.save(os.path.join(fpath, "spks.npy"), spks)
        np.save(os.path.join(fpath, "db.npy"), db)
        np.save(os.path.join(fpath, "stat.npy"), stat)
        np.save(os.path.join(fpath, "settings.npy"), settings)

        # save as matlab file
        if settings["io"]["save_mat"]:
            matpath = os.path.join(db["save_path"], "Fall.mat")
            save_mat(db, stat, F, Fneu, spks, iscell, redcell)

    return (stat, db, settings, F, Fneu, spks, 
            iscell[:,0], iscell[:,1], 
            redcell[:,0], redcell[:,1], hasred)

compute_dydx #

compute_dydx(db1)

Compute pixel offsets (dy, dx) for tiling multiple planes/ROIs into a combined view.

If the databases do not contain "dx" and "dy" fields, arranges planes in a grid that best tiles a square. If offsets are present (e.g. mesoscope ROIs), uses the physical offsets and further tiles across planes.

Parameters:

Name Type Description Default
db1 list of dict

List of per-plane database dictionaries. Each must contain "Ly" and "Lx". May contain "dx" and "dy" for physical ROI offsets.

 required

Returns:

Name Type Description
dy ndarray

Y-offsets (in pixels) for each plane, shape (len(db1),).
dx ndarray

X-offsets (in pixels) for each plane, shape (len(db1),).
Source code in suite2p/io/save.py 
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
def compute_dydx(db1):
    """
    Compute pixel offsets (dy, dx) for tiling multiple planes/ROIs into a combined view.

    If the databases do not contain "dx" and "dy" fields, arranges planes in a
    grid that best tiles a square. If offsets are present (e.g. mesoscope ROIs),
    uses the physical offsets and further tiles across planes.

    Parameters
    ----------
    db1 : list of dict
        List of per-plane database dictionaries. Each must contain "Ly" and "Lx".
        May contain "dx" and "dy" for physical ROI offsets.

    Returns
    -------
    dy : numpy.ndarray
        Y-offsets (in pixels) for each plane, shape (len(db1),).
    dx : numpy.ndarray
        X-offsets (in pixels) for each plane, shape (len(db1),).
    """
    db = db1[0].copy()
    dx = np.zeros(len(db1), np.int64)
    dy = np.zeros(len(db1), np.int64)
    if ("dx" not in db) or ("dy" not in db) or (db["dx"] is None) or (db["dy"] is None):
        Lx = db["Lx"]
        Ly = db["Ly"]
        nX = np.ceil(np.sqrt(db["Ly"] * db["Lx"] * len(db1)) / db["Lx"])
        nX = int(nX)
        for j in range(len(db1)):
            dx[j] = (j % nX) * Lx
            dy[j] = int(j / nX) * Ly
    else:
        dx = np.array([o["dx"] for o in db1])
        dy = np.array([o["dy"] for o in db1])
        unq = np.unique(np.vstack((dy, dx)), axis=1)
        nrois = unq.shape[1]
        if nrois < len(db1):
            nplanes = len(db1) // nrois
            Lx = np.array([o["Lx"] for o in db1])
            Ly = np.array([o["Ly"] for o in db1])
            ymax = (dy + Ly).max()
            xmax = (dx + Lx).max()
            nX = np.ceil(np.sqrt(ymax * xmax * nplanes) / xmax)
            nX = int(nX)
            nY = int(np.ceil(len(db1) / nX))
            for j in range(nplanes):
                for k in range(nrois):
                    dx[j * nrois + k] += (j % nX) * xmax
                    dy[j * nrois + k] += int(j / nX) * ymax
    return dy, dx

save_mat #

save_mat(ops, stat, F, Fneu, spks, iscell, redcell, F_chan2=None, Fneu_chan2=None)

Save Suite2p results to a MATLAB .mat file.

Converts pathlib paths to strings, replaces None values with empty arrays, and writes all results to "Fall.mat" in the save path specified in ops.

Parameters:

Name Type Description Default
ops dict

Suite2p options dictionary. Must contain "save_path". The "date_proc" field is converted to a string if present.

 required
stat ndarray

Array of ROI statistics dictionaries (one per detected cell).

 required
F ndarray

Fluorescence traces of shape (n_cells, n_frames).

 required
Fneu ndarray

Neuropil fluorescence traces of shape (n_cells, n_frames).

 required
spks ndarray

Deconvolved spike traces of shape (n_cells, n_frames).

 required
iscell ndarray

Cell classification array of shape (n_cells, 2), with columns for binary label and probability.

 required
redcell ndarray

Red channel cell classification array, or None.

 required
F_chan2 ndarray

Second channel fluorescence traces of shape (n_cells, n_frames).

 None
Fneu_chan2 ndarray

Second channel neuropil fluorescence traces of shape (n_cells, n_frames).

 None
Source code in suite2p/io/save.py 
 14
 15
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
def save_mat(ops, stat, F, Fneu, spks, iscell, redcell,
             F_chan2=None, Fneu_chan2=None):
    """
    Save Suite2p results to a MATLAB .mat file.

    Converts pathlib paths to strings, replaces None values with empty arrays,
    and writes all results to "Fall.mat" in the save path specified in `ops`.

    Parameters
    ----------
    ops : dict
        Suite2p options dictionary. Must contain "save_path". The "date_proc"
        field is converted to a string if present.
    stat : numpy.ndarray
        Array of ROI statistics dictionaries (one per detected cell).
    F : numpy.ndarray
        Fluorescence traces of shape (n_cells, n_frames).
    Fneu : numpy.ndarray
        Neuropil fluorescence traces of shape (n_cells, n_frames).
    spks : numpy.ndarray
        Deconvolved spike traces of shape (n_cells, n_frames).
    iscell : numpy.ndarray
        Cell classification array of shape (n_cells, 2), with columns for
        binary label and probability.
    redcell : numpy.ndarray
        Red channel cell classification array, or None.
    F_chan2 : numpy.ndarray, optional
        Second channel fluorescence traces of shape (n_cells, n_frames).
    Fneu_chan2 : numpy.ndarray, optional
        Second channel neuropil fluorescence traces of shape (n_cells, n_frames).
    """
    ops_matlab = ops.copy()
    if ops_matlab.get("date_proc"):
        try:
            ops_matlab["date_proc"] = str(
                datetime.strftime(ops_matlab["date_proc"], "%Y-%m-%d %H:%M:%S.%f"))
        except:
            pass
    for k in ops_matlab.keys():
        if isinstance(ops_matlab[k], (pathlib.WindowsPath, pathlib.PosixPath)):
            ops_matlab[k] = os.fspath(ops_matlab[k].absolute())
        elif isinstance(ops_matlab[k], list) and len(ops_matlab[k]) > 0:
            if isinstance(ops_matlab[k][0], (pathlib.WindowsPath, pathlib.PosixPath)):
                ops_matlab[k] = [os.fspath(p.absolute()) for p in ops_matlab[k]]
                logger.info(f"{k}: {ops_matlab[k]}")

    stat = np.array(stat, dtype=object)


    # Check for None values in ops_matlab and replace with empty arrays in recursive manner
    def replace_none(d):
        for k, v in d.items():
            if v is None:
                logger.warning(f"'{k}' is None, replacing with empty array")
                d[k] = np.array([])
            elif isinstance(v, dict):
                replace_none(v)
    replace_none(ops_matlab)

    # Handle None variables by replacing with empty arrays
    if redcell is None:
        logger.warning("redcell is None, replacing with empty array")
        redcell = np.array([])
    if iscell is None:
        logger.warning("iscell is None, replacing with empty array")
        iscell = np.array([])

    if F_chan2 is None:
        scipy.io.savemat(
            file_name=os.path.join(ops["save_path"], "Fall.mat"), mdict={
                "stat": stat,
                "ops": ops_matlab,
                "F": F,
                "Fneu": Fneu,
                "spks": spks,
                "iscell": iscell,
                "redcell": redcell
            })
    else:
        scipy.io.savemat(
            file_name=os.path.join(ops["save_path"], "Fall.mat"), mdict={
                "stat": stat,
                "ops": ops_matlab,
                "F": F,
                "Fneu": Fneu,
                "spks": spks,
                "iscell": iscell,
                "redcell": redcell,
                "F_chan2": F_chan2,
                "Fneu_chan2": Fneu_chan2
            })

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

sbx_to_binary #

sbx_to_binary(settings, ndeadcols=-1, ndeadrows=0)

finds scanbox files and writes them to binaries

Parameters:

Name Type Description Default
settings dictionary

"nplanes", "data_path", "save_path", "save_folder", "fast_disk", "nchannels", "keep_movie_raw", "look_one_level_down"

 required

Returns:

Type Description
settings : dictionary of first plane

"Ly", "Lx", settings["reg_file"] or settings["raw_file"] is created binary
Source code in suite2p/io/sbx.py 
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
def sbx_to_binary(settings, ndeadcols=-1, ndeadrows=0):
    """  finds scanbox files and writes them to binaries

    Parameters
    ----------
    settings : dictionary
        "nplanes", "data_path", "save_path", "save_folder", "fast_disk",
        "nchannels", "keep_movie_raw", "look_one_level_down"

    Returns
    -------
        settings : dictionary of first plane
            "Ly", "Lx", settings["reg_file"] or settings["raw_file"] is created binary

    """
    if not HAS_SBX:
        raise ImportError("sbxreader is required for this file type, please 'pip install sbxreader'")

    settings1 = init_settings(settings)
    # the following should be taken from the metadata and not needed but the files are initialized before...
    nplanes = settings1[0]["nplanes"]
    nchannels = settings1[0]["nchannels"]
    # open all binary files for writing
    settings1, sbxlist, reg_file, reg_file_chan2 = find_files_open_binaries(settings1)
    iall = 0
    for j in range(settings1[0]["nplanes"]):
        settings1[j]["nframes_per_folder"] = np.zeros(len(sbxlist), np.int32)
    ik = 0
    if "sbx_ndeadcols" in settings1[0].keys():
        ndeadcols = int(settings1[0]["sbx_ndeadcols"])
    if "sbx_ndeadrows" in settings1[0].keys():
        ndeadrows = int(settings1[0]["sbx_ndeadrows"])

    if ndeadcols == -1 or ndeadrows == -1:
        # compute dead rows and cols from the first file
        tmpsbx = sbx_memmap(sbxlist[0])
        # do not remove dead rows in non-multiplane mode
        # This number should be different for each plane since the artifact is larger
        # for larger ETL jumps.
        if nplanes > 1 and ndeadrows == -1:
            colprofile = np.array(np.mean(tmpsbx[0][0][0], axis=1))
            ndeadrows = np.argmax(np.diff(colprofile)) + 1
        else:
            ndeadrows = 0
        # do not remove dead columns in unidirectional scanning mode
        # do this only if ndeadcols is -1
        if tmpsbx.metadata["scanning_mode"] == "bidirectional" and ndeadcols == -1:
            ndeadcols = tmpsbx.ndeadcols
        else:
            ndeadcols = 0
        del tmpsbx
        logger.info("Removing {0} dead columns while loading sbx data.".format(ndeadcols))
        logger.info("Removing {0} dead rows while loading sbx data.".format(ndeadrows))

    settings1[0]["sbx_ndeadcols"] = ndeadcols
    settings1[0]["sbx_ndeadrows"] = ndeadrows

    for ifile, sbxfname in enumerate(sbxlist):
        f = sbx_memmap(sbxfname)
        nplanes = f.shape[1]
        nchannels = f.shape[2]
        nframes = f.shape[0]
        iblocks = np.arange(0, nframes, settings1[0]["batch_size"])
        if iblocks[-1] < nframes:
            iblocks = np.append(iblocks, nframes)

        # data = nframes x nplanes x nchannels x pixels x pixels
        if nchannels > 1:
            nfunc = settings1[0]["functional_chan"] - 1
        else:
            nfunc = 0
        # loop over all frames
        for ichunk, onset in enumerate(iblocks[:-1]):
            offset = iblocks[ichunk + 1]
            im = np.array(f[onset:offset, :, :, ndeadrows:, ndeadcols:]) // 2
            im = im.astype(np.int16)
            im2mean = im.mean(axis=0).astype(np.float32) / len(iblocks)
            for ichan in range(nchannels):
                nframes = im.shape[0]
                im2write = im[:, :, ichan, :, :]
                for j in range(0, nplanes):
                    if iall == 0:
                        settings1[j]["meanImg"] = np.zeros((im.shape[3], im.shape[4]),
                                                      np.float32)
                        if nchannels > 1:
                            settings1[j]["meanImg_chan2"] = np.zeros(
                                (im.shape[3], im.shape[4]), np.float32)
                        settings1[j]["nframes"] = 0
                    if ichan == nfunc:
                        settings1[j]["meanImg"] += np.squeeze(im2mean[j, ichan, :, :])
                        reg_file[j].write(
                            bytearray(im2write[:, j, :, :].astype("int16")))
                    else:
                        settings1[j]["meanImg_chan2"] += np.squeeze(im2mean[j, ichan, :, :])
                        reg_file_chan2[j].write(
                            bytearray(im2write[:, j, :, :].astype("int16")))

                    settings1[j]["nframes"] += im2write.shape[0]
                    settings1[j]["nframes_per_folder"][ifile] += im2write.shape[0]
            ik += nframes
            iall += nframes

    # write settings files
    do_registration = settings1[0]["do_registration"]
    do_nonrigid = settings1[0]["nonrigid"]
    for settings in settings1:
        settings["Ly"] = im.shape[3]
        settings["Lx"] = im.shape[4]
        if not do_registration:
            settings["yrange"] = np.array([0, settings["Ly"]])
            settings["xrange"] = np.array([0, settings["Lx"]])
        #settings["meanImg"] /= settings["nframes"]
        #if nchannels>1:
        #    settings["meanImg_chan2"] /= settings["nframes"]
        np.save(settings["settings_path"], settings)
    # close all binary files and write settings files
    for j in range(0, nplanes):
        reg_file[j].close()
        if nchannels > 1:
            reg_file_chan2[j].close()
    return settings1[0]

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

send_jobs #

send_jobs(save_folder, host=None, username=None, password=None, server_root=None, local_root=None, n_cores=8)

send each plane to compute on server separately

add your own host, username, password and path on server for where to save the data

Source code in suite2p/io/server.py 
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
def send_jobs(save_folder, host=None, username=None, password=None, server_root=None,
              local_root=None, n_cores=8):
    """ send each plane to compute on server separately

    add your own host, username, password and path on server 
    for where to save the data

    """
    if not HAS_PARAMIKO:
        raise ImportError("paramiko required, please 'pip install paramiko'")

    if host is None:
        raise Exception("No server specified, please edit suite2p/io/server.py")

    # server_root is different from where you created the binaries, which is local_root
    nparts = len(Path(local_root).parts)
    # e.g. if server is Z:/path on local computer, and server_root+path on remote, then nparts=1
    save_folder_server = Path(*Path(save_folder).parts[nparts:])
    save_folder_server = Path(server_root) / save_folder_server
    save_path0_server = Path(*Path(save_folder_server).parts[:-1])
    save_folder_name = Path(save_folder).parts[-1]
    logger.info("save path on server: ", unix_path(save_path0_server))
    ssh = ssh_connect(host, username, password)

    # create bash file in home directory to run
    run_script = Path.home().joinpath(".suite2p/run_script.sh")
    if run_script.exists():
        os.remove(run_script)
    with open(run_script, "x", newline="") as f:
        f.write("#!/bin/bash\n")
        # server specific commands to activate python
        f.write("source ~/add_anaconda.sh\n")
        f.write("eval $(~/anaconda4/bin/conda shell.bash hook)\n")
        # activate suite2p environment
        f.write("source activate suite2p\n")
        # run suite2p single plane command with settings as argument
        f.write('python -m suite2p --single_plane --settings "$@"')

    ssh.exec_command("rm ~/run_script.sh")
    ssh.exec_command("chmod 777 ~/")
    ftp_client = ssh.open_sftp()
    ftp_client.put(run_script, "run_script.sh")
    ssh.exec_command("chmod 777 run_script.sh")

    pdirs = natsorted(glob.glob(save_folder + "/*/"))
    for k, pdir in enumerate(pdirs):
        ipl = int(Path(pdir).parts[-1][5:])
        logger.info(">>>>>>>>>> PLANE %d <<<<<<<<<" % ipl)
        settings_path_orig = pdir + "settings.npy"
        op = np.load(settings_path_orig, allow_pickle=True).item()
        fast_disk_orig = Path(op["fast_disk"])

        ## change paths
        op["save_path0"] = unix_path(save_path0_server)
        op["save_folder"] = save_folder_name
        save_path = save_path0_server / save_folder_name / ("plane%d" % ipl)
        op["save_path"] = unix_path(save_path)
        op["fast_disk"] = unix_path(save_path)
        op["settings_path"] = unix_path(save_path / "settings.npy")
        logger.info(op["settings_path"])
        ## move binary files to server if needed
        # check if file structure needs to be created on remote server
        copy = False
        try:
            ftp_client.stat(op["save_path"])
        except IOError:
            logger.info("copying files")
            ftp_client.mkdir(op["save_path"])
            copy = True
        op["reg_file"] = unix_path(save_path / "data.bin")
        if "raw_file" in op:
            op["raw_file"] = unix_path(save_path / "data_raw.bin")
            if copy:
                ftp_client.put(fast_disk_orig / "data_raw.bin", op["raw_file"])
            if "raw_file_chan2" in op:
                op["raw_file_chan2"] = unix_path(save_path / "data_chan2_raw.bin")
                if copy:
                    ftp_client.put(fast_disk_orig / "data_raw_chan2.bin",
                                   op["raw_file_chan2"])
        else:
            if copy:
                ftp_client.put(fast_disk_orig / "data.bin", op["reg_file"])
            if "reg_file_chan2" in op:
                op["reg_file_chan2"] = unix_path(save_path / "data_chan2.bin")
                if copy:
                    ftp_client.put(fast_disk_orig / "data_chan2.bin",
                                   op["reg_file_chan2"])

        # save final version of settings and send to server
        np.save(settings_path_orig, op)
        if copy:
            logger.info("copying settings")
            ftp_client.put(settings_path_orig, op["settings_path"])

        # run plane (server-specific command)
        run_command = '''bsub -n %d -J test_s2p%d -R"select[avx512]" -o out%d.txt "~/run_script.sh "%s" > log%d.txt''' % (
            n_cores, ipl, ipl, op["settings_path"], ipl)
        stdin, stdout, stderr = ssh.exec_command(run_command)
        logger.info(stdout.readlines()[0])

    ftp_client.close()

    logger.info("Command done, closing SSH connection")
    ssh.close()

ssh_connect #

ssh_connect(host, username, password, verbose=True)

from paramiko example

Source code in suite2p/io/server.py 
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
def ssh_connect(host, username, password, verbose=True):
    """ from paramiko example """
    i = 0
    while True:
        if verbose:
            logger.info("Trying to connect to %s (attempt %i/30)" % (host, i + 1))
        try:
            ssh = paramiko.SSHClient()
            ssh.set_missing_host_key_policy(paramiko.AutoAddPolicy())
            ssh.connect(host, username=username, password=password)
            if verbose:
                logger.info("Connected to %s" % host)
            break
        except paramiko.AuthenticationException:
            logger.info("Authentication failed when connecting to %s" % host)
            sys.exit(1)
        except:
            logger.info("Could not SSH to %s, waiting for it to start" % host)
            i += 1
            time.sleep(2)
        # If we could not connect within time limit
        if i == 30:
            logger.info("Could not connect to %s. Giving up" % host)
            sys.exit(1)
    return ssh

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

ome_to_binary #

ome_to_binary(dbs, settings, reg_file, reg_file_chan2)

Convert non-interleaved OME-TIFF files to binary, splitting channels by filename.

Designed for recordings where channels are stored in separate TIFF files (distinguished by "Ch1"/"Ch2" in the filename) rather than interleaved pages. Supports both single-page and multi-page TIFFs, and handles bidirectional Bruker plane ordering.

Parameters:

Name Type Description Default
dbs list of dict

Database dictionaries for each plane. Must contain keys "file_list", "first_files", "batch_size", "nplanes", "nchannels", "functional_chan", and optionally "bruker_bidirectional". Updated in-place with "Ly", "Lx", "nframes", "frames_per_file", "frames_per_folder", "meanImg", and "meanImg_chan2".

 required
settings dict

Suite2p settings dictionary, saved alongside each plane's database.

 required
reg_file list of file objects

Opened binary files for writing each plane's functional channel data.

 required
reg_file_chan2 list of file objects

Opened binary files for writing each plane's second channel data (used only when nchannels > 1).

 required

Returns:

Name Type Description
dbs list of dict

Updated database dictionaries with image dimensions, frame counts, and mean images populated.
Source code in suite2p/io/tiff.py 
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
def ome_to_binary(dbs, settings, reg_file, reg_file_chan2):
    """
    Convert non-interleaved OME-TIFF files to binary, splitting channels by filename.

    Designed for recordings where channels are stored in separate TIFF files
    (distinguished by "Ch1"/"Ch2" in the filename) rather than interleaved pages.
    Supports both single-page and multi-page TIFFs, and handles bidirectional
    Bruker plane ordering.

    Parameters
    ----------
    dbs : list of dict
        Database dictionaries for each plane. Must contain keys "file_list",
        "first_files", "batch_size", "nplanes", "nchannels", "functional_chan",
        and optionally "bruker_bidirectional". Updated in-place with "Ly", "Lx",
        "nframes", "frames_per_file", "frames_per_folder", "meanImg", and
        "meanImg_chan2".
    settings : dict
        Suite2p settings dictionary, saved alongside each plane's database.
    reg_file : list of file objects
        Opened binary files for writing each plane's functional channel data.
    reg_file_chan2 : list of file objects
        Opened binary files for writing each plane's second channel data
        (used only when nchannels > 1).

    Returns
    -------
    dbs : list of dict
        Updated database dictionaries with image dimensions, frame counts, and
        mean images populated.
    """
    t0 = time.time()

    nplanes = dbs[0]["nplanes"]
    nchannels = dbs[0]["nchannels"]

    # get file list from dbs (already populated by run_s2p.py)
    fs = dbs[0]["file_list"]
    first_files = dbs[0]["first_files"]
    batch_size = dbs[0]["batch_size"]
    use_sktiff = not HAS_SCANIMAGE

    fs_Ch1, fs_Ch2 = [], []
    for f in fs:
        if f.find("Ch1") > -1:
            if dbs[0]["functional_chan"] == 1:
                fs_Ch1.append(f)
            else:
                fs_Ch2.append(f)
        else:
            if dbs[0]["functional_chan"] == 1:
                fs_Ch2.append(f)
            else:
                fs_Ch1.append(f)

    if len(fs_Ch2) == 0:
        dbs[0]["nchannels"] = 1
        nchannels = 1
    logger.info(f"nchannels = {nchannels}")

    # loop over all tiffs
    TiffReader = ScanImageTiffReader if HAS_SCANIMAGE else TiffFile
    with TiffReader(fs_Ch1[0]) as tif:
        if HAS_SCANIMAGE:
            n_pages = tif.shape()[0] if len(tif.shape()) > 2 else 1
            shape = tif.shape()[-2:]
        else:
            n_pages = len(tif.pages)
            im0 = tif.pages[0].asarray()
            shape = im0.shape

    for db in dbs:
        db["nframes"] = 0
        db["frames_per_folder"] = np.zeros(first_files.sum(), "int")
        db["frames_per_file"] = np.ones(len(fs_Ch1), "int") if n_pages==1 else np.zeros(len(fs_Ch1), "int")
        db["meanImg"] = np.zeros(shape, np.float32)
        if nchannels > 1:
            db["meanImg_chan2"] = np.zeros(shape, np.float32)

    bruker_bidirectional = dbs[0].get("bruker_bidirectional", False)
    iplanes = np.arange(0, nplanes)
    if not bruker_bidirectional:
        iplanes = np.tile(iplanes[np.newaxis, :],
                          int(np.ceil(len(fs_Ch1) / nplanes))).flatten()
        iplanes = iplanes[:len(fs_Ch1)]
    else:
        iplanes = np.hstack((iplanes, iplanes[::-1]))
        iplanes = np.tile(iplanes[np.newaxis, :],
                          int(np.ceil(len(fs_Ch1) / (2 * nplanes)))).flatten()
        iplanes = iplanes[:len(fs_Ch1)]

    itot = 0
    for ik, file in enumerate(fs_Ch1):
        ip = iplanes[ik]
        # read tiff
        if n_pages==1:    
            with TiffReader(file) as tif:
                im = tif.data()  if HAS_SCANIMAGE else tif.pages[0].asarray()
            if im.dtype.type == np.uint16:
                im = (im // 2)
            im = im.astype(np.int16)

            # write to binary
            dbs[ip]["nframes"] += 1
            dbs[ip]["frames_per_folder"][0] += 1
            dbs[ip]["meanImg"] += im.astype(np.float32)
            reg_file[ip].write(bytearray(im))
            #gc.collect()
        else:
            tif, Ltif = open_tiff(file, not HAS_SCANIMAGE)
            # keep track of the plane identity of the first frame (channel identity is assumed always 0)
            ix = 0
            while 1:
                im = read_tiff(file, tif, Ltif, ix, batch_size, use_sktiff)
                if im is None:
                    break
                nframes = im.shape[0]
                ix += nframes
                itot += nframes
                reg_file[ip].write(bytearray(im))
                dbs[ip]["meanImg"] += im.astype(np.float32).sum(axis=0)
                dbs[ip]["nframes"] += im.shape[0]
                dbs[ip]["frames_per_file"][ik] += nframes
                dbs[ip]["frames_per_folder"][0] += nframes
                if itot % 1000 == 0:
                    logger.info("%d frames of binary, time %0.2f sec." % (itot, time.time() - t0))
                gc.collect()            

    if nchannels > 1:
        itot = 0
        for ik, file in enumerate(fs_Ch2):
            ip = iplanes[ik]
            if n_pages==1:
                with TiffReader(file) as tif:
                    im = tif.data() if HAS_SCANIMAGE else tif.pages[0].asarray()
                if im.dtype.type == np.uint16:
                    im = (im // 2)
                im = im.astype(np.int16)
                dbs[ip]["meanImg_chan2"] += im.astype(np.float32)
                reg_file_chan2[ip].write(bytearray(im))
            else:
                tif, Ltif = open_tiff(file, not HAS_SCANIMAGE)
                ix = 0
                while 1:
                    im = read_tiff(file, tif, Ltif, ix, batch_size, use_sktiff)
                    if im is None:
                        break
                    nframes = im.shape[0]
                    ix += nframes
                    itot += nframes
                    dbs[ip]["meanImg_chan2"] += im.astype(np.float32).sum(axis=0)
                    reg_file_chan2[ip].write(bytearray(im))
                    if itot % 1000 == 0:
                        logger.info("%d frames of binary, time %0.2f sec." % (itot, time.time() - t0))
                    gc.collect()

    # update dbs with image dimensions and mean images
    for db in dbs:
        db["Ly"], db["Lx"] = shape
        db["meanImg"] /= db["nframes"]
        if nchannels > 1:
            db["meanImg_chan2"] /= db["nframes"]
        # Save db and settings to each plane folder
        np.save(db["db_path"], db)
        np.save(db["settings_path"], settings)

    return dbs

open_tiff #

open_tiff(file, sktiff)

Open a TIFF file and return the reader object and the number of pages.

Uses either tifffile or ScanImageTiffReader depending on sktiff.

Parameters:

Name Type Description Default
file str

Path to the TIFF file.

 required
sktiff bool

If True, use tifffile (TiffFile). If False, use ScanImageTiffReader.

 required

Returns:

Name Type Description
tif TiffFile or ScanImageTiffReader

Opened TIFF reader object.
Ltif int

Number of pages (frames) in the TIFF file.
Source code in suite2p/io/tiff.py 
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
def open_tiff(file, sktiff):
    """
    Open a TIFF file and return the reader object and the number of pages.

    Uses either tifffile or ScanImageTiffReader depending on `sktiff`.

    Parameters
    ----------
    file : str
        Path to the TIFF file.
    sktiff : bool
        If True, use tifffile (TiffFile). If False, use ScanImageTiffReader.

    Returns
    -------
    tif : TiffFile or ScanImageTiffReader
        Opened TIFF reader object.
    Ltif : int
        Number of pages (frames) in the TIFF file.
    """
    if sktiff:
        tif = TiffFile(file)
        Ltif = len(tif.pages)
    else:
        tif = ScanImageTiffReader(file)
        Ltif = 1 if len(tif.shape()) < 3 else tif.shape()[0]  # single page tiffs
    return tif, Ltif

read_tiff #

read_tiff(file, tif, Ltif, ix, batch_size, use_sktiff)

Read a batch of frames from an open TIFF file starting at index ix.

Reads up to batch_size frames and converts the data to int16. Returns None if ix is past the end of the file.

Parameters:

Name Type Description Default
file str

Path to the TIFF file (used by tifffile's imread when use_sktiff is True).

 required
tif TiffFile or ScanImageTiffReader

Already-opened TIFF reader object.

 required
Ltif int

Total number of pages (frames) in the TIFF file.

 required
ix int

Starting frame index to read from.

 required
batch_size int

Maximum number of frames to read in this batch.

 required
use_sktiff bool

If True, read with tifffile (imread). If False, read with ScanImageTiffReader.

 required

Returns:

Name Type Description
im ndarray or None

Frames as an int16 array of shape (nfr, Ly, Lx), or None if ix >= Ltif.
Source code in suite2p/io/tiff.py 
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
def read_tiff(file, tif, Ltif, ix, batch_size, use_sktiff):
    """
    Read a batch of frames from an open TIFF file starting at index `ix`.

    Reads up to `batch_size` frames and converts the data to int16. Returns None
    if `ix` is past the end of the file.

    Parameters
    ----------
    file : str
        Path to the TIFF file (used by tifffile's `imread` when `use_sktiff` is True).
    tif : TiffFile or ScanImageTiffReader
        Already-opened TIFF reader object.
    Ltif : int
        Total number of pages (frames) in the TIFF file.
    ix : int
        Starting frame index to read from.
    batch_size : int
        Maximum number of frames to read in this batch.
    use_sktiff : bool
        If True, read with tifffile (`imread`). If False, read with ScanImageTiffReader.

    Returns
    -------
    im : numpy.ndarray or None
        Frames as an int16 array of shape (nfr, Ly, Lx), or None if `ix >= Ltif`.
    """
    if ix >= Ltif:
        return None
    nfr = min(Ltif - ix, batch_size)
    if use_sktiff:
        im = imread(file, key=range(ix, ix + nfr))
    elif Ltif == 1:
        im = tif.data()
    else:
        im = tif.data(beg=ix, end=ix + nfr)
    # for single-page tiffs, add 1st dim
    if len(im.shape) < 3:
        im = np.expand_dims(im, axis=0)

    # check if uint16
    if im.dtype.type == np.uint16:
        im = (im // 2).astype(np.int16)
    elif im.dtype.type == np.int32:
        im = (im // 2).astype(np.int16)
    elif im.dtype.type != np.int16:
        im = im.astype(np.int16)

    if im.shape[0] > nfr:
        im = im[:nfr, :, :]

    return im

tiff_to_binary #

tiff_to_binary(dbs, settings, reg_file, reg_file_chan2)

Read TIFF files and write interleaved plane/channel data to binary files.

Iterates over all TIFF files listed in dbs[0]["file_list"], de-interleaves planes and channels, and writes each plane's frames to the corresponding binary file. Also computes per-plane mean images and frame counts.

Parameters:

Name Type Description Default
dbs list of dict

Database dictionaries for each plane/ROI. Must contain keys "file_list", "first_files", "batch_size", "force_sktiff", "nplanes", "nchannels", and optionally "nrois", "swap_order", "lines". Updated in-place with "Ly", "Lx", "nframes", "frames_per_file", "frames_per_folder", "meanImg", and "meanImg_chan2".

 required
settings dict

Suite2p settings dictionary, saved alongside each plane's database.

 required
reg_file list of file objects

Opened binary files for writing each plane's functional channel data.

 required
reg_file_chan2 list of file objects

Opened binary files for writing each plane's second channel data (used only when nchannels > 1).

 required

Returns:

Name Type Description
dbs list of dict

Updated database dictionaries with image dimensions, frame counts, and mean images populated.
Source code in suite2p/io/tiff.py 
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
def tiff_to_binary(dbs, settings, reg_file, reg_file_chan2):
    """
    Read TIFF files and write interleaved plane/channel data to binary files.

    Iterates over all TIFF files listed in `dbs[0]["file_list"]`, de-interleaves
    planes and channels, and writes each plane's frames to the corresponding binary
    file. Also computes per-plane mean images and frame counts.

    Parameters
    ----------
    dbs : list of dict
        Database dictionaries for each plane/ROI. Must contain keys "file_list",
        "first_files", "batch_size", "force_sktiff", "nplanes", "nchannels", and
        optionally "nrois", "swap_order", "lines". Updated in-place with "Ly", "Lx",
        "nframes", "frames_per_file", "frames_per_folder", "meanImg", and
        "meanImg_chan2".
    settings : dict
        Suite2p settings dictionary, saved alongside each plane's database.
    reg_file : list of file objects
        Opened binary files for writing each plane's functional channel data.
    reg_file_chan2 : list of file objects
        Opened binary files for writing each plane's second channel data
        (used only when nchannels > 1).

    Returns
    -------
    dbs : list of dict
        Updated database dictionaries with image dimensions, frame counts, and
        mean images populated.
    """

    t0 = time.time()

    fs = dbs[0]["file_list"]
    first_files = dbs[0]["first_files"]

    # try tiff readers
    batch_size = dbs[0]["batch_size"]
    use_sktiff = True if dbs[0]["force_sktiff"] else use_sktiff_reader(fs[0], batch_size=batch_size)

    nplanes, nchannels = dbs[0]["nplanes"], dbs[0]["nchannels"]
    # processing for multiple ROIs (nrois > 1 if mesoscope recording)
    nrois = dbs[0].get("nrois", 1)
    batch_size = nplanes * nchannels * math.ceil(batch_size / (nplanes * nchannels))

    # loop over all tiffs
    which_folder = -1
    ntotal = 0
    swap = dbs[0].get("swap_order", False)
    for ifile, file in enumerate(fs):
        # open tiff
        tif, Ltif = open_tiff(file, use_sktiff)
        # keep track of the plane identity of the first frame (channel identity is assumed always 0)
        if first_files[ifile]:
            which_folder += 1
            iplane = 0
        ix = 0
        while 1:
            im = read_tiff(file, tif, Ltif, ix, batch_size, use_sktiff)
            if im is None:
                break          
            nframes = im.shape[0]
            for j in range(0, nplanes):
                if ifile == 0 and ix == 0:
                    Ly, Lx = im.shape[1], im.shape[2]
                    for k in range(nrois):
                        jk = j*nrois + k
                        Ly = (dbs[jk]["lines"][-1] + 1 - dbs[jk]["lines"][0] 
                              if nrois > 1 else Ly)
                        dbs[jk]["Ly"], dbs[jk]["Lx"] = Ly, Lx
                        dbs[jk]["nframes"] = 0
                        dbs[jk]["frames_per_file"] = np.zeros(len(fs), "int")
                        dbs[jk]["frames_per_folder"] = np.zeros(first_files.sum(), "int")
                        dbs[jk]["meanImg"] = np.zeros((Ly, Lx), "float64")
                        if nchannels > 1:
                            dbs[jk]["meanImg_chan2"] = np.zeros((Ly, Lx), "float64")

                if not swap:
                    i0 = nchannels * ((iplane + j) % nplanes)
                else:
                    i0 = (iplane + j) % (nplanes*nchannels)

                if nchannels > 1:
                    nfunc = dbs[jk]["functional_chan"] - 1
                else:
                    nfunc = 0

                #print(i0, int(i0) + (swap+1)*nfunc)

                im2write = im[int(i0) + (swap+1)*nfunc:nframes:nplanes * nchannels]

                for k in range(nrois):
                    jk = j*nrois + k
                    if nrois > 1:
                        imk = im2write[:, dbs[jk]["lines"][0] : dbs[jk]["lines"][-1] + 1]
                    else:
                        imk = im2write
                    reg_file[jk].write(bytearray(imk))
                    dbs[jk]["meanImg"] += imk.sum(axis=0).astype("float64")
                    dbs[jk]["nframes"] += imk.shape[0]
                    dbs[jk]["frames_per_file"][ifile] += imk.shape[0]
                    dbs[jk]["frames_per_folder"][which_folder] += imk.shape[0]

                if nchannels > 1:
                    #print(int(i0) + (swap+1)*(1 - nfunc))
                    im2write = im[int(i0) + (swap+1)*(1 - nfunc):nframes:nplanes * nchannels]
                    for k in range(nrois):
                        jk = j*nrois + k
                        if nrois > 1:
                            imk = im2write[:, dbs[jk]["lines"][0] : dbs[jk]["lines"][-1] + 1]
                        else:
                            imk = im2write
                        reg_file_chan2[jk].write(bytearray(imk))
                        dbs[jk]["meanImg_chan2"] += imk.sum(axis=0).astype("float64")
            if not swap:
                iplane = (iplane - nframes / nchannels) % nplanes
            else:
                iplane = (iplane - nframes) % (nchannels * nplanes)
            ix += nframes
            ntotal += nframes
            if ntotal % (batch_size * 4) == 0:
                logger.info("%d frames of binary, time %0.2f sec." %
                      (ntotal, time.time() - t0))
        gc.collect()
    # write dbs and settings files
    for db in dbs:
        db["meanImg"] /= db["nframes"]
        if nchannels > 1:
            db["meanImg_chan2"] /= db["nframes"]
        np.save(db["db_path"], db)
        np.save(db["settings_path"], settings)

    # close all binary files and write settings files
    for jk in range(0, nplanes * nrois):
        reg_file[jk].close()
        if nchannels > 1:
            reg_file_chan2[jk].close()
    return dbs

use_sktiff_reader #

use_sktiff_reader(tiff_filename, batch_size)

Test whether ScanImageTiffReader can read a TIFF file.

Attempts to read a small batch from the file using ScanImageTiffReader. If it succeeds, returns False (use ScanImageTiffReader). If it fails or the package is not installed, returns True (fall back to tifffile).

Parameters:

Name Type Description Default
tiff_filename str

Path to the TIFF file to test.

 required
batch_size int

Number of frames to attempt reading as a test.

 required

Returns:

Name Type Description
use_sktiff bool

True if tifffile should be used, False if ScanImageTiffReader works.
Source code in suite2p/io/tiff.py 
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
def use_sktiff_reader(tiff_filename, batch_size):
    """
    Test whether ScanImageTiffReader can read a TIFF file.

    Attempts to read a small batch from the file using ScanImageTiffReader. If it
    succeeds, returns False (use ScanImageTiffReader). If it fails or the package is
    not installed, returns True (fall back to tifffile).

    Parameters
    ----------
    tiff_filename : str
        Path to the TIFF file to test.
    batch_size : int
        Number of frames to attempt reading as a test.

    Returns
    -------
    use_sktiff : bool
        True if tifffile should be used, False if ScanImageTiffReader works.
    """
    if HAS_SCANIMAGE:
        try:
            with ScanImageTiffReader(tiff_filename) as tif:
                tif.data() if len(tif.shape()) < 3 else tif.data(
                    beg=0, end=np.minimum(batch_size,
                                          tif.shape()[0] - 1))
            return False
        except:
            logger.info(
                "NOTE: ScanImageTiffReader not working for this tiff type, using tifffile"
            )
            return True
    else:
        logger.info("NOTE: ScanImageTiffReader not installed, using tifffile")
        return True

Copyright © 2023 Howard Hughes Medical Institute, Authored by Carsen Stringer and Marius Pachitariu.

find_files_open_binaries #

find_files_open_binaries(settings)

Find input files and open binary files for writing.

Parameters:

Name Type Description Default
settings dict

Suite2p settings dictionary.

 required

Returns:

Type Description
None

No longer used!!
Source code in suite2p/io/utils.py 
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
def find_files_open_binaries(settings):
    """
    Find input files and open binary files for writing.

    Parameters
    ----------
    settings : dict
        Suite2p settings dictionary.

    Returns
    -------
    None
        No longer used!!
    """
    return None

get_file_list #

get_file_list(db)

Build the list of input files to process from the database configuration.

Supports three modes: an explicit file list (db["file_list"]), subfolder-based discovery (db["subfolders"]), or recursive search with optional one-level-down lookup (db["look_one_level_down"]).

Parameters:

Name Type Description Default
db dict

Database dictionary. Must contain "data_path" (list of str). Optionally contains "file_list" (list of str), "subfolders" (list of str), "look_one_level_down" (bool), and "input_format" (str, default "tif").

 required

Returns:

Name Type Description
fsall list of str

List of all file paths to process.
first_files ndarray

Boolean array of length len(fsall), where True marks the first file from each folder.
Source code in suite2p/io/utils.py 
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
def get_file_list(db):
    """
    Build the list of input files to process from the database configuration.

    Supports three modes: an explicit file list (db["file_list"]), subfolder-based
    discovery (db["subfolders"]), or recursive search with optional one-level-down
    lookup (db["look_one_level_down"]).

    Parameters
    ----------
    db : dict
        Database dictionary. Must contain "data_path" (list of str). Optionally
        contains "file_list" (list of str), "subfolders" (list of str),
        "look_one_level_down" (bool), and "input_format" (str, default "tif").

    Returns
    -------
    fsall : list of str
        List of all file paths to process.
    first_files : numpy.ndarray
        Boolean array of length len(fsall), where True marks the first file from
        each folder.
    """
    data_path = db["data_path"]
    input_format = db.get("input_format", "tif")
    # use a user-specified list of tiffs
    if db.get("file_list", None) is not None:
        fsall = []
        for f in db["file_list"]:
            fsall.append(os.path.join(data_path[0], f))
        first_files = np.zeros(len(fsall), dtype="bool")
        first_files[0] = True
        logger.info(f"** Found {len(fsall)} files - converting to binary **")
    else:
        if len(data_path) == 1 and db.get("subfolders", None) is not None:
                fold_list = []
                for folder_down in db["subfolders"]:
                    fold = os.path.join(data_path[0], folder_down)
                    fold_list.append(fold)
        else:
            fold_list = data_path
        fsall = []
        first_files = []
        for k, fld in enumerate(fold_list):
            fs, firsts = list_files(fld, db["look_one_level_down"],
                                    EXTS[input_format])
            fsall.extend(fs)
            first_files.extend(list(firsts))
        if len(fsall) == 0:
            logger.info(f"Could not find any {EXTS[input_format]} files in {data_path}")
            raise Exception("no files found")
        else:
            first_files = np.array(first_files).astype("bool")
            logger.info(f"** Found {len(fsall)} files - converting to binary **")
    return fsall, first_files

get_suite2p_path #

get_suite2p_path(path)

Find the root suite2p folder within a given path.

Walks the path components backwards to locate the last occurrence of a folder named "suite2p" and returns the path up to and including that folder.

Parameters:

Name Type Description Default
path str or Path

File or directory path that should contain a "suite2p" folder component.

 required

Returns:

Name Type Description
new_path Path

Path truncated at the "suite2p" folder.

Raises:

Type Description
FileNotFoundError

If "suite2p" is not found in any component of path.
Source code in suite2p/io/utils.py 
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
def get_suite2p_path(path: Path) -> Path:
    """
    Find the root `suite2p` folder within a given path.

    Walks the path components backwards to locate the last occurrence of a folder
    named "suite2p" and returns the path up to and including that folder.

    Parameters
    ----------
    path : str or pathlib.Path
        File or directory path that should contain a "suite2p" folder component.

    Returns
    -------
    new_path : pathlib.Path
        Path truncated at the "suite2p" folder.

    Raises
    ------
    FileNotFoundError
        If "suite2p" is not found in any component of `path`.
    """

    path = Path(path)  # In case `path` is a string

    # Cheap sanity check
    if "suite2p" in str(path):
        # Walk the folders in path backwards
        for path_idx in range(len(path.parts) - 1, 0, -1):
            if path.parts[path_idx] == "suite2p":
                new_path = Path(path.parts[0])
                for path_part in path.parts[1:path_idx + 1]:
                    new_path = new_path.joinpath(path_part)
                break
    else:
        raise FileNotFoundError("The `suite2p` folder was not found in path")
    return new_path

init_dbs #

init_dbs(db0)

Initialize per-plane database dictionaries and create output directories.

Creates a deep copy of db0 for each plane (and each ROI for mesoscope recordings), setting up save paths, binary file paths, and fast-disk directories.

Parameters:

Name Type Description Default
db0 dict

Base database dictionary. Must contain "nplanes", "nchannels", "keep_movie_raw", "save_path0", and "save_folder". Optionally contains "fast_disk", "iplane", "lines", "dy", and "dx" (for mesoscope recordings with multiple ROIs).

 required

Returns:

Name Type Description
dbs list of dict

List of per-plane database dictionaries, each with added keys "save_path", "fast_disk", "settings_path", "db_path", "reg_file", and optionally "raw_file", "reg_file_chan2", "raw_file_chan2", "lines", "dy", "dx", "iroi", "iplane".
Source code in suite2p/io/utils.py 
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
def init_dbs(db0):
    """
    Initialize per-plane database dictionaries and create output directories.

    Creates a deep copy of `db0` for each plane (and each ROI for mesoscope recordings),
    setting up save paths, binary file paths, and fast-disk directories.

    Parameters
    ----------
    db0 : dict
        Base database dictionary. Must contain "nplanes", "nchannels",
        "keep_movie_raw", "save_path0", and "save_folder". Optionally contains
        "fast_disk", "iplane", "lines", "dy", and "dx" (for mesoscope recordings
        with multiple ROIs).

    Returns
    -------
    dbs : list of dict
        List of per-plane database dictionaries, each with added keys "save_path",
        "fast_disk", "settings_path", "db_path", "reg_file", and optionally
        "raw_file", "reg_file_chan2", "raw_file_chan2", "lines", "dy", "dx",
        "iroi", "iplane".
    """
    nplanes = db0["nplanes"]
    nchannels = db0["nchannels"]
    keep_movie_raw = db0["keep_movie_raw"]
    nfolders = nplanes
    iplane = db0.get("iplane", np.arange(0, nplanes))
    has_lines = False
    if "lines" in db0 and db0["lines"] is not None and len(db0["lines"]) > 0:
        nrois = len(db0["lines"])
        db0["nrois"] = nrois
        nfolders *= nrois
        logger.info(f"NOTE: nplanes={nplanes}, nrois={nrois} => nfolders = {nfolders}")
        # replicate lines across planes if nplanes > 1
        if nplanes > 1:
            lines0, dy0, dx0 = db0["lines"].copy(), db0["dy"].copy(), db0["dx"].copy()
            dy0, dx0 = np.array(dy0), np.array(dx0)
            dy = np.tile(dy0[np.newaxis, :], (nplanes, 1)).flatten()
            dx = np.tile(dx0[np.newaxis, :], (nplanes, 1)).flatten()
            lines = []
            [lines.extend(lines0) for _ in range(nplanes)]
            iroi = np.tile(np.arange(nrois)[np.newaxis,:], (nplanes, 1)).flatten()
            iplane = np.tile(np.arange(nplanes)[:, np.newaxis], (1, nrois)).flatten()
        else:
            lines, dy, dx = db0["lines"].copy(), db0["dy"].copy(), db0["dx"].copy()
            iroi = np.arange(nrois)
            iplane = np.zeros(nrois, "int")
        has_lines = True 

    dbs = []
    if db0.get("fast_disk", None) is None or len(db0["fast_disk"]) == 0:
        db0["fast_disk"] = db0["save_path0"]
    fast_disk = db0["fast_disk"]

    # compile dbs into list across planes
    for j in range(0, nfolders):
        db = copy.deepcopy(db0)
        db["save_path"] = os.path.join(db["save_path0"], db["save_folder"], f"plane{j}")
        fast_disk = os.path.join(db["fast_disk"], "suite2p", f"plane{j}")
        db["fast_disk"] = fast_disk
        db["settings_path"] = os.path.join(db["save_path"], "settings.npy")
        db["db_path"] = os.path.join(db["save_path"], "db.npy")
        db["reg_file"] = os.path.join(fast_disk, "data.bin")
        if keep_movie_raw:
            db["raw_file"] = os.path.join(fast_disk, "data_raw.bin")
        if has_lines:
            db["lines"], db["dy"], db["dx"] = lines[j], dy[j], dx[j]
            db["iroi"] = iroi[j]
        db["iplane"] = iplane[j]
        if nchannels > 1:
            db["reg_file_chan2"] = os.path.join(fast_disk, "data_chan2.bin")
            if keep_movie_raw:
                db["raw_file_chan2"] = os.path.join(fast_disk, "data_chan2_raw.bin")

        os.makedirs(db["fast_disk"], exist_ok=True)
        os.makedirs(db["save_path"], exist_ok=True)
        dbs.append(db)
    return dbs

list_files #

list_files(froot, look_one_level_down, exts)

Collect files matching the given extensions from a folder, optionally including subfolders.

Parameters:

Name Type Description Default
froot str

Root directory to search for files.

 required
look_one_level_down bool

If True, also search immediate subdirectories of froot.

 required
exts list of str

Glob patterns to match (e.g. [".tif", ".tiff"]).

 required

Returns:

Name Type Description
fs list of str

Naturally sorted list of matching file paths.
first_files ndarray

Boolean array of length len(fs), where True marks the first file from each folder (used to track folder boundaries).
Source code in suite2p/io/utils.py 
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
def list_files(froot, look_one_level_down, exts):
    """
    Collect files matching the given extensions from a folder, optionally including subfolders.

    Parameters
    ----------
    froot : str
        Root directory to search for files.
    look_one_level_down : bool
        If True, also search immediate subdirectories of `froot`.
    exts : list of str
        Glob patterns to match (e.g. ["*.tif", "*.tiff"]).

    Returns
    -------
    fs : list of str
        Naturally sorted list of matching file paths.
    first_files : numpy.ndarray
        Boolean array of length len(fs), where True marks the first file from each
        folder (used to track folder boundaries).
    """
    fs = []
    first_files = np.zeros(0, "bool")
    for e in exts:
        lpath = os.path.join(froot, e)
        fs.extend(glob.glob(lpath))
    fs = natsorted(set(fs))
    if len(fs) > 0:
        first_files = np.zeros(len(fs), "bool")
        first_files[0] = True
    lfs = len(fs)
    if look_one_level_down:
        fdir = natsorted(glob.glob(os.path.join(froot, "*/")))
        for folder_down in fdir:
            fsnew = []
            for e in exts:
                lpath = os.path.join(folder_down, e)
                fsnew.extend(glob.glob(lpath))
            fsnew = natsorted(set(fsnew))
            if len(fsnew) > 0:
                fs.extend(fsnew)
                first_files = np.append(first_files, np.zeros((len(fsnew),), "bool"))
                first_files[lfs] = True
                lfs = len(fs)
    return fs, first_files