more work related to graffy/lipase#3: for each pixel, the profile of variance along the circle is now computed. This image basically measures the non-roundedness for each pixel.

src/ij-plugins/Ipr/Lipase/Radial_Profile.py

@@ -1,8 +1,10 @@
 #@ ImagePlus (label="the input image") INPUT_IMAGE
-#@ Float (label="maximal radius", value=10.0, min=0.0, max=100.0, style="slider") MAX_RADIUS
+#@ Float (label="maximal radius", value=32.0, min=0.0, max=100.0, style="slider") MAX_RADIUS
-#@ Int (label="number of angular sectors", value=4, min=0, max=100, style="slider") NUM_ANGULAR_SECTORS
+#@ Integer (label="number of radial sectors", value=8, min=0, max=100, style="slider") NUM_RADIAL_SECTORS
+#@ Integer (label="number of angular sectors", value=4, min=0, max=100, style="slider") NUM_ANGULAR_SECTORS
 
-#@output ImagePlus RADIAL_PROFILE
+#@output ImagePlus RADIAL_PROFILES
+#@output ImagePlus ANGULAR_VARIANCE_AVG_IMAGE
 """This script is supposed to be launched from fiji's jython interpreter
 
 This imagej plugin computes the radial profile of each pixel of the input image. The resulting profiles are stored in a single hyperstack, where the profile data are stored along the channel axis
@@ -25,17 +27,20 @@ from jarray import zeros, array  # pylint: disable=import-error
 def run_script():
     global INPUT_IMAGE  # pylint:disable=global-variable-not-assigned
     global MAX_RADIUS  # pylint:disable=global-variable-not-assigned
-    global RADIAL_PROFILE  # pylint:disable=global-variable-not-assigned
+    global NUM_RADIAL_SECTORS  # pylint:disable=global-variable-not-assigned
     global NUM_ANGULAR_SECTORS  # pylint:disable=global-variable-not-assigned
+    global RADIAL_PROFILES  # pylint:disable=global-variable-not-assigned
+    global ANGULAR_VARIANCE_AVG_IMAGE  # pylint:disable=global-variable-not-assigned
     IImageEngine.set_instance(IJImageEngine())
 
     src_image = IImageEngine.get_instance().create_image(width=1, height=1, pixel_type=PixelType.U8)
     src_image.ij_image = INPUT_IMAGE  # pylint: disable=undefined-variable
 
-    detector = CircularSymmetryDetector(max_radius=MAX_RADIUS, num_angular_sectors=NUM_ANGULAR_SECTORS)  # pylint: disable=undefined-variable
+    detector = CircularSymmetryDetector(max_radius=MAX_RADIUS, num_angular_sectors=NUM_ANGULAR_SECTORS, num_radial_sectors=NUM_RADIAL_SECTORS)  # pylint: disable=undefined-variable
-    radial_profiles = detector.compute_radial_profiles(src_image)
+    radial_profiles, angular_variance_avg_image = detector.compute_radial_profiles(src_image)
 
-    RADIAL_PROFILE = radial_profiles.hyperstack
+    RADIAL_PROFILES = radial_profiles.hyperstack
+    ANGULAR_VARIANCE_AVG_IMAGE = angular_variance_avg_image.ij_image
 
 
 # note : when launched from fiji, __name__ doesn't have the value "__main__", as when launched from python

src/lipase/circsymdetector.py

@@ -168,12 +168,17 @@ class CircularSymmetryDetector:
         """ Computes for each pixel the radial profile (with this pixel as center)
 
         :param IImage src_image:
