idamnjanovic@40: function Dictionary = SMALL_rlsdla(X, params)
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: CODE_SPARSITY = 1;
idamnjanovic@40: CODE_ERROR = 2;
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: % Determine which method will be used for sparse representation step -
idamnjanovic@40: % Sparsity or Error mode
idamnjanovic@40: 
idamnjanovic@40: if (isfield(params,'codemode'))
idamnjanovic@40:     switch lower(params.codemode)
idamnjanovic@40:         case 'sparsity'
idamnjanovic@40:             codemode = CODE_SPARSITY;
idamnjanovic@40:             thresh = params.Tdata;
idamnjanovic@40:         case 'error'
idamnjanovic@40:             codemode = CODE_ERROR;
idamnjanovic@40:             thresh = params.Edata;
idamnjanovic@40:             
idamnjanovic@40:         otherwise
idamnjanovic@40:             error('Invalid coding mode specified');
idamnjanovic@40:     end
idamnjanovic@40: elseif (isfield(params,'Tdata'))
idamnjanovic@40:     codemode = CODE_SPARSITY;
idamnjanovic@40:     thresh = params.Tdata;
idamnjanovic@40: elseif (isfield(params,'Edata'))
idamnjanovic@40:     codemode = CODE_ERROR;
idamnjanovic@40:     thresh = params.Edata;
idamnjanovic@40:     
idamnjanovic@40: else
idamnjanovic@40:     error('Data sparse-coding target not specified');
idamnjanovic@40: end
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: % max number of atoms %
idamnjanovic@40: 
idamnjanovic@40: if (codemode==CODE_ERROR && isfield(params,'maxatoms'))
idamnjanovic@40:   maxatoms = params.maxatoms;
idamnjanovic@40: else
idamnjanovic@40:   maxatoms = -1;
idamnjanovic@40: end
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: % Forgetting factor
idamnjanovic@40: 
idamnjanovic@40: if (isfield(params,'forgettingMode'))
idamnjanovic@40:      switch lower(params.forgettingMode)
idamnjanovic@40:        case 'fix'
idamnjanovic@40:                 if (isfield(params,'forgettingFactor'))
idamnjanovic@40:                     lambda=params.forgettingFactor;
idamnjanovic@40:                 else
idamnjanovic@40:                     lambda=1;
idamnjanovic@40:                 end
idamnjanovic@40:        otherwise
idamnjanovic@40:               error('This mode is still not implemented');
idamnjanovic@40:      end
idamnjanovic@40: elseif (isfield(params,'forgettingFactor'))
idamnjanovic@40:      lambda=params.forgettingFactor;
idamnjanovic@40: else
idamnjanovic@40:      lambda=1;
idamnjanovic@40: end
idamnjanovic@40: 
idamnjanovic@40: % determine dictionary size %
idamnjanovic@40: 
idamnjanovic@40: if (isfield(params,'initdict'))
idamnjanovic@40:   if (any(size(params.initdict)==1) && all(iswhole(params.initdict(:))))
idamnjanovic@40:     dictsize = length(params.initdict);
idamnjanovic@40:   else
idamnjanovic@40:     dictsize = size(params.initdict,2);
idamnjanovic@40:   end
idamnjanovic@40: end
idamnjanovic@40: if (isfield(params,'dictsize'))    % this superceedes the size determined by initdict
idamnjanovic@40:   dictsize = params.dictsize;
idamnjanovic@40: end
idamnjanovic@40: 
idamnjanovic@40: if (size(X,2) < dictsize)
idamnjanovic@40:   error('Number of training signals is smaller than number of atoms to train');
idamnjanovic@40: end
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: % initialize the dictionary %
idamnjanovic@40: 
idamnjanovic@40: if (isfield(params,'initdict'))
idamnjanovic@40:   if (any(size(params.initdict)==1) && all(iswhole(params.initdict(:))))
idamnjanovic@40:     D = X(:,params.initdict(1:dictsize));
idamnjanovic@40:   else
idamnjanovic@40:     if (size(params.initdict,1)~=size(X,1) || size(params.initdict,2)<dictsize)
idamnjanovic@40:       error('Invalid initial dictionary');
idamnjanovic@40:     end
idamnjanovic@40:     D = params.initdict(:,1:dictsize);
idamnjanovic@40:   end
idamnjanovic@40: else
idamnjanovic@40:   data_ids = find(colnorms_squared(X) > 1e-6);   % ensure no zero data elements are chosen
idamnjanovic@40:   perm = randperm(length(data_ids));
idamnjanovic@40:   D = X(:,data_ids(perm(1:dictsize)));
idamnjanovic@40: end
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@40: % normalize the dictionary %
idamnjanovic@40: 
idamnjanovic@40: D = normcols(D);
idamnjanovic@40: 
idamnjanovic@40: % Training data
idamnjanovic@40: 
idamnjanovic@40: data=X;
idamnjanovic@65: cnt=size(data,2);
idamnjanovic@40: %
idamnjanovic@40: 
idamnjanovic@40: C=(100000*thresh)*eye(dictsize);
idamnjanovic@40: w=zeros(dictsize,1);
idamnjanovic@40: u=zeros(dictsize,1);
idamnjanovic@40: 
idamnjanovic@40: 
idamnjanovic@65: for i = 1:cnt
idamnjanovic@40: 
idamnjanovic@40:    if (codemode == CODE_SPARSITY)
idamnjanovic@40:         w = ompmex(D,data(:,i),[],[],thresh,1,-1,0);
idamnjanovic@40:    else
idamnjanovic@40:         w = omp2mex(D,data(:,i),[],[],[],thresh,0,-1,maxatoms,0);
idamnjanovic@40:    end
idamnjanovic@40:    
idamnjanovic@40:    spind=find(w);
idamnjanovic@40:    
idamnjanovic@40:    residual = data(:,i) - D * w;
idamnjanovic@40:    
idamnjanovic@40:    if (lambda~=1) 
idamnjanovic@40:        C = C *(1/ lambda);
idamnjanovic@40:    end
idamnjanovic@40:    
idamnjanovic@40:    u = C(:,spind) * w(spind);
idamnjanovic@40:   
idamnjanovic@40:     
idamnjanovic@40:    alfa = 1/(1 + w' * u);
idamnjanovic@40:     
idamnjanovic@40:    D = D + (alfa * residual) * u';
idamnjanovic@40:     
idamnjanovic@40:    
idamnjanovic@40:    C = C - (alfa * u)* u';
idamnjanovic@40:    
idamnjanovic@40: end
idamnjanovic@40: 
idamnjanovic@40: Dictionary = D;
idamnjanovic@40: 
idamnjanovic@40: end