compare_3d_algs

 Compare different 3D reconstructions
 
 algno=1     np_inv_solve            np_calc_image_prior
 algno=2     prior_laplace     np_inv_solve
 algno=3     prior_TV        np_inv_solve
 algno=4     prior_TV        inv_solve_TV_pdipm

0001 function imgr= compare_3d_algs( algno )
0002 % Compare different 3D reconstructions
0003 %
0004 % algno=1     np_inv_solve            np_calc_image_prior
0005 % algno=2     prior_laplace     np_inv_solve
0006 % algno=3     prior_TV        np_inv_solve
0007 % algno=4     prior_TV        inv_solve_TV_pdipm
0008 
0009 % (C) 2005 Andy Adler. License: GPL version 2 or version 3
0010 % $Id: compare_3d_algs.m 3428 2012-07-02 20:56:41Z bgrychtol $
0011 
0012 calc_colours('ref_level','auto');
0013 
0014 imb=  mk_common_model('n3r2',16);
0015 e= size(imb.fwd_model.elems,1);
0016 sigma= ones(e,1);
0017 img= eidors_obj('image','');
0018 img.elem_data= sigma;
0019 img.fwd_model= imb.fwd_model;
0020 vh= fwd_solve( img );
0021 
0022 load datacom A B;
0023 sigma(A)= 1.2; sigma(B)=0.8;
0024 clear A B;
0025 img.elem_data= sigma;
0026 vi= fwd_solve( img );
0027 
0028 sig= sqrt(norm(vi.meas - vh.meas));
0029 m= size(vi.meas,1);
0030 vi.meas = vi.meas + .001*sig*randn(m,1);
0031 
0032 %show_slices(img);
0033 inv3d= eidors_obj('inv_model', 'EIT inverse');
0034 inv3d.reconst_type= 'difference';
0035 inv3d.jacobian_bkgnd.value = 1;
0036 inv3d.fwd_model= imb.fwd_model;
0037 inv3d.fwd_model.np_fwd_solve.perm_sym= '{y}';
0038 
0039      iidx=1;
0040 switch algno
0041    case 1,
0042      inv3d.hyperparameter.value = 1e-4;
0043      inv3d.solve=       'np_inv_solve';
0044      inv3d.R_prior=     'np_calc_image_prior';
0045      inv3d.np_calc_image_prior.parameters= [3 1]; %  deg=1, w=1
0046 
0047    case 2,
0048      inv3d.hyperparameter.value = 1e-3;
0049      inv3d.RtR_prior=    'prior_laplace';
0050      inv3d.solve=        'np_inv_solve';
0051 
0052    case 3,
0053      inv3d.hyperparameter.value = 1e-2;
0054      inv3d.R_prior=      'prior_TV';
0055      inv3d.solve=        'np_inv_solve';
0056 
0057    case 4,
0058      inv3d.hyperparameter.value = 1e-4;
0059      inv3d.parameters.max_iterations= 5;
0060      inv3d.parameters.term_tolerance= 1e-3;
0061      inv3d.parameters.keep_iterations= 1;
0062      inv3d.R_prior=      'prior_TV';
0063      inv3d.solve=        'inv_solve_TV_pdipm';
0064      iidx=[1,2,5];
0065 
0066    case 5,
0067      inv3d.hyperparameter.value = 1e-3;
0068      inv3d.R_prior=      'prior_TV';
0069      inv3d.solve=        'aa_inv_total_var';
0070      inv3d.parameters.max_iterations= 5;
0071 
0072    otherwise,
0073      error('action unknown');
0074 end
0075 
0076 %
0077 % Step 3: Reconst and show image
0078 %
0079 imgr= inv_solve( inv3d, vh, vi);
0080 imgr.elem_data= imgr.elem_data(:,iidx);
0081 show_slices(imgr, [.5:1:2.5]'*[Inf,Inf,1]);
0082