function kappa=kappa2d(phi,h)

[m,n]=size(phi); 
phixx=(phi(3:m,2:n-1)-2*phi(2:m-1,2:n-1)+phi(1:m-2,2:n-1))/h^2;
phiyy=(phi(2:m-1,3:n)-2*phi(2:m-1,2:n-1)+phi(2:m-1,1:n-2))/h^2;
phixy=(phi(3:m,3:n)-phi(1:m-2,3:n)-phi(3:m,1:n-2)+phi(1:m-2,1:n-2))/(4*h*h);
phix=(phi(3:m,2:n-1)-phi(1:m-2,2:n-1))/(2*h);
phiy=(phi(2:m-1,3:n)-phi(2:m-1,1:n-2))/(2*h);
mag=sqrt(phix.^2+phiy.^2);
mag(mag==0)=1;
kappa=(phixx.*phiy.^2-2*phiy.*phix.*phixy+phiyy.*phix.^2)./mag.^3;

kappa1=zeros(m,n);
kappa1(2:m-1,2:n-1)=kappa;
kappa=kappa1;