function [F,S,P,R] = ncut_rw_demo(lambda)
%
%NCUT_RW_DEMO Ncut and random walk demo with a toy dataset
%
% [F,S,P,R] = ncut_rw_demo(lambda);
%
% lambda 2 x 1 "Weigths" for the distance and color cues
% F      2 x 1 Cell array of distance matrices
% F{1}.data 9 x 9 Distances
% F{2}.data 9 x 9 Color differences
% S         9 x 9 Similarity matrix
% P         9 x 9 Transition matrix
% R         3 x 3 Aggregated transition matrix
%
% Markus Herrgard
% $ Revision: 1.0 $ $ Date: 10/11/01 15:44:00 $

% Distance data

%ddata = [-1 1; -1.2 1.3; 1 1; 1.2 0.9; 0.9 1; 0 -1; 0.1 -1.1; -0.1 -0.9; 0 -0.8 ];
%ddata = [0.1 0; -0.1 0; 1 1; 1.2 0.9; 0.9 1; 0 -1; 0.1 -1.1; -0.1 -0.9; 0 -0.8 ];
ddata = [-0.5 -1; -0.6 -0.9; 1 1; 1.2 0.9; 0.9 1; 0 -1; 0.1 -1.1; -0.1 -0.9; 0 -0.8 ];
%ddata = [0.5 0.5; 0.6 0.6; 1 1; 1.2 0.9; 0.9 1; 0 -1; 0.1 -1.1; -0.1 -0.9; 0 -0.8 ];
%ddata = [-1 0; 0 0; 2 0];

% Color data 

%cdata = [0.1 1 1];
cdata = [0.1 0.2 0.1 0.1 0.2 0.8 0.9 0.9 1.0];
%cdata = [1:10];

%data = init_gmm(3,0.01,10);

[n,d] = size(ddata);

% Draw the original data
figure(1);
clf;
dc = demo_colormap;
colormap(dc);
for i =1:n
    h = rectangle('Curvature',[1 1],'Position',[ddata(i,1) ddata(i,2) 0.1 0.1]);
    set(h,'FaceColor',dc(cdata(i)*10,:));
end
set(gcf,'Color',[1 1 1]);
colorbar

% Calculate the similarity matrices F1 and F2
F1 = zeros(n);
F2 = zeros(n);
for i = 1:n
    for j = 1:n
        F1(i,j) = simfn(ddata(i,:),ddata(j,:));
        F2(i,j) = simfn(cdata(i),cdata(j));
    end
end

% Calculate similarity matrix S
S = exp(lambda(1)*F1).*exp(lambda(2)*F2);

% Plot S
figure(2)
clf;
colormap(dc);
imagesc(S);
colorbar;
%axis off;
set(gcf,'Color',[1 1 1]);

% Calculate transition matrix P
D = diag(sum(S,2));
P = inv(D)*S;

% Plot P
figure(3);
clf;
colormap(dc);
imagesc(P);
colorbar;
%axis off;
set(gcf,'Color',[1 1 1]);

% Figure out eigenvectors and eigenvalues of P
[X,L] = eig(P);
[L,I] = sort(diag(L));
X = X(:,flipud(I));
L = flipud(L);

% Plot eigenvectors and eigenvalues
figure(4);
clf;
plot_ev(X,L)
set(gcf,'Color',[1 1 1]);

% Calculate the aggregated transition matrix R
U = [1 1 0 0 0 0 0 0 0 ;
    0 0 1 1 1 0 0 0 0;
    0 0 0 0 0 1 1 1 1];
U = inv(diag(sum(U,2)))*U;
V = U';
V = inv(diag(sum(V,2)))*V;
R = U*P*V;

% Plot R
figure(5);
clf;
colormap(dc);
imagesc(R);
colorbar;
set(gcf,'Color',[1 1 1]);

% Return the similarities in a cell array
% (to be used with the supervised segmentation method)
F{1}.data = F1;
F{2}.data = F2;

function data = init_gmm(ng,sig,npts)

mix = gmm(2,ng,'spherical')
mix.centres = [-1 1;0 -1;1 1];
mix.covars = ones(1,ng)*sig;
data = gmmsamp(mix,npts);

function f = simfn(x,y)

f = sum((x-y).^2);

function cm = demo_colormap()

clear hot

cm = hot;
cm = cm(linspace(1,length(hot),10),:);

function plot_ev(V,L)

K = length(L);

Kx = floor(sqrt(K));
if (Kx^2 == K)
    Ky = Kx;
else
    Ky = Kx + 1;
end
K = 3;
Kx = 1;
Ky = 3;

for i = 1:K
    subplot(Kx,Ky,i);
    hold on
    plot(1:9,V(:,i),'ko');
    title(['\lambda_' num2str(i) '=' num2str(L(i))]);
    axis([0 9 -1 1]);
    hold off
end