clc
close all
clear all

% im = imread('leuk.jpg');
% im = rgb2gray(im);
% imshow(im);
% return
% im = double(im);
% [x, y] = meshgrid(1:size(im, 1), 1:size(im, 2));
% 
% figure;
% mesh(x, y, im);
% 
% [dx, dy] = gradient(im);

% figure
% grd = dx.*dx + dy.*dy;
% mesh(x, y, grd);

im100 = imread('leuk100.bmp');
im50 = imread('leuk50.bmp');

im100 = rgb2gray(im100);
im50 = rgb2gray(im50);

im100 = double(im100);
im50 = double(im50); 

[dx50, dy50] = gradient(im50);
[dx100, dy100] = gradient(im100);

col_start = 36; col_end = 58;
line = 90;

g50 = dx50.*dx50 + dy50.*dy50;
g100 = dx100.*dx100 + dy100.*dy100;

g100 = g100(line, col_start:col_end);
g50 = g50(line/2, col_start/2:col_end/2);

plot(1:length(g100), g100, 'r', 1:2:length(g100), g50, 'b');
legend('Full size', 'Half size');
title('Square of the gradient magnitude');
grid on;

gaus = [1 2 1]/4;
img = im100;
figure;

colors = 'bgry';
for s = 1:4    
    
    [dx, dy] = gradient(img);  
    grad = dx.*dx + dy.*dy;
    
    g_line{s} = grad(line, col_start:col_end);
    plot(1:length(g_line{s}), g_line{s}, colors(s));
    hold on;
    img = conv2(conv2(img, gaus,'same'), gaus','same');             
    
end

