From cce7ee89a8dd4138d1dfccf96950a1e0936e93b6 Mon Sep 17 00:00:00 2001
From: Guillaume Cerutti <guillaume.cerutti@inria.fr>
Date: Thu, 3 Aug 2023 09:57:03 +0200
Subject: [PATCH] add wall mesh extraction + basic quantification plugins

---
 setup.py                                      |   6 +-
 .../wallMeshImageSignal.py                    | 107 ++++++++++++++++++
 .../meshFromImage/segmentationWallMesh.py     |  69 +++++++++++
 3 files changed, 180 insertions(+), 2 deletions(-)
 create mode 100644 src/gnomon_package_tissueimagemesh/algorithm/cellImageQuantification/wallMeshImageSignal.py
 create mode 100644 src/gnomon_package_tissueimagemesh/algorithm/meshFromImage/segmentationWallMesh.py

diff --git a/setup.py b/setup.py
index 3624595..802e225 100644
--- a/setup.py
+++ b/setup.py
@@ -43,11 +43,13 @@ setup_kwds = dict(
         'cellImageQuantification': [
             'surfaceMeshCellCurvature = gnomon_package_tissueimagemesh.algorithm.cellImageQuantification.surfaceMeshCellCurvature',
             'surfaceMeshCellLayer = gnomon_package_tissueimagemesh.algorithm.cellImageQuantification.surfaceMeshCellLayer',
-            'surfaceMeshCellProperty = gnomon_package_tissueimagemesh.algorithm.cellImageQuantification.surfaceMeshCellProperty'
+            'surfaceMeshCellProperty = gnomon_package_tissueimagemesh.algorithm.cellImageQuantification.surfaceMeshCellProperty',
+            'wallMeshImageSignal = gnomon_package_tissueimagemesh.algorithm.cellImageQuantification.wallMeshImageSignal',
         ],
         'meshFromImage': [
             'imageSurfaceMesh = gnomon_package_tissueimagemesh.algorithm.meshFromImage.imageSurfaceMesh',
-            'segmentationSurfaceMesh = gnomon_package_tissueimagemesh.algorithm.meshFromImage.segmentationSurfaceMesh'
+            'segmentationSurfaceMesh = gnomon_package_tissueimagemesh.algorithm.meshFromImage.segmentationSurfaceMesh',
+            'segmentationWallMesh = gnomon_package_tissueimagemesh.algorithm.meshFromImage.segmentationWallMesh'
         ],
         'pointCloudQuantification': [
             'nucleiCurvatureNukem = gnomon_package_tissueimagemesh.algorithm.pointCloudQuantification.nucleiCurvatureNukem',
diff --git a/src/gnomon_package_tissueimagemesh/algorithm/cellImageQuantification/wallMeshImageSignal.py b/src/gnomon_package_tissueimagemesh/algorithm/cellImageQuantification/wallMeshImageSignal.py
new file mode 100644
index 0000000..0775221
--- /dev/null
+++ b/src/gnomon_package_tissueimagemesh/algorithm/cellImageQuantification/wallMeshImageSignal.py
@@ -0,0 +1,107 @@
+import logging
+from copy import deepcopy
+
+import numpy as np
+import pandas as pd
+
+from dtkcore import d_real, d_inliststring, d_inliststringlist
+
+from gnomon.utils import algorithmPlugin
+from gnomon.utils.decorators import cellImageInput, imageInput, meshInput, cellImageOutput, dataFrameOutput
+from gnomon.core import gnomonAbstractCellImageQuantification
+
+from timagetk import TissueImage3D
+from timagetk.features.pandas_tools import cells_to_dataframe
+
+from timagetk_geometry.signal_quantification.wall_mesh import quantify_wall_topomesh_signal_intensity
+
+
+@algorithmPlugin(version="0.3.0", coreversion="0.81.0", name="Wall Mesh Signal")
+@meshInput('wall_topomesh', data_plugin="gnomonMeshDataPropertyTopomesh")
+@imageInput('img', data_plugin='gnomonImageDataMultiChannelImage')
+@cellImageInput("tissue", data_plugin="gnomonCellImageDataTissueImage")
+@cellImageOutput("out_tissue", data_plugin="gnomonCellImageDataTissueImage")
+@dataFrameOutput('df', data_plugin="gnomonDataFrameDataPandas")
+class wallMeshImageSignal(gnomonAbstractCellImageQuantification):
+    """Measure cell-wall image signal intensity using a triangle mesh.
+
+    The method computes a value of signal for each cell wall represented in a
+    wall triangle mesh by averaging the image intensity at the level of the
+    mesh vertices.
+
+    """
+
+    def __init__(self):
+        super().__init__()
+
+        self.tissue = {}
+        self.img = {}
+        self.wall_topomesh = {}
+
+        self.out_tissue = {}
+        self.df = {}
+
+        self._parameters = {}
+        self._parameters['channels'] = d_inliststringlist("Channels", [""],[""], "Channels on which to quantify signal intensity at the cell walls")
+        self._parameters['gaussian_sigma'] = d_real("Sigma", 0.5, 0., 50., 2, "Standard deviation of the Gaussian kernel used to smooth signal")
+
+        self._parameter_groups = {}
+
+    def __del__(self):
+        pass
+
+    def refreshParameters(self):
+        if len(self.img)>0:
+            img = list(self.img.values())[0]
+
+            if len(img.channel_names) == 1:
+                if 'channels' in self._parameters.keys():
+                    del self._parameters['channels']
+            else:
+                if not 'channels' in self._parameters.keys():
+                    self._parameters['channels'] = d_inliststringlist("channels", [""], [""], "Channels on which to quantify signal intensity at the cell walls")
+                self._parameters['channels'].setValues(img.channel_names)
+                self._parameters['channels'].setValue(img.channel_names)
+
+    def run(self):
+        self.out_tissue = {}
+        self.df = {}
+
+        for time in self.tissue.keys():
+            if not time in self.img.keys():
+                logging.error("Impossible to quantify! No signal image provided")
+                self.df = {}
+                return
+            if not time in self.wall_topomesh.keys():
+                logging.error("Impossible to quantify! No wall mesh provided")
+                self.df = {}
+                return
+            tissue = self.tissue[time]
+
+            out_tissue = TissueImage3D(tissue, not_a_label=0, background=1)
+            out_tissue.cells = deepcopy(tissue.cells)
+            out_tissue.cells.image = out_tissue
+            out_tissue.walls = deepcopy(tissue.walls)
+            out_tissue.walls.image = out_tissue
+
+            img = self.img[time]
+            wall_topomesh = self.wall_topomesh[time]
+
+            wall_signals = quantify_wall_topomesh_signal_intensity(
+                signal_img=img,
+                wall_topomesh=wall_topomesh,
+                wall_sigma=self['gaussian_sigma']
+            )
+
+            for channel in wall_signals:
+                signal_name = channel if channel != "" else "image_signal"
+                if ('channels' not in self._parameters) or (signal_name in self['channels']):
+                    out_tissue.walls.set_feature(signal_name, wall_signals[channel])
+
+            self.out_tissue[time] = out_tissue
+
+            self.df[time] = cells_to_dataframe(out_tissue.cells, cell_ids=out_tissue.cell_ids())
+            self.df[time]['time'] = time
+
+        if len(self.df) > 0:
+            self.df = {0: pd.concat([self.df[time] for time in self.df.keys()])}
diff --git a/src/gnomon_package_tissueimagemesh/algorithm/meshFromImage/segmentationWallMesh.py b/src/gnomon_package_tissueimagemesh/algorithm/meshFromImage/segmentationWallMesh.py
new file mode 100644
index 0000000..16fece9
--- /dev/null
+++ b/src/gnomon_package_tissueimagemesh/algorithm/meshFromImage/segmentationWallMesh.py
@@ -0,0 +1,69 @@
+import numpy as np
+
+
+from dtkcore import d_bool, d_int, d_inliststringlist, d_real
+
+import gnomon.core
+from gnomon.core import gnomonAbstractMeshFromImage
+
+from gnomon.utils import corePlugin
+from gnomon.utils.decorators import cellImageInput, meshOutput
+
+from timagetk_geometry.image_surface.wall_mesh import tissue_image_wall_topomesh
+
+from cellcomplex.property_topomesh.analysis import compute_topomesh_property
+
+
+@corePlugin(version="0.3.0", coreversion="0.81.0", name="Cell Wall Meshes")
+@cellImageInput("seg_img", data_plugin='gnomonCellImageDataTissueImage')
+@meshOutput("topomesh", data_plugin="gnomonMeshDataPropertyTopomesh")
+class segmentationWallMesh(gnomonAbstractMeshFromImage):
+    """Extract the cell walls of a segmented image as triangular meshes.
+    """
+
+    def __init__(self):
+        super().__init__()
+
+        self.seg_img = {}
+        self.topomesh = {}
+
+        self._parameters = {}
+        self._parameters['resampling_voxelsize'] = d_real("Resampling voxelsize", 0, 0, 1, 2, "Cubic voxelsize in which to resample the image before mesh extraction")
+        self._parameters['smoothing'] = d_bool("Smoothing", True, "Whether to smooth and decimate the wall meshes")
+        self._parameters['wall_types'] = d_inliststringlist('Walls', ["anticlinal_L1"], ["anticlinal_L1", "epidermis_L1", "all"], "Which cell walls to extract in the resulting mesh")
+        self._parameters['edge_length'] = d_real('Edge length', 1, 0, 10, 1, "Target length of mesh edges")
+
+        self._parameter_groups = {}
+
+    def refreshParameters(self):
+        if len(self.seg_img)>0:
+            seg_img = list(self.seg_img.values())[0]
+
+            self._parameters['resampling_voxelsize'].setMin(0)
+            self._parameters['resampling_voxelsize'].setMax(10*np.max(seg_img.voxelsize))
+            self._parameters['resampling_voxelsize'].setValue(4*np.mean(seg_img.voxelsize))
+
+    def run(self):
+        self.topomesh = {}
+
+        for time in self.seg_img.keys():
+            seg_img = self.seg_img[time]
+
+            resampling_voxelsize = None if self['resampling_voxelsize'] == 0 else self['resampling_voxelsize']
+            topomesh = tissue_image_wall_topomesh(
+                tissue=seg_img,
+                resampling_voxelsize=resampling_voxelsize,
+                wall_types=self['wall_types'],
+                smoothing=self['smoothing'],
+                target_edge_length=self['edge_length']
+            )
+
+            for degree in [1, 2, 3]:
+                compute_topomesh_property(topomesh, 'barycenter', degree)
+
+            for degree in [0, 1, 2, 3]:
+                positions = topomesh.wisp_property('barycenter', degree)
+                for k, dim in enumerate(['x', 'y', 'z']):
+                    topomesh.update_wisp_property('barycenter_'+dim, degree, dict(zip(positions.keys(), positions.values()[:, k])))
+
+            self.topomesh[time] = topomesh
-- 
GitLab