vertex_area

Compute area around vertex, two variants available: "one_ring" and "mixed", see paper [1] for more details.

  1. Meyer, Mark, et al. "Discrete differential-geometry operators for triangulated 2-manifolds." Visualization and mathematics III. Springer Berlin Heidelberg, 2003. 35-57.

Contents

Syntax

va = vertex_area(face,vertex)
va = vertex_area(face,vertex,type)

Description

face  : double array, nf x 3, connectivity of mesh
vertex: double array, nv x 3, vertex of mesh
type  : string, area type, either "one_ring", or "mixed"
va: double array, nv x 1, area of all vertex.

Example

va = vertex_area(face,vertex)
va = vertex_area(face,vertex,'one_ring') % same as last
va = vertex_area(face,vertex,'mixed')

Contribution

Author : Wen Cheng Feng
Created: 2014/03/03
Revised: 2014/03/28 by Wen, add code and doc
Copyright 2014 Computational Geometry Group
Department of Mathematics, CUHK
http://www.math.cuhk.edu.hk/~lmlui
function va = vertex_area(face,vertex,type)
if ~exist('type','var')
    type = 'one_ring';
end
fa = face_area(face,vertex);

switch type
    case 'one_ring'
        [he,heif] = compute_halfedge(face);
        va = accumarray(he(:,1),fa(heif));
    case 'mixed'
        dvf12 = vertex(face(:,2),:)-vertex(face(:,1),:);
        dvf23 = vertex(face(:,3),:)-vertex(face(:,2),:);
        dvf31 = vertex(face(:,1),:)-vertex(face(:,3),:);
        c1 = vector_cot(dvf12,-dvf31);
        c2 = vector_cot(dvf23,-dvf12);
        c3 = vector_cot(dvf31,-dvf23);
        vaf1 = (dot(dvf12,dvf12,2).*c3+dot(dvf31,dvf31,2).*c2)/8;
        vaf2 = (dot(dvf23,dvf23,2).*c1+dot(dvf12,dvf12,2).*c3)/8;
        vaf3 = (dot(dvf31,dvf31,2).*c2+dot(dvf23,dvf23,2).*c1)/8;
        ind1 = c1<0;
        vaf1(ind1) = fa(ind1)/2;
        vaf2(ind1) = fa(ind1)/4;
        vaf3(ind1) = fa(ind1)/4;
        ind2 = c2<0;
        vaf1(ind2) = fa(ind2)/4;
        vaf2(ind2) = fa(ind2)/2;
        vaf3(ind2) = fa(ind2)/4;
        ind3 = c3<0;
        vaf1(ind3) = fa(ind3)/4;
        vaf2(ind3) = fa(ind3)/4;
        vaf3(ind3) = fa(ind3)/2;
        va = accumarray(face(:),[vaf1;vaf2;vaf3]);
    otherwise
        error('unknown area type.')
end

function vc = vector_cot(v1,v2)
% cot of angle between two vectors
cs = dot(v1,v2,2);
d1 = dot(v1,v1,2);
d2 = dot(v2,v2,2);
vc = cs./sqrt(d1.*d2-cs.*cs);