+# -*- coding: utf-8 -*-
#!/usr/bin/env python
# Copyright (C) 2010 Monash University
"""
from PIL import Image
-
import sys, numpy, os, random
from scipy import ndimage
#
# 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
# ==========================================
# 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)
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
# =======================
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):
"""
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()
# 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)))
# (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
)
if Show_infection:
a_output[infection,0] += 128
+ #a_output[stain, 1] += 128
+
# Hough transform
status('hough')
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')
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)
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:
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) ]]
--- /dev/null
+% Hough Transform for circles.
+% img : logical image, zeros are the circles edges
+% radius_range is a vector : [min, max]
+function [acc_votes, acc_radius] = HTCircle(img, radius_range)
+
+ [rows, columns] = size(img);
+
+ acc_votes = zeros(rows, columns);
+ acc_radius = zeros(rows, columns);
+
+ phi_range = linspace(0, 2 * pi, 100);
+ phi_range(end) = [];
+
+ indexes_non_zero = find(img)';
+ s = size(img);
+
+ for r = radius_range(1):radius_range(2)
+ acc = zeros(rows, columns);
+ for i = indexes_non_zero
+ [y, x] = ind2sub(s, i);
+ for phi = phi_range
+ x0 = round(x - r * cos(phi));
+ y0 = round(y - r * sin(phi));
+ if x0 < columns && x0 > 0 && y0 < rows && y0 > 0
+ acc(y0, x0) = acc(y0, x0) + 1;
+ end
+ end
+ end
+
+ for x = 1:columns
+ for y = 1:rows
+ if acc(y, x) > acc_votes(y, x)
+ acc_votes(y, x) = acc(y, x);
+ acc_radius(y, x) = r;
+ end
+ end
+ end
+ end
+
+ acc_votes = imgaussfilt(acc_votes, 1.0);
+end
+
--- /dev/null
+% Hough Transform for circles with gradient.
+% img : gradient image
+% radius_range is a vector : [min, max]
+function [acc_votes, acc_radius] = HTCircleGradient(img, radius_range)
+
+ [rows, columns] = size(img);
+
+ acc_votes = zeros(rows, columns);
+ acc_radius = zeros(rows, columns);
+
+ phi_range = linspace(0, 2 * pi, 100);
+ phi_range(end) = [];
+
+ indexes_non_zero = find(img)';
+ s = size(img);
+
+ for r = radius_range(1):radius_range(2)
+ acc = zeros(rows, columns);
+ for i = indexes_non_zero
+ [y, x] = ind2sub(s, i);
+ for phi = phi_range
+ x0 = round(x - r * cos(phi));
+ y0 = round(y - r * sin(phi));
+ if x0 < columns && x0 > 0 && y0 < rows && y0 > 0
+ acc(y0, x0) = acc(y0, x0) + img(y, x);
+ end
+ end
+ end
+
+ for x = 1:columns
+ for y = 1:rows
+ if acc(y, x) > acc_votes(y, x)
+ acc_votes(y, x) = acc(y, x);
+ acc_radius(y, x) = r;
+ end
+ end
+ end
+ end
+
+ acc_votes = imgaussfilt(acc_votes, 1.0);
+end
+
--- /dev/null
+% Hough Transform for circles.
+% x_dir : the x component of the normalized gradient
+% y_dir : the y component of the normalized gradient
+% edges : a binary image with edge = 1
+% radius_range is a vector : [min, max]
+function [acc_votes, acc_radius] = HTCircleMa(x_dir, y_dir, edges, radius_range)
+
+ [rows, columns] = size(edges);
+
+ acc_votes = zeros(rows, columns);
+ acc_radius = zeros(rows, columns);
+
+ indexes_edges = find(edges)';
+ s = size(edges);
+
+ for r = radius_range(1):radius_range(2)
+ acc = zeros(rows, columns);
+ for i = indexes_edges
+ [y, x] = ind2sub(s, i);
+
+ x0 = round(x + x_dir(i) * r);
+ y0 = round(y + y_dir(i) * r);
+
+ if x0 < columns && x0 > 0 && y0 < rows && y0 > 0
+ acc(y0, x0) = acc(y0, x0) + 1;
+ end
+ end
+
+ for x = 1:columns
+ for y = 1:rows
+ if acc(y, x) > acc_votes(y, x)
+ acc_votes(y, x) = acc(y, x);
+ acc_radius(y, x) = r;
+ end
+ end
+ end
+ end
+
+ acc_votes = imgaussfilt(acc_votes, 1.0);
+end
+
--- /dev/null
+% Hough Transform for circles.
+% img : logical image, zeros are the circles edges
+% radius_range is a vector : [min, max]
+function [acc_votes, acc_radius1, acc_radius2, acc_alpha] = HTEllipse(img, radius_range)
+
+ [rows, columns] = size(img);
+
+ acc_votes = zeros(rows, columns);
+ acc_radius1 = zeros(rows, columns);
+ acc_radius2 = zeros(rows, columns);
+ acc_alpha = zeros(rows, columns);
+
+ alpha_range = linspace(0, pi / 2, 8);
+ alpha_range(end) = [];
+
+ phi_range = linspace(0, 2 * pi, 32);
+ phi_range(end) = [];
+
+ indexes_non_zero = find(img)';
+ s = size(img);
+
+ for alpha = alpha_range
+ transform_mat = [cos(alpha) -sin(alpha); sin(alpha) cos(alpha)];
+ for r1 = radius_range(1):radius_range(2)
+ for r2 = radius_range(1):radius_range(2)
+ acc = zeros(rows, columns);
+ for i = indexes_non_zero
+ [y, x] = ind2sub(s, i);
+ for phi = phi_range
+ p0 = [x; y] - transform_mat * [r1 * cos(phi); r2 * sin(phi)];
+ x0 = round(p0(1));
+ y0 = round(p0(2));
+
+ if x0 < columns && x0 > 0 && y0 < rows && y0 > 0
+ acc(y0, x0) = acc(y0, x0) + 1;
+ end
+ end
+ end
+
+ for x = 1:columns
+ for y = 1:rows
+ if acc(y, x) > acc_votes(y, x)
+ acc_votes(y, x) = acc(y, x);
+ acc_radius1(y, x) = r1;
+ acc_radius2(y, x) = r2;
+ acc_alpha(y, x) = alpha;
+ end
+ end
+ end
+ end
+ end
+ end
+
+ acc_votes = imgaussfilt(acc_votes, 1.0);
+end
+
--- /dev/null
+% Hough Transform for circles.
+% img : logical image, zeros are the circles edges
+% radius_range is a vector : [min, max]
+function [acc_votes, acc_radius1, acc_radius2, acc_alpha] = HTEllipse2(x_dir, y_dir, edges, radius_range)
+
+ [rows, columns] = size(edges);
+
+ acc_votes = zeros(rows, columns);
+ acc_radius1 = zeros(rows, columns);
+ acc_radius2 = zeros(rows, columns);
+ acc_alpha = zeros(rows, columns);
+
+ alpha_range = linspace(0, pi / 2, 8);
+ alpha_range(end) = [];
+% alpha_range = pi/8;
+
+ indexes_edges = find(edges)';
+ s = size(edges);
+
+ for alpha = alpha_range
+ rot_mat = [cos(alpha) -sin(alpha); sin(alpha) cos(alpha)];
+ rot_mat_inv = [cos(-alpha) -sin(-alpha); sin(-alpha) cos(-alpha)];
+ for r1 = radius_range(1):radius_range(2)
+ for r2 = radius_range(1):radius_range(2)
+ acc = zeros(rows, columns);
+ for i = indexes_edges
+ [y, x] = ind2sub(s, i);
+ v = [-x_dir(i); -y_dir(i)];
+ u = rot_mat_inv * v;
+ if u(1) == 0
+ if u(2) > 0
+ theta = pi / 2.;
+ else
+ theta = 3. * pi / 2.;
+ end
+ elseif u(2) == 0
+ if u(1) > 0
+ theta = 0;
+ else
+ theta = pi;
+ end
+ elseif u(2) < 0
+ t = u(2) / u(1);
+ theta = 2*pi + 2 * atan(...
+ (-r1 - r2*t*sqrt((r1^2+r2^2*t^2)/(r2^2*t^2)))...
+ /(r2*t));
+ else
+ t = u(2) / u(1);
+ theta = 2 * atan(...
+ (-r1 + r2*t*sqrt((r1^2+r2^2*t^2)/(r2^2*t^2)))...
+ /(r2*t));
+ end
+ C = round([x; y] - rot_mat * [r1*cos(theta); r2*sin(theta)]);
+ if C(1) < columns && C(1) > 0 && C(2) < rows && C(2) > 0
+ acc(C(2),C(1)) = acc(C(2),C(1)) + 1;
+ end
+ end
+
+ for x = 1:columns
+ for y = 1:rows
+ if acc(y, x) > acc_votes(y, x)
+ acc_votes(y, x) = acc(y, x);
+ acc_radius1(y, x) = r1;
+ acc_radius2(y, x) = r2;
+ acc_alpha(y, x) = alpha;
+ end
+ end
+ end
+ end
+ end
+ end
+
+ acc_votes = imgaussfilt(acc_votes, 1.0);
+end
+
--- /dev/null
+% center : [x y]
+function draw_ellipses(axe, ellipses, plot_3d, plot_center)
+ if ~exist('plot_3d', 'var')
+ plot_3d = false;
+ end
+
+ if ~exist('plot_center', 'var')
+ plot_center = false;
+ end
+
+ hold on;
+
+ if plot_3d
+ z_space = [0 1];
+ else
+ z_space = [0];
+ end
+
+ for i = 1:length(ellipses)
+ e = ellipses{i};
+ transform_mat = [cos(e.alpha) -sin(e.alpha); sin(e.alpha) cos(e.alpha)];
+
+ p_previous = [];
+ for phi = linspace(0, 2 * pi, 36);
+ for z = z_space
+ p = [e.x0; e.y0] + transform_mat * [e.r1 * cos(phi); e.r2 * sin(phi)]; % p : [x y]
+ if ~isempty(p_previous)
+ line = plot3(axe, [p(1) p_previous(1)], [p(2) p_previous(2)], [z z], 'LineWidth', 1, 'Color', 'yellow');
+ line.Color(4) = 0.5;
+ end
+ if z == 1
+ line = plot3(axe, [p(1) p(1)], [p(2) p(2)], [0 1], 'LineWidth', 0.5, 'Color', 'yellow');
+ line.Color(4) = 0.5;
+ end
+ end
+ p_previous = p;
+ end
+
+ if plot_center
+ plot(e.x0, e.y0, 'r.', 'MarkerSize', 20);
+ end
+ end
+end
+
--- /dev/null
+function [m] = expand_matrix(m)
+ m = [m(:,1) m m(:,end)];
+ m = [m(1,:); m; m(end,:)];
+end
\ No newline at end of file
--- /dev/null
+function [ellipses] = find_ellipses(acc_votes, acc_radius_1, acc_radius_2, acc_alpha)
+ if ~exist('acc_radius_2', 'var')
+ acc_radius_2 = acc_radius_1;
+ end
+
+ if ~exist('acc_alpha', 'var')
+ acc_alpha = zeros(size(acc_votes));
+ end
+
+ ellipses = {}; % struct('x0', 'y0', 'r1', 'r2', 'alpha')
+ max_votes = max(acc_votes(:));
+ min_votes = min(acc_votes(:));
+ threshold_votes = (max_votes - min_votes) * 0.7 + min_votes;
+
+ acc_votes_suppressed = acc_votes;
+ while true
+ [votes, index_max] = max(acc_votes_suppressed(:));
+ if votes <= threshold_votes
+ break;
+ end
+
+ acc_votes_suppressed(index_max) = 0; % Suppress the vote.
+
+ [y, x] = ind2sub(size(acc_votes_suppressed), index_max);
+
+ % If the center ellipse isn't in a known ellipse then add it to 'ellipses'.
+ accept_ellipse = true;
+ for i = 1:length(ellipses)
+ ellipse = ellipses{i};
+% module = sqrt((x - ellipse.x)^2 + (y - ellipse.y)^2);
+% phi = asin((y - ellipse.y) / module);
+%
+% p = [cos(a) -sin(a); sin(a) cos(a)] * [ellipse.r1 * cos(phi); ellipse.r2 * sin(phi)];
+% module_ellipse = sqrt(p(1)^2 + p(2)^2);
+
+ a = ellipse.alpha;
+ x0 = ellipse.x0;
+ y0 = ellipse.y0;
+ r1 = ellipse.r1;
+ r2 = ellipse.r2;
+ n = ((x - x0) * cos(a) + (y - y0) * sin(a))^2 / r1^2 + ((x - x0) * sin(a) - (y - y0) * cos(a))^2 / r2^2;
+ if n <= 1
+ accept_ellipse = false;
+ break;
+ end
+ end
+
+ if accept_ellipse
+ e = struct(...
+ 'x0', x,...
+ 'y0', y,...
+ 'r1', acc_radius_1(index_max),...
+ 'r2', acc_radius_2(index_max),...
+ 'alpha', acc_alpha(index_max));
+ ellipses{end+1} = e;
+ fprintf('Ellipse : (x0: %.1f, y0: %.1f, r1: %.1f, r2: %.1f, alpha: %.2f)\n', e.x0, e.y0, e.r1, e.r2, e.alpha);
+ end
+ end
+end
+
--- /dev/null
+function [hlink] = show_figures(name, images, original, ellipses)
+ [hlink, axes] = show_votes_figure(name, images);
+ draw_ellipses(axes(4), ellipses, true);
+ savefig(['results/' name '.fig']);
+
+ show_image_ellipses(name, original, ellipses);
+ savefig(['results/' name ' - output.fig']);
+end
+
--- /dev/null
+function show_image_ellipses(name, img, ellipses)
+ f = figure;
+ f.Position = [10 10 700 700];
+ f.Name = name;
+ grayscale = linspace(0, 1, 255);
+ colormap([grayscale' grayscale' grayscale']);
+ image(img);
+ axis image;
+ draw_ellipses(f.CurrentAxes, ellipses, false, true);
+end
+
--- /dev/null
+function [hlink, axes] = show_votes_figure(name, images)
+ f = figure;
+ f.Name = name;
+ f.Position = [100 100 1300 900];
+
+ norm = @(A) (A - min(A(:)))/(max(A(:)) - min(A(:)));
+
+ for i = 1:size(images, 2)
+ axes(i) = subplot(2, 2, i);
+ %colormap(hsv);
+ surf(norm(double(images{i}{1})), 'FaceColor','interp', 'FaceLighting', 'gouraud', 'EdgeColor', 'none');
+ camlight left;
+ shading interp;
+ colormap(jet(10))
+ axis tight;
+ title(images{i}{2});
+ view(axes(i), [-39.5 70]);
+ end
+
+ hlink = linkprop(axes, {'CameraPosition', 'CameraUpVector'});
+
+ rotate3d on;
+end
+
+
+
--- /dev/null
+clear;
+
+images_path = 'D:\Master Thesis\src\Tests_hough\images\';
+image_names = {
+% 'heap_1.png'
+%
+% 'rbc_single.png'
+% 'rbc_single_blurred.png'
+
+% 'rbc_single_oblong_1.png'
+ 'rbc_single_oblong_2.png'
+% 'rbc_single_oblong_3.png'
+%
+% 'two_rbc_1.png'
+% 'two_rbc_2.png'
+% 'two_rbc_3.png'
+% 'wbc_1.png'
+% 'wbc_2.png'
+ };
+
+for k = 1:length(image_names)
+ image_name = image_names{k};
+
+% test_hough_simple(images_path, image_name);
+% close all;
+
+% test_hough_gradient(images_path, image_name, false)
+% close all;
+
+% test_hough_ma(images_path, image_name)
+% close all;
+
+% test_hough_ellipses(images_path, image_name);
+% close all;
+
+ test_hough_ellipses2(images_path, image_name);
+% close all;
+end
\ No newline at end of file
--- /dev/null
+function test_hough_ellipses(images_path, image_name)
+ i = imread([images_path image_name]);
+ green = i(:,:,2);
+ green_filtered = imgaussfilt(double(green), 1.0);
+ edges = edge(green_filtered);
+
+ disp('Applying Houg...');
+ tic
+ [acc_votes, acc_radius1, acc_radius2, acc_alpha] = HTEllipse(edges, [20, 40]);
+ toc
+
+ disp('Finding ellipses...');
+ tic
+ ellipses = find_ellipses(acc_votes, acc_radius1, acc_radius2, acc_alpha);
+ toc
+
+ show_figures(['Hough ellipses - ' image_name], {...
+ {green, 'original'},...
+ {green_filtered, 'original filtered with a gaussian, sigma=1'},...
+ {edges, 'edges (sobel)'},...
+ {expand_matrix(acc_votes), 'votes'}}, green, ellipses);
+end
\ No newline at end of file
--- /dev/null
+function test_hough_ellipses2(images_path, image_name)
+ i = imread([images_path image_name]);
+
+ green = i(:,:,2);
+ green_filtered = imgaussfilt(green, 1.0);
+ otsu = ~im2bw(green_filtered, graythresh(green_filtered));
+
+ se = strel('rectangle', [3 3]);
+ edges = imdilate(otsu, se) - imerode(otsu, se);
+
+ s = [-1 0 1; -2 0 2; -1 0 1];
+ x_edge = conv2(double(green_filtered), s, 'valid');
+ y_edge = conv2(double(green_filtered), s', 'valid');
+
+ magnitudes = sqrt(x_edge .^ 2 + y_edge .^ 2);
+
+ % Unit gradient.
+ x_dir = x_edge ./ magnitudes;
+ y_dir = y_edge ./ magnitudes;
+
+ disp('Applying Houg...');
+ tic
+ [acc_votes, acc_radius1, acc_radius2, acc_alpha] = HTEllipse2(x_dir, y_dir, edges(2:end-1, 2:end-1), [20, 40]);
+ toc
+
+ disp('Finding ellipses...');
+ tic
+ ellipses = find_ellipses(acc_votes, acc_radius1, acc_radius2, acc_alpha);
+ toc
+
+ name = ['Hough ellipses 2 - ' image_name];
+ show_figures(name, {...
+ {green_filtered, 'original filtered with a gaussian, sigma=1'},...
+ {edges, 'edges (otsu -> dilate - erode)'},...
+ {expand_matrix(x_edge), 'composante x du gradient'},...
+ {expand_matrix(acc_votes), 'votes'}}, green, ellipses);
+
+ f = figure;
+ img = green_filtered;
+ grayscale = linspace(0, 1, 255);
+ colormap([grayscale' grayscale' grayscale']);
+ img(~edges) = 0;
+ image(img);
+ axis image;
+ hold on;
+ [X, Y] = meshgrid(1:size(green, 2), 1:size(green, 1));
+ quiver(X, Y, expand_matrix(x_dir), expand_matrix(y_dir), 'Color', 'yellow');
+ savefig(['results/' name ' - gradient.fig']);
+end
\ No newline at end of file
--- /dev/null
+function test_hough_gradient(images_path, image_name, use_edges_mask)
+ i = imread([images_path image_name]);
+
+ green = i(:,:,2);
+ green_filtered = imgaussfilt(green, 1.0);
+
+ edges = zeros(size(green));
+ if use_edges_mask
+ otsu = ~im2bw(green_filtered, graythresh(green_filtered));
+ se = strel('rectangle', [3 3]);
+ edges = imdilate(otsu, se) - imerode(otsu, se);
+ end
+
+ s = [-1 0 1; -2 0 2; -1 0 1];
+ x_edge = conv2(cast(green_filtered, 'double'), s , 'valid');
+ y_edge = conv2(cast(green_filtered, 'double'), s', 'valid');
+
+ magnitudes = sqrt(x_edge .^ 2 + y_edge .^ 2);
+
+ if use_edges_mask
+ magnitudes(~edges(2:end-1, 2:end-1)) = 0;
+ end
+
+ disp('Applying Houg...');
+ tic
+ [acc_votes, acc_radius] = HTCircleGradient(magnitudes, [25, 35]);
+ toc
+
+ disp('Finding ellipses...');
+ tic
+ ellipses = find_ellipses(acc_votes, acc_radius);
+ toc
+
+ name = ['Hough gradient - ' image_name];
+ if use_edges_mask
+ name = [name ' (with edges)'];
+ end
+
+ show_figures(name, {...
+ {green_filtered, 'original filtered with a gaussian, sigma=1'},...
+ {edges, 'edges (otsu -> dilate - erode)'},...
+ {magnitudes, 'magnitudes'},...
+ {expand_matrix(acc_votes), 'votes'}}, green, ellipses);
+end
+
--- /dev/null
+function test_hough_ma(images_path, image_name)
+ i = imread([images_path image_name]);
+
+ green = i(:,:,2);
+ green_filtered = imgaussfilt(green, 1.0);
+ otsu = ~im2bw(green_filtered, graythresh(green_filtered));
+
+ se = strel('rectangle', [3 3]);
+ edges = imdilate(otsu, se) - imerode(otsu, se);
+
+ s = [-1 0 1; -2 0 2; -1 0 1];
+ x_edge = conv2(double(green_filtered), s, 'valid');
+ y_edge = conv2(double(green_filtered), s', 'valid');
+
+ magnitudes = sqrt(x_edge .^ 2 + y_edge .^ 2);
+
+ % Unit gradient.
+ x_dir = x_edge ./ magnitudes;
+ y_dir = y_edge ./ magnitudes;
+
+ disp('Applying Houg...');
+ tic
+ [acc_votes, acc_radius] = HTCircleMa(x_dir, y_dir, edges(2:end-1, 2:end-1), [25, 35]);
+ toc
+
+ disp('Finding ellipses...');
+ tic
+ ellipses = find_ellipses(acc_votes, acc_radius);
+ toc
+
+ name = ['Hough Ma - ' image_name];
+ show_figures(name,...
+ {{green_filtered, 'original filtered with a gaussian, sigma=1'},...
+ {edges, 'edges (otsu -> dilate - erode)'},...
+ {expand_matrix(x_edge), 'composante x du gradient'},...
+ {expand_matrix(acc_votes), 'votes'}}, green, ellipses);
+
+ f = figure;
+ img = green_filtered;
+ grayscale = linspace(0, 1, 255);
+ colormap([grayscale' grayscale' grayscale']);
+ img(~edges) = 0;
+ image(img);
+ axis image;
+ hold on;
+ [X, Y] = meshgrid(1:size(green, 2), 1:size(green, 1));
+ quiver(X, Y, expand_matrix(x_dir), expand_matrix(y_dir), 'Color', 'yellow');
+ savefig(['results/' name ' - gradient.fig']);
+end
--- /dev/null
+function test_hough_simple(images_path, image_name)
+ i = imread([images_path image_name]);
+ green = i(:, :, 2);
+ green_filtered = imgaussfilt(cast(i(:, :, 2), 'double'), 1.0);
+ edges = edge(green_filtered);
+
+ disp('Applying Houg...');
+ tic
+ [acc_votes, acc_radius] = HTCircle(edges, [20, 40]);
+ toc
+
+ disp('Finding ellipses...');
+ tic
+ ellipses = find_ellipses(acc_votes, acc_radius);
+ toc
+
+ show_figures(['Hough simple - ' image_name], {...
+ {green, 'original'},...
+ {green_filtered, 'original filtered with a gaussian, sigma=1'},...
+ {edges, 'edges (sobel)'},...
+ {acc_votes, 'votes'}}, green, ellipses);
+end
\ No newline at end of file
projects / master-thesis.git / commitdiff
