aa_e_move_image_prior

 AA_E_MOVE_IMAGE_PRIOR calculate image prior
 Reg= aa_e_move_image_prior( inv_model )
 Reg        => output regularization term
 inv_model  => inverse model struct
 Parameters:
   inv_model.image_prior.parameters(1) -> relative weighting
     of movement vs image fraction of hyperparameter
     => Default = 100
   inv_model.aa_e_move_image_prior.RegC.func = Cond Reg fcn

 For image portion, we use a Laplace prior, as 
 -1 for each adjacent element, and 3 (in 2D) or 4 (in 3D)
 for the element itself

 For the movmenent portion, we define a smoothness
 constraint, such that Rij = -1 for adjacent electrodes

 If used with a dual model (ie coarse2fine mapping), ensure
    imdl.prior_use_fwd_not_rec = 1;

0001 function Reg= aa_e_move_image_prior( inv_model );
0002 % AA_E_MOVE_IMAGE_PRIOR calculate image prior
0003 % Reg= aa_e_move_image_prior( inv_model )
0004 % Reg        => output regularization term
0005 % inv_model  => inverse model struct
0006 % Parameters:
0007 %   inv_model.image_prior.parameters(1) -> relative weighting
0008 %     of movement vs image fraction of hyperparameter
0009 %     => Default = 100
0010 %   inv_model.aa_e_move_image_prior.RegC.func = Cond Reg fcn
0011 %
0012 % For image portion, we use a Laplace prior, as
0013 % -1 for each adjacent element, and 3 (in 2D) or 4 (in 3D)
0014 % for the element itself
0015 %
0016 % For the movmenent portion, we define a smoothness
0017 % constraint, such that Rij = -1 for adjacent electrodes
0018 %
0019 % If used with a dual model (ie coarse2fine mapping), ensure
0020 %    imdl.prior_use_fwd_not_rec = 1;
0021 
0022 % (C) 2005 Andy Adler. License: GPL version 2 or version 3
0023 % $Id: aa_e_move_image_prior.html 2819 2011-09-07 16:43:11Z aadler $
0024 
0025 % relative strengths of conductivity and movement priors
0026 if isfield( inv_model,'aa_e_move_image_prior')
0027    hp_move= inv_model.aa_e_move_image_prior.parameters(1);
0028 else
0029    hp_move= 10;
0030 end
0031 
0032 pp= aa_fwd_parameters( inv_model.fwd_model );
0033 n_elec = pp.n_elec;
0034 
0035 % calc conductivity portion
0036 try
0037    RegCfcn = inv_model.aa_e_move_image_prior.RegC.func;
0038 catch
0039    RegCfcn = @laplace_image_prior;
0040 end
0041 try; inv_model = rmfield(inv_model, 'R_prior'); end
0042 try; inv_model = rmfield(inv_model, 'prior_use_fwd_not_rec'); end
0043 inv_model.RtR_prior = RegCfcn;
0044 
0045 pp= aa_fwd_parameters( inv_model.fwd_model );
0046 
0047 
0048 RegC= calc_RtR_prior( inv_model); 
0049 szC = size(RegC,1);
0050 
0051 % calc movement portion
0052 RegM = movement_image_prior( pp.n_dims, pp.n_elec );
0053 
0054 % zero blocks
0055 RegCM= sparse( szC, pp.n_dims*pp.n_elec );
0056 
0057 Reg= [RegC,           RegCM;
0058       RegCM', hp_move^2*RegM ];
0059 
0060 % For the movmenent portion, we define a smoothness
0061 % constraint, such that Rij = -1 for adjacent electrodes
0062 %
0063 % TODO: this approach assumes that electrodes close in number
0064 %   are close to each other. This isn't necessarily right.
0065 function RegM = movement_image_prior( dims, elecs );
0066 
0067    % movement constraint in each dimention
0068    idx =(0:elecs-1)';
0069    im1= rem(idx-1+elecs,elecs);
0070    ip1= rem(idx+1,elecs); 
0071    mv= sparse([idx,idx,idx]+1,[im1,idx,ip1]+1,ones(elecs,1)*[-1,2.1,-1]);
0072 
0073    RegM= spalloc(dims*elecs,dims*elecs, 3*dims*elecs);
0074 
0075    for i=0:dims-1;
0076      m_idx= idx + i*elecs + 1;
0077      RegM( m_idx, m_idx ) = mv;
0078    end
0079