wolffd@0: function base = som_clspread(sM,base,cldist,Ne,verbosity)
wolffd@0: 
wolffd@0: % SOM_CLSPREAD Partition the given data by flooding.
wolffd@0: %
wolffd@0: %  part = som_clspread(sM,part,cldist,[Ne],[verbos])
wolffd@0: %
wolffd@0: %  Input and output arguments ([]'s are optional):
wolffd@0: %   sM       (struct) map or data struct
wolffd@0: %            (matrix) size dlen x dim, the data set            
wolffd@0: %   base     (vector) initial partition, where if base(i) is
wolffd@0: %                      0         i should be assigned to some cluster
wolffd@0: %                      NaN       i should not be assigned to any cluster
wolffd@0: %                      otherwise i belongs to cluster base(i)
wolffd@0: %   cldist   (string) cluster distance measure: 'single', 'average',
wolffd@0: %                     'complete', 'neighf', 'ward', 'centroid', 'BMU'  
wolffd@0: %   [Ne]     (scalar) 0 = not constrined to neighborhood
wolffd@0: %                     1 = constrained   
wolffd@0: %            (matrix) size dlen x dlen, indicating possible connections
wolffd@0: %   [verbos] (scalar) 1 (default) = show status bar
wolffd@0: %                     0  = don't
wolffd@0: %
wolffd@0: % See also SOM_CLDIST. 
wolffd@0: 
wolffd@0: % Copyright (c) 2000 by Juha Vesanto
wolffd@0: % Contributed to SOM Toolbox on XXX by Juha Vesanto
wolffd@0: % http://www.cis.hut.fi/projects/somtoolbox/
wolffd@0:  
wolffd@0: % Version 2.0beta juuso 220800
wolffd@0: 
wolffd@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0: %% input arguments
wolffd@0: 
wolffd@0: q = 2; 
wolffd@0: 
wolffd@0: % map/data
wolffd@0: if isstruct(sM), 
wolffd@0:   switch sM.type, 
wolffd@0:    case 'som_map',  M = sM.codebook; mask = sM.mask; sT = sM.topol; 
wolffd@0:    case 'som_data', M = sM.data; mask = []; sT = []; 
wolffd@0:   end
wolffd@0: else M = sM; mask = []; sT = []; 
wolffd@0: end
wolffd@0: [dlen dim] = size(M); 
wolffd@0: if isempty(mask), mask = ones(dim,1); end
wolffd@0: 
wolffd@0: % simple option
wolffd@0: if any(strcmp(cldist,{'closest','BMU'})), 
wolffd@0:   i0 = find(base==0);
wolffd@0:   i1 = find(base>0);
wolffd@0:   bmus = som_bmus(M(i1,:),M(i0,:));
wolffd@0:   base(i0) = base(i1(bmus));
wolffd@0:   return; 
wolffd@0: end
wolffd@0: 
wolffd@0: % constrained clustering
wolffd@0: if nargin<4, Ne = []; end
wolffd@0: if prod(size(Ne))==1,  
wolffd@0:   if Ne & isempty(sT),
wolffd@0:     warning('Cannot use constrained clustering.'); Ne = 0; 
wolffd@0:   end
wolffd@0:   if Ne, Ne = som_unit_neighs(sT); else Ne = []; end
wolffd@0: end
wolffd@0: if ~isempty(Ne), 
wolffd@0:   Ne([0:dlen-1]*dlen+[1:dlen]) = 1; % set diagonal elements = 1
wolffd@0:   if all(Ne(:)>0), Ne = []; end
wolffd@0: end
wolffd@0: 
wolffd@0: if nargin<5, verbosity = 1; end
wolffd@0: 
wolffd@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0: %% initialize
wolffd@0: 
wolffd@0: if size(base,1)==1, base = base'; end
wolffd@0: 
wolffd@0: cid = unique(base(isfinite(base) & base~=0)); % cluster IDs
wolffd@0: nc = length(cid);    
wolffd@0: uind = find(base==0); % unclustered points
wolffd@0: nu = length(uind); 
wolffd@0: if nu==0, return; end
wolffd@0: 
wolffd@0: % initial clusters
wolffd@0: clinds = cell(nc,1); for i=1:nc, clinds{i} = find(base==i); end
wolffd@0: clinds2 = cell(nu,1); for i=1:nu, clinds2{i} = uind(i); end
wolffd@0: 
wolffd@0: % neighborhood function values
wolffd@0: if strcmp(cldist,'neighf')   
wolffd@0:   if isempty(sT), error('Cannot use neighf linkage.'); end
wolffd@0:   q = som_unit_dists(sT).^2; 
wolffd@0:   r = sM.trainhist(end).radius_fin^2; 
wolffd@0:   if isnan(r) | isempty(r), r = 1; end 
wolffd@0:   switch sM.neigh,
wolffd@0:    case 'bubble',   q = (q <= r);
wolffd@0:    case 'gaussian', q = exp(-q/(2*r));
wolffd@0:    case 'cutgauss', q = exp(-q/(2*r)) .* (q <= r);
wolffd@0:    case 'ep',       q = (1-q/r) .* (q <= r);
wolffd@0:   end
wolffd@0: end
wolffd@0: 
wolffd@0: % distance of each cluster to the unclustered points
wolffd@0: if any(strcmp(cldist,{'single','average','complete','neighf'})), 
wolffd@0:   M = som_mdist(M,2,mask,Ne); 
wolffd@0: end 
wolffd@0: Cd = som_cldist(M,clinds,clinds2,cldist,q,mask); 
wolffd@0: 			      
wolffd@0: % check out from Ne which of the clusters are not connected
wolffd@0: if ~isempty(Ne) & any(strcmp(cldist,{'centroid','ward'})),
wolffd@0:   Clconn = sparse(nc,nu);   
wolffd@0:   for i=1:nc, for j=1:nu, Clconn(i,j) = any(any(Ne(clinds{i},uind(j)))); end, end
wolffd@0:   Cd(Clconn==0) = Inf; 
wolffd@0: else
wolffd@0:   Clconn = []; 
wolffd@0: end
wolffd@0: 
wolffd@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0: %% action
wolffd@0: 
wolffd@0: if verbosity,
wolffd@0:   nu0 = nu; 
wolffd@0:   h = waitbar(1-nu/nu0,'Assigning unclustered points'); % tracking
wolffd@0: end
wolffd@0: 
wolffd@0: while 1, 
wolffd@0: 
wolffd@0:   % find closest unclustered point
wolffd@0:   [dk,k] = min(Cd,[],2);  % min distance from each unclustered point
wolffd@0:   [d,c]  = min(dk);       % cluster to which it is assigned  
wolffd@0:   k = k(c); 
wolffd@0: 
wolffd@0:   if ~isfinite(d), 
wolffd@0:     break; 
wolffd@0:   end
wolffd@0: 
wolffd@0:   % add k to cluster c
wolffd@0:   base(uind(k)) = cid(c);   
wolffd@0:   clinds{c} = [clinds{c}; uind(k)];
wolffd@0:   
wolffd@0:   % remove point k
wolffd@0:   notk = [1:k-1,k+1:nu]; 
wolffd@0:   nu = nu-1; if nu<=0, break; end  
wolffd@0:   Cd = Cd(:,notk); 
wolffd@0:   uind = uind(notk); 
wolffd@0:   clinds2 = clinds2(notk); 
wolffd@0:   if ~isempty(Clconn), Clconn = Clconn(:,notk); end
wolffd@0: 
wolffd@0:   % update cluster distances to c
wolffd@0:   Cd(c,:) = som_cldist(M,clinds(c),clinds2,cldist,q,mask); 
wolffd@0:   if ~isempty(Clconn), 
wolffd@0:     for j=1:nu, Clconn(c,j) = any(any(Ne(clinds{c},uind(j)))); end
wolffd@0:     Cd(c,find(Clconn(c,:)==0)) = Inf; 
wolffd@0:   end
wolffd@0:   
wolffd@0:   if verbosity, waitbar(1-nu/nu0,h); end % tracking
wolffd@0: 
wolffd@0: end
wolffd@0: if verbosity, close(h); end
wolffd@0: 
wolffd@0: return; 
wolffd@0: 
wolffd@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
wolffd@0: 
wolffd@0: