Auto Fiji ImageJ Cell Count

AutoFijiCellCount.py

#!/usr/bin/env python3
import sys
from pathlib import Path
from natsort import humansorted
import imagej
import scyjava
import tqdm
import openpyxl

# python3.10
# pip install tqdm natsort pyimagej scyjava openpyxl

ij = imagej.init( 'sc.fiji:fiji' , mode='headless' )

RANGE = "25-150"
CIRCULARITY = "0.30-1.00"

def get_files_in_directory( source_path ):
	ALLOWED_EXTENSIONS = [ ".czi" , ".tiff" ]
	_files = source_path.glob( "**/*" )
	_files = [x for x in _files if x.is_file()]
	_files = [x for x in _files if x.suffix.lower() in ALLOWED_EXTENSIONS]
	_files = humansorted(_files)
	return _files

from scyjava import jimport
ImagePlus = jimport("ij.ImagePlus")

def image_j_threshold_and_analyze_particles(input_image_path):

	print(f"\nProcessing: {input_image_path}")

	from scyjava import jimport

	Color       = jimport("java.awt.Color")
	Font        = jimport("java.awt.Font")
	Analyzer    = jimport("ij.plugin.filter.Analyzer")
	Overlay     = jimport("ij.gui.Overlay")
	Roi         = jimport("ij.gui.Roi")
	TextRoi     = jimport("ij.gui.TextRoi")

	# --- Load once for metadata/calibration ---
	base = ij.io().open(str(input_image_path))
	base = ij.py.to_imageplus(base)

	width      = base.getWidth()
	height     = base.getHeight()
	image_area = width * height
	n_channels = base.getNChannels()

	cal  = base.getCalibration()
	px_w = cal.pixelWidth if cal.pixelWidth else 1.0
	px_h = cal.pixelHeight if cal.pixelHeight else 1.0

	base.close()

	# AREA + CENTROID + BOUNDING_RECT
	Analyzer.setMeasurements(1 | 16 | 512)

	results = []

	for c in range(1, n_channels + 1):
		print(f"  Channel {c}/{n_channels}")

		# --- reload clean parent each loop ---
		imp = ij.io().open(str(input_image_path))
		imp = ij.py.to_imageplus(imp)

		# --- isolate channel exactly like your working flow ---
		ij.IJ.run(imp, "Duplicate...", f"title=ch{c} duplicate channels={c}")
		orig = ij.WindowManager.getImage(f"ch{c}")
		imp.close()

		# --- analysis copy only ---
		analysis = orig.duplicate()

		# --- save grayscale original ---
		# grey_path = input_image_path.parent / f"{input_image_path.stem}_channel_{c}_gray.png"
		# ij.IJ.save(orig, str(grey_path))

		# --- threshold/mask only analysis image ---
		ij.IJ.setAutoThreshold(analysis, "Moments dark")
		ij.IJ.run(analysis, "Convert to Mask", "method=Moments background=Dark")

		# --- save mask ---
		mask_path = input_image_path.parent / f"{input_image_path.stem}_channel_{c}_mask.png"
		ij.IJ.save(analysis, str(mask_path))

		# --- reset RT before analyze ---
		rt = ij.ResultsTable.getResultsTable()
		if rt is not None:
			rt.reset()

		Analyzer.setMeasurements(1 | 16 | 512)

		ij.IJ.run(
			analysis,
			"Analyze Particles...",
			f"size={RANGE} circularity={CIRCULARITY} show=Nothing clear include"
		)

		rt = ij.ResultsTable.getResultsTable()

		particles = []

		if rt is None or rt.size() == 0:
			count = 0
			total_area = 0
		else:
			cols = list(rt.getHeadings())
			count = rt.size()
			total_area = 0

			has_bbox     = all(k in cols for k in ("BX", "BY", "Width", "Height"))
			has_centroid = ("X" in cols and "Y" in cols)

			for i in range(count):
				area = rt.getValue("Area", i)
				total_area += area

				# ResultsTable values are calibrated in headless Fiji.
				# Convert back to pixel coords for drawing.
				if has_bbox:
					bx = int(round(rt.getValue("BX", i)     / px_w))
					by = int(round(rt.getValue("BY", i)     / px_h))
					bw = int(round(rt.getValue("Width", i)  / px_w))
					bh = int(round(rt.getValue("Height", i) / px_h))

					# clamp sane minimums
					if bw < 1: bw = 1
					if bh < 1: bh = 1

					cx = bx + bw // 2
					cy = by + bh // 2

				elif has_centroid:
					cx = int(round(rt.getValue("X", i) / px_w))
					cy = int(round(rt.getValue("Y", i) / px_h))
					bx, by, bw, bh = cx - 6, cy - 6, 12, 12

				else:
					cx, cy, bx, by, bw, bh = 0, 0, 0, 0, 1, 1

				# clip into image bounds
				if bx < 0: bx = 0
				if by < 0: by = 0
				if bx >= width:  bx = width - 1
				if by >= height: by = height - 1
				if bx + bw > width:   bw = max(1, width - bx)
				if by + bh > height:  bh = max(1, height - by)

				particles.append((i + 1, cx, cy, bx, by, bw, bh))

		avg_size     = total_area / count if count else 0
		percent_area = (total_area / image_area) * 100 if image_area else 0

		print(f"    Count={count}, TotalArea={total_area:.2f}")

		# --- draw traces on ORIGINAL using overlay + flatten ---
		# --- DRAW ON MASK (NOT GRAY) ---
		if particles:
			try:
				traced = analysis.duplicate()

				# 🔥 force RGB so colors show
				ij.IJ.run(traced, "RGB Color", "")

				ip = traced.getProcessor()

				Color = scyjava.jimport("java.awt.Color")
				Font  = scyjava.jimport("java.awt.Font")

				ip.setColor(Color(0, 255, 0))
				ip.setFont(Font("SansSerif", Font.BOLD, 12))

				for (num, cx, cy, bx, by, bw, bh) in particles:
					ip.drawRect(bx, by, bw, bh)
					ip.drawString(str(num), cx - 5, cy + 5)

				traced.updateAndDraw()

				traced_path = input_image_path.parent / f"{input_image_path.stem}_channel_{c}_traced.png"
				ij.IJ.save(traced, str(traced_path))

				print(f"    Saved traced PNG ({len(particles)} particles): {traced_path.name}")

				traced.close()

			except Exception as e:
				import traceback
				print(f"    Warning: tracing failed: {e}")
				traceback.print_exc()

		results.append({
			"image":        input_image_path.stem,
			"channel":      c,
			"count":        count,
			"total_area":   total_area,
			"average_size": avg_size,
			"percent_area": percent_area,
			"mean":         0
		})

		analysis.close()
		orig.close()

	return results
