diff --git a/1_identify_fibers.py b/1_identify_fibers.py
index f33076fda4f01cead2b932bec5c3b6dd76c9dfc6..186b49230ba18181cd03f60525ba4333e16e31b6 100755
--- a/1_identify_fibers.py
+++ b/1_identify_fibers.py
@@ -1,57 +1,78 @@
+# @ String (visibility=MESSAGE, value="<html><b> Welcome to Myosoft - identify fibers! </b></html>") msg1
+# @ File (label="Select folder with your images", description="select folder with your images", style="directory") src_dir
+# @ String(label="Extension for the images to look for", value="czi") filename_filter
+# @ File (label="Select directory for output", style="directory") output_dir
+# @ File(label="Cellpose environment folder", style="directory", description="Folder with the cellpose env") cellpose_dir
+# @ Boolean (label="close image after processing", description="tick this box when using batch mode", value=False) close_raw
+# @ String (visibility=MESSAGE, value="<html><b> Morphometric Gates </b></html>") msg2
+# @ Integer (label="Min Area [um²]", value=10) minAr
+# @ Integer (label="Max Area [um²]", value=6000) maxAr
+# @ Double (label="Min Circularity", value=0.5) minCir
+# @ Double (label="Max Circularity", value=1) maxCir
+# @ Integer (label="Min perimeter [um]", value=5) minPer
+# @ Integer (label="Max perimeter [um]", value=300) maxPer
+# @ String (visibility=MESSAGE, value="<html><b> Expand ROIS to match fibers </b></html>") msg3
+# @ Double (label="ROI expansion [microns]", value=1) enlarge_radius
+# @ String (visibility=MESSAGE, value="<html><b> channel positions in the hyperstack </b></html>") msg5
+# @ Integer (label="Membrane staining channel number", style="slider", min=1, max=5, value=1) membrane_channel
+# @ Integer (label="Fiber staining (MHC) channel number (0=skip)", style="slider", min=0, max=5, value=3) fiber_channel
+# @ Integer (label="minimum fiber intensity (0=auto)", description="0 = automatic threshold detection", value=0) min_fiber_intensity
+# @ CommandService command
+
+# @ RoiManager rm
+# @ ResultsTable rt
+
+
# this is a python rewrite of the original ijm published at
# https://github.com/Hyojung-Choo/Myosoft/blob/Myosoft-hub/Scripts/central%20nuclei%20counter.ijm
+# ─── Requirements ─────────────────────────────────────────────────────────────
+
+# List of update sites needed for the code
+# * TrackMate-Cellpose
+# * IMCF
+# * PTBIOP
+# * CLIJ-CLIJ2
+
+# ─── Imports ──────────────────────────────────────────────────────────────────
+
# IJ imports
# TODO: are the imports RoiManager and ResultsTable needed when using the services?
-from ij import IJ, WindowManager as wm
-from ij.plugin import Duplicator, RoiEnlarger, RoiScaler
-from trainableSegmentation import WekaSegmentation
-from de.biovoxxel.toolbox import Extended_Particle_Analyzer
-from ij.measure import ResultsTable
+import os
+import sys
# Bio-formats imports
-from loci.plugins import BF
-from loci.plugins.in import ImporterOptions
-
+# from loci.plugins import BF
+# from loci.plugins.in import ImporterOptions
# python imports
import time
-import os
-#@ String (visibility=MESSAGE, value="<html><b> Welcome to Myosoft - identify fibers! </b></html>") msg1
-#@ File (label="Select directory with classifiers", style="directory") classifiers_dir
-#@ File (label="Select directory for output", style="directory") output_dir
-#@ File (label="Select image file", description="select your image") path_to_image
-#@ Boolean (label="close image after processing", description="tick this box when using batch mode", value=False) close_raw
-#@ String (visibility=MESSAGE, value="<html><b> Morphometric Gates </b></html>") msg2
-#@ Integer (label="Min Area [um²]", value=10) minAr
-#@ Integer (label="Max Area [um²]", value=6000) maxAr
-#@ Float (label="Min Circularity", value=0.5) minCir
-#@ Float (label="Max Circularity", value=1) maxCir
-#@ Float (label="Min solidity", value=0.0) minSol
-#@ Float (label="Max solidity", value=1) maxSol
-#@ Integer (label="Min perimeter [um]", value=5) minPer
-#@ Integer (label="Max perimeter [um]", value=300) maxPer
-#@ Integer (label="Min min ferret [um]", value=0.1) minMinFer
-#@ Integer (label="Max min ferret [um]", value=100) maxMinFer
-#@ Integer (label="Min ferret AR", value=0) minFAR
-#@ Integer (label="Max ferret AR", value=8) maxFAR
-#@ Float (label="Min roundess", value=0.2) minRnd
-#@ Float (label="Max roundess", value=1) maxRnd
-#@ String (visibility=MESSAGE, value="<html><b> Expand ROIS to match fibers </b></html>") msg3
-#@ Float (label="ROI expansion [microns]", value=1) enlarge
-#@ String (visibility=MESSAGE, value="<html><b> channel positions in the hyperstack </b></html>") msg5
-#@ Integer (label="Membrane staining channel number", style="slider", min=1, max=5, value=1) membrane_channel
-#@ Integer (label="Fiber staining (MHC) channel number (0=skip)", style="slider", min=0, max=5, value=3) fiber_channel
-#@ Integer (label="minimum fiber intensity (0=auto)", description="0 = automatic threshold detection", value=0) min_fiber_intensity
-#@ Integer (label="sub-tiling to economize RAM", style="slider", min=1, max=8, value=4) tiling_factor
-
-#@ RoiManager rm
-#@ ResultsTable rt
+from ch.epfl.biop.ij2command import Labels2CompositeRois
+
+# TrackMate imports
+from fiji.plugin.trackmate import Logger, Model, Settings, TrackMate
+from fiji.plugin.trackmate.action import LabelImgExporter
+from fiji.plugin.trackmate.cellpose import CellposeDetectorFactory
+from fiji.plugin.trackmate.cellpose.CellposeSettings import PretrainedModel
+from fiji.plugin.trackmate.features import FeatureFilter
+from fiji.plugin.trackmate.providers import (
+ SpotAnalyzerProvider,
+ SpotMorphologyAnalyzerProvider,
+)
+from fiji.plugin.trackmate.tracking.jaqaman import SparseLAPTrackerFactory
+from ij import IJ
+from ij import WindowManager as wm
+from ij.measure import ResultsTable
+from ij.plugin import Duplicator, ImageCalculator, RoiEnlarger
+from imcflibs import pathtools
+from imcflibs.imagej import bioformats as bf
+from imcflibs.imagej import misc
+
+# ─── Functions ────────────────────────────────────────────────────────────────
def fix_ij_options():
- """put IJ into a defined state
- """
+ """Put IJ into a defined state."""
# disable inverting LUT
IJ.run("Appearance...", " menu=0 16-bit=Automatic")
# set foreground color to be white, background black
@@ -68,7 +89,7 @@ def fix_ij_options():
def fix_ij_dirs(path):
- """use forward slashes in directory paths
+ """use forward slashes in directory paths.
Parameters
----------
@@ -87,61 +108,72 @@ def fix_ij_dirs(path):
return fixed_path
-def open_image_with_BF(path_to_file):
- """ use Bio-Formats to opens the first image from an image file path
+def fix_BF_czi_imagetitle(imp):
+ """Fix the title of an image read using the bio-formats importer.
+
+ The title is modified to remove the ".czi" extension and replace
+ spaces with underscores.
Parameters
----------
- path_to_file : string
- path to the image file
+ imp : ij.ImagePlus
+ The image to be processed.
Returns
-------
- ImagePlus
- the first imp stored in a give file
+ string
+ The modified title of the image.
"""
- options = ImporterOptions()
- options.setColorMode(ImporterOptions.COLOR_MODE_GRAYSCALE)
- options.setAutoscale(True)
- options.setId(path_to_file)
- imps = BF.openImagePlus(options) # is an array of ImagePlus
-
- return imps[0]
-
-
-def fix_BF_czi_imagetitle(imp):
- image_title = os.path.basename( imp.getShortTitle() )
+ image_title = os.path.basename(imp.getShortTitle())
+ # remove the ".czi" extension
image_title = image_title.replace(".czi", "")
+ # replace spaces with underscores
image_title = image_title.replace(" ", "_")
+ # remove any double underscores
image_title = image_title.replace("_-_", "")
+ # remove any double underscores
image_title = image_title.replace("__", "_")
+ # remove any "#" characters
image_title = image_title.replace("#", "Series")
return image_title
-def preprocess_membrane_channel(imp):
- """apply myosoft pre-processing steps for the membrane channel
+def do_background_correction(imp, gaussian_radius=20):
+ """Perform background correction on an image.
+
+ This is done by applying a Gaussian blur to the image and then dividing the
+ original image by the blurred image.
Parameters
----------
- imp : ImagePlus
- a single channel image of the membrane staining
+ imp : ij.ImagePlus
+ The image to be corrected.
+ gaussian_radius : int
+ The radius of the Gaussian filter to be used. Default value is 20.
+
+ Returns
+ -------
+ ij.ImagePlus
+ The background-corrected image.
"""
- IJ.run(imp, "Enhance Contrast", "saturated=0.35")
- IJ.run(imp, "Apply LUT", "")
- IJ.run(imp, "Enhance Contrast", "saturated=1")
- IJ.run(imp, "8-bit", "")
- IJ.run(imp, "Invert", "")
- IJ.run(imp, "Convolve...", "text1=[-1.0 -1.0 -1.0 -1.0 -1.0\n-1.0 -1.0 -1.0 -1.0 0\n-1.0 -1.0 24.0 -1.0 -1.0\n-1.0 -1.0 -1.0 -1.0 -1.0\n-1.0 -1.0 -1.0 -1.0 0] normalize")
+ imp_bgd = imp.duplicate()
+ IJ.run(
+ imp_bgd,
+ "Gaussian Blur...",
+ "sigma=" + str(gaussian_radius) + " scaled",
+ )
+ return ImageCalculator.run(imp, imp_bgd, "Divide create 32-bit")
def get_threshold_from_method(imp, channel, method):
- """returns the threshold value of chosen IJ AutoThreshold method in desired channel
+ """Get the value of automated threshold method.
+
+ Returns the threshold value of chosen IJ AutoThreshold method in desired channel.
Parameters
----------
- imp : ImagePlus
+ imp : ij.ImagePlus
the imp from which to get the threshold value
channel : integer
the channel in which to get the treshold
@@ -153,7 +185,7 @@ def get_threshold_from_method(imp, channel, method):
list
the upper and the lower threshold (integer values)
"""
- imp.setC(channel) # starts at 1
+ imp.setC(channel) # starts at 1
ip = imp.getProcessor()
ip.setAutoThreshold(method + " dark")
lower_thr = ip.getMinThreshold()
@@ -163,52 +195,194 @@ def get_threshold_from_method(imp, channel, method):
return lower_thr, upper_thr
-def apply_weka_model(model_path, imp, tiles_per_dim):
- """apply a pretrained WEKA model to an ImagePlus
+def run_tm(
+ implus,
+ channel_seg,
+ cellpose_env,
+ seg_model,
+ diam_seg,
+ channel_sec=0,
+ quality_thresh=[0, 0],
+ intensity_thresh=[0, 0],
+ circularity_thresh=[0, 0],
+ perimeter_thresh=[0, 0],
+ area_thresh=[0, 0],
+ crop_roi=None,
+ use_gpu=True,
+):
+ """
+ Function to run TrackMate on open data, applying filters to spots.
Parameters
----------
- model_path : string
- path to the model file
- imp : ImagePlus
- ImagePlus to apply the model to
- tiles_per_dim : integer
- tiles the imp to save RAM
+ implus : ij.ImagePlus
+ ImagePlus on which to run the function
+ channel_seg : int
+ Channel of interest
+ cellpose_env : str
+ Path to the cellpose environment
+ seg_model : PretrainedModel
+ Model to use for the segmentation
+ diam_seg : float
+ Diameter to use for segmentation
+ channel_sec : int, optional
+ Secondary channel to use for segmentation, by default 0
+ quality_thresh : float, optional
+ Threshold for quality filtering, by default None
+ intensity_thresh : float, optional
+ Threshold for intensity filtering, by default None
+ circularity_thresh : float, optional
+ Threshold for circularity filtering, by default None
+ perimeter_thresh : float, optional
+ Threshold for perimeter filtering, by default None
+ area_thresh : float, optional
+ Threshold for area filtering, by default None
+ crop_roi : ROI, optional
+ ROI to crop on the image, by default None
+ use_gpu : bool, optional
+ Boolean to use GPU or not, by default True
Returns
-------
- ImagePlus
- the result of the WEKA segmentation. One channel per class.
+ ij.ImagePlus
+ Label image with the segmented objects
"""
- segmentator = WekaSegmentation()
- segmentator.loadClassifier( model_path )
- result = segmentator.applyClassifier( imp, [tiles_per_dim, tiles_per_dim], 0, True ) #ImagePlus imp, int[x,y,z] tilesPerDim, int numThreads (0=all), boolean probabilityMaps
-
- return result
-
-def process_weka_result(imp):
- """apply myosoft pre-processing steps for the imp after WEKA classification to prepare it
- for ROI detection with the extended particle analyzer
+ # Get image dimensions and calibration
+ dims = implus.getDimensions()
+ cal = implus.getCalibration()
+
+ # If the image has more than one slice, adjust the dimensions
+ if implus.getNSlices() > 1:
+ implus.setDimensions(dims[2], dims[4], dims[3])
+
+ # Set ROI if provided
+ if crop_roi is not None:
+ implus.setRoi(crop_roi)
+
+ # Initialize TrackMate model
+ model = Model()
+ model.setLogger(Logger.IJTOOLBAR_LOGGER)
+
+ # Prepare settings for TrackMate
+ settings = Settings(implus)
+ settings.detectorFactory = CellposeDetectorFactory()
+
+ # Configure detector settings
+ settings.detectorSettings["TARGET_CHANNEL"] = channel_seg
+ settings.detectorSettings["OPTIONAL_CHANNEL_2"] = channel_sec
+ settings.detectorSettings["CELLPOSE_PYTHON_FILEPATH"] = os.path.join(
+ cellpose_env, "python.exe"
+ )
+ settings.detectorSettings["CELLPOSE_MODEL_FILEPATH"] = os.path.join(
+ os.environ["USERPROFILE"], ".cellpose", "models"
+ )
+ settings.detectorSettings["CELLPOSE_MODEL"] = seg_model
+ settings.detectorSettings["CELL_DIAMETER"] = diam_seg
+ settings.detectorSettings["USE_GPU"] = use_gpu
+ settings.detectorSettings["SIMPLIFY_CONTOURS"] = True
+
+ settings.initialSpotFilterValue = -1.0
+
+ # Add spot analyzers
+ spotAnalyzerProvider = SpotAnalyzerProvider(1)
+ spotMorphologyProvider = SpotMorphologyAnalyzerProvider(1)
+
+ for key in spotAnalyzerProvider.getKeys():
+ settings.addSpotAnalyzerFactory(spotAnalyzerProvider.getFactory(key))
+
+ for key in spotMorphologyProvider.getKeys():
+ settings.addSpotAnalyzerFactory(spotMorphologyProvider.getFactory(key))
+
+ # Apply spot filters based on thresholds
+ if any(quality_thresh):
+ settings = set_trackmate_filter(settings, "QUALITY", quality_thresh)
+ if any(intensity_thresh):
+ settings = set_trackmate_filter(
+ settings, "MEAN_INTENSITY_CH" + str(channel_seg), intensity_thresh
+ )
+ if any(circularity_thresh):
+ settings = set_trackmate_filter(settings, "CIRCULARITY", circularity_thresh)
+ if any(area_thresh):
+ settings = set_trackmate_filter(settings, "AREA", area_thresh)
+ if any(perimeter_thresh):
+ settings = set_trackmate_filter(settings, "PERIMETER", perimeter_thresh)
+
+ # print(settings)
+
+ # Configure tracker
+ settings.trackerFactory = SparseLAPTrackerFactory()
+ settings.trackerSettings = settings.trackerFactory.getDefaultSettings()
+ # settings.addTrackAnalyzer(TrackDurationAnalyzer())
+ settings.trackerSettings["LINKING_MAX_DISTANCE"] = 3.0
+ settings.trackerSettings["GAP_CLOSING_MAX_DISTANCE"] = 3.0
+ settings.trackerSettings["MAX_FRAME_GAP"] = 2
+
+ # Initialize TrackMate with model and settings
+ trackmate = TrackMate(model, settings)
+ trackmate.computeSpotFeatures(True)
+ trackmate.computeTrackFeatures(False)
+
+ # Check input validity
+ if not trackmate.checkInput():
+ sys.exit(str(trackmate.getErrorMessage()))
+ return
+
+ # Process the data
+ if not trackmate.process():
+ if "[SparseLAPTracker] The spot collection is empty." in str(
+ trackmate.getErrorMessage()
+ ):
+ return IJ.createImage(
+ "Untitled",
+ "8-bit black",
+ implus.getWidth(),
+ implus.getHeight(),
+ implus.getNFrames(),
+ )
+ else:
+ sys.exit(str(trackmate.getErrorMessage()))
+ return
+
+ # Export the label image
+ # sm = SelectionModel(model)
+ exportSpotsAsDots = False
+ exportTracksOnly = False
+ label_imp = LabelImgExporter.createLabelImagePlus(
+ trackmate, exportSpotsAsDots, exportTracksOnly, False
+ )
+ label_imp.setDimensions(1, dims[3], dims[4])
+ label_imp.setCalibration(cal)
+ implus.setDimensions(dims[2], dims[3], dims[4])
+ return label_imp
+
+
+def set_trackmate_filter(settings, filter_name, filter_value):
+ """Sets a TrackMate spot filter with specified filter name and values.
Parameters
----------
- imp : ImagePlus
- a single channel (= desired class) of the WEKA classification result imp
+ settings : Settings
+ TrackMate settings object to which the filter will be added.
+ filter_name : str
+ The name of the filter to be applied.
+ filter_value : list
+ A list containing two values for the filter. The first value is
+ applied as an above-threshold filter, and the second as a below-threshold filter.
"""
- IJ.run(imp, "8-bit", "")
- IJ.run(imp, "Median...", "radius=3")
- IJ.run(imp, "Gaussian Blur...", "sigma=2")
- IJ.run(imp, "Auto Threshold", "method=MaxEntropy")
- IJ.run(imp, "Invert", "")
+ filter = FeatureFilter(filter_name, filter_value[0], True)
+ settings.addSpotFilter(filter)
+ filter = FeatureFilter(filter_name, filter_value[1], False)
+ settings.addSpotFilter(filter)
+ return settings
def delete_channel(imp, channel_number):
- """delete a channel from target imp
+ """Delete a channel from target imp.
Parameters
----------
- imp : ImagePlus
+ imp : ij.ImagePlus
the imp from which to delete target channel
channel_number : integer
the channel number to be deleted. starts at 0.
@@ -217,59 +391,14 @@ def delete_channel(imp, channel_number):
IJ.run(imp, "Delete Slice", "delete=channel")
-def run_extended_particle_analyzer( imp, eda_parameters ):
- """identifies ROIs in target imp using the extended particle analyzer of the BioVoxxel toolbox
- with given parameters
+def measure_in_all_rois(imp, channel, rm):
+ """Gives measurements for all ROIs in ROIManager.
- Parameters
- ----------
- imp : ImagePlus
- the image on which to run the EPA on. Should be 8-bit thresholded
- eda_parameters : array
- all user defined parameters to restrict ROI identification
- """
- epa = Extended_Particle_Analyzer()
- epa.readInputImageParameters(imp)
- epa.setDefaultParameterFields()
-
- # expose all parameters explicitly
- epa.usePixel = False
- epa.usePixelForOutput = False
- epa.Area = str(eda_parameters[0]) + "-" + str(eda_parameters[1])
- epa.Extent = "0.00-1.00"
- epa.Perimeter = str(eda_parameters[2]) + "-" + str(eda_parameters[3])
- epa.Circularity = str(eda_parameters[4]) + "-" + str(eda_parameters[5])
- epa.Roundness = str(eda_parameters[6]) + "-" + str(eda_parameters[7])
- epa.Solidity = str(eda_parameters[8]) + "-" + str(eda_parameters[9])
- epa.Compactness = "0.00-1.00"
- epa.AR = "0-Infinity"
- epa.FeretAR = str(eda_parameters[10]) + "-" + str(eda_parameters[11])
- epa.EllipsoidAngle = "0-180"
- epa.MaxFeret = "0-Infinity"
- epa.MinFeret = str(eda_parameters[12]) + "-" + str(eda_parameters[13])
- epa.FeretAngle = "0-180"
- epa.COV = "0.00-1.00"
- epa.Output = "Nothing"
- epa.Redirect = "None"
- epa.Correction = "None"
- epa.Reset = False
- epa.DisplayResults = False
- epa.ClearResults = False
- epa.Summarize = False
- epa.AddToManager = True
- epa.ExcludeEdges = False
- epa.IncludeHoles = False
-
- epa.defineParticleAnalyzers()
- epa.particleAnalysis( imp.getProcessor(), imp, imp.getTitle() )
-
-
-def measure_in_all_rois( imp, channel, rm ):
- """measures in all ROIS on a given channel of imp all parameters that are set in IJ "Set Measurements"
+ Measures in all ROIS on a given channel of imp all parameters that are set in IJ "Set Measurements".
Parameters
----------
- imp : ImagePlus
+ imp : ij.ImagePlus
the imp to measure on
channel : integer
the channel to measure in. starts at 1.
@@ -277,12 +406,12 @@ def measure_in_all_rois( imp, channel, rm ):
a reference of the IJ-RoiManager
"""
imp.setC(channel)
- rm.runCommand(imp,"Deselect")
- rm.runCommand(imp,"Measure")
+ rm.runCommand(imp, "Deselect")
+ rm.runCommand(imp, "Measure")
-def change_all_roi_color( rm, color ):
- """change the color of all ROIs in the RoiManager
+def change_all_roi_color(rm, color):
+ """Cchange the color of all ROIs in the RoiManager.
Parameters
----------
@@ -292,13 +421,13 @@ def change_all_roi_color( rm, color ):
the desired color. e.g. "green", "red", "yellow", "magenta" ...
"""
number_of_rois = rm.getCount()
- for roi in range( number_of_rois ):
+ for roi in range(number_of_rois):
rm.select(roi)
rm.runCommand("Set Color", color)
-def change_subset_roi_color( rm, selected_rois, color ):
- """change the color of selected ROIs in the RoiManager
+def change_subset_roi_color(rm, selected_rois, color):
+ """Change the color of selected ROIs in the RoiManager.
Parameters
----------
@@ -316,21 +445,21 @@ def change_subset_roi_color( rm, selected_rois, color ):
def show_all_rois_on_image(rm, imp):
- """shows all ROIs in the ROiManager on imp
+ """Shows all ROIs in the ROiManager on imp.
Parameters
----------
rm : RoiManager
a reference of the IJ-RoiManager
- imp : ImagePlus
+ imp : ij.ImagePlus
the imp on which to show the ROIs
"""
imp.show()
- rm.runCommand(imp,"Show All")
+ rm.runCommand(imp, "Show All")
def save_all_rois(rm, target):
- """save all ROIs in the RoiManager as zip to target path
+ """Save all ROIs in the RoiManager as zip to target path.
Parameters
----------
@@ -342,8 +471,8 @@ def save_all_rois(rm, target):
rm.runCommand("Save", target)
-def save_selected_rois( rm, selected_rois, target ):
- """save selected ROIs in the RoiManager as zip to target path
+def save_selected_rois(rm, selected_rois, target):
+ """Save selected ROIs in the RoiManager as zip to target path.
Parameters
----------
@@ -360,8 +489,8 @@ def save_selected_rois( rm, selected_rois, target ):
rm.runCommand("Deselect")
-def enlarge_all_rois( amount_in_um, rm, pixel_size_in_um ):
- """enlarges all ROIs in the RoiManager by x scaled units
+def enlarge_all_rois(amount_in_um, rm, pixel_size_in_um):
+ """Enlarges all ROIs in the RoiManager by x scaled units.
Parameters
----------
@@ -380,13 +509,15 @@ def enlarge_all_rois( amount_in_um, rm, pixel_size_in_um ):
rm.addRoi(enlarged_roi)
-def select_positive_fibers( imp, channel, rm, min_intensity ):
- """For all ROIs in the RoiManager, select ROIs based on intensity measurement in given channel of imp.
+def select_positive_fibers(imp, channel, rm, min_intensity):
+ """Select ROIs in ROIManager based on intensity in specific channel.
+
+ For all ROIs in the RoiManager, select ROIs based on intensity measurement in given channel of imp.
See https://imagej.nih.gov/ij/developer/api/ij/process/ImageStatistics.html
Parameters
----------
- imp : ImagePlus
+ imp : ij.ImagePlus
the imp on which to measure
channel : integer
the channel on which to measure. starts at 1
@@ -412,8 +543,10 @@ def select_positive_fibers( imp, channel, rm, min_intensity ):
return selected_rois
-def preset_results_column( rt, column, value):
- """pre-set all rows in given column of the IJ-ResultsTable with desired value
+def preset_results_column(rt, column, value):
+ """Pre-set values in selected column from the ResultsTable.
+
+ Pre-set all rows in given column of the IJ-ResultsTable with desired value.
Parameters
----------
@@ -424,14 +557,14 @@ def preset_results_column( rt, column, value):
value : string or float or integer
the value to be set
"""
- for i in range( rt.size() ):
+ for i in range(rt.size()):
rt.setValue(column, i, value)
rt.show("Results")
-def add_results( rt, column, row, value ):
- """adds a value in desired rows of a given column
+def add_results(rt, column, row, value):
+ """Adds a value in desired rows of a given column.
Parameters
----------
@@ -444,27 +577,27 @@ def add_results( rt, column, row, value ):
value : string or float or integer
the value to be set
"""
- for i in range( len( row ) ):
+ for i in range(len(row)):
rt.setValue(column, row[i], value)
rt.show("Results")
-def enhance_contrast( imp ):
- """use "Auto" Contrast & Brightness settings in each channel of imp
+def enhance_contrast(imp):
+ """Use "Auto" Contrast & Brightness settings in each channel of imp.
Parameters
----------
- imp : ImagePlus
+ imp : ij.ImagePlus
the imp on which to change C&B
"""
- for channel in range( imp.getDimensions()[2] ):
- imp.setC(channel + 1) # IJ channels start at 1
+ for channel in range(imp.getDimensions()[2]):
+ imp.setC(channel + 1) # IJ channels start at 1
IJ.run(imp, "Enhance Contrast", "saturated=0.35")
def renumber_rois(rm):
- """rename all ROIs in the RoiManager according to their number
+ """Rename all ROIs in the RoiManager according to their number.
Parameters
----------
@@ -472,12 +605,12 @@ def renumber_rois(rm):
a reference of the IJ-RoiManager
"""
number_of_rois = rm.getCount()
- for roi in range( number_of_rois ):
- rm.rename( roi, str(roi + 1) )
+ for roi in range(number_of_rois):
+ rm.rename(roi, str(roi + 1))
def setup_defined_ij(rm, rt):
- """set up a clean and defined Fiji user environment
+ """Set up a clean and defined Fiji user environment.
Parameters
----------
@@ -487,124 +620,165 @@ def setup_defined_ij(rm, rt):
a reference of the IJ-ResultsTable
"""
fix_ij_options()
- rm.runCommand('reset')
+ rm.runCommand("reset")
rt.reset()
IJ.log("\\Clear")
-execution_start_time = time.time()
-
-setup_defined_ij(rm, rt)
-
-print rt.size()
-
-# open image using Bio-Formats
-path_to_image = fix_ij_dirs(path_to_image)
-raw = open_image_with_BF(path_to_image)
-
-# get image info
-raw_image_calibration = raw.getCalibration()
-raw_image_title = fix_BF_czi_imagetitle(raw)
-print("raw image title: ", str(raw_image_title))
-
-# take care of paths and directories
-output_dir = fix_ij_dirs(output_dir) + "/" + str(raw_image_title) + "/1_identify_fibers"
-print("output_dir: ", str(output_dir))
-
-if not os.path.exists( str(output_dir) ):
- os.makedirs( str(output_dir) )
-
-classifiers_dir = fix_ij_dirs(classifiers_dir)
-primary_model = classifiers_dir + "/" + "primary.model"
-secondary_model = classifiers_dir + "/" + "secondary_central_nuclei.model"
-
-# update the log for the user
-IJ.log( "Now working on " + str(raw_image_title) )
-if raw_image_calibration.scaled() == False:
- IJ.log("Your image is not spatially calibrated! Size measurements are only possible in [px].")
-IJ.log( " -- settings used -- ")
-IJ.log( "area = " + str(minAr) + "-" + str(maxAr) )
-IJ.log( "perimeter = " + str(minPer) + "-" + str(maxPer) )
-IJ.log( "circularity = " + str(minCir) + "-" + str(maxCir) )
-IJ.log( "roundness = " + str(minRnd) + "-" + str(maxRnd) )
-IJ.log( "solidity = " + str(minSol) + "-" + str(maxSol) )
-IJ.log( "feret_ar = " + str(minFAR) + "-" + str(maxFAR) )
-IJ.log( "min_feret = " + str(minMinFer) + "-" + str(maxMinFer) )
-IJ.log( "ROI expansion [microns] = " + str(enlarge) )
-IJ.log( "Membrane channel = " + str(membrane_channel) )
-IJ.log( "MHC positive fiber channel = " + str(fiber_channel) )
-IJ.log( "sub-tiling = " + str(tiling_factor) )
-IJ.log( " -- settings used -- ")
-
-# image (pre)processing and segmentation (-> ROIs)
-membrane = Duplicator().run(raw, membrane_channel, membrane_channel, 1, 1, 1, 1) # imp, firstC, lastC, firstZ, lastZ, firstT, lastT
-preprocess_membrane_channel(membrane)
-weka_result1 = apply_weka_model(primary_model, membrane, tiling_factor )
-delete_channel(weka_result1, 1)
-weka_result2 = apply_weka_model(secondary_model, weka_result1, tiling_factor )
-delete_channel(weka_result2, 1)
-weka_result2.setCalibration(raw_image_calibration)
-process_weka_result(weka_result2)
-IJ.saveAs(weka_result2, "Tiff", output_dir + "/" + raw_image_title + "_all_fibers_binary")
-eda_parameters = [minAr, maxAr, minPer, maxPer, minCir, maxCir, minRnd, maxRnd, minSol, maxSol, minFAR, maxFAR, minMinFer, maxMinFer]
-raw.show() # EPA will not work if no image is shown
-run_extended_particle_analyzer(weka_result2, eda_parameters)
-
-# modify rois
-rm.hide()
-raw.hide()
-enlarge_all_rois( enlarge, rm, raw_image_calibration.pixelWidth )
-renumber_rois(rm)
-save_all_rois( rm, output_dir + "/" + raw_image_title + "_all_fiber_rois.zip" )
-
-# check for positive fibers
-if fiber_channel > 0:
- if min_fiber_intensity == 0:
- min_fiber_intensity = get_threshold_from_method(raw, fiber_channel, "Mean")[0]
- IJ.log( "automatic intensity threshold detection: True" )
-
- IJ.log( "fiber intensity threshold: " + str(min_fiber_intensity) )
- change_all_roi_color(rm, "blue")
- positive_fibers = select_positive_fibers( raw, fiber_channel, rm, min_fiber_intensity )
- change_subset_roi_color(rm, positive_fibers, "magenta")
- save_selected_rois( rm, positive_fibers, output_dir + "/" + raw_image_title + "_mhc_positive_fiber_rois.zip")
-
-# measure size & shape, save
-IJ.run("Set Measurements...", "area perimeter shape feret's redirect=None decimal=4")
-IJ.run("Clear Results", "")
-measure_in_all_rois( raw, membrane_channel, rm )
-
-rt = ResultsTable.getResultsTable("Results")
-
-print rt.size()
-
-if fiber_channel > 0:
- print rt.size()
- preset_results_column( rt, "MHC Positive Fibers (magenta)", "NO" )
- print rt.size()
- add_results( rt, "MHC Positive Fibers (magenta)", positive_fibers, "YES")
- print rt.size()
-
-rt.save(output_dir + "/" + raw_image_title + "_all_fibers_results.csv")
-print "saved the all_fibers_results.csv"
-# dress up the original image, save a overlay-png, present original to the user
-rm.show()
-raw.show()
-show_all_rois_on_image( rm, raw )
-raw.setDisplayMode(IJ.COMPOSITE)
-enhance_contrast( raw )
-IJ.run("From ROI Manager", "") # ROIs -> overlays so they show up in the saved png
-qc_duplicate = raw.duplicate()
-IJ.saveAs(qc_duplicate, "PNG", output_dir + "/" + raw_image_title + "_all_fibers")
-qc_duplicate.close()
-wm.toFront( raw.getWindow() )
-IJ.run("Remove Overlay", "")
-raw.setDisplayMode(IJ.GRAYSCALE)
-show_all_rois_on_image( rm, raw )
-total_execution_time_min = (time.time() - execution_start_time) / 60.0
-IJ.log("total time in minutes: " + str(total_execution_time_min))
-IJ.log( "~~ all done ~~" )
-IJ.selectWindow("Log")
-IJ.saveAs("Text", str(output_dir + "/" + raw_image_title + "_all_fibers_Log"))
-if close_raw == True:
- raw.close()
\ No newline at end of file
+# ─── Main Code ────────────────────────────────────────────────────────────────
+
+if __name__ == "__main__":
+ execution_start_time = time.time()
+ IJ.log("\\Clear")
+ misc.timed_log("Script starting")
+ setup_defined_ij(rm, rt)
+
+ file_list = pathtools.listdir_matching(
+ src_dir.getPath(), filename_filter, fullpath=True
+ )
+
+ out_dir_info = pathtools.parse_path(output_dir)
+
+ for index, file in enumerate(file_list):
+ # open image using Bio-Formats
+ file_info = pathtools.parse_path(file)
+ misc.progressbar(index + 1, len(file_list), 1, "Opening : ")
+ raw = bf.import_image(file_info["full"])[0]
+
+ # get image info
+ raw_image_calibration = raw.getCalibration()
+ raw_image_title = fix_BF_czi_imagetitle(raw)
+ print("raw image title: ", str(raw_image_title))
+
+ # take care of paths and directories
+ output_dir = os.path.join(
+ out_dir_info["full"], str(raw_image_title), "1_identify_fibers"
+ )
+ print("output_dir: ", str(output_dir))
+
+ if not os.path.exists(str(output_dir)):
+ os.makedirs(str(output_dir))
+
+ # update the log for the user
+ misc.timed_log("Now working on " + str(raw_image_title))
+ if raw_image_calibration.scaled() is False:
+ IJ.log(
+ "Your image is not spatially calibrated! Size measurements are only possible in [px]."
+ )
+ # Only print it once since we'll use the same settings everytime
+ if index == 0:
+ IJ.log(" -- settings used -- ")
+ IJ.log("area = " + str(minAr) + "-" + str(maxAr))
+ IJ.log("perimeter = " + str(minPer) + "-" + str(maxPer))
+ IJ.log("circularity = " + str(minCir) + "-" + str(maxCir))
+ IJ.log("ROI expansion [microns] = " + str(enlarge_radius))
+ IJ.log("Membrane channel = " + str(membrane_channel))
+ IJ.log("MHC positive fiber channel = " + str(fiber_channel))
+ # IJ.log("sub-tiling = " + str(tiling_factor))
+ IJ.log(" -- settings used -- ")
+
+ # image (pre)processing and segmentation (-> ROIs)# imp, firstC, lastC, firstZ,
+ # lastZ, firstT, lastT
+ membrane = Duplicator().run(raw, membrane_channel, membrane_channel, 1, 1, 1, 1)
+ imp_bgd_corrected = do_background_correction(membrane)
+ IJ.run("Conversions...", "scale")
+ IJ.run(imp_bgd_corrected, "16-bit", "")
+
+ imp_result = run_tm(
+ imp_bgd_corrected,
+ 1,
+ cellpose_dir.getPath(),
+ PretrainedModel.CYTO2,
+ 30.0,
+ area_thresh=[minAr, maxAr],
+ circularity_thresh=[minCir, maxCir],
+ perimeter_thresh=[minPer, maxPer],
+ )
+ IJ.saveAs(
+ imp_result,
+ "Tiff",
+ os.path.join(output_dir, raw_image_title + "_all_fibers_binary"),
+ )
+
+ command.run(Labels2CompositeRois, True, "rm", rm, "imp", imp_result).get()
+
+ enlarge_all_rois(enlarge_radius, rm, raw_image_calibration.pixelWidth)
+ renumber_rois(rm)
+ save_all_rois(
+ rm, os.path.join(output_dir, raw_image_title + "_all_fiber_rois.zip")
+ )
+
+ # check for positive fibers
+ if fiber_channel > 0:
+ if min_fiber_intensity == 0:
+ min_fiber_intensity = get_threshold_from_method(
+ raw, fiber_channel, "Mean"
+ )[0]
+ IJ.log("automatic intensity threshold detection: True")
+
+ IJ.log("fiber intensity threshold: " + str(min_fiber_intensity))
+ change_all_roi_color(rm, "blue")
+ positive_fibers = select_positive_fibers(
+ raw, fiber_channel, rm, min_fiber_intensity
+ )
+ change_subset_roi_color(rm, positive_fibers, "magenta")
+ save_selected_rois(
+ rm,
+ positive_fibers,
+ os.path.join(
+ output_dir, raw_image_title + "_mhc_positive_fiber_rois.zip"
+ ),
+ )
+
+ # measure size & shape, save
+ IJ.run(
+ "Set Measurements...",
+ "area perimeter shape feret's redirect=None decimal=4",
+ )
+ IJ.run("Clear Results", "")
+ measure_in_all_rois(raw, membrane_channel, rm)
+
+ rt = ResultsTable.getResultsTable("Results")
+
+ # print(rt.size())
+
+ if fiber_channel > 0:
+ # print(rt.size())
+ preset_results_column(rt, "MHC Positive Fibers (magenta)", "NO")
+ # print(rt.size())
+ add_results(rt, "MHC Positive Fibers (magenta)", positive_fibers, "YES")
+ # print(rt.size())
+
+ rt.save(os.path.join(output_dir, raw_image_title + "_all_fibers_results.csv"))
+ # print("saved the all_fibers_results.csv")
+ # dress up the original image, save a overlay-png, present original to the user
+ rm.show()
+ raw.show()
+ show_all_rois_on_image(rm, raw)
+ raw.setDisplayMode(IJ.COMPOSITE)
+ enhance_contrast(raw)
+ IJ.run(
+ "From ROI Manager", ""
+ ) # ROIs -> overlays so they show up in the saved png
+ qc_duplicate = raw.duplicate()
+ IJ.saveAs(
+ qc_duplicate, "PNG", output_dir + "/" + raw_image_title + "_all_fibers"
+ )
+ qc_duplicate.close()
+ wm.toFront(raw.getWindow())
+ IJ.run("Remove Overlay", "")
+ raw.setDisplayMode(IJ.GRAYSCALE)
+ show_all_rois_on_image(rm, raw)
+
+ IJ.selectWindow("Log")
+ IJ.saveAs("Text", str(output_dir + "/" + raw_image_title + "_all_fibers_Log"))
+ membrane.close()
+ imp_bgd_corrected.close()
+ imp_result.close()
+
+ if close_raw == True:
+ raw.close()
+
+ total_execution_time_min = (time.time() - execution_start_time) / 60.0
+ IJ.log("total time in minutes: " + str(total_execution_time_min))
+ IJ.log("~~ all done ~~")
diff --git a/README.md b/README.md
index 81a4b80af5f90a525212f294b2eb855e6be3a186..5b106e1f781e1988516ad1adad9ca4a199793ee0 100644
--- a/README.md
+++ b/README.md
@@ -9,22 +9,22 @@ Original publication: <https://doi.org/10.1371/journal.pone.0229041>
Original code: <https://github.com/Hyojung-Choo/Myosoft/tree/Myosoft-hub>
-## `1_identify_fibers.py`
+## [`1_identify_fibers.py`](1_identify_fibers.py)
-- Will identify all fibers based on the membrane staining using WEKA pixel
- classification, filter them according to the morphometric gates and save the
+- Will identify all fibers based on the membrane staining using [Cellpose](https://github.com/MouseLand/cellpose) segmentation, filter them according to the morphometric gates and save the
corresponding ROIs.
-- Will now also save the WEKA segmentation as a binary so it can be edited
+ - Need to be installed ont the machine where the script is run. Follow [this guide](https://wiki.biozentrum.unibas.ch/display/IMCF/Cellpose+python+environment) to create the environment.
+- Will now also save the Cellpose segmentation as a binary so it can be edited
manually. If you do so, you need to run the "extended particle analyzer"
manually as well to choose & apply the morphometric gates.
- Can be run in batch.
-## `2a_identify_MHC_positive_fibers.py`
+## [`2a_identify_MHC_positive_fibers.py`](2a_identify_MHC_positive_fibers.py)
- Allows to manual re-run the MHC positive fiber detection. Useful in case you
would like to re-run detection with a manual threshold for an image.
-## `2b_central_nuclei_counter.py`
+## [`2b_central_nuclei_counter.py`](2b_central_nuclei_counter.py)
- Will identify centralized nuclei given a ROI-zip together with its
corresponding image.
@@ -32,14 +32,14 @@ Original code: <https://github.com/Hyojung-Choo/Myosoft/tree/Myosoft-hub>
information of a MHC staining channel.
- The ROI color code is annotated in the results table.
-## `2c_fibertyping.py`
+## [`2c_fibertyping.py`](2c_fibertyping.py)
- Identifies positive fibers in up to 3 channels given a ROI-zip together with
its corresponding image.
- Includes identification of double and triple positive combinations.
- The ROI color code is annotated in the results table.
-## `3_manual_rerun.py`
+## [`3_manual_rerun.py`](3_manual_rerun.py)
- Requires an already open image with an already populated ROI manager.
- Allows to manually select measurement parameters and the measurement channel.