inv_solve_TV_pdipm

PURPOSE ^

INV_SOLVE_TV_PDIPM inverse solver for Andrea Borsic's

SYNOPSIS ^

function img= inv_solve_TV_pdipm( inv_model, data1, data2)

DESCRIPTION ^

 INV_SOLVE_TV_PDIPM inverse solver for Andrea Borsic's
   Total Variation solver for use with difference EIT
 img= inv_solve_TV_pdipm( inv_model, data1, data2)
 img        => output image
 inv_model  => inverse model struct
 data1      => differential data at earlier time
 data2      => differential data at later time
 Parameters
   alpha1
   alpha2
   want_dual_variable  (set to 1 if you want access to dual)
 Termination parameters
  max_iters =  inv_model.parameters.max_iteration (default 10)
      Max number of iterations before stopping
  min change = inv_model.parameters.min_change   (default 0)
      Min Change in objective fcn (norm(y-Jx)^2 + hp*TV(x)) before stopping

CROSS-REFERENCE INFORMATION ^

This function calls: This function is called by:

SUBFUNCTIONS ^

SOURCE CODE ^

0001 function img= inv_solve_TV_pdipm( inv_model, data1, data2)
0002 % INV_SOLVE_TV_PDIPM inverse solver for Andrea Borsic's
0003 %   Total Variation solver for use with difference EIT
0004 % img= inv_solve_TV_pdipm( inv_model, data1, data2)
0005 % img        => output image
0006 % inv_model  => inverse model struct
0007 % data1      => differential data at earlier time
0008 % data2      => differential data at later time
0009 % Parameters
0010 %   alpha1
0011 %   alpha2
0012 %   want_dual_variable  (set to 1 if you want access to dual)
0013 % Termination parameters
0014 %  max_iters =  inv_model.parameters.max_iteration (default 10)
0015 %      Max number of iterations before stopping
0016 %  min change = inv_model.parameters.min_change   (default 0)
0017 %      Min Change in objective fcn (norm(y-Jx)^2 + hp*TV(x)) before stopping
0018 
0019 % (C) 2002-2009 Andrea Borsic and Andy Adler. License: GPL version 2 or version 3
0020 % $Id: inv_solve_TV_pdipm.m 5289 2016-11-02 08:43:12Z fab-b $
0021 
0022 
0023 p= get_params(inv_model);
0024 
0025 dva = calc_difference_data( data1, data2, inv_model.fwd_model);
0026 % TEST CODE -> Put elsewhere
0027 back_val = get_good_background(inv_model, data1);
0028 inv_model.jacobian_bkgnd.value= back_val;
0029 
0030 sol= [];
0031 for i=1:size(dva,2)
0032    [soln,dual_x]=primaldual_tvrecon_lsearch(inv_model, dva(:,i), ...
0033         p.maxiter,p.alpha1,p.alpha2, p.tol, p.beta, p.min_change);
0034 
0035    if ~p.keepiters
0036       soln=soln(:,end);
0037    end
0038 
0039    sol=[sol, soln];
0040 end
0041 
0042 img.name= 'solved by inv_solve_TV_pdipm';
0043 img.elem_data = sol;
0044 img.fwd_model= inv_model.fwd_model;
0045 try if inv_model.inv_solve_TV_pdipm.want_dual_variable
0046    img.dual_data = dual_x;
0047 end; end
0048 
0049 function p = get_params(inv_model);
0050    
0051    try   
0052        p.alpha1= inv_model.inv_solve_TV_pdipm.alpha1;
0053    catch
0054        p.alpha1= 1e-2;
0055    end
0056 
0057    try   
0058        p.beta= inv_model.inv_solve_TV_pdipm.beta;
0059    catch
0060        p.beta= 1e-4;
0061    end
0062 
0063    p.alpha2= calc_hyperparameter( inv_model);
0064 
0065    try   
0066        p.min_change = inv_model.parameters.min_change;
0067    catch
0068        p.min_change = 0;
0069    end
0070 
0071    try   
0072        p.maxiter = inv_model.parameters.max_iterations;
0073    catch
0074        p.maxiter= 10;
0075    end
0076 
0077    try   
0078        p.keepiters = inv_model.parameters.keep_iterations;
0079    catch
0080        p.keepiters= 0;
0081    end
0082 
0083    p.tol = 0; % TODO
0084 
0085 function back_val = get_good_background(inv_mdl, data1);
0086 
0087    % Create homogeneous model
0088    IM= eidors_obj('image','');
0089    IM.fwd_model= inv_mdl.fwd_model;
0090    s= ones(size(IM.fwd_model.elems,1),1);
0091    IM.elem_data= s;
0092 
0093    vsim= fwd_solve( IM);
0094    back_val=abs( data1\vsim.meas );
0095    back_val=1;
0096

Generated on Wed 21-Jun-2017 09:29:07 by m2html © 2005