choose_noise_figure

 CHOOSE_NOISE_FIGURE: choose hyperparameter based on NF calculation
 HP= choose_noise_figure( inv_model );
 inv_model  => inverse model struct

 In order to use this function, it is necessary to specify
 inv_model.hyperparameter. has the following fields
 hpara.func         = @choose_noise_figure;
 hpara.noise_figure = NF Value requested
 hpara.tgt_elems    = vector of element numbers of contrast in centre

 The NF parameter is from the definition in Adler & Guardo (1996).

 SNR_z = sumsq(z0) / var(z) = sum(z0.^2) / trace(Rn)
 SNR_x = sumsq(A*x0) / var(A*x) = sum((A*x0).^2) / trace(ABRnB'A')
   where Rn = noise covariance and A_ii = area of element i
 NF = SNR_z / SNR_x

0001 function HP= choose_noise_figure( inv_model );
0002 % CHOOSE_NOISE_FIGURE: choose hyperparameter based on NF calculation
0003 % HP= choose_noise_figure( inv_model );
0004 % inv_model  => inverse model struct
0005 %
0006 % In order to use this function, it is necessary to specify
0007 % inv_model.hyperparameter. has the following fields
0008 % hpara.func         = @choose_noise_figure;
0009 % hpara.noise_figure = NF Value requested
0010 % hpara.tgt_elems    = vector of element numbers of contrast in centre
0011 %
0012 % The NF parameter is from the definition in Adler & Guardo (1996).
0013 %
0014 % SNR_z = sumsq(z0) / var(z) = sum(z0.^2) / trace(Rn)
0015 % SNR_x = sumsq(A*x0) / var(A*x) = sum((A*x0).^2) / trace(ABRnB'A')
0016 %   where Rn = noise covariance and A_ii = area of element i
0017 % NF = SNR_z / SNR_x
0018 
0019 % (C) 2006 Andy Adler. License: GPL version 2 or version 3
0020 % $Id: choose_noise_figure.html 2819 2011-09-07 16:43:11Z aadler $
0021 
0022 reqNF= inv_model.hyperparameter.noise_figure;
0023 try % remove the value field
0024    inv_model.hyperparameter = rmfield(inv_model.hyperparameter,'value');
0025 end
0026 
0027 HP = eidors_obj('get-cache', inv_model, 'HP_for_NF');
0028 if ~isempty(HP)
0029    eidors_msg('choose_noise_figure: using cached value', 2);
0030    return
0031 end
0032 
0033 HP = HP_for_NF_search(reqNF,inv_model);
0034 
0035    
0036 eidors_cache('boost_priority',3);
0037 eidors_obj('set-cache', inv_model, 'HP_for_NF', HP);
0038 eidors_msg('choose_noise_figure: setting cached value', 2);
0039 eidors_cache('boost_priority',-3);
0040 
0041 
0042 function [HP,NF,SE] = HP_for_NF_search(dNF,imdl);
0043    hp= search1(dNF, imdl, 1);
0044 
0045    dx= hp-linspace(-0.7,0.7,5);
0046    hp= search2(dNF, imdl, dx);
0047 
0048    dx= hp-linspace(-0.2,0.2,5);
0049    hp= search2(dNF,imdl, dx);
0050 
0051    dx= hp-linspace(-0.1,0.1,5); 
0052    hp= search2(dNF,imdl, dx);
0053 
0054    dx= hp-linspace(-0.05,0.05,21);
0055    hp= search2(dNF,imdl, dx);
0056 
0057    HP= 10^-hp;
0058   
0059 function hp= search1(dNF, imdl, hp)
0060   [NF,SE]=calc_noise_figure( imdl, 10^(-hp));
0061    if     NF+3*SE < dNF; dir = 1;
0062    elseif NF-3*SE > dNF; dir = -1;
0063    else   dir = 0; end
0064    while  dir*NF+3*SE < dir*dNF %>
0065      hp= hp+0.5*dir;
0066      [NF,SE]=calc_noise_figure( imdl, 10^(-hp));
0067    end
0068 
0069 function hp=search2(dNF, imdl, dx)
0070    hp = [];
0071    it = 0;
0072    nf = zeros(1,length(dx));
0073    while isempty(hp) && it < 2
0074        it = it+1;
0075        %if it > 1 keyboard, end
0076        for k=1:length(dx)
0077            nf(k)=nf(k)+calc_noise_figure( imdl, 10^-dx(k), 10 );
0078        end
0079        log_nf = log10(nf/it);
0080        p= polyfit( dx, log_nf-log10(dNF), 1);
0081        hp = roots(p);
0082        hp = hp( hp<max(dx) & hp>min(dx) );  %USE if poly>1
0083    end
0084    if isempty(hp)
0085        %fallback
0086        [jnk,idx] = min(abs(log_nf-log10(dNF)));
0087        hp = dx(idx);
0088    end