AVERAGEMESH Compute average mesh from a list of meshes. AVGMESH = averageMesh(MESHLIST) Where MESHLIST is a cell array of meshes, computes an average mesh that minimizes the sum of squared distances between average mesh vertices and all other mesh vertices. The method is to choose an arbitrary reference mesh, and to iterate a series of smoothin of the reference mesh, computation of nearest vertex neighbors, and computing average position of nearest neighbors. [AVGMESH, DISTS] = averageMesh(MESHLIST) Returns also a cell array containing for each iteration, the standard deviation of distances to other individual meshes. [AVGMESH, DISTS, VERT_ITERS] = averageMesh(MESHLIST) Also returns for each iteration, the positions of the average vertices. [AVGMESH, DISTS] = averageMesh(..., PNAME, PVALUE) Provides addition options are parameter name-value pairs. Available options are: * verbose: (logical, default false) display or not information about process * nIters: number of smooth-projection iterations to perform. Default value is 10. * smoothingIteration: the number of smoothing operations to apply on average mesh at each iteration. Default value is 3. Example averageMesh See also meshes3d, smoothMesh
0001 function [avgMesh, distsIters, verticesIters] = averageMesh(meshList, varargin) 0002 %AVERAGEMESH Compute average mesh from a list of meshes. 0003 % 0004 % AVGMESH = averageMesh(MESHLIST) 0005 % Where MESHLIST is a cell array of meshes, computes an average mesh that 0006 % minimizes the sum of squared distances between average mesh vertices 0007 % and all other mesh vertices. 0008 % The method is to choose an arbitrary reference mesh, and to iterate a 0009 % series of smoothin of the reference mesh, computation of nearest vertex 0010 % neighbors, and computing average position of nearest neighbors. 0011 % 0012 % [AVGMESH, DISTS] = averageMesh(MESHLIST) 0013 % Returns also a cell array containing for each iteration, the standard 0014 % deviation of distances to other individual meshes. 0015 % 0016 % [AVGMESH, DISTS, VERT_ITERS] = averageMesh(MESHLIST) 0017 % Also returns for each iteration, the positions of the average vertices. 0018 % 0019 % [AVGMESH, DISTS] = averageMesh(..., PNAME, PVALUE) 0020 % Provides addition options are parameter name-value pairs. Available 0021 % options are: 0022 % * verbose: (logical, default false) display or not information about 0023 % process 0024 % * nIters: number of smooth-projection iterations to perform. Default 0025 % value is 10. 0026 % * smoothingIteration: the number of smoothing operations to apply on 0027 % average mesh at each iteration. Default value is 3. 0028 % 0029 % 0030 % Example 0031 % averageMesh 0032 % 0033 % See also 0034 % meshes3d, smoothMesh 0035 0036 % ------ 0037 % Author: David Legland 0038 % E-mail: david.legland@inrae.fr 0039 % Created: 2020-01-31, using Matlab 9.7.0.1247435 (R2019b) Update 2 0040 % Copyright 2020-2024 INRAE - BIA Research Unit - BIBS Platform (Nantes) 0041 0042 %% Parse input values 0043 0044 % default values 0045 nIters = 10; 0046 verbose = false; 0047 smoothingIterations = 3; 0048 0049 % parse input arguments 0050 while length(varargin) > 1 0051 name = varargin{1}; 0052 if ~ischar(name) 0053 error('require parameter name-value pairs'); 0054 end 0055 0056 if strcmpi(name, 'verbose') 0057 verbose = varargin{2}; 0058 elseif strcmpi(name, 'nIters') 0059 nIters = varargin{2}; 0060 elseif strcmpi(name, 'smoothingIterations') 0061 smoothingIterations = varargin{2}; 0062 else 0063 error(['Unknown parameter name: ' name]); 0064 end 0065 varargin(1:2) = []; 0066 end 0067 0068 0069 %% Initialisations 0070 0071 nMeshes = length(meshList); 0072 0073 % initialize kd-trees to accelerate nearest-neighbor searches 0074 treeList = cell(nMeshes, 1); 0075 for iMesh = 1:nMeshes 0076 treeList{iMesh} = KDTreeSearcher(meshList{iMesh}.vertices); 0077 end 0078 0079 % choose arbitrary initial mesh 0080 avgMesh = struct('vertices', meshList{1}.vertices, 'faces', meshList{1}.faces); 0081 0082 verticesIters = cell(1, nIters); 0083 distsIters = cell(1, nIters); 0084 0085 0086 %% Main iteration 0087 0088 % iterates smoothing + computation of average projections 0089 for iIter = 1:nIters 0090 if verbose 0091 fprintf('iter %d/%d\n', iIter, nIters); 0092 end 0093 % apply smoothing to current average mesh 0094 avgMesh = smoothMesh(avgMesh, smoothingIterations); 0095 0096 % create new array for average vertices 0097 newVerts = zeros(size(avgMesh.vertices)); 0098 dists = zeros(size(avgMesh.vertices, 1), 1); 0099 0100 % iterate over all meshes 0101 for iMesh = 1:nMeshes 0102 if verbose 0103 fprintf(' mesh %d/%d\n', iMesh, nMeshes); 0104 end 0105 0106 % for each vertex of reference mesh, find index of closest vertex 0107 % in current mesh 0108 inds = knnsearch(treeList{iMesh}, avgMesh.vertices); 0109 0110 % keep position of closest vertex to update new position 0111 closest = treeList{iMesh}.X(inds,:); 0112 newVerts = newVerts + closest; 0113 0114 % keep distance to closest index to build variability map 0115 dists = dists + sum((closest - avgMesh.vertices).^2, 2); 0116 end 0117 0118 % update position of new vertices 0119 newVerts = newVerts / nMeshes; 0120 verticesIters{iIter} = newVerts; 0121 avgMesh.vertices = newVerts; 0122 0123 % keep list of distances 0124 dists = sqrt(dists / nMeshes); 0125 distsIters{iIter} = dists; 0126 end 0127 0128 0129 % figure; drawMesh(refMesh, 'lineStyle', 'none', 'faceColor', [.5 .5 .5]) 0130 % axis equal; view(3); hold on; axis([-2.5 2.5 -2 2 -3.5 3.5]); light; 0131 % lighting gouraud 0132 % title('Average mesh'); 0133 % print(gcf, 'averageMesh_initial.png', '-dpng'); 0134