www.pudn.com > image_mva_0.rar > mcr.m

function [c,s] = mcr(data,c0,ccon,scon,ittol); 
%MCR Multivariate curve resolution with constraints 
%  Input: 
%  data (m by n)- is the matrix to be analyzed; 
%  c0 (m by k) is the initial guess for c 
%  k - number of components to resolve 
% 
% Output: 
%  c (m by k)- concentration profiles 
%  s  (k by n)- pure spectra  
% 
%  (ccon) is an option describing constraints on (c):  
%    ccon = 0, no contraints 
%    ccon = 2, non-negative 
% 
%  (scon) is an option describing constraints on (s):  
%    scon = 0, no contraints  
%    scon = 2, non-negative  
% 
%  (ittol) is an optional convergance criteria, 
%    ittol < 1, convergance tolerance (default 1e-4), and 
%    ittol >=1, max number of iterations (default 100) 
% 
%   
% Example: [c,s] = mcr(data,c0,0,0); 
% 
 
 
if nargin<3 
  ccon  = 0; 
elseif isempty(ccon) 
  ccon  = 0; 
end 
 
if nargin<4 
  scon  = 0; 
elseif isempty(scon) 
  scon  = 0; 
end 
 
if nargin<5 
  itmax = 100; itmin = 1e-8; ittol = itmax; 
elseif isempty(ittol) 
  itmax = 100; itmin = 1e-8; 
elseif ittol<1 
  itmax = 1e6; itmin = ittol; 
  if ittol < 1e-8 
    itmin = 1e-8; 
  end 
else 
  itmax = ittol; itmin = 1e-8; 
  if itmax > 1e6 
    itmax = 1e6; 
  end 
end 
 
ka      = size(c0,2); 
kac     = ones(1,ka);  
s0      = zeros(ka,size(data,2)); 
it      = 0; 
 
if nargin<6 
  cc    = 0; cc1 = 0; 
elseif isempty(cc) 
  cc    = 0; cc1 = 0; 
elseif sum(sum(isfinite(cc))') 
  if size(cc,1)~=size(c0,1)|size(cc,2)~=size(c0,2) 
    error('cc must be the same size as c0') 
  end 
  cc1   = 1;  
  ii    = find(isfinite(cc)); 
  c0(ii) = cc(ii); 
  i1    = find(sum(isfinite(cc))); 
  kac(1,i1) = zeros(1,length(i1)); 
else 
  cc    = 0; cc1 = 0; 
end 
c       = c0; 
 
if nargin<8 
  sclc  = [1:size(data,1)]; 
elseif isempty(sclc)|(length(sclc)~=size(data,1)) 
  sclc  = [1:size(data,1)]; 
end 
 
if nargin<9 
  scls  = [1:size(data,2)]; 
elseif isempty(scls)|(length(scls)~=size(data,2)) 
  scls  = [1:size(data,2)]; 
end 
 
if nargin<10 
  nnlstol = max(size(data))*norm(data,1)*eps; 
elseif isempty(nnlstol) 
  nnlstol = max(size(data))*norm(data,1)*eps; 
end 
 
if nargin<7 
  sc    = 0; sc1 = 0; 
elseif isempty(sc) 
  sc    = 0; sc1 = 0; 
elseif sum(sum(isfinite(sc))') 
  if size(sc,1)~=size(s0,1)|size(sc,2)~=size(s0,2) 
    error('size of sc must be size(c0,2) by size(data,2)') 
  end 
  sc1   = 1;  
  ii    = find(isfinite(sc)); 
  s0(ii) = sc(ii); 
  i1    = find(sum(isfinite(sc'))); 
  kac(1,i1) = zeros(1,length(i1)); 
else 
  sc    = 0; sc1 = 0; 
end 
  kac   = find(kac); 
 
switch scon  
case 0 
  s0    = c\data; 
case 1 
  for ii=1:length(scls) 
    s0(:,ii) = nnls(c,data(:,ii),nnlstol); 
  end 
case 2 
  for ii=1:length(scls) 
    s0(:,ii) = fastnnls(c,data(:,ii),nnlstol); 
  end 
end 
if sc1 
  i2  = find(isfinite(sc));     
  s0(i2) = sc(i2); 
end 
s       = s0; 
if ~isempty(kac) 
  for i1=1:length(kac) 
    s(kac(i1),:) = s(kac(i1),:)/norm(s(kac(i1),:)'); 
  end 
end 
 
while it