function yi=pccp(x,y,xi)
n=max(size(xi));
m=max(size(x));
dx=(xi(n)-xi(1))/(n-1);
xv=xi(1);
j=1;
for i=1:n
	if (xv<=x(2))
		yv=linear(x(1),y(1),x(2),y(2),xv);
	end
	if (xv>=x(m-1))
		yv=linear(x(m-1),y(m-1),x(m),y(m),xv);
	end
	if (xv>x(2)&xv<x(m-1))
		% use pccp interpolation
		x0=x(j-1);		y0=y(j-1);
		x1=x(j);		y1=y(j);
		x2=x(j+1);		y2=y(j+1);
		x3=x(j+2);		y3=y(j+2);
		yv=pccp1(x0,y0,x1,y1,x2,y2,x3,y3,xv);
	end
xv=xv+dx;
yi(i)=yv;
	if (xv>x(j+1))
		j=j+1;
	end
end