Roughed in chaetoceros project

extensions/chaeo/products.py

+from PIL import Image, ImageDraw, ImageFont
+
+def draw_boxes_on_2d_image(img, boxes, **kwargs):
+    pilimg = Image.fromarray(np.copy(img))  # drawing modifies array in-place
+    draw = ImageDraw.Draw(pilimg)
+    font_size = kwargs.get('font_size', 18)
+    linewidth = kwargs.get('linewidth', 4)
+
+    draw.font = ImageFont.truetype(font="arial.ttf", size=font_size)
+
+    for box in boxes:
+        y0 = box['info'].y0
+        y1 = box['info'].y1
+        x0 = box['info'].x0
+        x1 = box['info'].x1
+        xm = round((x0 + x1) / 2)
+
+        la = box['info'].label
+        zi = box['info'].zi
+
+        draw.rectangle([(x0, y0), (x1, y1)], outline='white', width=linewidth)
+
+        if kwargs.get('add_label') is True:
+            draw.text((xm, y0), f'Z{zi:04d}-L{la:04d}', fill='white', anchor='mb')
+
+    return pilimg
+
+def generate_patches(
+        desc, stack, boxes, rescale_clip=0.0,
+        pad_to=256,
+        proj=lambda x: x.max(axis=0),
+        prefix='patch',
+        **kwargs
+):
+    patch_dir = root / 'output' / 'patches' / desc
+    patch_dir.mkdir(parents=True, exist_ok=True)
+
+    for box in boxes:
+        obj = box['info']
+        sl = box['slice']
+        rbb = box['relative_bounding_box']
+
+        patch = proj(stack[sl])
+        patch_fname = f'{prefix}-la{obj.label:04d}-zi{obj.zi:04d}'
+
+        if rescale_clip is not None:
+            patch = rescale(patch, rescale_clip)
+
+        if kwargs.get('draw_bounding_box') is True:
+            x0 = rbb['x0']
+            y0 = rbb['y0']
+            x1 = rbb['x1']
+            y1 = rbb['y1']
+
+            pilimg = Image.fromarray(patch)  # drawing modifies array in-place
+            draw = ImageDraw.Draw(pilimg)
+            draw.rectangle([(x0, y0), (x1, y1)], outline='white', width=kwargs.get('linewidth', 1))
+            patch = np.array(pilimg)
+
+        if pad_to:
+            patch = pad(patch, pad_to)
+
+        imsave(
+            patch_dir / (patch_fname + '.png'),
+            resample(patch),
+            check_contrast=False,
+        )
+    print(f'Successfully wrote {len(boxes)} patches to:\n{patch_dir}')
+
+
+def generate_3d_patches(  # in extensions.chaeo.products
+        desc, stack, boxes, rescale_clip=0.0,
+        pad_to=256,
+        prefix='patch',
+        proj=lambda x: x,
+):
+    patch_dir = root / 'output' / '3d_patches' / desc
+    patch_dir.mkdir(parents=True, exist_ok=True)
+
+    for box in boxes:
+        obj = box['info']
+        sl = box['slice']
+        patch = proj(stack[sl])
+        patch_fname = f'{prefix}-la{obj.label:04d}-zi{obj.zi:04d}'
+
+        if rescale_clip is not None:
+            patch = rescale(patch, rescale_clip)
+
+        if pad_to:
+            patch = pad_3d(patch, pad_to)
+
+        imwrite(
+            patch_dir / (patch_fname + '.tif'),
+            patch,
+            imagej=True
+        )
+    print(f'Successfully wrote {len(boxes)} patches to:\n{patch_dir}')
\ No newline at end of file

extensions/chaeo/test_zstack.py

+import unittest
+
+from extensions.chaeo.zstack import build_stack_mask
+from model_server.accessors import InMemoryDataAccessor
+
+class TestZStackDerivedDataProducts(unittest.TestCase):
+
+    def setUp(self) -> None:
+        # need test data incl obj map
+        self.obmap = None
+        self.stack = None
+
+    def test_zmask_makes_correct_boxes(self):
+        zmask, meta = build_stack_mask(
+            'test_zmask_with boxes',
+            self.stack,
+            self.obmap,
+            mask_type='boxes',
+        )
+        zmask_acc = InMemoryDataAccessor(zmask)
+        self.assertTrue(zmask_acc.is_object_map())
+
+        # assert dimensionality of zmask
+        self.assertEqual(zmask.shape_dict['Z'] > 1)
+        self.assertEqual(zmask.shape_dict['C'] == 1)
+
+        # assert non-trivial meta info in boxes
+        pass
+
+    def test_zmask_makes_correct_contours(self):
+        pass
\ No newline at end of file

