THESE/codes/snake/basic_code/PatchNormalsDouble.m

   1 function [Nx,Ny,Nz]=PatchNormalsDouble(Fa,Fb,Fc,Vx,Vy,Vz)\r
   2 %\r
   3 %  [Nx,Ny,Nz]=PatchNormalsDouble(Fa,Fb,Fc,Vx,Vy,Vz)\r
   4 %\r
   5 \r
   6 FV.vertices=zeros(length(Vx),3);\r
   7 FV.vertices(:,1)=Vx;\r
   8 FV.vertices(:,2)=Vy;\r
   9 FV.vertices(:,3)=Vz;\r
  10 \r
  11 % Get all edge vectors\r
  12 e1=FV.vertices(Fa,:)-FV.vertices(Fb,:);\r
  13 e2=FV.vertices(Fb,:)-FV.vertices(Fc,:);\r
  14 e3=FV.vertices(Fc,:)-FV.vertices(Fa,:);\r
  15 \r
  16 % Normalize edge vectors\r
  17 e1_norm=e1./repmat(sqrt(e1(:,1).^2+e1(:,2).^2+e1(:,3).^2),1,3); \r
  18 e2_norm=e2./repmat(sqrt(e2(:,1).^2+e2(:,2).^2+e2(:,3).^2),1,3); \r
  19 e3_norm=e3./repmat(sqrt(e3(:,1).^2+e3(:,2).^2+e3(:,3).^2),1,3);\r
  20 \r
  21 % Calculate Angle of face seen from vertices\r
  22 Angle =  [acos(dot(e1_norm',-e3_norm'));acos(dot(e2_norm',-e1_norm'));acos(dot(e3_norm',-e2_norm'))]';\r
  23 \r
  24 % Calculate normal of face\r
  25  Normal=cross(e1,e3);\r
  26 \r
  27 % Calculate Vertice Normals \r
  28 VerticeNormals=zeros([size(FV.vertices,1) 3]);\r
  29 for i=1:size(Fa,1),\r
  30     VerticeNormals(Fa(i),:)=VerticeNormals(Fa(i),:)+Normal(i,:)*Angle(i,1);\r
  31     VerticeNormals(Fb(i),:)=VerticeNormals(Fb(i),:)+Normal(i,:)*Angle(i,2);\r
  32     VerticeNormals(Fc(i),:)=VerticeNormals(Fc(i),:)+Normal(i,:)*Angle(i,3);\r
  33 end\r
  34 \r
  35 V_norm=sqrt(VerticeNormals(:,1).^2+VerticeNormals(:,2).^2+VerticeNormals(:,3).^2)+eps;\r
  36 VerticeNormals=VerticeNormals./repmat(V_norm,1,3);\r
  37 Nx=VerticeNormals(:,1);\r
  38 Ny=VerticeNormals(:,2);\r
  39 Nz=VerticeNormals(:,3);\r