def main():

	p = Path(sys.argv[1])
	files   = []
	results = []
	if p.is_file():
		files.append(p)
	elif p.is_dir():
		files = get_files_in_directory(p)
	for f in tqdm.tqdm(files):
		res = image_j_threshold_and_analyze_particles(f)
		if res:
			results.append(res)
	if results:
		wb = openpyxl.Workbook()
		ws = wb.active
		ws.title = "Results"

		headers = [
			"image", "fraction",
			"ch1_count", "ch1_total_area", "ch1_avg_size", "ch1_percent_area",
			"ch2_count", "ch2_total_area", "ch2_avg_size", "ch2_percent_area",
		]
		ws.append(headers)

		for image_results in results:
			ch1 = image_results[0] if len(image_results) > 0 else None
			ch2 = image_results[1] if len(image_results) > 1 else None

			fraction = ch1["count"] / ch2["count"] if (ch1 and ch2 and ch2["count"]) else 0

			row = [
				ch1["image"]        if ch1 else "unknown",
				fraction,
				ch1["count"]        if ch1 else 0,
				ch1["total_area"]   if ch1 else 0,
				ch1["average_size"] if ch1 else 0,
				ch1["percent_area"] if ch1 else 0,
				ch2["count"]        if ch2 else 0,
				ch2["total_area"]   if ch2 else 0,
				ch2["average_size"] if ch2 else 0,
				ch2["percent_area"] if ch2 else 0,
			]
			ws.append(row)

		output_path = Path("results.xlsx")
		wb.save(output_path)
		print(f"\nSaved results to {output_path}")

if __name__ == "__main__":
	main()

https://docs.google.com/presentation/d/1TXpvDIMVnWfA9GeINPQamX7Nh9AK9OwLeO6Q9yZg1Dg/edit?usp=sharing