extensions/chaeo/zstack.py

+import numpy as np
+import pandas as pd
+
+from skimage.measure import find_contours, regionprops_table
+
+# build a single boolean 3d mask (objects v. bboxes) and return bounding boxes
+def build_stack_mask(desc, obmap, stack, filters=None, mask_type='contour', expand_box_by=(0, 0)): # TODO: specify boxes data type
+    """
+
+    filters: dict of (min, max) tuples
+    expand_box_by: (xy, z) pixelsf
+    """
+
+    # validate inputs
+    assert len(stack.shape) == 3, stack.shape
+    assert mask_type in ('contour', 'box'), mask_type # TODO: replace with call to validator
+
+    for k in filters.keys():
+        assert k in ('area', 'solidity')
+        vmin, vmax = filters[k]
+        assert vmin >= 0
+
+    # build object query
+    query_str = 'label > 0'  # always true
+    if filters:
+        for k in filters.keys():
+            assert k in ('area', 'solidity')
+            vmin, vmax = filters[k]
+            assert vmin >= 0
+            query_str = query_str + f' & {k} > {vmin} & {k} < {vmax}'
+
+    # build dataframe of objects, assign z index to each object
+    argmax = stack.argmax(axis=0)
+    df = (
+        pd.DataFrame(
+            regionprops_table(
+                obmap,
+                intensity_image=argmax,
+                properties=('label', 'area', 'intensity_mean', 'solidity', 'bbox')
+            )
+        )
+        .query(query_str)
+        .rename(
+            columns={
+                'bbox-0': 'y0',
+                'bbox-1': 'x0',
+                'bbox-2': 'y1',
+                'bbox-3': 'x1',
+            }
+        )
+    )
+    df['zi'] = df['intensity_mean'].round().astype('int')
+
+    # make an object map where label is replaced by focus position in stack and background is -1
+    lut = np.zeros(obmap.max() + 1) - 1
+    lut[df.label] = df.zi
+
+    # convert bounding boxes to slices
+    ebxy, ebz = expand_box_by
+    nz, h, w = stack.shape
+
+    boxes = []
+    for ob in df.itertuples(name='LabeledObject'):
+        y0 = max(ob.y0 - ebxy, 0)
+        y1 = min(ob.y1 + ebxy, h - 1)
+        x0 = max(ob.x0 - ebxy, 0)
+        x1 = min(ob.x1 + ebxy, w - 1)
+        z0 = max(ob.zi - ebz, 0)
+        z1 = min(ob.zi + ebz, nz)
+
+        # relative bounding box positions
+        rbb = {
+            'y0': ob.y0 - y0,
+            'y1': ob.y1 - y0,
+            'x0': ob.x0 - x0,
+            'x1': ob.x1 - x0,
+        }
+
+        sl = np.s_[z0: z1 + 1, y0: y1, x0: x1]
+
+        # compute contours
+        obmask = (obmap == ob.label)
+        contour = find_contours(obmask)
+        mask = obmask[ob.y0: ob.y1, ob.x0: ob.x1]
+
+        boxes.append({
+            'info': ob,
+            'slice': sl,
+            'relative_bounding_box': rbb,
+            'contour': contour,
+            'mask': mask
+        })
+
+    # build mask z-stack
+    zi_st = np.zeros(stack.shape, dtype='bool')
+    if mask_type == 'contour':
+        zi_map = (lut[obmap] + 1.0).astype('int')
+        idxs = np.array([zi_map]) - 1
+        np.put_along_axis(zi_st, idxs, 1, axis=0)
+
+        # change background level from to 0 in final frame
+        zi_st[-1, :, :][obmap == 0] = 0
+
+    elif mask_type == 'box':
+        for bb in boxes:
+            sl = bb['slice']
+            zi_st[sl] = 1
+
+    return zi_st, boxes
\ No newline at end of file

