Daniel@0: function Md = som_mdist(D,q,mask,Ne)
Daniel@0: 
Daniel@0: % SOM_MDIST Mutual (or pairwise) distance matrix for the given data.
Daniel@0: % 
Daniel@0: %   Md = som_mdist(D,[q],[mask],[Ne])
Daniel@0: %
Daniel@0: %    Md = som_mdist(D); 
Daniel@0: %    Md = som_mdist(D,Inf); 
Daniel@0: %    Md = som_mdist(D,2,Ne); 
Daniel@0: %
Daniel@0: %  Input and output arguments ([]'s are optional):
Daniel@0: %   D        (matrix) size dlen x dim, the data set
Daniel@0: %            (struct) map or data struct
Daniel@0: %   [q]      (scalar) distance norm, default = 2
Daniel@0: %   [mask]   (vector) size dim x 1, the weighting mask
Daniel@0: %   [Ne]     (matrix) size dlen x dlen, sparse matrix 
Daniel@0: %                     indicating which distances should be 
Daniel@0: %                     calculated (ie. less than Infinite) 
Daniel@0: %
Daniel@0: % See also PDIST. 
Daniel@0: 
Daniel@0: % Copyright (c) 2000 by Juha Vesanto
Daniel@0: % Contributed to SOM Toolbox on XXX by Juha Vesanto
Daniel@0: % http://www.cis.hut.fi/projects/somtoolbox/
Daniel@0:  
Daniel@0: % Version 2.0beta juuso 220800
Daniel@0: 
Daniel@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Daniel@0: 
Daniel@0: % mask
Daniel@0: if nargin<3, mask = []; end
Daniel@0: 
Daniel@0: % the data
Daniel@0: if isstruct(D), 
Daniel@0:   switch D.type, 
Daniel@0:    case 'som_map', if isempty(mask), mask = D.mask; end, D = D.codebook; 
Daniel@0:    case 'som_data', D = D.data; 
Daniel@0:    otherwise, error('Bad first argument');
Daniel@0:   end
Daniel@0: end
Daniel@0: nans = sum(isnan(D),2);
Daniel@0: if any(nans>0), 
Daniel@0:   D(find(nans>0),:) = 0; 
Daniel@0:   warning('Distances of vectors with NaNs are not calculated.'); 
Daniel@0: end
Daniel@0: [dlen dim] = size(D);
Daniel@0: 
Daniel@0: % distance norm
Daniel@0: if nargin<2 | isempty(q) | isnan(q), q = 2; end
Daniel@0: 
Daniel@0: % mask
Daniel@0: if isempty(mask), mask = ones(dim,1); end
Daniel@0: 
Daniel@0: % connections 
Daniel@0: if nargin<4, Ne = []; end
Daniel@0: if ~isempty(Ne), 
Daniel@0:   l = size(Ne,1); Ne([0:l-1]*l+[1:l]) = 1; % set diagonal elements = 1
Daniel@0: end
Daniel@0: 
Daniel@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Daniel@0: 
Daniel@0: m = mask; 
Daniel@0: o = ones(dlen,1); 
Daniel@0: l = dlen; 
Daniel@0: Md = zeros(dlen);
Daniel@0: calculate_all = isempty(Ne); 
Daniel@0: 
Daniel@0: if ~calculate_all, Md(Ne==0) = Inf; end
Daniel@0: 
Daniel@0: for i=1:l-1,
Daniel@0:   j=(i+1):l; 
Daniel@0:   if ~calculate_all, j=find(Ne(i,j))+i; end
Daniel@0:   C=D(j,:)-D(i*o(1:length(j)),:);
Daniel@0:   switch q, 
Daniel@0:    case 1,    Md(j,i)=abs(C)*m;
Daniel@0:    case 2,    Md(j,i)=sqrt((C.^2)*m);  
Daniel@0:    case Inf,  Md(j,i)=max(diag(m)*abs(C),[],2);
Daniel@0:    otherwise, Md(j,i)=((abs(C).^q)*m).^(1/q);
Daniel@0:   end   
Daniel@0:   Md(i,j) = Md(j,i)';
Daniel@0: end
Daniel@0: 
Daniel@0: Md(find(nans>0),:) = NaN; 
Daniel@0: Md(:,find(nans>0)) = NaN; 
Daniel@0: 
Daniel@0: return;
Daniel@0: 
Daniel@0: %%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
Daniel@0: