Add some tests on different variations of Hough algorithm.
[master-thesis.git] / Documentation / ma_autocount / plasmodium-autocount-0.1.py
index b2fa3d2..d97e059 100644 (file)
@@ -1,3 +1,4 @@
+# -*- coding: utf-8 -*-
 #!/usr/bin/env python
 
 #    Copyright (C) 2010 Monash University
@@ -23,7 +24,6 @@ Identify red blood cells and red blood cells infected with plasmodium in a set o
 """
 
 from PIL import Image
-
 import sys, numpy, os, random
 from scipy import ndimage
 
@@ -58,20 +58,32 @@ N_processors = 1
 #
 # What to show in diagnostic output images
 
-Show_edges = False
+Show_edges = True
 Show_stain = False
 Show_infection = False
 
-Show_hand_marks = True
+Show_hand_marks = False
 Show_computer_marks = True
-Show_cell_shapes = False
+Show_cell_shapes = True
 
 
 # =================
 # Image scaling
 
-Resize_width = 680
-Resize_height = 512
+Resize = False
+
+# Tailles pour les images dans '08.10.2015'
+Resize_width = 1296 #680
+Resize_height = 972 #512
+
+# Tailles pour les images dans '22.05.2015' (50%)
+#Resize_width = 1296
+#Resize_height = 972
+
+
+# [GB] Hemato
+# Resize_width = 2864
+# Resize_height = 2885
 
 
 # ==========================================
@@ -109,7 +121,7 @@ Infection_pixels_required = 1
 
 # Size of cells, for Hough transform
 Hough_min_radius = 5  
-Hough_max_radius = 20 
+Hough_max_radius = 32 # 20
 
 # Portion of cell used for classification as infected/weird
 # (may be abutting an infected cell, so be cautious)
@@ -117,7 +129,7 @@ Hough_max_radius = 20
 Cell_inside_radius = 0.9
 
 # Do not allow another cell center to be detected within this radius
-Cell_suppression_radius = 1.25
+Cell_suppression_radius = 1.05 # 1.25
 
 
 # =======================
@@ -127,17 +139,13 @@ Cell_suppression_radius = 1.25
 Max_dark_stain_proportion = 0.728
 
 # Minimum/maximum radius (as determined by Hough transform)
-Min_cell_radius = 9
-Max_cell_radius = 20
+Min_cell_radius = 11 # 9
+Max_cell_radius = 32 # 20
 
 # Proportion of radius that center of mass can be away from Hough determined center
 Max_offcenter = 0.5
 
 
-
-
-
-
 _offset_cache = { }
 def offsets(inner_radius, outer_radius):
     """
@@ -244,14 +252,14 @@ def process(filename_in, filename_out, filename_coords, name):
     print 'Processing', filename_in, '->', filename_out
 
     # Resize (mostly done for speed increase)
-    i = i.resize((Resize_width,Resize_height), Image.ANTIALIAS)
+    i = i.resize((Resize_width,Resize_height), Image.ANTIALIAS) if Resize else i
 
     # Convert to numpy array
-    a = numpy.asarray(i).astype('float64')
+    a = numpy.asarray(i).astype('float64') # [GB] Shape: (512, 680, 3)
     height = a.shape[0]
     width = a.shape[1]
 
-    a_gray = numpy.sum(a, 2)
+    a_gray = numpy.sum(a, 2) # [GB] Somme de chaque composante RGB (non-normalisé)
     
     # Make a copy of the array to doodle various diagnostic markings on
     a_output = a.copy()    
@@ -262,33 +270,39 @@ def process(filename_in, filename_out, filename_coords, name):
 
 
     # Denoise/unsharp mask
-    a_enhanced = gaussian_blur(a, Enhance_a) - Enhance_b*gaussian_blur(a,Enhance_c)
-    
+    # [GB] (Concerne les trois channels RGB).
+    a_enhanced = gaussian_blur(a, Enhance_a) - Enhance_b * gaussian_blur(a, Enhance_c)
     
     # Split into foreground and background using k-medians
     status('fg/bg')
     
+    # [GB] Chaque canal de couleur est linéarisé dans un tableau.
     a_raveled = numpy.reshape(a_enhanced,(width*height,a_enhanced.shape[2]))
-    average = numpy.average(a_raveled, axis=0)
-    # Initial guess
-    color_bg = average*1.5
-    color_fg = average*0.5
     
+    # [GB] Moyenne de chaque canal.
+    average = numpy.average(a_raveled, axis=0) 
+    
+    # Initial guess of the medians
+    color_bg = average * 1.5 # [GB] Plus clair
+    color_fg = average * 0.5 # [GB] Plus sombre
+    
+    # Fifty times: (0 to 4).
     for i in xrange(5):
-        d_bg = a_raveled-color_bg[None,:]
+        d_bg = a_raveled - color_bg[None,:]
         d_bg *= d_bg
-        d_bg = numpy.sum(d_bg,1)
-        d_fg = a_raveled-color_fg[None,:]
+        d_bg = numpy.sum(d_bg, 1)  # [GB] Somme des carrés des couleurs
+        d_fg = a_raveled - color_fg[None,:]
         d_fg *= d_fg
-        d_fg = numpy.sum(d_fg,1)
-        fg = d_fg*Background_weight < d_bg
-        color_fg = numpy.median(a_raveled[fg,:],axis=0)
+        d_fg = numpy.sum(d_fg, 1)
+        fg = d_fg * Background_weight < d_bg
+        color_fg = numpy.median(a_raveled[fg,:],axis=0) # [GB] axis = 0 -> par colonne (ordonnées). Le résultat est alors 3 valeurs, une pour chaque couleur
         color_bg = numpy.median(a_raveled[~fg,:],axis=0)
