shapeworks::Mesh
Module: Mesh Classes
#include <Mesh.h>
Public Types
Name | |
---|---|
enum | FieldType |
enum | AlignmentType |
enum | DistanceMethod |
enum | CurvatureType |
enum | SubdivisionType |
Public Functions
Name | |
---|---|
Mesh(const std::string & pathname) | |
void | set_id(int id) |
int | get_id() const |
Mesh(vtkSmartPointer< vtkPolyData > meshPtr) | |
Mesh(const Mesh & orig) | |
Mesh(Mesh && orig) | |
Mesh & | operator=(const Mesh & orig) |
Mesh(const Eigen::MatrixXd & points, const Eigen::MatrixXi & faces) | |
Mesh & | operator=(Mesh && orig) |
Mesh & | operator+=(const Mesh & otherMesh) append two meshes |
vtkSmartPointer< vtkPolyData > | getVTKMesh() const return the current mesh |
Mesh & | write(const std::string & pathname, bool binaryFile =false) writes mesh, format specified by filename extension |
Mesh & | coverage(const Mesh & otherMesh, bool allowBackIntersections =true, double angleThreshold =0, double backSearchRadius =0) determines coverage between current mesh and another mesh (e.g. acetabular cup / femoral head) |
Mesh & | smooth(int iterations =0, double relaxation =0.0) applies laplacian smoothing |
Mesh & | smoothSinc(int iterations =0, double passband =0.0) applies vtk windowed sinc smoothing |
Mesh & | remesh(int numVertices, double adaptivity =1.0) applies remeshing using approximated centroidal voronoi diagrams for a given number of vertices and adaptivity |
Mesh & | remeshPercent(double percentage, double adaptivity =1.0) applies remeshing using approximated centroidal voronoi diagrams for a given percentage of vertices and adaptivity |
Mesh & | invertNormals() handle flipping normals |
Mesh & | reflect(const Axis & axis, const Vector3 & origin =makeVector({0.0, 0.0, 0.0})) reflect meshes with respect to a specified center and specific axis |
MeshTransform | createTransform(const Mesh & target, AlignmentType align =Similarity, unsigned iterations =10) |
Mesh & | applyTransform(const MeshTransform transform) applies the given transformation to the mesh |
Mesh & | rotate(const double angle, const Axis axis) applies the given rotation to the given axis |
Mesh & | fillHoles(double hole_size =1000.0) finds holes in a mesh and closes them |
Mesh & | clean() clean mesh |
Mesh & | probeVolume(const Image & image) samples image data values at point locations specified by image |
Mesh & | clip(const Plane plane) clips a mesh using a cutting plane |
Mesh & | translate(const Vector3 & v) helper to translate mesh |
Mesh & | scale(const Vector3 & v) helper to scale mesh |
PhysicalRegion | boundingBox() const computes bounding box of current mesh |
Mesh & | fixElement() fix element winding of mesh |
Mesh & | fixNonManifold() Attempt to fix non-manifold edges. |
bool | detectNonManifold() Detect if mesh contain non-manifold edges. |
bool | detectTriangular() Detect if mesh is triangular;. |
std::vector< Field > | distance(const Mesh & target, const DistanceMethod method =PointToCell) const |
Mesh & | clipClosedSurface(const Plane plane) clips a mesh using a cutting plane resulting in a closed surface |
Mesh & | computeNormals() computes and adds oriented point and cell normals |
Point3 | closestPoint(const Point3 point, double & distance, vtkIdType & face_id) const |
int | closestPointId(const Point3 point) const returns closest point id in this mesh to the given point in space |
bool | isPointInside(const Point3 point) const returns if the given point is inside the mesh |
double | geodesicDistance(int source, int target) const computes geodesic distance between two vertices (specified by their indices) on mesh |
Field | geodesicDistance(const Point3 landmark) const computes geodesic distance between a point (landmark) and each vertex on mesh |
Field | geodesicDistance(const std::vector< Point3 > curve) const computes geodesic distance between a set of points (curve) and each vertex on mesh |
Field | curvature(const CurvatureType type =Principal) const computes curvature using principal (default) or gaussian or mean algorithms |
void | computeFieldGradient(const std::string & field) const compute the gradient of a scalar field for all vertices |
Eigen::Vector3d | computeFieldGradientAtPoint(const std::string & field, const Point3 & query) const compute the gradient of a scalar field at a point |
double | interpolateFieldAtPoint(const std::string & field, const Point3 & query) const interpolate a scalar field at a given point |
Mesh & | applySubdivisionFilter(const SubdivisionType type =Butterfly, int subdivision =1) applies subdivision filter (butterfly (default) or loop) |
Image | toImage(PhysicalRegion region =PhysicalRegion(), Point3 spacing =Point3({1., 1., 1.})) const rasterizes specified region to create binary image of desired dims (default: unit spacing) |
Image | toDistanceTransform(PhysicalRegion region =PhysicalRegion(), const Point3 spacing =Point3({1., 1., 1.}), const Dims padding =Dims({1, 1, 1})) const converts specified region to distance transform image (default: unit spacing) with (logical) padding |
Mesh & | computeThickness(Image & image, Image * dt =nullptr, double max_dist =10000, double median_radius =5.0, std::string distance_mesh ="") assign cortical thickness values from mesh points |
Mesh & | computeLandmarkGeodesics(const std::vector< Point3 > & landmarks) compute geodesic distances to landmarks and assign as fields |
Point3 | center() const center of mesh |
Point3 | centerOfMass() const center of mass of mesh |
int | numPoints() const number of points |
int | numFaces() const number of faces |
Eigen::MatrixXd | points() const matrix with number of points with (x,y,z) coordinates of each point |
Eigen::MatrixXi | faces() const matrix with number of faces with indices of the three points from which each face is composed |
Point3 | getPoint(int id) const (x,y,z) coordinates of vertex at given index |
IPoint3 | getFace(int id) const return indices of the three points with which the face at the given index is composed |
std::vector< std::string > | getFieldNames() const print all field names in mesh |
Mesh & | setField(const std::string name, Array array, const FieldType type) sets the given field for points or faces with array (*does not copy array's values) |
Field | getField(const std::string & name, const FieldType type) const gets a pointer to the requested field of points or faces, null if field doesn't exist |
void | setFieldValue(const std::string & name, int idx, double value) sets the given index of field to value |
double | getFieldValue(const std::string & name, int idx) const gets the value at the given index of field (NOTE: returns first component of vector fields) |
Eigen::VectorXd | getMultiFieldValue(const std::string & name, int idx) const gets the multi value at the given index of [vertex] field |
bool | compareAllPoints(const Mesh & other_mesh) const compare if values of the points in two (corresponding) meshes are (eps)equal |
bool | compareAllFaces(const Mesh & other_mesh) const compare if face indices in two (corresponding) meshes are equal |
bool | compareAllFields(const Mesh & other_mesh, const double eps =-1.0) const compare if all fields in two meshes are (eps)equal |
bool | compareField(const Mesh & other_mesh, const std::string & name1, const std::string & name2 ="", const double eps =-1.0) const compare field of meshes to be (eps)equal (same field for both if only one specified) |
bool | compare(const Mesh & other_mesh, const double eps =-1.0) const compare meshes |
bool | operator==(const Mesh & other) const compare meshes |
double | getFFCValue(Eigen::Vector3d query) const Gets values for FFCs. |
Eigen::Vector3d | getFFCGradient(Eigen::Vector3d query) const Gets gradients for FFCs. |
vtkSmartPointer< vtkPoints > | getIGLMesh(Eigen::MatrixXd & V, Eigen::MatrixXi & F) const Formats mesh into an IGL format. |
vtkSmartPointer< vtkPolyData > | clipByField(const std::string & name, double value) Clips the mesh according to a field value. |
vtkSmartPointer< vtkStaticCellLocator > | getCellLocator() const Returns the cell locator. |
int | getClosestFace(const Point3 & point) const |
Eigen::Vector3d | computeBarycentricCoordinates(const Eigen::Vector3d & pt, int face) const Computes baricentric coordinates given a query point and a face number. |
std::vector< std::string > | getSupportedTypes() Return supported file types. |
Friends
Name | |
---|---|
struct | SharedCommandData |
Detailed Description
class shapeworks::Mesh;
This class encapsulates a Mesh and operations that can be performed on meshes
Public Types Documentation
enum FieldType
Enumerator | Value | Description |
---|---|---|
Point | ||
Face |
enum AlignmentType
Enumerator | Value | Description |
---|---|---|
Rigid | ||
Similarity | ||
Affine |
enum DistanceMethod
Enumerator | Value | Description |
---|---|---|
PointToPoint | ||
PointToCell |
enum CurvatureType
Enumerator | Value | Description |
---|---|---|
Principal | ||
Gaussian | ||
Mean |
enum SubdivisionType
Enumerator | Value | Description |
---|---|---|
Butterfly | ||
Loop |
Public Functions Documentation
function Mesh
Mesh(
const std::string & pathname
)
function set_id
inline void set_id(
int id
)
function get_id
inline int get_id() const
function Mesh
inline Mesh(
vtkSmartPointer< vtkPolyData > meshPtr
)
function Mesh
inline Mesh(
const Mesh & orig
)
function Mesh
inline Mesh(
Mesh && orig
)
function operator=
inline Mesh & operator=(
const Mesh & orig
)
function Mesh
Mesh(
const Eigen::MatrixXd & points,
const Eigen::MatrixXi & faces
)
function operator=
inline Mesh & operator=(
Mesh && orig
)
function operator+=
Mesh & operator+=(
const Mesh & otherMesh
)
append two meshes
function getVTKMesh
inline vtkSmartPointer< vtkPolyData > getVTKMesh() const
return the current mesh
function write
Mesh & write(
const std::string & pathname,
bool binaryFile =false
)
writes mesh, format specified by filename extension
function coverage
Mesh & coverage(
const Mesh & otherMesh,
bool allowBackIntersections =true,
double angleThreshold =0,
double backSearchRadius =0
)
determines coverage between current mesh and another mesh (e.g. acetabular cup / femoral head)
function smooth
Mesh & smooth(
int iterations =0,
double relaxation =0.0
)
applies laplacian smoothing
function smoothSinc
Mesh & smoothSinc(
int iterations =0,
double passband =0.0
)
applies vtk windowed sinc smoothing
function remesh
Mesh & remesh(
int numVertices,
double adaptivity =1.0
)
applies remeshing using approximated centroidal voronoi diagrams for a given number of vertices and adaptivity
function remeshPercent
Mesh & remeshPercent(
double percentage,
double adaptivity =1.0
)
applies remeshing using approximated centroidal voronoi diagrams for a given percentage of vertices and adaptivity
function invertNormals
Mesh & invertNormals()
handle flipping normals
function reflect
Mesh & reflect(
const Axis & axis,
const Vector3 & origin =makeVector({0.0, 0.0, 0.0})
)
reflect meshes with respect to a specified center and specific axis
function createTransform
MeshTransform createTransform(
const Mesh & target,
AlignmentType align =Similarity,
unsigned iterations =10
)
creates transform to target mesh using specified AlignmentType (Mesh::Rigid, Mesh::Similarity, Mesh::Affine) for specified number of iterations
function applyTransform
Mesh & applyTransform(
const MeshTransform transform
)
applies the given transformation to the mesh
function rotate
Mesh & rotate(
const double angle,
const Axis axis
)
applies the given rotation to the given axis
function fillHoles
Mesh & fillHoles(
double hole_size =1000.0
)
finds holes in a mesh and closes them
function clean
Mesh & clean()
clean mesh
function probeVolume
Mesh & probeVolume(
const Image & image
)
samples image data values at point locations specified by image
function clip
Mesh & clip(
const Plane plane
)
clips a mesh using a cutting plane
function translate
Mesh & translate(
const Vector3 & v
)
helper to translate mesh
function scale
Mesh & scale(
const Vector3 & v
)
helper to scale mesh
function boundingBox
PhysicalRegion boundingBox() const
computes bounding box of current mesh
function fixElement
Mesh & fixElement()
fix element winding of mesh
function fixNonManifold
Mesh & fixNonManifold()
Attempt to fix non-manifold edges.
function detectNonManifold
bool detectNonManifold()
Detect if mesh contain non-manifold edges.
function detectTriangular
bool detectTriangular()
Detect if mesh is triangular;.
function distance
std::vector< Field > distance(
const Mesh & target,
const DistanceMethod method =PointToCell
) const
Computes distance from each vertex to closest cell or point in target mesh, specified as PointToCell (default) or PointToPoint. Returns Fields containing distance to target and ids of the associated cells or points.
function clipClosedSurface
Mesh & clipClosedSurface(
const Plane plane
)
clips a mesh using a cutting plane resulting in a closed surface
function computeNormals
Mesh & computeNormals()
computes and adds oriented point and cell normals
function closestPoint
Point3 closestPoint(
const Point3 point,
double & distance,
vtkIdType & face_id
) const
Returns closest point on this mesh to the given point in space. In addition, returns by reference:
- the distance of the point in space from this mesh
- the face_id containing the closest point
function closestPointId
int closestPointId(
const Point3 point
) const
returns closest point id in this mesh to the given point in space
function isPointInside
bool isPointInside(
const Point3 point
) const
returns if the given point is inside the mesh
function geodesicDistance
double geodesicDistance(
int source,
int target
) const
computes geodesic distance between two vertices (specified by their indices) on mesh
function geodesicDistance
Field geodesicDistance(
const Point3 landmark
) const
computes geodesic distance between a point (landmark) and each vertex on mesh
function geodesicDistance
Field geodesicDistance(
const std::vector< Point3 > curve
) const
computes geodesic distance between a set of points (curve) and each vertex on mesh
function curvature
Field curvature(
const CurvatureType type =Principal
) const
computes curvature using principal (default) or gaussian or mean algorithms
function computeFieldGradient
void computeFieldGradient(
const std::string & field
) const
compute the gradient of a scalar field for all vertices
function computeFieldGradientAtPoint
Eigen::Vector3d computeFieldGradientAtPoint(
const std::string & field,
const Point3 & query
) const
compute the gradient of a scalar field at a point
function interpolateFieldAtPoint
double interpolateFieldAtPoint(
const std::string & field,
const Point3 & query
) const
interpolate a scalar field at a given point
function applySubdivisionFilter
Mesh & applySubdivisionFilter(
const SubdivisionType type =Butterfly,
int subdivision =1
)
applies subdivision filter (butterfly (default) or loop)
function toImage
Image toImage(
PhysicalRegion region =PhysicalRegion(),
Point3 spacing =Point3({1., 1., 1.})
) const
rasterizes specified region to create binary image of desired dims (default: unit spacing)
function toDistanceTransform
Image toDistanceTransform(
PhysicalRegion region =PhysicalRegion(),
const Point3 spacing =Point3({1., 1., 1.}),
const Dims padding =Dims({1, 1, 1})
) const
converts specified region to distance transform image (default: unit spacing) with (logical) padding
function computeThickness
Mesh & computeThickness(
Image & image,
Image * dt =nullptr,
double max_dist =10000,
double median_radius =5.0,
std::string distance_mesh =""
)
assign cortical thickness values from mesh points
function computeLandmarkGeodesics
Mesh & computeLandmarkGeodesics(
const std::vector< Point3 > & landmarks
)
compute geodesic distances to landmarks and assign as fields
function center
Point3 center() const
center of mesh
function centerOfMass
Point3 centerOfMass() const
center of mass of mesh
function numPoints
inline int numPoints() const
number of points
function numFaces
inline int numFaces() const
number of faces
function points
Eigen::MatrixXd points() const
matrix with number of points with (x,y,z) coordinates of each point
function faces
Eigen::MatrixXi faces() const
matrix with number of faces with indices of the three points from which each face is composed
function getPoint
Point3 getPoint(
int id
) const
(x,y,z) coordinates of vertex at given index
function getFace
IPoint3 getFace(
int id
) const
return indices of the three points with which the face at the given index is composed
function getFieldNames
std::vector< std::string > getFieldNames() const
print all field names in mesh
function setField
Mesh & setField(
const std::string name,
Array array,
const FieldType type
)
sets the given field for points or faces with array (*does not copy array's values)
function getField
Field getField(
const std::string & name,
const FieldType type
) const
gets a pointer to the requested field of points or faces, null if field doesn't exist
function setFieldValue
void setFieldValue(
const std::string & name,
int idx,
double value
)
sets the given index of field to value
function getFieldValue
double getFieldValue(
const std::string & name,
int idx
) const
gets the value at the given index of field (NOTE: returns first component of vector fields)
function getMultiFieldValue
Eigen::VectorXd getMultiFieldValue(
const std::string & name,
int idx
) const
gets the multi value at the given index of [vertex] field
function compareAllPoints
bool compareAllPoints(
const Mesh & other_mesh
) const
compare if values of the points in two (corresponding) meshes are (eps)equal
function compareAllFaces
bool compareAllFaces(
const Mesh & other_mesh
) const
compare if face indices in two (corresponding) meshes are equal
function compareAllFields
bool compareAllFields(
const Mesh & other_mesh,
const double eps =-1.0
) const
compare if all fields in two meshes are (eps)equal
function compareField
bool compareField(
const Mesh & other_mesh,
const std::string & name1,
const std::string & name2 ="",
const double eps =-1.0
) const
compare field of meshes to be (eps)equal (same field for both if only one specified)
function compare
bool compare(
const Mesh & other_mesh,
const double eps =-1.0
) const
compare meshes
function operator==
inline bool operator==(
const Mesh & other
) const
compare meshes
function getFFCValue
double getFFCValue(
Eigen::Vector3d query
) const
Gets values for FFCs.
function getFFCGradient
Eigen::Vector3d getFFCGradient(
Eigen::Vector3d query
) const
Gets gradients for FFCs.
function getIGLMesh
vtkSmartPointer< vtkPoints > getIGLMesh(
Eigen::MatrixXd & V,
Eigen::MatrixXi & F
) const
Formats mesh into an IGL format.
function clipByField
vtkSmartPointer< vtkPolyData > clipByField(
const std::string & name,
double value
)
Clips the mesh according to a field value.
function getCellLocator
inline vtkSmartPointer< vtkStaticCellLocator > getCellLocator() const
Returns the cell locator.
function getClosestFace
int getClosestFace(
const Point3 & point
) const
function computeBarycentricCoordinates
Eigen::Vector3d computeBarycentricCoordinates(
const Eigen::Vector3d & pt,
int face
) const
Computes baricentric coordinates given a query point and a face number.
function getSupportedTypes
static inline std::vector< std::string > getSupportedTypes()
Return supported file types.
Friends
friend SharedCommandData
friend struct SharedCommandData(
SharedCommandData
);
Updated on 2024-03-17 at 12:58:44 -0600