-        :rtype IHyperstack: each radial profile is stored in the channel coordinate of the hyperstack
+        :rtype IHyperstack, Image : 
+        :returns:
+            - radial_profiles (:py:class:`IHyperstack`) - each radial profile is stored in the channel coordinate of the hyperstack
+            - angular_variance_avg (:py:class:`IImage`) - each pixel stores the average variance of signal along the circles of different diameters in the neighborhood
         """
+        # https://stackoverflow.com/questions/39759503/how-to-document-multiple-return-values-using-restructuredtext-in-python-2
         ie = IImageEngine.get_instance()
         ie.debugger.on_image(src_image, 'src_image')
         projector_base = CircularSymmetryProjectorBase(self.max_radius, num_angular_sectors=self.num_angular_sectors, num_radial_sectors=self.num_radial_sectors, oversampling_scale=1)
         radial_profile_image = ie.create_hyperstack(width=src_image.get_width(), height=src_image.get_height(), num_channels=projector_base.num_radial_sectors, num_slices=1, num_frames=1, pixel_type=PixelType.F32)
+        angular_variances_image = ie.create_hyperstack(width=src_image.get_width(), height=src_image.get_height(), num_channels=projector_base.num_radial_sectors, num_slices=1, num_frames=1, pixel_type=PixelType.F32)
         for radius_index in range(projector_base.num_radial_sectors):
             circle_stack = ie.create_hyperstack(width=src_image.get_width(), height=src_image.get_height(), num_channels=projector_base.num_angular_sectors, num_slices=1, num_frames=1, pixel_type=PixelType.F32)
             for angular_sector_index in range(projector_base.num_angular_sectors):
@@ -184,6 +189,15 @@ class CircularSymmetryDetector:
                 circle_stack.set_image(projection, frame_index=0, slice_index=0, channel_index=angular_sector_index)
             # compute the image in which each pixel contains the mean value of the source image in the given circular region around the pixel
             circle_mean_image = ie.compute_mean(StackImageFeeder(circle_stack))
-
             radial_profile_image.set_image(circle_mean_image, frame_index=0, slice_index=0, channel_index=radius_index)
-        return radial_profile_image
+
+            square_diffs_stack = circle_stack  # note : circle_stack is reused to store square_diffs_stack
+            for angular_sector_index in range(projector_base.num_angular_sectors):
+                angular_sector_image = circle_stack.get_image(frame_index=0, slice_index=0, channel_index=angular_sector_index)
+                diff_image = ie.subtract(angular_sector_image, circle_mean_image)
+                square_diff_image = ie.multiply(diff_image, diff_image)
+                square_diffs_stack.set_image(square_diff_image, frame_index=0, slice_index=0, channel_index=angular_sector_index)
+            variance_image = ie.compute_mean(StackImageFeeder(square_diffs_stack))
+            angular_variances_image.set_image(variance_image, frame_index=0, slice_index=0, channel_index=radius_index)
+        angular_variance_avg_image = ie.compute_mean(StackImageFeeder(angular_variances_image))
+        return radial_profile_image, angular_variance_avg_image

src/lipase/imageengine.py

@@ -170,7 +170,6 @@ class IImageFeeder(ABC):
         creates an hyperstack from this image feeder
         :rtype IHyperStack
         """
-
         it = iter(self)
 
         image = it.next()
@@ -235,6 +234,9 @@ class StackImageFeeder(IImageFeeder):
             raise StopIteration
         else:
             image = self.hyperstack.get_image(frame_index=self.next_frame_index, slice_index=self.next_slice_index, channel_index=self.next_channel_index)
+            # mean_value = image.get_mean_value()
+            # print("mean_value", mean_value)
+            # assert mean_value > 0.001
         
         # compute next image index
         self.next_slice_index += 1
@@ -380,7 +382,17 @@ class IImageEngine(ABC):
         :param IImage image2:
         :rtype IImage: image1-image2
         """
-        
+
+    @abc.abstractmethod
+    def multiply(self, image1, image2):
+        """
+        computes the difference image1-image2
+
+        :param IImage image1:
+        :param IImage image2:
+        :rtype IImage: image1*image2
+        """
+
     @abc.abstractmethod
     def divide(self, numerator_image, denominator_image):
         """

src/lipase/imagej/ijimageengine.py

@@ -290,6 +290,11 @@ class IJImageEngine(IImageEngine):
         result_image = ic.run("Subtract create float", image1.ij_image, image2.ij_image)
         return IJImage(self, result_image)
 
+    def multiply(self, image1, image2):
+        ic = ImageCalculator()
+        result_image = ic.run("Multiply create float", image1.ij_image, image2.ij_image)
+        return IJImage(self, result_image)
+
     def divide(self, numerator_image, denominator_image):
         ic = ImageCalculator()
         result_image = ic.run("Divide create float", numerator_image.ij_image, denominator_image.ij_image)