iso_s_smooth.m

用来实现三维阻抗及光学断层成像重建的matlab程序

function [sm2] = iso_s_smooth(simp,w);
%function [sm2] = iso_f_smooth(simp,w);
%
%Calculates a second order discrete Gaussian smoothing operator
%
%
%
%simp = The simplices matrix.
%w    = smoothing weight, w=1...k, default value = 1
%sm2  = The second order smoothing regulariser. 


if nargin<2 
   w=1;
end

if w<0
   error('Weight must be possitive');
end


ndsrch = [simp,(1:size(simp,1))'];
sm2 =spalloc(size(simp,1),size(simp,1),80*size(simp,1));

   
for i=1:size(ndsrch,1)
      
    t_id = ndsrch(i,1:size(simp,2));
      
    X2simp = [];
    Xsimp = [];
    
 for j =1:length(t_id)
      
      t_nd = t_id(j);
      
      X1 = find(simp(:,1)==t_nd);
      X2 = find(simp(:,2)==t_nd);
      X3 = find(simp(:,3)==t_nd);
      X4 = find(simp(:,3)==t_nd);
      
      XX = unique([X1;X2;X3;X4]);
      
      Xsimp = [Xsimp;XX];
   end
   
   Joined = simp(Xsimp,:);
   
   Jr = unique(Joined);
   
   for k=1:length(Jr)
      
      t_2nd = Jr(k);
      
      XX1 = find(simp(:,1)==t_2nd);
      XX2 = find(simp(:,2)==t_2nd);
      XX3 = find(simp(:,3)==t_2nd);
      XX4 = find(simp(:,4)==t_2nd);
      
      XXX = unique([XX1;XX2;XX3;XX4]);
      
      X2simp = [X2simp;XXX];
      
   end
   
   Xdif = setdiff(Xsimp,i);
   Xdif2 = setdiff(X2simp,i);
   
   sm2(i,Xdif) = -0.5*w;
   sm2(i,Xdif2) = -1/w;
   sm2(i,i) = abs(sum(sm2(i,:)));
end



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
% This is part of the EIDORS suite.
% Copyright (c) N. Polydorides 2001
% Copying permitted under terms of GNU GPL
% See enclosed file gpl.html for details.
% EIDORS 3D version 1.0
% MATLAB version 5.3 R11
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%