now a clean trap image is computed as the median value of all traps

src/lipase/imageengine.py

@@ -113,11 +113,35 @@ class IImageFeeder(ABC):
 
     @abc.abstractmethod
     def next(self):
-        """returns the nex image in the collection
+        """returns the next image in the collection
         for iterator
         """
         pass
 
+    @abc.abstractmethod
+    def get_num_images(self):
+        pass
+
+    def create_hyperstack(self):
+        """
+        creates an hyperstack from this image feeder
+        :rtype IHyperStack
+        """
+
+        it = iter(self)
+
+        image = it.next()
+        ie = IImageEngine.get_instance()
+        hyperstack = ie.create_hyperstack(width=image.width(), height=image.height(), num_channels=1, num_slices=1, num_frames=self.get_num_images(), pixel_type=image.get_pixel_type())
+        frame_index = 0
+        for image in it:
+            hyperstack.set_image(image, frame_index=frame_index)
+            frame_index += 1
+        print(frame_index)
+        print(self.get_num_images())
+        
+        assert frame_index == self.get_num_images()
+        return hyperstack
 
 class FileImageFeeder(IImageFeeder):
 
@@ -137,6 +161,9 @@ class FileImageFeeder(IImageFeeder):
             return image
         raise StopIteration
 
+    def get_num_images(self):
+        return len(self.image_filepaths)
+
     def add_image(self, image_filepath):
         self.image_filepaths.append(image_filepath)
 
@@ -160,29 +187,42 @@ class StackImageFeeder(IImageFeeder):
         return self
 
     def next(self):
-        print(self.next_channel_index, self.next_frame_index, self.next_slice_index)
+        print("channel %d/%d frame %d/%d slice %d/%d" % (self.next_channel_index, self.hyperstack.num_channels(), self.next_frame_index, self.hyperstack.num_frames(), self.next_slice_index, self.hyperstack.num_slices()))
-        print(self.hyperstack.num_channels, self.hyperstack.num_frames, self.hyperstack.num_slices)
         if self.end_is_reached:
             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)
         
         # compute next image index
-        if self.next_slice_index < self.hyperstack.num_slices():
         self.next_slice_index += 1
-        else:
+        if self.next_slice_index == self.hyperstack.num_slices():
             self.next_slice_index = 0
-            if self.next_frame_index < self.hyperstack.num_frames():
             self.next_frame_index += 1
-            else:
+            if self.next_frame_index == self.hyperstack.num_frames():
                 self.next_frame_index = 0
-                if self.next_channel_index < self.hyperstack.num_channels():
                 self.next_channel_index += 1
-                else:
+                if self.next_channel_index == self.hyperstack.num_channels():
                     self.end_is_reached = True
 
+
+        # if self.next_slice_index < self.hyperstack.num_slices() - 1:
+        #     self.next_slice_index += 1
+        # else:
+        #     self.next_slice_index = 0
+        #     if self.next_frame_index < self.hyperstack.num_frames() - 1:
+        #         self.next_frame_index += 1
+        #     else:
+        #         self.next_frame_index = 0
+        #         if self.next_channel_index < self.hyperstack.num_channels():
+        #             self.next_channel_index += 1
+        #         else:
+        #             self.end_is_reached = True
+        print("after : channel %d/%d frame %d/%d slice %d/%d" % (self.next_channel_index, self.hyperstack.num_channels(), self.next_frame_index, self.hyperstack.num_frames(), self.next_slice_index, self.hyperstack.num_slices()))
         return image
 
+    def get_num_images(self):
+        return self.hyperstack.num_frames()
+
 
 class IImageEngine(ABC):
     """
@@ -245,6 +285,14 @@ class IImageEngine(ABC):
         :rtype IImage:
         """
 
+    @abc.abstractmethod
+    def compute_median(self, image_feeder):
+        """Compute for each pixel position the median value at this position in all input images.
+
+        :param IImageFeeder image_feeder:
+        :rtype IImage:
+        """
+
     @abc.abstractmethod
     def mean_filter(self, image, radius):
         """Each pixel becomes an average of its neighbours within the given radius

src/lipase/imagej/ijimageengine.py

