Skip to content

Commit 7f27146

Browse files
adding plant prep for cpsam for #1437
1 parent 746f0d2 commit 7f27146

2 files changed

Lines changed: 157 additions & 3 deletions

File tree

Lines changed: 154 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,154 @@
1+
import os
2+
import numpy as np
3+
import matplotlib.pyplot as plt
4+
from scipy.ndimage import gaussian_filter, gaussian_filter1d
5+
from cellpose import io, transforms, utils, models, dynamics, metrics
6+
import h5py
7+
import gc
8+
from glob import glob
9+
import cv2
10+
from natsort import natsorted
11+
import shutil
12+
from pathlib import Path
13+
import pandas as pd
14+
from tqdm import trange
15+
import fastremap
16+
from skimage.measure import label
17+
18+
19+
def create_plant_train(lateral_root = True):
20+
""" Data from Wolny et al 2020 https://osf.io/uzq3w/overview
21+
22+
For lateral root we used Movie2 only for testing since it was fully segmented
23+
24+
"""
25+
if lateral_root:
26+
root = Path("/media/carsen/ssd3/plants/lateral_root/")
27+
# 0.1625 × 0.1625 × 0.250 um/px
28+
anisotropy = 0.25 / 0.1625
29+
rsz = 0.25
30+
else:
31+
root = Path("/media/carsen/ssd3/plants/ovules/")
32+
# 0.075 × 0.075 × 0.235 um/px
33+
anisotropy = 0.235 / 0.075
34+
rsz = 0.75
35+
print(anisotropy)
36+
37+
ly, lx = 256, 256
38+
pm = [(0,1,2), (1,0,2), (2,0,1)]
39+
pstr = ["YX", "ZX", "ZY"]
40+
41+
train_files = (root / "train").glob("*.h5")
42+
train_files = natsorted([tf for tf in train_files])
43+
test_files = (root / "test").glob("*.h5")
44+
test_files = natsorted([tf for tf in test_files])
45+
print(len(train_files), len(test_files))
46+
47+
np.random.seed(0)
48+
for k, files in enumerate([train_files, test_files]):
49+
for i, tf in enumerate(files):
50+
print(tf.stem)
51+
f = h5py.File(tf, "r")
52+
print(f.keys())
53+
img = np.array(f["raw"])
54+
masks = (np.array(f["label"])).astype("uint16")
55+
if lateral_root:
56+
if k==0: masks -= 1
57+
else:
58+
if i==0: masks[masks==1] = 0
59+
else: masks[masks==411] = 0
60+
else:
61+
ignore = np.array(f["label_with_ignore"])==-1
62+
zignore = ignore.mean(axis=(1,2))
63+
zmask = (masks>0).mean(axis=(1,2))
64+
zmin = np.nonzero(zignore < zmask*0.2)[0]
65+
if len(zmin) == 0: continue
66+
zmin = zmin[0]
67+
zmax = np.nonzero(zignore[zmin:] > zmask[zmin:]*0.2)[0]
68+
zmax = zmax[0] + zmin if len(zmax) > 0 else len(zmask)
69+
print(zmin, zmax)
70+
if zmax - zmin < 50: continue
71+
masks = masks[zmin:zmax]
72+
img = img[zmin:zmax]
73+
74+
fastremap.renumber(masks, in_place=True)
75+
print(f"ncells = {masks.max()}")
76+
### resize
77+
print(img.shape)
78+
if rsz!=1:
79+
Lyr = int(masks.shape[-2] * rsz)
80+
Lxr = int(masks.shape[-1] * rsz)
81+
masks_rsz = transforms.resize_image(masks, Ly=Lyr, Lx=Lxr, no_channels=True,
82+
interpolation=cv2.INTER_NEAREST).astype(masks.dtype)
83+
img_rsz = transforms.resize_image(img, Ly=Lyr, Lx=Lxr,
84+
no_channels=True).astype(img.dtype)
85+
else:
86+
masks_rsz = masks.copy()
87+
img_rsz = img.copy()
88+
# make isotropic
89+
Lyr = int(masks_rsz.shape[0] * anisotropy * rsz)
90+
Lxr = int(masks_rsz.shape[-1])
91+
masks_rsz = transforms.resize_image(masks_rsz.transpose(1,0,2), Ly=Lyr, Lx=Lxr,
92+
no_channels=True, interpolation=cv2.INTER_NEAREST).astype(masks.dtype).transpose(1,0,2)
93+
img_rsz = transforms.resize_image(img_rsz.transpose(1,0,2), Ly=Lyr, Lx=Lxr,
94+
no_channels=True).astype(img.dtype).transpose(1,0,2)
95+
print(img_rsz.shape)
96+
if lateral_root and k==1:
97+
masks_rsz = masks_rsz[:130]
98+
img_rsz = img_rsz[:130]
99+
else:
100+
th = 0.05 # fraction of mask pixels required for cropping
101+
for d, inds in enumerate([(1, 2), (0, 2), (0, 1)]):
102+
# compute fraction of mask pixels in each slice
103+
m = (masks_rsz > 0).mean(axis=inds)
104+
imin = max(0, np.nonzero(m>th)[0][0] - 10)
105+
imax = min(len(m), len(m) - np.nonzero(m[::-1]>th)[0][0] + 10)
106+
# slice d dimension
107+
masks_rsz = masks_rsz.take(range(imin, imax), axis=d)
108+
img_rsz = img_rsz.take(range(imin, imax), axis=d)
109+
110+
img = img_rsz.copy()
111+
masks = masks_rsz.copy()
112+
print(img.shape, masks.shape)
113+
114+
if k==1:
115+
folder = root / "test_3D"
116+
folder.mkdir(exist_ok=True)
117+
io.imsave(folder / f"{tf.stem}.tif", img_rsz)
118+
io.imsave(folder / f"{tf.stem}_masks.tif", masks_rsz)
119+
120+
for p in range(3):
121+
n_slices = 20 if p==0 else 10
122+
if p > 0:
123+
masks_rsz = masks.transpose(pm[p])
124+
img_rsz = img.transpose(pm[p])
125+
else:
126+
masks_rsz = masks.copy()
127+
img_rsz = img.copy()
128+
129+
Lz, Ly, Lx = masks_rsz.shape
130+
131+
# random z/y/x slices
132+
#nr = 400
133+
#iz = np.random.randint(Lz - 2*zpad, size=(nr,)) + zpad
134+
zpad = 10
135+
iz = np.linspace(zpad, Lz - zpad, n_slices*2, dtype=int)
136+
mr = masks_rsz[iz]
137+
mp = (mr>0).mean(axis=(-2,-1))
138+
igood = np.sort(mp.argsort()[::-1][:n_slices])
139+
mr = mr[igood]
140+
iz = iz[igood]
141+
imr = img_rsz[iz]
142+
mr = [label(mr[i]) for i in range(len(mr))]
143+
144+
# save training crops
145+
(root / "slices").mkdir(exist_ok=True)
146+
folder = root / "slices" / "train" if k == 0 else root / "slices" / "test"
147+
folder.mkdir(exist_ok=True)
148+
for i in range(len(iz)):
149+
M0 = mr[i].copy()
150+
Im0 = imr[i]
151+
M0 = utils.fill_holes_and_remove_small_masks(M0, 100)
152+
fastremap.renumber(M0, in_place=True)
153+
io.imsave(folder / f"{tf.stem}_{pstr[p]}_{i:03d}.tif", Im0)
154+
io.imsave(folder / f"{tf.stem}_{pstr[p]}_{i:03d}_masks.tif", M0)

