vb_divide_patch

 divide triangle by 4 :
 make new vertex at midpoint of triangle edge
  [V, F, Iparent] = vb_divide_patch(V, F, xxT, Vindex)
 
 xxT{n} : 頂点 n の 隣接３角面番号
 Iparent(n,1)   : vertex index for new vertex
 Iparent(n,2:3) : Parents vertex index for new vertex
   Vmid = ( V(Iparent(:,2),:) + V(Iparent(:,3),:) )/2
   V(Iparent(:,1),:) = Vmid./repmat( sqrt(sum(Vmid.^2,2)) ,1 , 3);

0001 function    [V, F, Iparent] = vb_divide_patch(V, F, xxT, Vindex)
0002 % divide triangle by 4 :
0003 % make new vertex at midpoint of triangle edge
0004 %  [V, F, Iparent] = vb_divide_patch(V, F, xxT, Vindex)
0005 %
0006 % xxT{n} : 頂点 n の 隣接３角面番号
0007 % Iparent(n,1)   : vertex index for new vertex
0008 % Iparent(n,2:3) : Parents vertex index for new vertex
0009 %   Vmid = ( V(Iparent(:,2),:) + V(Iparent(:,3),:) )/2
0010 %   V(Iparent(:,1),:) = Vmid./repmat( sqrt(sum(Vmid.^2,2)) ,1 , 3);
0011 
0012 NV = size(V,1);
0013 NF = size(F,1);
0014 
0015 % target patch index for divide
0016 DFindex = vb_cell_merge(xxT,Vindex);
0017 DFindex = unique(DFindex );
0018 Iremain = vb_setdiff2([1:NF],DFindex);
0019 
0020 % patch without change
0021 Fremain = F(Iremain, : );
0022 
0023 % target patch for divide
0024 divF = F(DFindex, : );
0025 % target patch number
0026 NDF = length(DFindex);
0027 
0028 % Vertex coordinate for each Triangle
0029 v1 = V(divF(:,1),:);
0030 v2 = V(divF(:,2),:);
0031 v3 = V(divF(:,3),:);
0032 
0033 % Midpoint for 3 edge
0034 f1 = (v2+v3)/2;
0035 f2 = (v3+v1)/2;
0036 f3 = (v1+v2)/2;
0037 
0038 
0039 % Midpoint vertex coordinate list
0040 Vf = [f1;
0041       f2;
0042       f3];
0043 
0044 % Edge index for midpoint
0045 Jedge  = [ divF(:,2),divF(:,3); 
0046            divF(:,3),divF(:,1); 
0047            divF(:,1),divF(:,2)];
0048 
0049 Jedge  = sort(Jedge,2);
0050 
0051 % [B,I,J] = unique(A)        : B = A(I)   and  A = B(J)
0052 % [B,I,J] = unique(A,'rows') : B = A(I,:) and  A = B(J,:))
0053 
0054 % Extract unique edge index for midpoint
0055 [Juniq, II, JJ ] = unique( Jedge, 'rows');
0056 
0057 % Unique edge index
0058 % Juniq = Jedge(II,:)
0059 % Jedge = Juniq(JJ,:)
0060 
0061 % Unique vertex coordinate for midpoint
0062 Vnew = Vf(II,:);
0063 Nnew = size(Vnew,1);
0064 
0065 % Append midpoint after original V
0066 V    = [ V ; Vnew]; 
0067 Inew = (1:Nnew) + NV;
0068 
0069 % Midpoint index from Vnew index
0070 % Vnew = Vf(II,:);
0071 % Vf   = Vnew(JJ,:);
0072 Jnew = JJ + NV;
0073 
0074 ix   = 1:NDF;
0075 
0076 % Index of midpoint for [ divF(:,2),divF(:,3)];
0077 J1   = Jnew(ix);
0078 % Index of midpoint for [ divF(:,3),divF(:,1)];
0079 J2   = Jnew(ix+NDF);
0080 % Index of midpoint for [ divF(:,1),divF(:,2) ];
0081 J3   = Jnew(ix+2*NDF);
0082 
0083 % Patch index
0084 Fadd = [ divF(:,1), J2, J3; ...
0085          divF(:,2), J3, J1; ...
0086          divF(:,3), J1, J2; ...
0087          J1,J2,J3 ];
0088 
0089 F = [Fremain; Fadd];
0090 
0091 
0092 % Parent index for new midpoint
0093 Iparent = [J1, divF(:,2), divF(:,3);
0094            J2, divF(:,3), divF(:,1);
0095            J3, divF(:,1), divF(:,2)];
0096 
0097 % [B,I,J] = unique(A)  : B = A(I)   and  A = B(J)
0098 [Iuniq , I, J]= unique(Iparent(:,1));
0099 
0100 Iparent = Iparent(I,:);
0101 
0102 return
0103 
0104 % check Iparent
0105 Nerr = Nnew - size(Iparent,1)
0106 
0107 vmid = V(Iparent(:,2),:) + V(Iparent(:,3),:) ;
0108 % Normalization
0109 vmid = vmid./repmat( sqrt(sum(vmid.^2,2)) ,1 , 3);
0110 
0111 err = sum((V(Iparent(:,1),:) - vmid).^2,2);
0112 
0113 errsum = sum(err)
0114 
0115 %   % Make 4 triangle patch for each triangle
0116 %   vv1  = [ v1; v2; v3; f1];
0117 %   vv2  = [ f2; f3; f1; f2];
0118 %   vv3  = [ f3; f1; f2; f3];
0119 %
0120 %   v1   = vv1;
0121 %   v2   = vv2;
0122 %   v3   = vv3;
0123 %