@@ -7,6 +7,7 @@ from ij import IJ, ImagePlus  # pylint: disable=import-error
 from ij.measure import ResultsTable # pylint: disable=import-error
 from ij.plugin import ImageCalculator  # pylint: disable=import-error
 from ij.plugin.filter import MaximumFinder # pylint: disable=import-error
+from ij.plugin import ZProjector # pylint: disable=import-error
 from ijopencv.ij import ImagePlusMatConverter  # pylint: disable=import-error
 from ijopencv.opencv import MatImagePlusConverter  # pylint: disable=import-error
 import org.bytedeco.javacpp.opencv_core as opencv_core  # pylint: disable=import-error
@@ -239,6 +240,21 @@ class IJImageEngine(IImageEngine):
         print('max_image', max_image)
         return IJImage(self, max_image)
 
+    def compute_median(self, image_feeder):
+        """Computes for each pixel position the median value at this position in all input images.
+
+        :param IImageFeeder image_feeder:
+        :rtype IJmage:
+        """
+        hyperstack = image_feeder.create_hyperstack()
+        # https://imagej.nih.gov/ij/developer/api/ij/plugin/ZProjector.html
+        projector = ZProjector()
+        median_image = projector.run(hyperstack.hyperstack, 'median')
+        # imagej_run_image_command(image=hyperstack.hyperstack, command="Z Project...", options="projection=Median")
+        # # after median computation, the resulting image is expected to be the selected one
+        # median_image = IJ.getImage()  # get the currently selected image
+        return IJImage(self, median_image)
+
     def mean_filter(self, image, radius):
         """Implement interface method."""
         IJ.run(image.ij_image, "Mean...", "radius=%d" % radius)

src/lipase/traps_detector.py

@@ -1,10 +1,10 @@
 
 from catalog import ImageCatalog, Sequence
 # from imageengine import Toto
-from imageengine import IImageEngine
+from imageengine import IImageEngine, Aabb, StackImageFeeder
 from imageengine import PixelType
 from preprocessing import WhiteEstimator, correct_non_uniform_lighting
-
+from maxima_finder import Match
 from template_matcher import TemplateMatcher
 
 class TrapsDetector(object):
@@ -15,7 +15,7 @@ class TrapsDetector(object):
         """
         self.template_matcher = template_matcher
 
-    def compute_traps_mask(self, sequence, channel_id, trap_aabb):
+    def compute_traps_mask(self, sequence, channel_id, template_trap_aabb):
         """Remove the traps in the input sequence
 
         :param Sequence sequence:
@@ -33,10 +33,31 @@ class TrapsDetector(object):
         uniform_stack = correct_non_uniform_lighting(sequence, channel_id, white_estimator=WhiteEstimator(open_size=75, close_size=75, average_size=75))
 
         first_image = uniform_stack.get_image(frame_index=0)
-        trap_image = first_image.get_subimage(trap_aabb)
+        template_trap_image = first_image.get_subimage(template_trap_aabb)
         # ie.save_as_tiff(trap_image, '/home/graffy/Desktop/template.tiff')
 
-        matches = self.template_matcher.match_template(first_image, trap_image)
+        matches = self.template_matcher.match_template(first_image, template_trap_image)
+        num_traps_per_frame = len(matches)
+        num_traps = uniform_stack.num_frames() * num_traps_per_frame
+        
+        ie = IImageEngine.get_instance()
+        traps_stack = ie.create_hyperstack(width=template_trap_aabb.width, height=template_trap_aabb.height, num_slices=1, num_frames=num_traps, num_channels=1, pixel_type=PixelType.F32)
+
+        trap_index = 0
+        for frame_index in range(uniform_stack.num_frames()):
+            uniform_frame = uniform_stack.get_image(frame_index=frame_index)
+            for frame_trap_index in range(num_traps_per_frame):
+                match = matches[frame_trap_index]
+                trap_aabb = Aabb(x_min=match.x, y_min=match.y, x_max=match.x + template_trap_aabb.width - 1, y_max=match.y + template_trap_aabb.height - 1)
+                assert trap_aabb.width == template_trap_aabb.width
+                trap_image = uniform_frame.get_subimage(trap_aabb)
+                traps_stack.set_image(frame_index=trap_index, image=trap_image)
+                trap_index += 1
+
+        image_feeder = StackImageFeeder(traps_stack)
+
+        clean_trap_image = ie.compute_median(image_feeder)
+        ie.save_as_tiff(clean_trap_image, './clean_trap_image.tiff')
 
         #non_uniform_stack = sequence.as_stack()
         #uniform_stack = IImageEngine.get_instance().divide(non_uniform_stack, white_estimate)

