% 080408 gabriel@pads.ufrj.br

close all;

% Gerando Dados

clear all; K = 8; P = 100;
randn('state',0); a = randn(1,1); rand('state',0);
X0 = randn(2,K);
X=[];
for k=1:size(X0,2),   
    X = [X repmat(X0(:,k),1,P)+randn(2,P)*0.1];    
end;

figure; plot(X(1,:),X(2,:),'k.',X0(1,:),X0(2,:),'b.'); grid on; hold on; voronoi(X0(1,:),X0(2,:)); [Cn,p,Dmin] = ExemploFuncaoJ(X,X0); Dmin = Dmin/size(X,2) % Dmin = 0.0195

% Metodo 1: Linde-Buzo-Gray (Generalized Lloyd Algorithm)

C0 = mean(X,2);
D0 = mean(sum((X-repmat(C0,1,size(X,2))).^2,1));
C = C0; Di = D0; e = 1e-3; delta = 1e-5; 
n = 0; fim = 0; split = 1; 

while not(fim),
    n = n + 1; BKD(n,:) = [size(C,2) Di];
    if split,
        Phi = e*randn(size(C)); C = [C-Phi C+Phi];
    end;
    split = 0;
    [Cn,p,D] = ExemploFuncaoJ(X,C);
    for k = 1:length(p),
        C(:,k) = Cn(:,k)/p(k); p(k) = p(k)/size(X,2);
    end;
    D = D/size(X,2);
    if abs((D-Di)/Di) < delta,
        if size(C,2) < K, split = 1; else fim = 1; end;
    end;
    Di = D;
end;

figure; plot(X(1,:),X(2,:),'k.',C(1,:),C(2,:),'b.'); grid on; hold on; voronoi(C(1,:),C(2,:)); plot(C0(1,:),C0(2,:),'r.'); D % D = 0.0672

% Gerando Dados Novamente

clear all; K = 8; P = 100;
randn('state',0); a = randn(1,1); rand('state',0);
X0 = randn(2,K);
X=[];
for k=1:size(X0,2),   
    X = [X repmat(X0(:,k),1,P)+randn(2,P)*0.1];    
end;

% Metodo 2: Simulated Annealing

C0 = randn(2,K); [Cn,p,D0] = ExemploFuncaoJ(X,C0); Ci = C0; Di = D0/size(X,2);
e = 5e-2; t = 1; T0 = 5e-3; T = T0/log2(1+t); N = 10000; L = 8;
n = 0; fim = 0; Df = Di; Cf = Ci; BKD = zeros(N*(L-1),2);

while not(fim),
    n = n + 1; BKD(n,:) = [T Di];
    C = Ci + e*randn(size(Ci));
    [Cn,p,D] = ExemploFuncaoJ(X,C); D = D/size(X,2); % Notacao: D = D(C) ao inves de J = J(X)!
    if D<Df, Df = D; Cf = C; end;
    if D<Di, 
        Ci = C; Di = D;
    else
        a = rand(1,1); 
        if a < exp((Di-D)/T), Ci = C; Di = D; end;
    end;
    if (rem(n,N)==0),
        % [T Df]
        t = t+1;
        T = T0/log2(1+t);
        if t==L, fim = 1; end;
    end;    
end;
C = Cf; 

figure; plot(X(1,:),X(2,:),'k.',C(1,:),C(2,:),'b.'); grid on; hold on; voronoi(C(1,:),C(2,:)); plot(C0(1,:),C0(2,:),'r.'); Df % D = 0.0223

