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