AA_INV_SOLVE inverse solver using approach of Adler&Guardo 1996 img= aa_inv_solve( inv_model, data1, data2) img => output image (or vector of images) inv_model => inverse model struct data1 => differential data at earlier time data2 => differential data at later time both data1 and data2 may be matrices (MxT) each of M measurements at T times if either data1 or data2 is a vector, then it is expanded to be the same size matrix
0001 function img= aa_inv_solve( inv_model, data1, data2) 0002 % AA_INV_SOLVE inverse solver using approach of Adler&Guardo 1996 0003 % img= aa_inv_solve( inv_model, data1, data2) 0004 % img => output image (or vector of images) 0005 % inv_model => inverse model struct 0006 % data1 => differential data at earlier time 0007 % data2 => differential data at later time 0008 % 0009 % both data1 and data2 may be matrices (MxT) each of 0010 % M measurements at T times 0011 % if either data1 or data2 is a vector, then it is expanded 0012 % to be the same size matrix 0013 0014 % (C) 2005 Andy Adler. License: GPL version 2 or version 3 0015 % $Id: aa_inv_solve.html 2819 2011-09-07 16:43:11Z aadler $ 0016 0017 dv = calc_difference_data( data1, data2, inv_model.fwd_model); 0018 sol = get_RM( inv_model ) * dv; 0019 0020 img.name= 'solved by aa_inv_solve'; 0021 img.elem_data = sol; 0022 img.fwd_model= inv_model.fwd_model; 0023 0024 function RM = get_RM( inv_model ); 0025 % The one_step reconstruction matrix is cached 0026 RM = eidors_obj('get-cache', inv_model, 'aa_inv_solve'); 0027 if ~isempty(RM) 0028 eidors_msg('aa_inv_solve: using cached value', 3); 0029 return; 0030 end 0031 0032 img_bkgnd= calc_jacobian_bkgnd( inv_model ); 0033 J = calc_jacobian( img_bkgnd); 0034 0035 RtR = calc_RtR_prior( inv_model ); 0036 W = calc_meas_icov( inv_model ); 0037 hp = calc_hyperparameter( inv_model ); 0038 0039 RM= (J'*W*J + hp^2*RtR)\J'*W; 0040 0041 eidors_obj('set-cache', inv_model, 'aa_inv_solve', RM); 0042 eidors_msg('aa_inv_solve: setting cached value', 3);