clc;
close all;
clear all;

N = 100;

m1 = [0; 0];
m2 = [5; 5];

var1 = 0.5;
var2 = 3;

S1 = sqrt(var1)*eye(2);
invS1 = inv(S1);
S2 = sqrt(var2)*eye(2);
invS2 = inv(S2);

x = linspace(m1(1), m2(1), N);
y = linspace(m1(2), m2(2), N);

class = zeros(N, 1);
d1 = zeros(N, 1);
d2 = d1;

for k = 1:N
    
    p = [x(k); y(k)];
    
    d1(k) = sqrt((p - m1)'*invS1*(p - m1));
    d2(k) = sqrt((p - m2)'*invS2*(p - m2));
    if (d1(k) < d2(k))
        class(k) = 1;
    else
        class(k) = 2;
    end
    
end

figure;

plot(x(class==1), y(class==1), 'b', x(class==2), y(class==2), 'r');
hold on;
circle(m1, var1, N, 'b');
circle(m2, var2, N, 'r');

x1 = sqrt(var1)*randn(N, 1) + m1(1);
y1 = sqrt(var1)*randn(N, 1) + m1(2);

x2 = sqrt(var2)*randn(N, 1) + m2(1);
y2 = sqrt(var2)*randn(N, 1) + m2(2);

plot(x1, y1, 'b.', x2, y2, 'r.');

figure;
plot(1:N, d1, 'b', 1:N, d2, 'r');
legend('d1', 'd2');