current.m

MATLAB二维电阻抗断层成像算法!用于医学成像,里面包括有限元剖分正问题,及反问题的算法.并且附有网络剖分数据表!

function [I,T]=Current(L,lg,style,rms,numpat);

%Current Supplies some current patterns for 2D EIT
% Function [I,T]=Current(L,lg,style,rms,numpat) calculates 
% different current patterns. The 'style' can be 'tri' (三角形注入), 
% 'adj' (相邻注入), 'opp' (相对注入) or 'ref' (一个参考点注入). 
%
% INPUT
% L = 电极数
% lg =剖分节点数
% style = 电流注入模式
% rms =注入电流的均方根值. Default=1 (optional).
% numpat =电流注入次数(optional, in certain cases forced to be L-1)
%
% OUTPUT
% I = 全电流矩阵包含电极处和内部节点的电流 
% T = 电极电流矩阵
%****************************
if nargin < 4 & style == 'tri'
 rms=1;
 numpat=L-1;
end

if nargin < 4 & style ~= 'tri'
 rms=1;
 numpat=L;
end

if nargin == 4 & style == 'tri'
 numpat=L-1;
end

if nargin == 4 & style ~= 'tri'
 numpat=L;
end
%****************************

switch style
case 'tri'   %三角形注入的电流矩阵计算
if numpat > L-1,
  str=['Too many current patterns. Automatically set to ', num2str(L-1),'.' ];
  disp(str)
  numpat =L-1;
 end
II=zeros(L,numpat);
l=(1:L)';
th=2*pi*l/L;
 for k=1:(numpat+1)/2,
  II(:,k)=rms*cos(k*th);
 end
 for k=(numpat+1)/2+1:numpat
  II(:,k)=rms*sin((k-L/2)*th);
 end
I=[zeros(lg,numpat);II];
T=II;

case 'adj'   %相邻电极注入的电流矩阵计算
II=zeros(L,numpat);
l=(1:L)';
II1=diag(ones(L,1));
II2=diag(ones(L-1,1),-1);
if numpat < L
 II=II1(:,1:numpat)-II2(:,1:numpat);
else
 II=II1(:,1:numpat-1)-II2(:,1:numpat-1);
 II=[II,[-1;zeros(L-2,1);1]];% Optional!!
end
I=[zeros(lg,numpat);rms*II];
T=II;

case 'ref' %参考电极注入的电流矩阵计算
 if numpat > L-1, 
  str=['Too many current patterns. Automatically set to ', num2str(L-1),'.' ];
  disp(str)
  numpat =L-1;
 end
II=zeros(L,numpat);
II1=diag(ones(L-1,1),-1);
II(1,:)=ones(1,numpat);
T=II-II1(:,1:numpat);
I=[zeros(lg,numpat);rms*T];


case 'opp'    %相对电极注入的电流矩阵计算
if numpat > L/2, 
str=['Too many current patterns. Automatically set to ', num2str(L/2),'.' ];
  disp(str)
numpat=L/2;end
II=zeros(L,numpat);
l=(1:L)';
II1=diag(ones(L,1));
II2=diag(ones(L,1),-8);
II=II1(:,1:numpat)-II2(1:L,1:numpat);
%II=[II(:,1:ceil(numpat/2)),-II(:,1:floor(numpat/2))];% Optional!!
I=[zeros(lg,numpat);rms*II];
T=II;
end