Author: Qianqian Fang <fangq at nmr.mgh.harvard.edu>
"iso2mesh" is a Matlab/Octave-based mesh generation toolbox. It contains a number of optimized mesh processing scripts and interfaces to several free meshing tools, such as surface mesh simplifier and extraction utilities (based on CGAL), tetgen, and mesh validation&repairing utility (based on JMeshLib). Iso2mesh toolbox provides a simple yet powerful way to create quality tetrahedral volumetric mesh directly from 3D binary, segmented or gray-scale images such as MRI or CT scans.
To be able to create volumetric finite-element mesh from 3D volumetric images is of great interest for researchers in the field of image-based modeling and analysis. Unfortunately, very limited software packages for this purpose can be found in publications and market; they are either limited in functionalities or very expensive (such as Amira from Mercury Computer Systems). This toolbox was developed as a free alternative to the expensive commertial tools and provides researchers a complete workflow from 3D binary image preprocessing (hole-filling and island-removing), to surface mesh modeling (extraction, remeshing, validation, repairing, and smoothing) and volumetric mesh creation. Users can either choose the streamlined wrapper scripts to create 3D mesh with a fully automated process, or to call invidual subroutines to perform specific meshing tasks, such as closed surface extraction or boundary extraction.
The details of this toolbox can be found in the following reference:
[node,elem,face]=v2m(img,isovalues,opt,maxvol,method) volumetric mesh generation from binary or gray-scale volumetric images shortcut for vol2mesh inputs and outputs are similar to those defined in vol2mesh
[no,el,regions,holes]=v2s(img,isovalues,opt,method) surface mesh generation from binary or gray-scale volumetric images shortcut for vol2surf inputs and outputs are similar to those defined in vol2surf
[node,elem,face]=s2m(v,f,keepratio,maxvol) volumetric mesh generation from a closed surface, shortcut for surf2mesh inputs and outputs are similar to those defined in surf2mesh
newnode=sms(node,face,iter,useralpha) simplified version of surface mesh smoothing input: node: node coordinates of a surface mesh face: face element list of the surface mesh iter: smoothing iteration number alpha: scaler, smoothing parameter, v(k+1)=alpha*v(k)+(1-alpha)*mean(neighbors) output: newnode: output, the smoothed node coordinates
[node,elem,face]=vol2mesh(img,ix,iy,iz,opt,maxvol,dofix,method,isovalues) convert a binary (or multi-valued) volume to tetrahedral mesh input: img: a volumetric binary image ix,iy,iz: subvolume selection indices in x,y,z directions opt: as defined in vol2surf.m maxvol: target maximum tetrahedral elem volume dofix: 1: perform mesh validation&repair, 0: skip repairing method: 'cgalsurf' or omit: use CGAL surface mesher 'simplify': use binsurface and then simplify 'cgalmesh': use CGAL 3.5 3D mesher for direct mesh generation [new] generally speaking, 'cgalmesh' is the most robust path if you want to product meshes from binary or multi-region volumes, however, its limitations include 1) only accept uint8 volume, and 2) can not extract meshes from gray-scale volumes. If ones goal is to process a gray-scale volume, he/she should use the 'cgalsurf' option. 'simplify' approach is not recommended unless other options failed. isovalues: a list of isovalues where the levelset is defined output: node: output, node coordinates of the tetrahedral mesh elem: output, element list of the tetrahedral mesh, the last column is the region ID face: output, mesh surface element list of the tetrahedral mesh the last column denotes the boundary ID
[no,el,regions,holes]=vol2surf(img,ix,iy,iz,opt,dofix,method,isovalues) converting a 3D volumetric image to surfaces input: img: a volumetric binary image; if img is empty, vol2surf will return user defined surfaces via opt.surf if it exists ix,iy,iz: subvolume selection indices in x,y,z directions opt: function parameters if method is 'cgalsurf': opt=a float number>1: max radius of the Delaunay sphere(element size) opt.radbound: same as above, max radius of the Delaunay sphere opt.distbound: maximum deviation from the specified isosurfaces opt(1,2,...).radbound: same as above, for each levelset if method is 'simplify': opt=a float number<1: compression rate for surf. simplification opt.keeyratio=a float less than 1: same as above, same for all surf. opt(1,2,..).keeyratio: setting compression rate for each levelset opt(1,2,..).maxsurf: 1 - only use the largest disjointed surface 0 - use all surfaces for that levelset opt(1,2,..).side: - 'upper': threshold at upper interface 'lower': threshold at lower interface opt(1,2,..).maxnode: - the maximum number of surface node per levelset opt(1,2,..).holes: user specified holes interior pt list opt(1,2,..).regions: user specified regions interior pt list opt(1,2,..).surf.{node,elem}: add additional surfaces opt(1,2,..).{A,B}: linear transformation for each surface dofix: 1: perform mesh validation&repair, 0: skip repairing method: - if method is 'simplify', iso2mesh will first call binsurface to generate a voxel-based surface mesh and then use meshresample/meshcheckrepair to create a coarser mesh; - if method is 'cgalsurf', iso2mesh will call the surface extraction program from CGAL to make surface mesh - if method is not specified, 'cgalsurf' is assumed by default isovalues: a list of isovalues where the levelset is defined output: no: list of nodes on the resulting suface mesh, 3 columns for x,y,z el: list of trianglular elements on the surface, [n1,n2,n3,region_id] regions: list of interior points for all sub-region, [x,y,z] holes: list of interior points for all holes, [x,y,z]
[node,elem,face]=surf2mesh(v,f,p0,p1,keepratio,maxvol,regions,holes,forcebox) create quality volumetric mesh from isosurface patches input parameters: v: input, isosurface node list, dimension (nn,3) if v has 4 columns, the last column specifies mesh density near each node f: input, isosurface face element list, dimension (be,3) p0: input, coordinates of one corner of the bounding box, p0=[x0 y0 z0] p1: input, coordinates of the other corner of the bounding box, p1=[x1 y1 z1] keepratio: input, percentage of elements being kept after the simplification maxvol: input, maximum tetrahedra element volume regions: list of regions, specifying by an internal point for each region holes: list of holes, similar to regions forcebox: 1: add bounding box, 0: automatic outputs: node: output, node coordinates of the tetrahedral mesh elem: output, element list of the tetrahedral mesh face: output, mesh surface element list of the tetrahedral mesh the last column denotes the boundary ID
img=surf2vol(node,face,xi,yi,zi) convert a triangular surface to a shell of voxels in a 3D image input: node: node list of the triangular surface, 3 columns for x/y/z face: triangle node indices, each row is a triangle xi,yi,zi: x/y/z grid for the resulting volume output: img: a volumetric binary image at position of ndgrid(xi,yi,zi)
[node,elem]=binsurface(img,nface) fast isosurface extraction from 3D binary images input: img: a 3D binary image nface: nface=3 or ignored - for triangular faces, nface=4 - square faces nface=0 - return a boundary mask image via node output elem: integer array with dimensions of NE x nface, each row represents a surface mesh face element node: node coordinates, 3 columns for x, y and z respectively the outputs of this subroutine can be easily plotted using patch('Vertices',node,'faces',elem,'FaceVertexCData',node(:,3), 'FaceColor','interp'); if the surface mesh has triangular faces, one can plot it with trisurf(elem,node(:,1),node(:,2),node(:,3))
[node,elem,face]=cgalv2m(vol,opt,maxvol) wrapper for CGAL 3D mesher (CGAL 3.5) convert a binary (or multi-valued) volume to tetrahedral mesh http://www.cgal.org/Manual/3.5/doc_html/cgal_manual/Mesh_3/Chapter_main.html input: vol: a volumetric binary image ix,iy,iz: subvolume selection indices in x,y,z directions opt: parameters for CGAL mesher, if opt is a structure, then opt.radbound: defines the maximum surface element size opt.angbound: defines the miminum angle of a surface triangle opt.surfaceapprox: defines the maximum distance between the center of the surface bounding circle and center of the element bounding sphere opt.reratio: maximum radius-edge ratio if opt is a scalar, it only specifies radbound. maxvol: target maximum tetrahedral elem volume output: node: output, node coordinates of the tetrahedral mesh elem: output, element list of the tetrahedral mesh, the last column is the region id face: output, mesh surface element list of the tetrahedral mesh the last column denotes the boundary ID note: each triangle will appear twice in the face list with each one attaches to each side of the interface. one can remove the redundant triangles by unique(face(:,1:3),'rows')
[node,elem,face]=cgals2m(v,f,opt,maxvol) wrapper for CGAL 3D mesher (CGAL 3.5) convert a binary (or multi-valued) volume to tetrahedral mesh http://www.cgal.org/Manual/3.5/doc_html/cgal_manual/Mesh_3/Chapter_main.html input: v: the node coordinate list of a surface mesh (nn x 3) f: the face element list of a surface mesh (be x 3) opt: parameters for CGAL mesher, if opt is a structure, then opt.radbound: defines the maximum surface element size opt.angbound: defines the miminum angle of a surface triangle opt.surfaceapprox: defines the maximum distance between the center of the surface bounding circle and center of the element bounding sphere opt.reratio: maximum radius-edge ratio if opt is a scalar, it only specifies radbound. maxvol: target maximum tetrahedral elem volume output: node: output, node coordinates of the tetrahedral mesh elem: output, element list of the tetrahedral mesh, the last column is the region id face: output, mesh surface element list of the tetrahedral mesh the last column denotes the boundary ID
[node,elem]=vol2restrictedtri(vol,thres,cent,brad,ang,radbound,distbound,maxnode) surface mesh extraction using CGAL mesher input: vol: a 3D volumetric image thres: a scalar as the threshold of of the extraction cent: a 3d position (x,y,z) which locates inside the resulting mesh, this is automatically computed from vol2surf brad: maximum bounding sphere squared of the resulting mesh ang: minimum angular constrains of the resulting tranglar elements (in degrees) radbound: maximum triangle delaunay circle radius distbound: maximum delaunay sphere distances maxnode: maximum number of surface nodes (even radbound is not reached) output: node: the list of 3d nodes in the resulting surface (x,y,z) elem: the element list of the resulting mesh (3 columns of integers)
img=surf2volz(node,face,xi,yi,zi) convert a triangular surface to a shell of voxels in a 3D image along the z-axis input: node: node list of the triangular surface, 3 columns for x/y/z face: triangle node indices, each row is a triangle xi,yi,zi: x/y/z grid for the resulting volume output: img: a volumetric binary image at position of ndgrid(xi,yi,zi)
facecell=finddisconnsurf(f) subroutine to extract disconnected surfaces from a cluster of surfaces Date: 2008/03/06 input: f: faces defined by node indices for all surface triangles output: facecell: separated disconnected surface node indices
openedge=surfedge(f) find the edge of an open surface or surface of a volume input: f: input, surface face element list, dimension (be,3) output: openedge: list of edges of the specified surface
openface=volface(t) find the surface patches of a volume input: t: input, volumetric element list, dimension (ne,4) output: openface: list of faces of the specified volume
loops=extractloops(edges) extract individual loops from an edge table of a loop collection input: edges: two column matrix recording the starting/ending points of all edge segments output: loops: output, a single vector separated by NaN, each segment is a close-polygon consisted by node IDs
[conn,connnum,count]=meshconn(elem,nn) create node neighbor list from a mesh input: elem: element table of a mesh nn : total node number of the mesh output: conn: output, a cell structure of length nn, conn{n} contains a list of all neighboring node ID for node n connnum: vector of length nn, denotes the neighbor number of each node count: total neighbor numbers
centroid=meshcentroid(v,f) compute the centroids of a mesh defined by nodes and elements (surface or tetrahedra) in R^n space input: v: surface node list, dimension (nn,3) f: surface face element list, dimension (be,3) output: centroid: centroid positions, one row for each element
nodevol=nodevolume(node,elem) calculate the Voronoi volume of each node in a simplex mesh input: node: node coordinates elem: element table of a mesh output: nodevol: volume values for all nodes
vol=elemvolume(node,elem,option) calculate the volume for a list of simplexes input: node: node coordinates elem: element table of a mesh option: if option='signed', the volume is the raw determinant, else, the results will be the absolute values output vol: volume values for all elements
[conn,connnum,count]=neighborelem(elem,nn) create node neighbor list from a mesh input: elem: element table of a mesh nn : total node number of the mesh output: conn: output, a cell structure of length nn, conn{n} contains a list of all neighboring elem ID for node n connnum: vector of length nn, denotes the neighbor number of each node count: total neighbor numbers
facenb=faceneighbors(t,opt) to find 4 face-neighboring elements of a tetrahedron input: t: tetrahedron element list, 4 columns of integers opt: if opt='surface', return boundary triangle list (should be the same as the face output from v2m) otherwise, return the element list for each element: each row contains 4 numbers, representing the element indices sharing triangular faces [1 2 3],[1 2 4],[1 3 4] and [2 3 4] in order, where 1~4 is the node local index. if the index is 0, indicating the face has no neighbor (i.e. a boundary face) output: facenb: see opt
f=maxsurf(facecell) return the surface with the maximum element number (not necessarily in area) from a cell arry of surfaces input: facecell: a cell array, each element is a face array output: f: the surface data (node indices) for the surface with most elements
mask=flatsegment(node,edge) decompose edge loops into flat segments alone arbitrary planes of the bounding box this code is fragile: it can not handle curves with many co-linear nodes near the corner point input: node: x,y,z coordinates of each node of the mesh edge: input, a single vector separated by NaN, each segment is a close-polygon consisted by node IDs output: mask: output, a cell, each element is a close-polygon on x/y/z plane
[newedge]=orderloopedge(edge) order the node list of a simple loop based on connection sequence input: edge: a loop consisted by a sequence of edges, each row is an edge with two integers: start/end node index output: newedge: reordered edge node list
[X,V,E,F]=mesheuler(face) Euler's charastistics of a mesh input: face: a closed surface mesh output: X: Euler's charastistics V: number of vertices E: number of edges F: number of faces
seg=bbxflatsegment(node,loop) decompose edge loops into flat segments alone x/y/z planes of the bounding box input: node: x,y,z coordinates of each node of the mesh loop: input, a single vector separated by NaN, each segment is a close-polygon consisted by node IDs output: seg: output, a single vector separated by NaN, each segment is a close-polygon on x/y/z plane
[node,elem]=meshcheckrepair(node,elem,opt) check and repair a surface mesh input/output: node: input/output, surface node list, dimension (nn,3) elem: input/output, surface face element list, dimension (be,3) opt: options, including 'duplicated': remove duplicated elements 'isolated': remove isolated nodes 'deep': call external jmeshlib to remove non-manifold vertices
newelem=meshreorient(node,elem) reorder nodes in a surface or tetrahedral mesh to ensure all elements are oriented consistently input: node: list of nodes elem: list of elements (each row are indices of nodes of each element) output: newelem: the element list with consistent ordering
elem=removedupelem(elem) remove doubly duplicated elements input: elem: list of elements (node indices) output: elem: element list after removing the duplicated elements
[newnode,newelem]=removedupnodes(node,elem) removing the duplicated nodes from a mesh input: elem: integer array with dimensions of NE x 4, each row contains the indices of all the nodes for each tetrahedron node: node coordinates, 3 columns for x, y and z respectively output: newnode: nodes without duplicates newelem: elements with only the unique nodes
[no,el]=removeisolatednode(node,elem) remove isolated nodes: nodes that are not included in any element input: node: list of node coordinates elem: list of elements of the mesh output: no: node coordinates after removing the isolated nodes el: element list of the resulting mesh
fnew=removeisolatedsurf(v,f,maxdiameter) remove disjointed surface fragment by using mesh diameter input: v: list of nodes of the input surface f: list of triangles of the input surface maxdiameter: maximum bounding box size for surface removal ouput: fnew: new face list after removing the components smaller than maxdiameter
f=surfaceclean(f,v) remove surface patches that are located inside the bounding box faces input: v: surface node list, dimension (nn,3) f: surface face element list, dimension (be,3) output: f: faces free of those on the bounding box
eid=getintersecttri(tmppath) get the IDs of self-intersecting elements from tetgen call this when tetgen complains about self-intersection input: tmppath: working dir, use mwpath('') in most cases output: eid: an array of all intersecting surface elements, one can read the corresponding node/elem by [no,el]=readoff(mwpath('post_vmesh.off'));
elem=delendelem(elem,mask) delete elements whose nodes are all edge nodes input/output: elem: input/output, surface/volumetric element list mask: of length of node number, =0 for internal nodes, =1 for edge nodes
[node,elem]=meshresample(v,f,keepratio) resample mesh using CGAL mesh simplification utility input: v: list of nodes f: list of surface elements (each row for each triangle) keepratio: decimation rate, a number less than 1, as the percentage of the elements after the sampling output: node: the node coordinates of the sampled surface mesh elem: the element list of the sampled surface mesh
[newno,newfc]=remeshsurf(node,face,opt) remesh a triangular surface and the output is guaranteed to be free of self-intersecting element. This function is similar to meshresample, but it can both downsample or upsample a mesh input: node: list of nodes on the input suface mesh, 3 columns for x,y,z face: list of trianglular elements on the surface, [n1,n2,n3,region_id] opt: function parameters opt.gridsize: resolution for the voxelization of the mesh opt.closesize: if there are openings, set the closing diameter opt.elemsize: the size of the element of the output surface if opt is a scalar, it defines the elemsize and gridsize=opt/4 output: newno: list of nodes on the resulting suface mesh, 3 columns for x,y,z newfc: list of trianglular elements on the surface, [n1,n2,n3,region_id]
p=smoothsurf(node,mask,conn,iter,useralpha,usermethod,userbeta) smoothing a surface mesh input: node: node coordinates of a surface mesh mask: flag whether a node is movable: 0 movable, 1 non-movable if mask=[], it assumes all nodes are movable conn: input, a cell structure of length size(node), conn{n} contains a list of all neighboring node ID for node n, this can be computed from meshconn function iter: smoothing iteration number useralpha: scaler, smoothing parameter, v(k+1)=alpha*v(k)+(1-alpha)*mean(neighbors) usermethod: smoothing method, including 'laplacian','laplacianhc' and 'lowpass' userbeta: scaler, smoothing parameter, for 'laplacianhc' output: p: output, the smoothed node coordinates recommendations Based on [Bade2006], 'Lowpass' method outperforms 'Laplacian-HC' in volume preserving and both are significantly better than the standard Laplacian method [Bade2006] R. Bade, H. Haase, B. Preim, "Comparison of Fundamental Mesh Smoothing Algorithms for Medical Surface Models," Simulation and Visualization, pp. 289-304, 2006.
[no,el,fc]=sortmesh(origin,node,elem,face) sort nodes and elements in a mesh so that the indexed nodes and elements are closer to each order (this may reduce cache-miss in a calculation) input: origin: sorting all nodes and elements with the distance and angles wrt this location, if origin=[], it will be node(1,:) node: list of nodes elem: list of elements (each row are indices of nodes of each element) ecol: list of columns in elem to participate sorting face: list of surface triangles (this can be omitted) fcol: list of columns in face to participate sorting output: no: node coordinates in the sorted order el: the element list in the sorted order fc: the surface triangle list in the sorted order (can be ignored) nodemap: the new node mapping order, no=node(nodemap,:)
saveasc(v,f,fname) save a surface mesh to FreeSurfer ASC mesh format input: v: input, surface node list, dimension (nn,3) f: input, surface face element list, dimension (be,3) fname: output file name
savedxf(node,elem,face,fname) save a surface mesh to DXF format input: node: input, surface node list, dimension (nn,3) elem: input, tetrahedral element list, dimension (ne,4) face: input, surface face element list, dimension (be,3) fname: output file name
saveinr(vol,fname) save a surface mesh to INR Format input: vol: input, a binary volume fname: output file name
saveoff(v,f,fname) save a surface mesh to Geomview Object File Format (OFF) input: v: input, surface node list, dimension (nn,3) f: input, surface face element list, dimension (be,3) fname: output file name
savesmf(v,f,fname) save a surface mesh to smf format input: v: input, surface node list, dimension (nn,3) f: input, surface face element list, dimension (be,3) fname: output file name
savesurfpoly(v,f,holelist,regionlist,p0,p1,fname) save a set of surfaces into poly format (for tetgen) input: v: input, surface node list, dimension (nn,3) if v has 4 columns, the last column specifies mesh density near each node f: input, surface face element list, dimension (be,3) holelist: list of holes, each hole is represented by an internal point regionlist: list of regions, similar to holelist p0: coordinate of one of the end of the bounding box p1: coordinate for the other end of the bounding box fname: output file name forcebox: non-empty: add bounding box, []: automatic if forcebox is a 8x1 vector, it will be used to specify max-edge size near the bounding box corners
savevrml(node,elem,face,fname) save a surface mesh to VRML 1.0 format input: node: input, surface node list, dimension (nn,3) elem: input, tetrahedral element list, dimension (ne,4) face: input, surface face element list, dimension (be,3) fname: output file name
[node,elem]=readasc(fname) read FreeSurfer ASC mesh format input: fname: name of the asc file output: node: node positions of the mesh elem: element list of the mesh
vol=readinr(fname) load a volume from an INR file input: fname: input file name output: dat: output, data read from the inr file
[node,elem,face]=readmedit(filename) read Medit mesh format input: fname: name of the medit data file output: node: node coordinates of the mesh elem: list of elements of the mesh face: list of surface triangles of the mesh
[node,elem]=readoff(fname) read Geomview Object File Format (OFF) input: fname: name of the OFF data file output: node: node coordinates of the mesh elem: list of elements of the mesh
[node,elem]=readsmf(fname) read simple model format (SMF) input: fname: name of the SMF data file output: node: node coordinates of the mesh elem: list of elements of the mesh
[node,elem,face]=readtetgen(fstub) read tetgen output files input: fstub: file name stub output: node: node coordinates of the tetgen mesh elem: tetrahedra element list of the tetgen mesh face: surface triangles of the tetgen mesh
flag=deletemeshfile(fname) delete a given work mesh file under the working directory input: fname: specified file name (without path) output: flag: not used
binname=mcpath(fname) get full executable path by prepending a command directory path parameters: input: fname: input, a file name string output: binname: output, full file name located in the bin directory if global variable ISO2MESH_BIN is set in 'base', it will use [ISO2MESH_BIN filesep cmdname] as the command full path, otherwise, let matlab pass the cmdname to the shell, which will search command in the directories listed in system $PATH variable.
tempname=meshtemppath(fname) get full temp-file name by prepend working-directory and current session name input: fname: input, a file name string output: tempname: output, full file name located in the working directory if global variable ISO2MESH_TEMP is set in 'base', it will use it as the working directory; otherwise, will use matlab function tempdir to return a working directory. if global variable ISO2MESH_SESSION is set in 'base', it will be prepended for each file name, otherwise, use supplied file name.
islands=bwislands(img) return the indices of non-zero elements in a 2D or 3D image grouped by connected regions in a cell array input: img: a 2D or 3D array output: islands: a cell array, each cell records the indices of the non-zero elements in img for a connected region (or an island)
resimg=fillholes3d(img,ballsize) close a 3D image with the speicified gap size and then fill the holes input: img: a 3D binary image ballsize: maximum gap size for image closing output: resimg: the image free of holes this function requires the image processing toolbox for matlab/octave
cleanimg=deislands2d(img,sizelim) remove isolated islands on a 2D image below speicified size limit input: img: a 2D binary image sizelim: a integer as the maximum pixel size of a isolated region output: cleanimg: a binary image after removing islands below sizelim
cleanimg=deislands3d(img,sizelim) remove isolated islands for 3D image (for each slice) input: img: a 3D volumetric image sizelim: maximum island size (in pixels) for each x/y/z slice output: cleanimg: 3D image after removing the islands
imgdiff=imedge3d(binimg,isdiff) Extract the boundary voxels from a binary image author: Aslak Grinsted <ag at glaciology.net> input: binimg: a 3D binary image isdiff: if isdiff=1, output will be all voxels which is different from its neighbors; if isdiff=0 or ignored, output will be the edge voxels of the non-zero regions in binimg output: imgdiff: a 3D logical array with the same size as binimg with 1 for voxels on the boundary and 0 otherwise
p=internalpoint(v,aloop) imperical function to find an internal point of a planar polygon input: v: x,y,z coordinates of each node of the mesh aloop: input, a single vector separated by NaN, each segment is a close-polygon consisted by node IDs output: p: output, [x y z] of an internal point of aloop
vol=smoothbinvol(vol,layer) convolve a 3x3 gaussian kernel to a binary image multiple times input: vol: a 3D volumetric image to be smoothed layer: number of iterations for the smoothing output: vol: the volumetric image after smoothing
vol=thickenbinvol(vol,layer) thickening a binary volume by a given pixel width this is similar to bwmorph(vol,'thicken',3) except this does it in 3d and only does thickening for non-zero elements (and hopefully faster) input: vol: a volumetric binary image layer: number of iterations for the thickenining output: vol: the volume image after the thickening
vol=thickenbinvol(vol,layer) thining a binary volume by a given pixel width this is similar to bwmorph(vol,'thin',3) except this does it in 3d and only does thickening for non-zero elements (and hopefully faster) input: vol: a volumetric binary image layer: number of iterations for the thickenining output: vol: the volume image after the thickening
hm=plotmesh(node,face,elem,opt) plot surface and volumetric meshes input: node: a node coordinate list, 3 columns for x/y/z face: a triangular surface face list elem: a tetrahedral element list opt: additional options for the plotting for simple point plotting, opt can be markers or color options, such as 'r.', or opt can be a logic statement to select a subset of the mesh, such as 'x>0 & y+z<1'; opt can have more than one items to combine these options, for example: plotmesh(...,'x>0','r.'); the range selector must appear before the color/marker specifier output: hm: handle or handles (vector) to the plotted surfaces example: h=plotmesh(node,'r.'); h=plotmesh(node,'x<20','r.'); h=plotmesh(node,face); h=plotmesh(node,face,'y>10'); h=plotmesh(node,face,'facecolor','r'); h=plotmesh(node,elem,'x<20'); h=plotmesh(node,elem,'x<20 & y>0'); h=plotmesh(node,face,elem); h=plotmesh(node,face,elem,'linestyle','--');
hm=plotsurf(node,face,opt) plot 3D surface meshes input: node: node coordinates, dimension (nn,3) face: triangular surface face list opt: additional options for the plotting, see trisurf output: hm: handle or handles (vector) to the plotted surfaces example: h=plotsurf(node,face); h=plotsurf(node,face,'facecolor','r');
hm=plottetra(node,elem,opt) plot 3D surface meshes input: node: node coordinates, dimension (nn,3) elem: tetrahedral element list opt: additional options for a patch object output: hm: handle or handles (vector) to the plotted surfaces example: h=plottetra(node,elem); h=plottetra(node,elem,'facealpha',0.5);
[cutpos,cutvalue,facedata]=qmeshcut(elem,node,value,plane) fast tetrahedral mesh cross-section plot input: elem: integer array with dimensions of NE x 4, each row contains the indices of all the nodes for each tetrahedron node: node coordinates, 3 columns for x, y and z respectively value: a scalar array with the length of node numbers, can have multiple columns plane: defines a plane by 3 points: plane=[x1 y1 z1;x2 y2 z2;x3 y3 z3] output: cutpos: all the intersections of mesh edges by the plane cutpos is similar to node, containing 3 columns for x/y/z cutvalue: interpolated values at the intersections, with row number being the num. of the intersections, column number being the same as "value". facedata: define the intersection polygons in the form of patch "Faces" the outputs of this subroutine can be easily plotted using patch('Vertices',cutpos,'Faces',facedata,'FaceVertexCData',cutvalue,... 'FaceColor','interp');
valnew=surfdiffuse(node,tri,val,ddt,iter,type1,opt) apply a smoothing/diffusion process on a surface input: node: list of nodes of the surface mesh tri: triangular element list of the surface val: vector, scalar value for each node ddt: diffusion coefficient multiplied by delta t iter: iterations for applying the smoothing type1: indices of the nodes which will not be updated opt: method, 'grad' for gradient based, and 'simple' for simple average output: valnew: nodal value vector after the smoothing
[node,elem,face]=vol2mesh(img,ix,iy,iz,thickness,elemnum,maxvol,A,B) convert a binary volume to tetrahedral mesh input: img, ix,iy,iz, elemnum and maxvol: see vol2mesh.m Amat: a 3x3 transformation matrix Bvec: a 3x1 vector Amat and Bvec maps the image index space to real world coordnate system by [x,y,z]_new=Amat*[x,y,z]_old+Bvec thickness: scale z-dimension of the mesh to specified thickness, if thickness==0, scaling is bypassed
isoctave=isoctavemesh determine whether the code is running in octave output: isoctave: 1 if in octave, otherwise 0
p=getvarfrom(ws,name) get variable value by name from specified work-space input: ws: name of the work-space, for example, 'base' name: name string of the variable output: p: the value of the specified variable, if the variable does not exist, return empty array
[a,b,c,d]=getplanefrom3pt(plane) define a plane equation ax+by+cz+d=0 from three 3D points input: plane: a 3x3 matrix with each row specifying a 3D point (x,y,z) output: a,b,c,d: the coefficient for plane equation ax+by+cz+d=0
exesuff=getexeext() get meshing external tool extension names for the current platform output: exesuff: file extension for iso2mesh tool binaries
This toolbox interacts with a number external meshing tools to perform the essential functionalities. These tools are listed below:
Note: iso2mesh and the above meshing utilities are considered as an "aggregate" rather than "derived work", based on the definitions in GPL FAQ (http://www.gnu.org/licenses/gpl-faq.html#MereAggregation) Therefore, the license of iso2mesh and these utilities are independent. The iso2mesh license only applies to the scripts and documents/data in this package and exclude those programs stored in the bin/ directory. The source codes of the modified meshing utilities are available separatly at iso2mesh's website and retain their upstream licenses.
Your acknowledgement of iso2mesh in your publications or presentations would be greatly appreciated by the author of this toolbox. The citation information can be found in the Introduction section.