Add chan2 manual labeling

suite2p/__main__.py

@@ -51,16 +51,7 @@ def main():
             logging.exception(f'fatal error in {"run_plane" if args.single_plane else "run_s2p"}:')
             raise
 
-    else:
-        # Check if the OS is macOS and the machine is Apple Silicon (ARM-based)
-        if platform.system() == "Darwin" and 'arm' in platform.processor().lower():
-            # Set the number of threads for OpenMP and OpenBLAS
-            os.environ["OMP_NUM_THREADS"] = "1"
-            os.environ["OPENBLAS_NUM_THREADS"] = "1"
-            print("Environment set to use 1 thread for OpenMP and OpenBLAS (Apple Silicon macOS).")
-        else:
-            print("Not macOS on Apple Silicon, proceeding without limiting threads.")
-
+    else:   
         from suite2p import gui
         gui.run()#statfile="C:/DATA/exs2p/suite2p/plane0/stat.npy")
 

suite2p/detection/detect.py

@@ -166,8 +166,7 @@ def detection_wrapper(f_reg, diameter=[12., 12.], tau=1., fs=30, meanImg_chan2=N
 
     if mov is None:
         nbins = settings["nbins"]
-        bin_size = int(max(1, n_frames // nbins, np.round(tau * fs)))
-        #bin_size = int(max(1, np.round(tau * fs)))
+        bin_size = settings.get("bin_size") or int(max(1, n_frames // nbins, np.round(tau * fs)))
         logger.info("Binning movie in chunks of %2.2d frames" % bin_size)
         mov = bin_movie(f_reg, bin_size, yrange=yrange, xrange=xrange,
                         badframes=badframes, nbins=nbins)

suite2p/detection/sparsedetect.py

@@ -565,6 +565,8 @@ def sparsery(mov, sdmov, highpass_neuropil,
     scale, estimate_mode = find_best_scale(I=I, spatial_scale=spatial_scale)
 
     spatscale_pix = 3 * 2**scale
+    if isinstance(spatscale_pix, np.ndarray):
+        spatscale_pix = spatscale_pix.item()
     mask_window = int(((spatscale_pix * 1.5) // 2) * 2)
     Th2 = threshold_scaling * 5 * max(
         1, scale)  # threshold for accepted peaks (scale it by spatial scale)

suite2p/extraction/extract.py

@@ -57,21 +57,22 @@ def extract_traces(f_in, cell_masks, neuropil_masks, batch_size=500,
     if device.type == 'mps':
         device = torch.device('cpu')
 
-    npix_neuropil = torch.Tensor([len(nm) for nm in neuropil_masks]).to(device)
-    # create coo tensor of neuropil and cell masks 
-    ccol_indices = [m for nm in neuropil_masks for m in nm]
-    row_indices = [k for k in range(len(neuropil_masks)) for m in neuropil_masks[k]]
-    inds = torch.Tensor([ccol_indices, row_indices]).to(device)
-    # convert to csc (tried creating csc directly but it was slow)
-    nmasks = torch.sparse_coo_tensor(inds, torch.ones(len(row_indices), device=device), 
-                                     size=(Ly*Lx, ncells))
-    nmasks = nmasks.to_sparse_csc()
+    if neuropil_masks is not None:
+        npix_neuropil = torch.Tensor([len(nm) for nm in neuropil_masks]).to(device)
+        # create coo tensor of neuropil masks
+        ccol_indices = [m for nm in neuropil_masks for m in nm]
+        row_indices = [k for k in range(len(neuropil_masks)) for m in neuropil_masks[k]]
+        inds = torch.Tensor([ccol_indices, row_indices]).to(device)
+        # convert to csc (tried creating csc directly but it was slow)
+        nmasks = torch.sparse_coo_tensor(inds, torch.ones(len(row_indices), device=device),
+                                         size=(Ly*Lx, ncells))
+        nmasks = nmasks.to_sparse_csc()
 
     ccol_indices = [m for cm in cell_masks for m in cm[0]]
     row_indices = [k for k in range(len(cell_masks)) for m in cell_masks[k][0]]
     cell_lam = torch.Tensor([l for cm in cell_masks for l in cm[1]]).to(device)
     inds = torch.Tensor([ccol_indices, row_indices]).to(device)
-    cmasks = torch.sparse_coo_tensor(inds, cell_lam, 
+    cmasks = torch.sparse_coo_tensor(inds, cell_lam,
                                      size=(Ly*Lx, ncells))
     cmasks = cmasks.to_sparse_csc()
 
@@ -87,11 +88,12 @@ def extract_traces(f_in, cell_masks, neuropil_masks, batch_size=500,
         tstart, tend = n * batch_size, min((n+1) * batch_size, n_frames)
         data = torch.from_numpy(f_in[tstart : tend]).to(device)
         data = data.reshape(-1, Ly*Lx).float()
-
-        Fneu_batch = (data @ nmasks) / npix_neuropil
-        Fneu[:, tstart : tend] = Fneu_batch.T.cpu().numpy()
-
-        F_batch = data @ cmasks 
+
+        if neuropil_masks is not None:
+            Fneu_batch = (data @ nmasks) / npix_neuropil
+            Fneu[:, tstart : tend] = Fneu_batch.T.cpu().numpy()
+
+        F_batch = data @ cmasks
         F[:, tstart : tend] = F_batch.T.cpu().numpy()
 
     return F, Fneu

suite2p/gui/drawroi.py

@@ -8,7 +8,7 @@
 import numpy as np
 import pyqtgraph as pg
 from qtpy import QtGui, QtCore
-from qtpy.QtWidgets import QPushButton, QLabel, QLineEdit, QMainWindow, QGridLayout, QButtonGroup, QMessageBox, QWidget
+from qtpy.QtWidgets import QPushButton, QLabel, QLineEdit, QMainWindow, QGridLayout, QButtonGroup, QMessageBox, QWidget, QVBoxLayout
 from matplotlib.colors import hsv_to_rgb
 from scipy import stats
 from scipy.ndimage import rotate
@@ -18,6 +18,10 @@
 from ..detection.stats import roi_stats
 from ..extraction import preprocess
 from ..extraction.dcnv import oasis
+from ..extraction.extract import extract_traces
+from ..io.binary import BinaryFile
+from ..parameters import default_settings
+from ..run_s2p import _assign_torch_device
 
 
 def masks_and_traces(settings, stat_manual, stat_orig):
@@ -27,6 +31,8 @@ def masks_and_traces(settings, stat_manual, stat_orig):
         returns: F (ROIs x time), Fneu (ROIs x time), F_chan2, Fneu_chan2, settings, stat
         F_chan2 and Fneu_chan2 will be empty if no second channel
     """
+    # Merge with defaults to ensure all required keys are present
+    settings = {**default_settings(), **settings}
 
     t0 = time.time()
 
@@ -35,11 +41,12 @@ def masks_and_traces(settings, stat_manual, stat_orig):
     for n in range(len(stat_orig)):
         stat_all.append(stat_orig[n])
 
-    stat_all = roi_stats(stat_all, settings["Ly"], settings["Lx"], aspect=settings.get("aspect", None),
+    stat_all = np.array(stat_all)
+    stat_all = roi_stats(stat_all, settings["Ly"], settings["Lx"],
                          diameter=settings["diameter"])
     cell_masks = [
         masks.create_cell_mask(stat, Ly=settings["Ly"], Lx=settings["Lx"],
-                               allow_overlap=settings["allow_overlap"]) for stat in stat_all
+                               allow_overlap=settings["extraction"]["allow_overlap"]) for stat in stat_all
     ]
     cell_pix = masks.create_cell_pix(stat_all, Ly=settings["Ly"], Lx=settings["Lx"])
     manual_roi_stats = stat_all[:len(stat_manual)]
@@ -48,13 +55,26 @@ def masks_and_traces(settings, stat_manual, stat_orig):
         ypixs=[stat["ypix"] for stat in manual_roi_stats],
         xpixs=[stat["xpix"] for stat in manual_roi_stats],
         cell_pix=cell_pix,
-        inner_neuropil_radius=settings["inner_neuropil_radius"],
-        min_neuropil_pixels=settings["min_neuropil_pixels"],
+        inner_neuropil_radius=settings["extraction"]["inner_neuropil_radius"],
+        min_neuropil_pixels=settings["extraction"]["min_neuropil_pixels"],
     )
     print("Masks made in %0.2f sec." % (time.time() - t0))
 
-    F, Fneu, F_chan2, Fneu_chan2 = extract_traces_from_masks(settings, manual_cell_masks,
-                                                             manual_neuropil_masks)
+    # Extract traces from binary file
+    Ly, Lx = settings["Ly"], settings["Lx"]
+    batch_size = settings["extraction"]["batch_size"]
+    device = _assign_torch_device(settings["torch_device"])
+    f_reg = BinaryFile(Ly, Lx, settings["reg_file"])
+    F, Fneu = extract_traces(f_reg, manual_cell_masks, manual_neuropil_masks, batch_size=batch_size, device=device)
+    f_reg.close()
+
+    # Handle chan2 if present
+    if "reg_file_chan2" in settings and settings["reg_file_chan2"]:
+        f_reg_chan2 = BinaryFile(Ly, Lx, settings["reg_file_chan2"])
+        F_chan2, Fneu_chan2 = extract_traces(f_reg_chan2, manual_cell_masks, manual_neuropil_masks, batch_size=batch_size, device=device)
+        f_reg_chan2.close()
+    else:
+        F_chan2, Fneu_chan2 = None, None
 
     # compute activity statistics for classifier
     npix = np.array([stat_orig[n]["npix"] for n in range(len(stat_orig))
@@ -69,7 +89,7 @@ def masks_and_traces(settings, stat_manual, stat_orig):
             manual_roi_stats[n]["iplane"] = stat_orig[0]["iplane"]
 
     # subtract neuropil and compute skew, std from F
-    dF = F - settings["neucoeff"] * Fneu
+    dF = F - settings["extraction"]["neuropil_coefficient"] * Fneu
     sk = stats.skew(dF, axis=1)
     sd = np.std(dF, axis=1)
 
@@ -81,10 +101,10 @@ def masks_and_traces(settings, stat_manual, stat_orig):
             np.mean(manual_roi_stats[n]["xpix"])
         ]
 
-    dF = preprocess(F=dF, baseline=settings["baseline"], win_baseline=settings["win_baseline"],
-                    sig_baseline=settings["sig_baseline"], fs=settings["fs"],
-                    prctile_baseline=settings["prctile_baseline"])
-    spks = oasis(F=dF, batch_size=settings["batch_size"], tau=settings["tau"], fs=settings["fs"])
+    dF = preprocess(F=dF, baseline=settings["dcnv_preprocess"]["baseline"], win_baseline=settings["dcnv_preprocess"]["win_baseline"],
+                    sig_baseline=settings["dcnv_preprocess"]["sig_baseline"], fs=settings["fs"],
+                    prctile_baseline=settings["dcnv_preprocess"]["prctile_baseline"], device=device)
+    spks = oasis(F=dF, batch_size=settings["extraction"]["batch_size"], tau=settings["tau"], fs=settings["fs"])
 
     return F, Fneu, F_chan2, Fneu_chan2, spks, settings, manual_roi_stats
 
@@ -101,8 +121,10 @@ def __init__(self, bid, Text, parent=None):
         self.show()
 
     def press(self, parent, bid):
+        for b in range(len(parent.views)):
+            parent.viewbtns.button(b).setStyleSheet(parent.styleUnpressed)
+        parent.viewbtns.button(bid).setStyleSheet(parent.stylePressed)
         parent.img0.setImage(parent.masked_images[:, :, :, bid])
-
         parent.win.show()
         parent.show()
 
@@ -121,7 +143,7 @@ def __init__(self, parent):
         # layout = QtGui.QFormLayout()
         self.cwidget.setLayout(self.l0)
         self.stylePressed = ("QPushButton {Text-align: left; "
-                             "background-color: rgb(100,50,100); "
+                             "background-color: rgb(100,100,100); "
                              "color:white;}")
         self.styleUnpressed = ("QPushButton {Text-align: left; "
                                "background-color: rgb(50,50,50); "
@@ -187,7 +209,7 @@ def __init__(self, parent):
         self.saveGUI = False
         self.closeGUI = QPushButton("Save and Quit")
         self.closeGUI.setFont(QtGui.QFont("Arial", 8, QtGui.QFont.Bold))
-        self.closeGUI.clicked.connect(self.close_GUI)
+        self.closeGUI.clicked.connect(lambda: self.close_GUI())
         self.closeGUI.setEnabled(False)
         self.closeGUI.setFixedWidth(100)
         self.closeGUI.setStyleSheet(self.styleUnpressed)
@@ -198,17 +220,28 @@ def __init__(self, parent):
             "W: mean img", "E: mean img (enhanced)", "R: correlation map",
             "T: max projection"
         ]
+        self.has_chan2 = "meanImg_chan2" in parent.ops
+        if self.has_chan2:
+            self.views.append("Y: mean img chan2")
         b = 0
         self.viewbtns = QButtonGroup(self)
+        view_container = QWidget()
+        view_vbox = QVBoxLayout()
+        view_vbox.setContentsMargins(0, 0, 0, 0)
+        view_vbox.setSpacing(4)
+        view_container.setLayout(view_vbox)
         for names in self.views:
             btn = ViewButton(b, "&" + names, self)
             self.viewbtns.addButton(btn, b)
-            self.l0.addWidget(btn, b, 4, 1, 1)
+            view_vbox.addWidget(btn)
             btn.setEnabled(True)
             b += 1
-        b = 0
-        self.viewbtns.button(b).setChecked(True)
-        self.viewbtns.button(b).setStyleSheet(self.stylePressed)
+        view_vbox.addStretch()
+        self.l0.addWidget(view_container, 0, 4, 3, 1)
+        for b in range(len(self.views)):
+            self.viewbtns.button(b).setStyleSheet(self.styleUnpressed)
+        self.viewbtns.button(0).setChecked(True)
+        self.viewbtns.button(0).setStyleSheet(self.stylePressed)
 
         self.l0.addWidget(QLabel("neuropil"), 13, 13, 1, 1)
 
@@ -247,9 +280,7 @@ def close_GUI(self):
 
         # Append new stat file with old and save
         print("Saving new stat")
-        stat_all = self.new_stat.copy()
-        for n in range(len(self.parent.stat)):
-            stat_all.append(self.parent.stat[n])
+        stat_all = np.concatenate((self.new_stat, self.parent.stat))
         np.save(os.path.join(self.parent.basename, "stat.npy"), stat_all)
         iscell_prob = np.concatenate(
             (self.parent.iscell[:, np.newaxis], self.parent.probcell[:, np.newaxis]),
@@ -291,41 +322,38 @@ def close_GUI(self):
         self.close()
 
     def normalize_img_add_masks(self):
-        masked_image = np.zeros(
-            ((self.Ly, self.Lx, 3, 4)))  # 3 for RGB and 4 for buttons
-        for i in np.arange(4):  # 4 because 4 buttons
+        nviews = len(self.views)
+        masked_image = np.zeros((self.Ly, self.Lx, 3, nviews))
+        yr = slice(self.parent.ops["yrange"][0], self.parent.ops["yrange"][1])
+        xr = slice(self.parent.ops["xrange"][0], self.parent.ops["xrange"][1])
+        for i in np.arange(nviews):
             if i == 0:
-                mimg = np.zeros((self.Ly, self.Lx), np.float32)
-                mimg[self.parent.ops["yrange"][0]:self.parent.ops["yrange"][1],
-                     self.parent.ops["xrange"][0]:self.parent.
-                     settings["xrange"][1]] = self.parent.ops["meanImg"][
-                         self.parent.ops["yrange"][0]:self.parent.ops["yrange"][1],
-                         self.parent.ops["xrange"][0]:self.parent.ops["xrange"][1]]
-
+                src = self.parent.ops["meanImg"]
             elif i == 1:
-                mimg = np.zeros((self.Ly, self.Lx), np.float32)
-                mimg[self.parent.ops["yrange"][0]:self.parent.ops["yrange"][1],
-                     self.parent.ops["xrange"][0]:self.parent.
-                     settings["xrange"][1]] = self.parent.ops["meanImgE"][
-                         self.parent.ops["yrange"][0]:self.parent.ops["yrange"][1],
-                         self.parent.ops["xrange"][0]:self.parent.ops["xrange"][1]]
+                src = self.parent.ops["meanImgE"]
             elif i == 2:
-                mimg = np.zeros((self.Ly, self.Lx), np.float32)
-                mimg[self.parent.ops["yrange"][0]:self.parent.ops["yrange"][1],
-                     self.parent.ops["xrange"][0]:self.parent.
-                     settings["xrange"][1]] = self.parent.ops["Vcorr"]
-
+                src = self.parent.ops["Vcorr"]
+            elif i == 3:
+                src = self.parent.ops.get("max_proj", None)
+            elif i == 4 and self.has_chan2:
+                src = self.parent.ops["meanImg_chan2"]
             else:
-                mimg = np.zeros((self.Ly, self.Lx), np.float32)
-                if "max_proj" in self.parent.ops:
-                    mimg[self.parent.ops["yrange"][0]:self.parent.ops["yrange"][1],
-                         self.parent.ops["xrange"][0]:self.parent.
-                         settings["xrange"][1]] = self.parent.ops["max_proj"]
-
-            mimg1 = np.percentile(mimg, 1)
-            mimg99 = np.percentile(mimg, 99)
-            mimg = (mimg - mimg1) / (mimg99 - mimg1)
-            mimg = np.maximum(0, np.minimum(1, mimg))
+                src = None
+
+            mimg = np.zeros((self.Ly, self.Lx), np.float32)
+            if src is not None:
+                mimg1 = np.percentile(src, 1)
+                mimg99 = np.percentile(src, 99)
+                if mimg99 > mimg1:
+                    src = (src - mimg1) / (mimg99 - mimg1)
+                else:
+                    src = np.zeros_like(src)
+                src = np.clip(src, 0, 1).astype(np.float32)
+                if src.shape[0] == self.Ly and src.shape[1] == self.Lx:
+                    mimg = src
+                else:
+                    mimg[yr, xr] = src
+
             masked_image[:, :, :, i] = self.create_masks_of_cells(mimg)
 
         return masked_image
@@ -373,6 +401,9 @@ def keyPressEvent(self, event):
             elif event.key() == QtCore.Qt.Key_T:
                 self.viewbtns.button(3).setChecked(True)
                 self.viewbtns.button(3).press(self, 3)
+            elif event.key() == QtCore.Qt.Key_Y and self.has_chan2:
+                self.viewbtns.button(4).setChecked(True)
+                self.viewbtns.button(4).press(self, 4)
 
     def add_ROI(self, pos=None):
         self.iROI = len(self.ROIs)
@@ -543,13 +574,13 @@ def remove(self, parent):
         parent.win.show()
         parent.show()
 
-    def rotate_ROI(self, parent, ellipse, xrange, yrange, posx, posy):
+    def rotate_ROI(self, parent, ellipse, xrange, yrange, center_x, center_y):
         #Rotates ROI depending on Rotatehandle degree
         ellipse = rotate(ellipse, angle=math.floor(self.ROI.angle()), order=0)
         ellipse = np.flip(ellipse, axis=0)
-        xrange = (np.arange(-1 * int(ellipse.shape[1] - 1), 1) + int(posx)).astype(np.int32)
-        yrange = (np.arange(-1 * int(ellipse.shape[0] - 1), 1) + int(posy)).astype(np.int32)
-        yrange += int(np.floor(ellipse.shape[0] / 2)) + 1
+        w, h = ellipse.shape[1], ellipse.shape[0]
+        xrange = (np.arange(-(w // 2), w - w // 2) + int(center_x)).astype(np.int32)
+        yrange = (np.arange(-(h // 2), h - h // 2) + int(center_y)).astype(np.int32)
         return ellipse, xrange, yrange
 
     def position(self, parent):
@@ -565,13 +596,16 @@ def position(self, parent):
         yrange = (np.arange(-1 * int(sizey), 1) + int(posy)).astype(np.int32)
         yrange += int(np.floor(sizey / 2)) + 1
         # what is ellipse circling?
-        br = self.ROI.boundingRect()
         ellipse = np.zeros((yrange.size, xrange.size), "bool")
         x, y = np.meshgrid(np.arange(0, xrange.size, 1), np.arange(0, yrange.size, 1))
         ellipse = ((y - br.center().y())**2 / (br.height() / 2)**2 +
                    (x - br.center().x())**2 / (br.width() / 2)**2) <= 1
+        center_scene = self.ROI.mapToScene(br.center())
+        center_view = parent.p0.mapSceneToView(center_scene)
+        center_x = center_view.x()
+        center_y = center_view.y()
         if self.ROI.angle() not in (0, 180, -180):
-            ellipse, xrange, yrange = self.rotate_ROI(parent, ellipse, xrange, yrange, posx, posy)
+            ellipse, xrange, yrange = self.rotate_ROI(parent, ellipse, xrange, yrange, center_x, center_y)
         #ensures that ROI is not placed outside of movie coordinates
         ellipse = ellipse[:, np.logical_and(xrange >= 0, xrange < parent.Lx)]
         xrange = xrange[np.logical_and(xrange >= 0, xrange < parent.Lx)]