zisejinghua
2011-12-23, 20:16
function [newx newy] = FourPoint(x,y,m)
n=numel(x);
oldx=x;
oldy=y;
w=1/16;
for i=1:1:m
for j=1:1:n
newx(2*j)=oldx(j);
newy(2*j)=oldy(j);
if(j<1.5)
newx(2*j+1)=(1/2+w)*(oldx(j)+oldx(j+1))-w*(oldx(n)+oldx(j+2));
newy(2*j+1)=(1/2+w)*(oldy(j)+oldy(j+1))-w*(oldy(n)+oldy(j+2));
end
if(j>n-1.5)
if(j>n-0.5)
newx(2*j+1)=(1/2+w)*(oldx(j)+oldx(1))-w*(oldx(j-1)+oldx(2));
newy(2*j+1)=(1/2+w)*(oldy(j)+oldy(1))-w*(oldy(j-1)+oldy(2));
else
newx(2*j+1)=(1/2+w)*(oldx(j)+oldx(j+1))-w*(oldx(j-1)+oldx(1));
newy(2*j+1)=(1/2+w)*(oldy(j)+oldy(j+1))-w*(oldy(j-1)+oldy(1));
end
end
end
newx(1)=newx(2*n+1);
newy(1)=newy(2*n+1);
n=2*n;
oldx=newx;
oldy=newy;
end
n=numel(x);
oldx=x;
oldy=y;
w=1/16;
for i=1:1:m
for j=1:1:n
newx(2*j)=oldx(j);
newy(2*j)=oldy(j);
if(j<1.5)
newx(2*j+1)=(1/2+w)*(oldx(j)+oldx(j+1))-w*(oldx(n)+oldx(j+2));
newy(2*j+1)=(1/2+w)*(oldy(j)+oldy(j+1))-w*(oldy(n)+oldy(j+2));
end
if(j>n-1.5)
if(j>n-0.5)
newx(2*j+1)=(1/2+w)*(oldx(j)+oldx(1))-w*(oldx(j-1)+oldx(2));
newy(2*j+1)=(1/2+w)*(oldy(j)+oldy(1))-w*(oldy(j-1)+oldy(2));
else
newx(2*j+1)=(1/2+w)*(oldx(j)+oldx(j+1))-w*(oldx(j-1)+oldx(1));
newy(2*j+1)=(1/2+w)*(oldy(j)+oldy(j+1))-w*(oldy(j-1)+oldy(1));
end
end
end
newx(1)=newx(2*n+1);
newy(1)=newy(2*n+1);
n=2*n;
oldx=newx;
oldy=newy;
end