extensions/ilastik/validators.py

+from model_server.accessors import GenericImageDataAccessor
+
+def is_object_map(img: GenericImageDataAccessor):
+    # TODO: implement
+    pass
\ No newline at end of file

model_server/batch.py

+def get_tif_seq(where, prefix=None, verbose=True):
+    files = [ff for ff in os.listdir(where) if ff.startswith(prefix)]
+    files.sort()
+    paths = [where / fn for fn in files]
+    ndas = [imread(pa.__str__()) for pa in paths]
+    assert all([n.shape == ndas[0].shape and n.dtype == ndas[0].dtype for n in ndas])
+    if verbose:
+        print(f'Successfully read {len(ndas)} x {ndas[0].shape} px images of type {ndas[0].dtype} from:\n{where}')
+    return ndas
+
+def write_tif_zstack(data, desc, dtype='uint16', subdir='output'):
+    where = root / subdir
+    filename = desc + '.tif'
+    if not os.path.exists(where):
+        os.makedirs(where)
+    abspath = where / filename
+    imwrite(
+        abspath,
+        np.array(data).astype(dtype),
+        imagej=True
+    )
+    print(f'Successfully wrote {len(data)} x {data[0].shape} px images of type {dtype} to:\n{abspath}')
\ No newline at end of file

model_server/process.py

+"""
+Image processing utility functions
+"""
+from math import ceil, floor
+
+import numpy as np
+from skimage.exposure import rescale_intensity
+
+
+def pad(im, mpx):  # now in model_server.batch
+    '''Pads and crops image width edge values to specified dimension'''
+    dh = 0.5 * (mpx - im.shape[0])
+    dw = 0.5 * (mpx - im.shape[1])
+
+    if dw < 0:
+        x0 = floor(-dw)
+        x1 = x0 + mpx
+        im = im[:, x0:x1]
+        dw = 0
+    if dh < 0:
+        y0 = floor(-dh)
+        y1 = y0 + mpx
+        im = im[y0:y1, :]
+        dh = 0
+
+    border = ((floor(dh), ceil(dh)), (floor(dw), ceil(dw)))
+    padded = np.pad(im, border, mode='constant')
+    if padded.shape != (mpx, mpx):
+        raise Exception(f'Incorrect image shape: {padded.shape} v. {(mpx, mpx)}')
+    return padded
+
+def pad_3d(im, mpx): # im: [z x h x w]
+    assert(len(im.shape) == 3)
+    nz, h, w = im.shape
+    padded = np.zeros((nz, mpx, mpx), dtype=im.dtype)
+    for zi in range(nz):
+        padded[zi, :, :] = pad(im[zi, :, :], mpx)
+    return padded
+
+def resample(nda, cmin=0, cmax=2**16): # now in model_server.batch
+    return rescale_intensity(
+        np.clip(nda, cmin, cmax),
+        in_range=(cmin, cmax + 1),
+        out_range=(0, 2**8)
+    ).astype('uint8')
+
+
+def rescale(nda, clip=0.0): # now in model_server.batch
+    clip_pct = (100.0 * clip, 100.0 * (1.0 - clip))
+    cmin, cmax = np.percentile(nda, clip_pct)
+    rescaled = rescale_intensity(nda, in_range=(cmin, cmax))
+    return rescaled
\ No newline at end of file

readme.md

+# model_server
+
+# How to extend service
+Add sub-package to extensions
+Add models that inherit from model_server.Model 
+In workflows, implement pipelines with File I/O via accessors.GenericImageDataAccessor
+(to decouple model logic from image data source)
+Set extensions-specific folders, etc. in conf relative to overall package root (set by user)
+As much as possible, set pipeline and model parameters with defaults and support overrides by optional API arguments; 
+this helps non-coding users control their jobs
+Set up API endpoints in router, following as much as possible existing conventions with load, infer, etc. keyword
+
+decouple data access from processing
+
+control either via batch runners or API (serial)
+
+workflow: combines data access with processing via models, produces primary outputs
\ No newline at end of file