![[Home]](/images/iso2mesh_2015_banner.png)
Showing revision 1.0
A registration toolbox for FEM surface meshes
Author: Qianqian Fang <fangq at nmr.mgh.harvard.edu>
- Martinos Center for Biomedical Imaging
- Massachusetts General Hospital (Harvard Medical School)
- Bldg. 149, 13th St., Charlestown, MA 02148
Table of Content
- 1. Introduction
- 2. List of functions
- 3. Acknowledgement
1. Introduction
"Metch", coined from "mesh" and "match", is a Matlab/Octave- based mesh registration toolbox. It provides straightforward functions to register point clouds (or surfaces) to a triangular/cubic surface mesh by calculating an optimal affine transformation (in terms of matrix A for scaling and rotation, and b for translation). It also allows one to project a point cloud onto the surface using surface norms and guarantee the conformness of the points to the surface.
2. List of functions
regpt2surf.m
[A,b,newpos]=regpt2surf(node,elem,p,pmask,A0,b0,cmask,maxiter)Perform point cloud registration to a triangular surface (surface can be either triangular or cubic), Gauss-Newton method is used for the calculation
parameters:
node: node coordinates of the surface mesh (nn x 3)
elem: element list of the surface mesh (3 columns for
triangular mesh, 4 columns for cubic surface mesh)
p: points to be registered, 3 columns for x,y and z respectively
pmask: a mask vector with the same length as p, determines the
method to handle the point, if pmask(i)=-1, the point is a free
node and can be move by the optimization, if pmask(i)=0, the
point is fixed; if pmask(i)=n>0, the distance between p(i,:)
and node(n,:) will be part of the object function and be optimized
A0: a 3x3 matrix, as the initial guess for the affine A matrix (rotation&scaling)
b0: a 3x1 vector, as the initial guess for the affine b vector (translation)
cmask: a binary 12x1 vector, determines which element of [A(:);b] will be optimized
maxiter: a integer, specifying the optimization iterations
outputs:
newpos: the registered positions for p, newpos=A*p'+b
proj2mesh.m
[newpt elemid weight]=proj2mesh(v,f,pt,nv,cn)
project a point cloud on to the surface mesh (surface can only be triangular)
parameters:
v: node coordinates of the surface mesh (nn x 3)
f: element list of the surface mesh (3 columns for
triangular mesh, 4 columns for cubic surface mesh)
pt: points to be projected, 3 columns for x,y and z respectively
nv: nodal norms (vector) calculated from nodesurfnorm.m
with dimensions of (size(v,1),3)
cn: a integer vector with the length of p, denoting the closest
surface nodes (indices of v) for each point in p. this
value can be calculated from dist2surf.m
if nv and cn are not supplied, proj2mesh will project the point
cloud onto the surface by the direction pointing to the centroid
of the mesh
outputs:
newpt: the projected points from p
elemid: a vector of length of p, denotes which surface trangle (in elem)
contains the projected point
weight: the barycentric coordinates for each projected points, these are
the weights
dist2surf.m
[d2surf,cn]=dist2surf(node,nv,p,cn)
calculate the distances from a point cloud to a surface, and return the indices of the closest surface node
parameters:
node: node coordinates of the surface mesh (nn x 3)
nv: nodal norms (vector) calculated from nodesurfnorm.m
with dimensions of (size(node,1),3), this can be
calculated from nodesurfnorm.m
pt: points to be calculated, 3 columns for x,y and z respectively
outputs:
d2surf: a vector of length of p, the distances from p(i) to the surface
cn: a integer vector with the length of p, the indicies of the closest surface node
getplanefrom3pt.m
[a,b,c,d]=getplanefrom3pt(plane)
calculate the plane equation coefficients for a plane (determined by 3 points), the plane equation is a*x+b*y+c*z+d=0
parameters:
plane: a 3x3 matrix, each row is a 3d point in form of (x,y,z)
this is used to define a plane
outputs:
a,b,c,d: the coefficients of the plane equation
linextriangle.m
[isinside,pt,coord]=linextriangle(p0,p1,plane)
calculate the intersection of a 3d line (passing two points) with a plane (determined by 3 points)
parameters:
p0: a 3d point in form of (x,y,z)
p1: another 3d point in form of (x,y,z), p0 and p1 determins the line
plane: a 3x3 matrix, each row is a 3d point in form of (x,y,z)
this is used to define a plane
outputs:
isinside: a boolean variable, 1 for the intersection is within the
3d triangle determined by the 3 points in plane; 0 is outside
pt: the coordinates of the intersection pint
coord: 1x3 vector, if isinside=1, coord will record the barycentric
coordinates for the intersection point within the triangle;
otherwise it will be all zeros.
mapcoord4.m
[A,b]=mapcoord4(pfrom,pto)
calculate an affine transform (A matrix and b vector) to map 4 vertices from one space to the other
parameters:
pfrom: 4x3 matrix, each row is a 3d point in original space
pto: 4x3 matrix, each row is a 3d point in the mapped space
triangular mesh, 4 columns for cubic surface mesh)
outputs:
the solution will satisfy the following equation: A*pfrom'+b=pto
nodesurfnorm.m
nv=nodesurfnorm(node,elem)
calculate a nodal norm for each vertix on a surface mesh (surface can only be triangular or cubic)
parameters:
node: node coordinates of the surface mesh (nn x 3)
elem: element list of the surface mesh (3 columns for
triangular mesh, 4 columns for cubic surface mesh)
pt: points to be projected, 3 columns for x,y and z respectively
outputs:
nv: nodal norms (vector) calculated from nodesurfnorm.m
with dimensions of (size(v,1),3)
trisurfnorm.m
ev=trisurfnorm(node,elem)calculate the surface norms for each element (surface can be either triangular or cubic)
parameters:
node: node coordinates of the surface mesh (nn x 3)
elem: element list of the surface mesh (3 columns for
triangular mesh, 4 columns for cubic surface mesh)
outputs:
ev: norm vector for each surface element
3. Acknowledgement
This toolbox was developed with the support from NIH grant titled "Dynamic Inverse Solutions for Multimodal Imaging" (R01EB006385)