nt_sns_cluster

 y=nt_sns_cluster(x,nneigbors,cluster_size) - sensor noise suppression within clusters

   y: denoised matrix

   x: matrix  to denoise (time * channels * repeats * ...)
   nneighbors: number of channels to use in projection
   cluster_size: target cluster size [default: size(x,1)/2]

 NoiseTools

0001  function x=nt_sns_cluster(x,nneighbors,cluster_size)
0002 % y=nt_sns_cluster(x,nneigbors,cluster_size) - sensor noise suppression within clusters
0003 %
0004 %   y: denoised matrix
0005 %
0006 %   x: matrix  to denoise (time * channels * repeats * ...)
0007 %   nneighbors: number of channels to use in projection
0008 %   cluster_size: target cluster size [default: size(x,1)/2]
0009 %
0010 % NoiseTools
0011 nt_greetings;
0012 
0013 
0014 if nargin<3; cluster_size=[]; end
0015 if nargin<2; error('!'); end
0016 
0017 sz=size(x);
0018 x=x(:,:,:); % merge higher dimensions into repeats
0019 x=nt_unfold(x);
0020 x0=x;
0021 
0022 if isempty(cluster_size); cluster_size=round(size(x,1)/2); end
0023 
0024 if size(x,2)<=cluster_size
0025     x=nt_sns(x,nneighbors);
0026 else
0027     NCLUSTERS=2;
0028     [C,A]=vl_kmeans(x,NCLUSTERS,'algorithm', 'elkan','initialization','plusplus','numrepetitions', 100);
0029     if numel(find(A==1)) && numel(find(A==2))
0030         xA=nt_sns_cluster(x(:,find(A==1)),nneighbors,cluster_size);
0031         xB=nt_sns_cluster(x(:,find(A==2)),nneighbors,cluster_size);
0032         x(:,find(A==1))=xA;
0033         x(:,find(A==2))=xB;
0034     end % else no split, return
0035 end
0036 x1=x;
0037 
0038 x=nt_fold(x,sz(1));
0039 x=reshape(x,sz);
0040 
0041 verbose=0;
0042 if nargout==0 || verbose;
0043     disp(['cluster size: ', num2str(sz(2)), ',  power ratio: ',num2str(nt_wpwr(x1)/nt_wpwr(x0))]);
0044 end
0045 
0046