fwd_solve_apparent_resistivity

 fwd_solve_apparent_resistivity: fwd_solve output as apparent resistivity
  This function is a wrapper to the fwd_solve; however, the output
  is converted into apparent resistivity units, rather than in
  voltage units

0001 function data =fwd_solve_apparent_resistivity(fwd_model, img)
0002 % fwd_solve_apparent_resistivity: fwd_solve output as apparent resistivity
0003 %  This function is a wrapper to the fwd_solve; however, the output
0004 %  is converted into apparent resistivity units, rather than in
0005 %  voltage units
0006 
0007 % $Id: fwd_solve_apparent_resistivity.m 5877 2018-12-21 20:41:28Z aadler $
0008 
0009 % correct input paralemeters if function was called with only img
0010 
0011 if ischar(fwd_model) && strcmp(fwd_model,'UNIT_TEST'); do_unit_test; return; end
0012 
0013 if nargin == 1
0014    img= fwd_model;
0015 elseif  strcmp(getfield(warning('query','EIDORS:DeprecatedInterface'),'state'),'on')
0016    warning('EIDORS:DeprecatedInterface', ...
0017       ['Calling FWD_SOLVE_APPARENT_RESISTIVITY with two arguments is deprecated and will cause' ...
0018        ' an error in a future version. First argument ignored.']);
0019 end
0020 fwd_model= img.fwd_model;
0021 
0022 solver = @eidors_default;
0023 try % replace solver if this one exists
0024    solver = fwd_model.fwd_solve_apparent_resistivity.solve;
0025 end
0026 
0027 img.fwd_model.solve = solver;
0028 
0029 data = fwd_solve(img);
0030 
0031 fctr = get_factor(img);
0032 
0033 data.meas = fctr * data.meas;
0034 data.name = ['apparent resistivity ' data.name];
0035 data.quantity = 'apparent resistivity';
0036 data.apparent_resistivity_factor= fctr;
0037 
0038 function fctr = get_factor(img)
0039 fctr = NaN;
0040 try
0041    fctr = img.fwd_model.apparent_resistivity_factor;
0042 end
0043 
0044 if ischar(fctr), fctr = str2func(fcstr); end;
0045 
0046 if isa(fctr, 'function_handle')
0047    fctr = feval(fctr, img);
0048 end
0049 
0050 if isnan(fctr)
0051    
0052    vh = fwd_solve(mk_image(img.fwd_model,1));
0053    n = length(vh.meas);
0054    fctr = spdiags(1./vh.meas,0,n,n);
0055 end
0056 
0057 
0058 function do_unit_test
0059    img = mk_image(mk_common_model('a2c2',8)); 
0060    img.elem_data(:) = 1.0; vrh = fwd_solve( img );
0061    img.elem_data(1) = 1.1; vri = fwd_solve( img );
0062 
0063    img.fwd_model.solve = @fwd_solve_apparent_resistivity;
0064    img.elem_data(:) = 1.0; vah = fwd_solve( img );
0065    img.elem_data(1) = 1.1; vai = fwd_solve( img );
0066 
0067    unit_test_cmp('homog is ones', vah.meas, ones(size(vah.meas)), 1e-10);
0068    unit_test_cmp('ratio is same', vai.meas./vah.meas, vri.meas./vrh.meas, 1e-10);