% -------------------------------------------------------------------------
% Function: 
% [RD,CD,order]=optics(x,k)
% -------------------------------------------------------------------------
% Aim: 
% Ordering objects of a data set to obtain the clustering structure 
% -------------------------------------------------------------------------
% Input: 
% x - data set (m,n); m-objects, n-variables
% k - number of objects in a neighborhood of the selected object
% (minimal number of objects considered as a cluster)
% -------------------------------------------------------------------------
% Output: 
% RD - vector with reachability distances (m,1)
% CD - vector with core distances (m,1)
% order - vector specifying the order of objects (1,m)
% -------------------------------------------------------------------------
% Example of use:
% x=[randn(30,2)*.4;randn(40,2)*.5+ones(40,1)*[4 4]];
% [RD,CD,order]=optics(x,4)
% -------------------------------------------------------------------------
% 

function [RD,CD,order]=optics(x,k)

[m,n]=size(x);
CD=zeros(1,m);%核心距离
RD=ones(1,m)*10^10;%可达距离

% Calculate Core Distances 得到每个点的核心距离 计算每个点和其他节点的距离，然后进行排序，得到前
% k个点 故而得到核心距离 
for i=1:m	
    D=sort(dist(x(i,:),x))
    CD(i)=D(k);  
end
order=[];
seeds=[1:m];
ind=1;
while ~isempty(seeds)
    ob=seeds(ind);      
    seeds(ind)=[] ;
    order=[order ob];
    mm=max([ones(1,length(seeds))*CD(ob);dist(x(ob,:),x(seeds,:))]);
    ii=(RD(seeds))>mm;
    RD(seeds(ii))=mm(ii);
    [i1 ind]=min(RD(seeds));
end   

%RD(1)=max(RD(2:m))+.1*max(RD(2:m));
kk=RD>k;
RD(kk)=10;
plot(RD(order))
function [D]=dist(i,x)

% function: [D]=dist(i,x)
%
% Aim: 
% Calculates the Euclidean distances between the i-th object and all objects in x	 
% Input: 
% i - an object (1,n)
% x - data matrix (m,n); m-objects, n-variables	    
%                                                                 
% Output: 
% D - Euclidean distance (m,1)

[m,n]=size(x);
D=(sum((((ones(m,1)*i)-x).^2)'));
D=sqrt(D);
if n==1
   D=abs((ones(m,1)*i-x))';
end