-        
-    fg = numpy.reshape(fg, (height,width))
+
+    fg = numpy.reshape(fg, (height,width)) # [GB] Image binaire du premier plan.
     
-    edge = ndimage.maximum_filter(fg,size=3) != ndimage.minimum_filter(fg,size=3)
+    edge = ndimage.maximum_filter(fg,size=3) != ndimage.minimum_filter(fg,size=3) # [GB] image binaire des contours du premier plan.
 
+    # [GB] délinéarisation des tableaux résultant du k-median (niveau de gris).
     d_fg = numpy.sqrt(numpy.reshape(d_fg, (height,width)))
     d_bg = numpy.sqrt(numpy.reshape(d_bg, (height,width)))
     
@@ -298,10 +312,15 @@ def process(filename_in, filename_out, filename_coords, name):
     # (should be invariant under changes in brightness) 
     mag_bg = numpy.average(color_bg)
     
+    # [GB] Les "dark stain" :
+    # * sont d'une intensité plus élevée que la moyenne du bg multiplié par un facteur (Dark_stain_level)
+    # * les valeurs rouges sont plus petites que la mediane du premier-plan pour les rouges
+    # * les valeurs vertes sont plus petites que la médiane du premier-plan pour les verts
+    # * appartiennent au foreground
     dark_stain = (
         (d_fg > mag_bg*Dark_stain_level) &
-        (a_enhanced[:,:,0] < color_fg[0]) & 
-        (a_enhanced[:,:,1] < color_fg[1]) &
+        (a_enhanced[:,:,0] < color_fg[0]) & # [GB] Red channel
+        (a_enhanced[:,:,1] < color_fg[1]) & # [GB] Green channel
         fg
     )
     
@@ -332,6 +351,8 @@ def process(filename_in, filename_out, filename_coords, name):
     if Show_infection:
         a_output[infection,0] += 128
         
+    #a_output[stain, 1] += 128
+        
         
     # Hough transform
     status('hough')
@@ -409,9 +430,11 @@ def process(filename_in, filename_out, filename_coords, name):
     candidates = [ ]
     for y in xrange(height):
         for x in xrange(width):
-            if hough[y,x] > thresh:
-                candidates.append((hough[y,x],y,x))
-    candidates.sort(reverse=True)
+            if hough[y, x] > thresh:
+                candidates.append((hough[y,x], y, x))
+    candidates.sort(reverse = True)
+    hough
+    
     
     status('classify')
     
@@ -420,21 +443,20 @@ def process(filename_in, filename_out, filename_coords, name):
     
     suppress = numpy.zeros((height,width), 'bool')
     claim_strength = numpy.zeros((height,width), 'float64') + 1e30
-    cells = [ ] # [(y,x)]
+    cells = [] # [(y,x)]
     n = 0
     n_infected = 0
-    n_late = 0
     for _,y,x in candidates:
         radius = hough_radius[y,x]        
-
+        
         # Make sure candidate is not near a previously detected cell,
         # and does not abut the edge of the image        
         if not suppress[y,x] and \
            y >= radius and x >= radius and y < height-radius and x < width-radius:
-            coords = clipped_coordinates(0,radius*Cell_inside_radius, y,x,height,width)
+            coords = clipped_coordinates(0, radius * Cell_inside_radius, y, x, height, width)
             
             #Only include foreground
-            this_fg = fg[coords[:,0],coords[:,1]]
+            this_fg = fg[coords[:,0], coords[:,1]] # [GB] Booleans
             coords = coords[this_fg]
             
             this_strength = ((coords[:,0]-y)**2 + (coords[:,1]-x)**2)
@@ -443,13 +465,12 @@ def process(filename_in, filename_out, filename_coords, name):
                 this_strength
             )
             
-            cells.append((y,x, coords,this_strength))
+            cells.append((y, x, coords, this_strength))
 
-           
+       
         # Suppress detection of cells too close to here
-        # (note: we do this irrespective of whether we treated this location as a hit)
-        
-        coords = clipped_coordinates(0,radius*Cell_suppression_radius, y,x,height,width)
+        # (note: we do this irrespective of whether we treated this location as a hit)        
+        coords = clipped_coordinates(0,radius * Cell_suppression_radius, y,x,height,width)
         suppress[coords[:,0],coords[:,1]] = True
         
         #for oy, ox in offsets(0, radius*Cell_suppression_radius): #suppression_offsets:
@@ -525,17 +546,18 @@ def process(filename_in, filename_out, filename_coords, name):
             color = (0,0,255)
             file_coords.write('%d,%d,Uninfected\n' % (x,y))
         
+        if Show_cell_shapes:            
+            a_output[coords[:,0],coords[:,1]] += [[ random.random()*64+32 for i in (0,1,2) ]]
+            
         # Mark detected cell on output image
         if Show_computer_marks:
-            set(a_output,y,x,color)
+            set(a_output, y, x, color)
             for i in xrange(1,3):
                 set(a_output,y-i,x,color)  
                 set(a_output,y+i,x,color)  
                 set(a_output,y,x-i,color)  
                 set(a_output,y,x+i,color)
 
-        if Show_cell_shapes:            
-            a_output[coords[:,0],coords[:,1]] += [[ random.random()*64+32 for i in (0,1,2) ]]