tests/test0001.py

@@ -7,13 +7,13 @@
 
 import unittest  # unittest2 doesn't exist in fiji
 import sys
-from lipase.imageengine import IImageEngine, PixelType, Aabb
+from lipase.imageengine import IImageEngine, PixelType, Aabb  # pylint: disable=import-error
-from lipase.imagej.ijimageengine import IJImageEngine, IJImage
+from lipase.imagej.ijimageengine import IJImageEngine, IJImage  # pylint: disable=import-error
-from lipase.preprocessing import WhiteEstimator, correct_non_uniform_lighting
+from lipase.preprocessing import WhiteEstimator, correct_non_uniform_lighting  # pylint: disable=import-error
-from lipase.maxima_finder import MaximaFinder
+from lipase.maxima_finder import MaximaFinder  # pylint: disable=import-error
-from lipase.template_matcher import TemplateMatcher
+from lipase.template_matcher import TemplateMatcher  # pylint: disable=import-error
-from lipase.traps_detector import TrapsDetector
+from lipase.traps_detector import TrapsDetector  # pylint: disable=import-error
-from lipase.catalog import ImageCatalog, Sequence
+from lipase.catalog import ImageCatalog, Sequence  # pylint: disable=import-error
 
 class TestLipase(unittest.TestCase):
 
@@ -31,19 +31,19 @@ class TestLipase(unittest.TestCase):
         print("uninitializing TestLipase instance")
         self.catalog = None
 
-    def test_estimate_white(self):
+    # def test_estimate_white(self):
-        sequence = self.catalog.sequences['res_soleil2018/GGH/GGH_2018_cin2_phiG_I_327_vis_-40_1/Pos2']
+    #     sequence = self.catalog.sequences['res_soleil2018/GGH/GGH_2018_cin2_phiG_I_327_vis_-40_1/Pos2']
-        white_estimator = WhiteEstimator(open_size=75, close_size=75, average_size=75)
+    #     white_estimator = WhiteEstimator(open_size=75, close_size=75, average_size=75)
-        white_estimate = white_estimator.estimate_white([sequence], ['DM300_327-353_fluo'])
+    #     white_estimate = white_estimator.estimate_white([sequence], ['DM300_327-353_fluo'])
-        # find_white_reference_image(white_estimate, sequence.get_white())
+    #     # find_white_reference_image(white_estimate, sequence.get_white())
-        print(white_estimate)
+    #     print(white_estimate)
-        IImageEngine.get_instance().save_as_tiff(white_estimate, './white_estimate.tiff')
+    #     IImageEngine.get_instance().save_as_tiff(white_estimate, './white_estimate.tiff')
-        # assert False, "hellooooo"
+    #     # assert False, "hellooooo"
-        print('end of test_estimate_white')
+    #     print('end of test_estimate_white')
 
-    def test_uniform_lighting_correction(self):
+    # def test_uniform_lighting_correction(self):
-        non_uniform_sequence = self.catalog.sequences['res_soleil2018/GGH/GGH_2018_cin2_phiG_I_327_vis_-40_1/Pos0']
+    #     non_uniform_sequence = self.catalog.sequences['res_soleil2018/GGH/GGH_2018_cin2_phiG_I_327_vis_-40_1/Pos0']
-        uniform_sequence = correct_non_uniform_lighting(non_uniform_sequence, 'DM300_nofilter_vis', white_estimator=WhiteEstimator(open_size=75, close_size=75, average_size=75))
+    #     uniform_sequence = correct_non_uniform_lighting(non_uniform_sequence, 'DM300_nofilter_vis', white_estimator=WhiteEstimator(open_size=75, close_size=75, average_size=75))
 
     def test_traps_detector(self):
         # the typical value of peaks is -500 and the value between peaks is below -2500