Code covered by the BSD License  

Highlights from
NLCSmoothReg

image thumbnail
from NLCSmoothReg by Michael Wendlandt
Computes a smooth regularized solution of an ill-posed discrete linear inverse problem.

[Lgfit,rho,eta]=LCurve(s,K,lambdavec,method,guess,tols,tolg,maxiter,plotLC,as,lb,ub,A,b,Aeq,beq,nonlcon);
function [Lgfit,rho,eta]=LCurve(s,K,lambdavec,method,guess,tols,tolg,maxiter,plotLC,as,lb,ub,A,b,Aeq,beq,nonlcon);

% determines a set of solutions $g$ according to a user-defined set of regularization parameters. 
% Two plots are generated for the detection of the optimal degree of
% regularization, the so-called L-curve and the discrepancy plot.

if nargin<10
    fprintf('\n\tToo less input arguments ! \n \n')  
else
    
    % get dimensions of g
    [gridpointsrow,gridpointscol]=size(guess);
    
    %start fitting
    time=0;
    fprintf('\n\tStart determination of L-curve !\n \n')
    for tel=1:length(lambdavec);
        tic;
        fprintf('\tlambda(%2d/%2d) =  %8.4e',tel,length(lambdavec),lambdavec(tel))
        if nargin==10
            [gfit,rho(tel),eta(tel)]=NLCSmoothReg(s,K,lambdavec(tel),method,guess,tols,tolg,maxiter,as);
            Lgfit(:,tel)=reshape(gfit,gridpointsrow*gridpointscol,1);
        elseif nargin==11
            [gfit,rho(tel),eta(tel)]=NLCSmoothReg(s,K,lambdavec(tel),method,guess,tols,tolg,maxiter,as,lb);
            Lgfit(:,tel)=reshape(gfit,gridpointsrow*gridpointscol,1);
        elseif nargin==12
            [gfit,rho(tel),eta(tel)]=NLCSmoothReg(s,K,lambdavec(tel),method,guess,tols,tolg,maxiter,as,lb,ub);
            Lgfit(:,tel)=reshape(gfit,gridpointsrow*gridpointscol,1);
        elseif nargin==14
            [gfit,rho(tel),eta(tel)]=NLCSmoothReg(s,K,lambdavec(tel),method,guess,tols,tolg,maxiter,as,lb,ub,A,b);
            Lgfit(:,tel)=reshape(gfit,gridpointsrow*gridpointscol,1);
        elseif nargin==16
            [gfit,rho(tel),eta(tel)]=NLCSmoothReg(s,K,lambdavec(tel),method,guess,tols,tolg,maxiter,as,lb,ub,A,b,Aeq,beq);
            Lgfit(:,tel)=reshape(gfit,gridpointsrow*gridpointscol,1);
        elseif nargin==17
            [gfit,rho(tel),eta(tel)]=NLCSmoothReg(s,K,lambdavec(tel),method,guess,tols,tolg,maxiter,as,lb,ub,A,b,Aeq,beq,nonlcon);
            Lgfit(:,tel)=reshape(gfit,gridpointsrow*gridpointscol,1);  
        end
        
        %determine elapsed and left total computing time
        t=toc;
        time=time + t;
        tottime=round(time * length(lambdavec) / tel);
        percentdone=round((time/tottime) * 100);
        
        timeelapsed=round(time);
        edays=round(timeelapsed/86400);
        ehours=round(mod(timeelapsed,86400)/3600);
        eminute=round(mod(mod(timeelapsed,86400),3600)/60);
        esec=mod(mod(mod(timeelapsed,86400),3600),60);
        
        timeleft=round(tottime-time);
        ldays=round(timeleft/86400);
        lhours=round(mod(timeleft,86400)/3600);
        lminute=round(mod(mod(timeleft,86400),3600)/60);
        lsec=mod(mod(mod(timeleft,86400),3600),60);
        
        fprintf('\n\ttime elapsed  = %2dd %2dh %2dm %2ds\n\ttime left     = %2dd %2dh %2dm %2ds, %3d%% done\n\n',edays,ehours,eminute,esec,ldays,lhours,lminute,lsec,percentdone)
    end
    
    rho=rho';
    eta=eta';
    
    if (length(lambdavec)>1) & (strcmp(plotLC,'on'))
        Lplot(lambdavec,rho,eta);
    end
    
end

Contact us at files@mathworks.com