www.pudn.com > balloonSnake(matlab).rar > snakedeform2.m


function [x,y] = snakedeform(x,y,alpha,beta,gamma,kappa, kappap, fx,fy,ITER)
> SNAKEDEFORM2 Deform snake in the given external force field with
> pressure force added
> [x,y] = snakedeform(x,y,alpha,beta,gamma,kappa,kappap,fx,fy,ITER)
>
> alpha: elasticity parameter
> beta: rigidity parameter
> gamma: viscosity parameter
> kappa: external force weight
> kappap: pressure force weight
> fx,fy: external force field

> Chenyang Xu and Jerry L. Prince, 4/1/95, 6/17/97
> Copyright (c) 1995-97 by Chenyang Xu and Jerry L. Prince
> Image Analysis and Communications Lab, Johns Hopkins University


> generates the parameters for snake

N = length(x);

alpha = alpha* ones(1,N);
beta = beta*ones(1,N);

> produce the five diagnal vectors
alpham1 = [alpha(2:N) alpha(1)];
alphap1 = [alpha(N) alpha(1:N-1)];
betam1 = [beta(2:N) beta(1)];
betap1 = [beta(N) beta(1:N-1)];

a = betam1;
b = -alpha - 2*beta - 2*betam1;
c = alpha + alphap1 +betam1 + 4*beta + betap1;
d = -alphap1 - 2*beta - 2*betap1;
e = betap1;

> generate the parameters matrix
A = diag(a(1:N-2),-2) + diag(a(N-1:N),N-2);
A = A + diag(b(1:N-1),-1) + diag(b(N), N-1);
A = A + diag(c);
A = A + diag(d(1:N-1),1) + diag(d(N),-(N-1));
A = A + diag(e(1:N-2),2) + diag(e(N-1:N),-(N-2));

invAI = inv(A + gamma * diag(ones(1,N)));

for count = 1:ITER,
vfx = interp2(fx,x,y,'*linear');
vfy = interp2(fy,x,y,'*linear');

> adding the pressure force
xp = [x(2:N);x(1)]; yp = [y(2:N);y(1)];
xm = [x(N);x(1:N-1)]; ym = [y(N);y(1:N-1)];

qx = xp-xm; qy = yp-ym;
pmag = sqrt(qx.*qx+qy.*qy);
px = qy./pmag; py = -qx./pmag;

> deform snake
x = invAI * (gamma* x + kappa*vfx + kappap.*px);
y = invAI * (gamma* y + kappa*vfy + kappap.*py);
end