MATLAB对三角网格进行线性细分

这是除去三角网格的LOOP细分外， Jesus Mena的另一篇代码，是对三角网格进行线性细分，先贴代码

function [newVertices, newFaces] =  linearSubdivision(vertices, faces)
% Linear subdivision for triangle meshes
%
%  Dimensions:
%    vertices: 3xnVertices
%    faces:    3xnFaces
%  
%  Author: Jesus Mena

	global edgeVertex;
    global newIndexOfVertices;
	newFaces = [];
	newVertices = vertices;

	nVertices = size(vertices,2);
	nFaces    = size(faces,2);
	edgeVertex= zeros(nVertices, nVertices);
	newIndexOfVertices = nVertices;

    % ------------------------------------------------------------------------ %
	% create a matrix of edge-vertices and a new triangulation (newFaces).
    % 
    % * edgeVertex(x,y): index of the new vertex between (x,y)
    %
    %  0riginal vertices: va, vb, vc.
    %  New vertices: vp, vq, vr.
    %
    %      vb                   vb             
    %     / \                  /  \ 
    %    /   \                vp--vq
    %   /     \              / \  / \
    % va ----- vc   ->     va-- vr --vc 
	%
    
	for i=1:nFaces
		[vaIndex, vbIndex, vcIndex] = deal(faces(1,i), faces(2,i), faces(3,i));
		
		vpIndex = addEdgeVertex(vaIndex, vbIndex);
		vqIndex = addEdgeVertex(vbIndex, vcIndex);
		vrIndex = addEdgeVertex(vaIndex, vcIndex);
		
		fourFaces = [vaIndex,vpIndex,vrIndex; vpIndex,vbIndex,vqIndex; vrIndex,vqIndex,vcIndex; vrIndex,vpIndex,vqIndex]';
		newFaces  = [newFaces, fourFaces]; 
    end;
    	
    % ------------------------------------------------------------------------ %
	% positions of the new vertices
	for v1=1:nVertices-1
		for v2=v1:nVertices
			vNIndex = edgeVertex(v1,v2);
            if (vNIndex~=0)
 				newVertices(:,vNIndex) = 1/2*(vertices(:,v1)+vertices(:,v2));
            end;
        end;
    end;
 	
end

% ---------------------------------------------------------------------------- %
function vNIndex = addEdgeVertex(v1Index, v2Index)
	global edgeVertex;
	global newIndexOfVertices;

	if (v1Index>v2Index) % setting: v1 <= v2
		vTmp = v1Index;
		v1Index = v2Index;
		v2Index = vTmp;
	end;
	
	if (edgeVertex(v1Index, v2Index)==0)  % new vertex
		newIndexOfVertices = newIndexOfVertices+1;
		edgeVertex(v1Index, v2Index) = newIndexOfVertices;
	end;

	vNIndex = edgeVertex(v1Index, v2Index);

    return;
end

命名为linearSubdivision.m文件。其代码的大致思想和之前LOOP细分是一样的，但是代码显然简单了不少，因为线性细分的规则是只需要在两点之间线性插值其中点，因此所用的代码简单非常多，相信不难理解，下面是用来绘图的plotMesh.m文件

function plotMesh(vertices, faces)
    hold on;
    trimesh(faces', vertices(1,:), vertices(2,:), vertices(3,:));
    colormap hot;
    axis tight;
    axis square; 
    axis off;
    view(3);
end

最后是用来测试的test.m文件

 vertices = [10 10 10; -100 10 -10; -100 -10 10; 10 -10 -10]';
faces = [1 2 3; 1 3 4; 1 4 2; 4 3 2]';

figure(2);
subplot(1,4,1);
plotMesh(vertices, faces);
for i=2:4
 subplot(1,4,i);
 [vertices, faces] = linearSubdivision(vertices, faces); 
 plotMesh(vertices, faces);
end

下面附上效果图

线性细分四面体	LOOP细分正四面体

对比一下之前的LOOP细分，可以看到二者有明显的不同。