paper/cpsam/figures.py

Lines changed: 3 additions & 3 deletions
Original file line numberDiff line numberDiff line change
@@ -218,7 +218,7 @@ def fig2(root, save_fig=False):
218218
axin.text(x0, 0, "false negatives (FN)", color="r", va="center")
219219

220220
axin.scatter(2, -1.5, marker="o", color=colors_ia[1], s=100, facecolor="none", lw=1.5)
221-
axin.text(1.8, -1.5, "human consensus", color=colors_ia[1], va="center", ha="right")
221+
axin.text(1.8, -1.5, "ground-truth", color=colors_ia[1], va="center", ha="right")
222222

223223
for k in range(2):
224224
if k==0:
@@ -245,7 +245,7 @@ def fig2(root, save_fig=False):
245245
ax.scatter(np.median(outline[:,0]), np.median(outline[:,1]), color="r",
246246
s=300, marker="x", lw=2, zorder=30)
247247
ax.axis("off")
248-
ax.set_title(f"Annotator {k+1}\n to human consensus", loc="center")
248+
ax.set_title(f"Annotator {k+1}\n to ground-truth", loc="center")
249249

250250
iplot_ia = np.array([0, 1])
251251
iplot_errors = np.arange(len(aps))
@@ -284,7 +284,7 @@ def fig2(root, save_fig=False):
284284
axin.set_xlim([-0.5, (len(aps)+1)+.5])
285285
axin.axis("off")
286286
colors_all = [*colors_ia, *[alg_dict[alg]["color"] for alg in algs]]
287-
labels_all = ["Annotator 2", "human consensus\n(estimate)"] + [alg_dict[alg]["name"] for alg in algs]
287+
labels_all = ["Annotator 2", "absolute bound\n(estimate)"] + [alg_dict[alg]["name"] for alg in algs]
288288
labels_all = [label.replace("MicroSAM", "MicroSAM (2026)") for label in labels_all]
289289
ax.set_xticks(np.arange(0, len(aps)+2))
290290
ax.set_xticklabels(labels_all, rotation=45, ha="right")

0 commit comments

Comments
 (0)