public class DistanceGrid extends ScalarGridBase
The field is implemented using a regular 3D grid composed of
nx
X ny
X nz
vertices along the x, y
and z directions, giving the grid a resolution of rx
,
ry
, rz
, cells along each of these directions,
where rx = nx-1
, ry = ny-1
, and rz =
nz-1
. The grid has widths w = (wx, wy, wz)
along each of
these directions, and cell widths given by wx/rx
,
wy/ry
, and wz/rz
.
Several coordinate frames are associated with the grid: the local
frame L, the world frame W, and the grid frame G. Details on
these are given in the documentation for InterpolatingGridBase
.
Support is provided to compute the distance field from sets of underlying mesh features (such as faces or vertices). These computations are done in the local frame.
Distances for any point within the grid can be obtained by interpolation
of the distances at surrounding vertices. Interpolation can be either
linear or quadratic. Linearly interpolated distances are obtained using
getLocalDistance(maspack.matrix.Point3d)
and are computed using trilinear interpolation
within the grid cell containing the point. If the query point is outside the
grid, then the distance is assigned the special value ScalarGridBase.OUTSIDE_GRID
. Normals can also be obtained, using the
getLocalDistanceAndNormal
methods. Normals are also
interpolated using trilinear interpolation, with the normals at each vertex
being computed on demand by mid-point differencing with adjacent
vertices. Note that this means that the normals do not correspond to
the gradient of the distance function; however, their interpolated values
are smoother than those of the gradient. Gradient values can be obtained
instead, using the getLocalDistanceAndGradient
methods. The
returned gradient is the true gradient of the distance function within the
interpolating cell and a linear function within that cell. All inputs and
outputs for the above described distance methods are assumed to be in local
coordinates. The method getWorldDistance(maspack.matrix.Point3d)
and associated
getWorldDistanceAndNormal
and
getWorldDistanceAndGradient
methods perform the same operations
in world coordinates.
Quadratically interpolated distances are obtained using getQuadDistance(maspack.matrix.Point3d)
. Quadratic interpolation is done within a composite
quadratic cell composed of 2 x 2 x 2 regular cells. To ensure that
all points within the grid can be assigned a unique quadratic cell, the grid
resolution is restricted so that rx
, ry
and
rz
are always even. Moreover, within each quadratic cell,
interpolation is done by partitioning the cell into six tetrahedra, each
with 10 nodes, and then performing quadratic interpolation within the
appropriate tetrahedron. This is done to ensure that the interpolation
function is exactly quadratic instead of tri-quadratic. Because the
interpolation is exactly quadratic, certain operations, such as intersecting
level set surfaces, can be performed by solving a simple quadratic equation.
The gradient of the quadratically interpolated distance function can also be
obtained using getQuadDistanceAndGradient(maspack.matrix.Vector3d, maspack.matrix.Matrix3d, maspack.matrix.Point3d)
; this will return
the gradient of the distance function within the quadratic cell containing
the query point. These gradient values are smoother than those associated
with linear interpolation and so (if normalized) can be used to provide
normal vectors. Consequently, quadratically interpolated normals, analagous
to those returned by getLocalDistanceAndNormal(maspack.matrix.Vector3d, double, double, double)
, are not currently
supported. All inputs and outputs for the above described quadratic distance
methods are assumed to be in local coordinates. The methods
getWorldQuadDistance(maspack.matrix.Point3d)
and
getWorldQuadDistanceAndGradient(maspack.matrix.Vector3d, maspack.matrix.Matrix3d, maspack.matrix.Point3d)
perform the same operations
in world coordinates.
Distances at the grid vertices can either be assigned directly, using
setVertexDistances(double[], boolean)
, or can be computed based on nearest distance
calculations to a set of point-based geometric features such as vertices,
edges, or faces, using computeDistances(java.util.List<? extends maspack.geometry.Feature>, boolean)
. The convenience methods
computeFromFeatures(java.util.List<? extends maspack.geometry.Feature>, double, maspack.matrix.RigidTransform3d, int, boolean)
and computeFromFeaturesOBB(java.util.List<? extends maspack.geometry.Feature>, double, int, boolean)
will both
compute distances and also fit the grid to the features. When computed from
features, each grid vertex is associated with a nearest
feature.
Feature based distance computation assumes that the feature points are
represented in local coordinates. For grid vertices close to the mesh,
nearby features are examined to determine the distance from the mesh to the
vertex. A sweep method is then used to propagate distance values throughout
the volume and determine whether vertices are inside or outside. If the
features are faces, it is possible to request that the grid be
signed. In this case, ray casts are used to determine whether
vertices are inside or outside the feature set, with interior and exterior
vertices being given negative and positive distance values,
respectively. The algorithm is based on C++ code provided by Robert Bridson
(UBC) in the makelevelset3
class of his
common code set.
Modifier and Type | Class and Description |
---|---|
static class |
DistanceGrid.DistanceMethod
Describes what method to use when generating distance values from
polygonal meshes.
|
class |
DistanceGrid.GridEdge |
static class |
DistanceGrid.TetDesc |
static class |
DistanceGrid.TetID
Identifies a sub tet within a hex cell.
|
Modifier and Type | Field and Description |
---|---|
static int |
culledQueries |
static DistanceGrid.DistanceMethod |
DEFAULT_DISTANCE_METHOD |
static int |
totalQueries |
boolean |
useNormalsForQuadGrad |
OUTSIDE_GRID
TRANSPARENT, TWO_DIMENSIONAL
Constructor and Description |
---|
DistanceGrid()
Default constructor.
|
DistanceGrid(DistanceGrid grid)
Creates a new distance grid that is a copy of an existing grid.
|
DistanceGrid(java.util.List<? extends Feature> features,
double marginFrac,
int maxRes,
boolean signed)
Creates a new distance grid for a specified list of features,
axis-aligned and centered on the features, with a uniform cell width
along all axes defined so that the resolution along the maximum width
axis is
maxRes . |
DistanceGrid(java.util.List<? extends Feature> features,
double marginFrac,
Vector3i resolution,
boolean signed)
Creates a new distance grid for a specified list of features,
axis-aligned and centered on the features, with a the cell resolution
along each axis given by
resolution . |
DistanceGrid(Vector3d widths,
Vector3i resolution,
RigidTransform3d TCL)
Creates a new distance grid with specified widths, resolution, and
position and orientation of the center given by
TCL . |
DistanceGrid(Vector3i resolution)
Creates a new distance grid, axis-aligned and centered on the origin,
with the specified resolution and x, y, z widths set to 1.
|
Modifier and Type | Method and Description |
---|---|
void |
checkGradient() |
void |
clearFeatures()
Clears the features, if any, associated with this distance
grid.
|
void |
computeDifference01(PolygonalMesh mesh) |
void |
computeDifference10(PolygonalMesh mesh) |
void |
computeDistances(java.util.List<? extends Feature> features,
boolean signed)
Computes the distance field for this grid, based on the nearest distances
to a supplied set of features.
|
void |
computeDistances(PolygonalMesh mesh,
boolean signed) |
void |
computeFromFeatures(java.util.List<? extends Feature> features,
double marginFrac,
RigidTransform3d TCL,
int maxRes,
boolean signed)
Fits this grid to a set of features and computes the corresponding
distance field.
|
void |
computeFromFeaturesOBB(java.util.List<? extends Feature> features,
double marginFrac,
int maxRes,
boolean signed)
Fits this grid in an oriented manner to a set of features and computes
the corresponding distance field.
|
void |
computeFromMesh(PolygonalMesh mesh,
double marginFrac,
RigidTransform3d TCL,
int maxRes,
boolean signed) |
void |
computeIntersection(PolygonalMesh mesh) |
void |
computeQuadCoefs(double[] a,
DistanceGrid.TetDesc tdesc) |
void |
computeQuadCoefsStub(double[] a,
DistanceGrid.TetDesc tdesc) |
void |
computeUnion(PolygonalMesh mesh) |
PolygonalMesh |
createQuadDistanceSurface(double val,
int res)
Creates a triangular mesh approximating the surface on which the
quadratically interpolated distance function equals
val . |
boolean |
epsilonEquals(DistanceGrid grid,
double tol)
Returns
true if this grid equals another within a prescribed
tolerance. |
boolean |
findQuadSurfaceIntersection(Point3d pi,
Point3d p0,
Point3d pa,
Vector3d nrm)
Experimental method.
|
boolean |
findQuadSurfaceTangent(Point3d pt,
Point3d p0,
Point3d pa,
Vector3d nrm)
Find a point
pt such that the line segment
pa-pt is tangent to the quadratic zero distance surface and
also lies in the plane defined by point p0 and normal
nrm . |
java.util.ArrayList<DistanceGrid.TetDesc> |
findSurfaceEdgeTets(java.util.ArrayList<DistanceGrid.GridEdge> edges) |
java.util.ArrayList<DistanceGrid.GridEdge> |
findSurfaceIntersectingEdges() |
void |
fitToFeatures(java.util.Collection<java.util.List<? extends Feature>> featureSets,
double marginFrac,
RigidTransform3d TCL,
int maxRes)
Fits the widths and center of this grid to a several sets of
features.
|
void |
fitToFeatures(java.util.List<? extends Feature> features,
double marginFrac,
RigidTransform3d TCL,
int maxRes)
Fits the widths and center of this grid to a set of features.
|
void |
fitToFeaturesOBB(java.util.Collection<java.util.List<? extends Feature>> featureSets,
double marginFrac,
int maxRes)
Fits the widths, center and local orientation of this grid to several
sets of features.
|
void |
fitToFeaturesOBB(java.util.List<? extends Feature> features,
double marginFrac,
int maxRes)
Fits the widths, center and local orientation of this grid to a set of
features.
|
Feature |
getClosestFeature(int idx)
Returns the closest Feature to the vertex indexed by
idx . |
DistanceGrid.DistanceMethod |
getDistanceMethod()
Returns the method used to compute distance values from polygonal meshes.
|
Feature[] |
getFeatures()
Returns the features, if any, associated with this distance
grid Features will be associated with the field if they were used to
compute it, via
computeFromFeatures(java.util.List<? extends maspack.geometry.Feature>, double, maspack.matrix.RigidTransform3d, int, boolean) , computeFromFeaturesOBB(java.util.List<? extends maspack.geometry.Feature>, double, int, boolean) , or one of the associated constructors. |
double |
getLocalDistance(Point3d point)
Calculates the distance at an arbitrary point in local coordinates using
multilinear interpolation of the vertex values for the grid cell
containing the point.
|
double |
getLocalDistanceAndGradient(Vector3d grad,
double x,
double y,
double z)
Calculates the distance and gradient at an arbitrary point in
local coordinates using multilinear interpolation, as described
for
getLocalDistanceAndGradient(Vector3d,Point3d) . |
double |
getLocalDistanceAndGradient(Vector3d grad,
Point3d point)
Calculates the distance and gradient at an arbitrary point in local
coordinates using multilinear interpolation of the vertex values for the
grid cell containing the point.
|
double |
getLocalDistanceAndNormal(Vector3d norm,
double x,
double y,
double z)
Calculates the distance and normal at an arbitrary point in local
coordinates using multilinear interpolation, as described
for
getLocalDistanceAndNormal(Vector3d,Point3d) . |
double |
getLocalDistanceAndNormal(Vector3d norm,
Matrix3d Dnrm,
Point3d point)
Calculates the distance and normal at an arbitrary point in local
coordinates using multilinear interpolation, as described for
getLocalDistanceAndNormal(Vector3d,Point3d) . |
double |
getLocalDistanceAndNormal(Vector3d norm,
Point3d point)
Calculates the distance and normal at an arbitrary point in local
coordinates using multilinear interpolation of the vertex values for the
grid cell containing the point.
|
Feature |
getNearestLocalFeature(Point3d nearest,
Point3d point)
Determines nearest feature to an arbitray point in local coordinates.
|
Feature |
getNearestWorldFeature(Point3d nearest,
Point3d point)
Determines nearest feature to an arbitray point in world coordinates.
|
java.lang.String |
getQuadCellInfo(Point3d point) |
double |
getQuadDistance(Point3d point)
Calculates the distance at an arbitrary point in local coordinates using
quadratic interpolation, as described for
getQuadDistanceAndGradient(maspack.matrix.Vector3d, maspack.matrix.Matrix3d, maspack.matrix.Point3d) . |
double |
getQuadDistanceAndGradient(Vector3d grad,
Matrix3d dgrad,
Point3d point)
Calculates the distance and gradient at an arbitrary point in local
coordinates using quadratic interpolation.
|
DistanceGrid.TetDesc |
getQuadTet(Point3d point)
Used for debugging
|
Vector3d |
getRenderVector(int xi,
int yj,
int zk) |
double |
getRenderVectorScale()
For grid subclasses wishing to render 3D vector information, this method
should be overridden to return a scaling factor to be applied to vectors
returned by
InterpolatingGridBase.getRenderVector(maspack.matrix.Vector3d, int, int, int) to determine the length of the
actual rendered vector. |
double |
getVertexDistance(Vector3i vxyz)
Queries the distance value at a specified vertex, as specified by x, y, z
indices.
|
double[] |
getVertexDistances()
Returns the full array of distances at each vertex.
|
double |
getWorldDistance(Point3d point)
Calculates the distance at an arbitrary point in world coordinates using
multilinear interpolation of the vertex values for the grid cell
containing the point.
|
double |
getWorldDistanceAndGradient(Vector3d grad,
double x,
double y,
double z)
Calculates the distance and gradient at an arbitrary point in
world coordinates using multilinear interpolation, as described
for
getLocalDistanceAndGradient(Vector3d,Point3d) . |
double |
getWorldDistanceAndGradient(Vector3d grad,
Point3d point)
Calculates the distance and gradient at an arbitrary point.
|
double |
getWorldDistanceAndNormal(Vector3d norm,
double x,
double y,
double z)
Calculates the distance and normal at an arbitrary point in world
coordinates using multilinear interpolation, as described
for
getLocalDistanceAndNormal(Vector3d,Point3d) . |
double |
getWorldDistanceAndNormal(Vector3d norm,
Point3d point)
Calculates the distance and normal at an arbitrary point in world
coordinates using multilinear interpolation, as described for
getLocalDistanceAndNormal(Vector3d,Point3d) . |
double |
getWorldQuadDistance(Point3d point)
Calculates the distance at an arbitrary point in world coordinates using
quadratic interpolation, as described for
getQuadDistanceAndGradient(maspack.matrix.Vector3d, maspack.matrix.Matrix3d, maspack.matrix.Point3d) . |
double |
getWorldQuadDistanceAndGradient(Vector3d grad,
Matrix3d dgrad,
Point3d point)
Calculates the distance and gradient at an arbitrary point in
world coordinates using quadratic interpolation, as described
for
getLocalDistanceAndGradient(Vector3d,Point3d) . |
boolean |
isSigned()
Queries whether or not this grid is signed.
|
void |
markOutsideQuadtets(double mindist)
Marks all tets whose nodes all have distance values
> mindist as
being "outside". |
void |
scaleDistance(double s)
Scales the distance units of the grid, which entails scaling the widths,
distances and coordinate transform information.
|
void |
set(DistanceGrid grid)
Sets this distance grid to be a copy of an existing grid.
|
void |
setDistanceMethod(DistanceGrid.DistanceMethod method)
Sets the method used to compute distance values from polygonal meshes.
|
void |
setDistancesAndFeatures(double[] distances,
java.util.List<? extends Feature> features,
int[] closestFeatures,
boolean signed)
Explicitly sets the distance field and features for this grid.
|
void |
setLocalToWorld(RigidTransform3d TLW)
Sets a transform that maps from local to world coordinates.
|
void |
setResolution(Vector3i resolution)
Sets the resolution for this grid along the x, y, and z axes.
|
void |
setVertexDistances(double[] distances,
boolean signed)
Explicitly sets the distance field for this grid, by setting the distance
values at each vertex.
|
void |
smooth()
Applies a simple Laplacian smoothing operation to the distance grid
|
void |
smooth(double lambda,
double mu)
Taubin Smoothing
|
void |
smooth(double lambda,
double mu,
int iters)
Applies Taubin smoothing
|
void |
smooth(int n)
Apply the smoothing operation n times
|
void |
zeroVertexDistances()
Explicitly zeros the distance field for this grid.
|
createDistanceSurface, createDistanceSurface, createDistanceSurface, epsilonEquals, getParameterType, setVertexValue
clearColors, createRenderProps, epsilonEquals, getCellVertex, getCellWidths, getCenter, getCenterAndOrientation, getClosestVertex, getDebug, getDefaultVertexColor, getGridToLocalTransformer, getLocalToWorld, getLocalToWorldTransformer, getLocalVertexCoords, getLocalVertexCoords, getOrientation, getRadius, getRenderHints, getRenderProps, getRenderRanges, getRenderVector, getResolution, getSelection, getVertexColor, getWidths, getWidths, getWorldCenter, getWorldOrientation, getWorldVertexCoords, getWorldVertexCoords, hasLocalToWorld, isSelectable, isVertexRenderingEnabled, isWritable, numSelectionQueriesNeeded, numVertices, parseRenderRanges, prerender, prerender, render, render, scan, set, setCenter, setCenterAndOrientation, setDebug, setDefaultVertexColor, setOrientation, setRenderProps, setRenderRanges, setVertexColor, setVertexRenderingEnabled, updateBounds, vertexToXyzIndices, write, xyzIndicesToVertex
public static DistanceGrid.DistanceMethod DEFAULT_DISTANCE_METHOD
public boolean useNormalsForQuadGrad
public static int totalQueries
public static int culledQueries
public DistanceGrid()
InterpolatingGridBase.scan(maspack.util.ReaderTokenizer, java.lang.Object)
is called immediately after.public DistanceGrid(java.util.List<? extends Feature> features, double marginFrac, int maxRes, boolean signed)
maxRes
. The grid is
created by calling
computeFromFeatures (features, marginFrac, null, maxRes, signed);
features
- features used to compute the distance fieldmarginFrac
- specifies the fractional amount that the
grid should be grown in each direction to better contain the featuresmaxRes
- specfies the resolution along the longest
width, with resolutions along other widths set to ensure uniform
cell size. Must be >
0.signed
- if true
, indicates that the field should be
signed. At present, signed fields can only be computed if all features
are faces (i.e., Face
).public DistanceGrid(java.util.List<? extends Feature> features, double marginFrac, Vector3i resolution, boolean signed)
resolution
. The grid is
created by calling
setResolution (resolution); computeFromFeatures (features, marginFrac, null, 0, signed);
features
- features used to compute the distance fieldmarginFrac
- specifies the fractional amount that the
grid should be grown in each direction to better contain the featuresresolution
- specfies the resolution along each of the
x, y, and z axes.signed
- if true
, indicates that the field should be
signed. At present, signed fields can only be computed if all features
are faces (i.e., Face
).public DistanceGrid(Vector3d widths, Vector3i resolution, RigidTransform3d TCL)
TCL
.
The grid distances are initialized to zero.widths
- widths along the x, y, and z axesresolution
- cell resolution along the x, y, and z axesTCL
- transformation giving the position and orientation
of the grid centerpublic DistanceGrid(Vector3i resolution)
resolution
- cell resolution along the x, y, and z axespublic DistanceGrid(DistanceGrid grid)
grid
- distance grid to copypublic void setDistanceMethod(DistanceGrid.DistanceMethod method)
null
sets the method to the default value,
which is DistanceGrid.DistanceMethod.BRIDSON
.method
- distance method for polygonal meshespublic DistanceGrid.DistanceMethod getDistanceMethod()
public void setLocalToWorld(RigidTransform3d TLW)
hasLocalToWorld
will subsequently return false
.setLocalToWorld
in class InterpolatingGridBase
TLW
- transform from local to world coordinatespublic void setResolution(Vector3i resolution)
resolution
argument is adjusted
accordingly). If features are present for this grid, as described for
getFeatures()
, then these are used to recompute the
distances. Otherwise, the distances are set to zero. If a color map is
present, it is cleared.setResolution
in class InterpolatingGridBase
resolution
- cell resolution along the x, y and z axes. Must be at
least 1 along each axis.public void set(DistanceGrid grid)
grid
- distance grid to copypublic void computeDistances(java.util.List<? extends Feature> features, boolean signed)
getFeatures()
. The grid resolution,
widths, center and orientation remain unchanged.features
- features used to compute the distance fieldsigned
- if true
, indicates that the field should be
signed. At present, signed fields can only be computed if all features
are faces (i.e., Face
).public void computeDistances(PolygonalMesh mesh, boolean signed)
public void setVertexDistances(double[] distances, boolean signed)
The input array is indexed such that for vertex indices xi, yj, zk, the corresponding index into this array is
idx = xi + nx*yj + (nx*ny)*zkwhere
nx
and ny
are the number
of vertices along x and y axes.distances
- distance for each vertex. Must have a length
>=
InterpolatingGridBase.numVertices()
.signed
- if true
, indicates that the field should be
considered signed.public double[] getVertexDistances()
setVertexDistances(double[], boolean)
for a description of how vertices are indexed with
respect to this array.null
if distances have
not yet been set.public double getVertexDistance(Vector3i vxyz)
vxyz
- x, y, z vertex indicespublic void zeroVertexDistances()
false
.public void setDistancesAndFeatures(double[] distances, java.util.List<? extends Feature> features, int[] closestFeatures, boolean signed)
distances
- distance for each vertex. Must have a length
>=
InterpolatingGridBase.numVertices()
.features
- list of features to be associated with this
grid.closestFeatures
- index (with respect to features
(
of the nearest feature to each vertex. Must have a length
>=
InterpolatingGridBase.numVertices()
.signed
- if true
, indicates that the field should be
considered signed.public void computeFromFeatures(java.util.List<? extends Feature> features, double marginFrac, RigidTransform3d TCL, int maxRes, boolean signed)
getFeatures()
. The existing grid resolution is
unchanged, unless maxRes
exceeds 0, in which
case the resolution along the longest width is set to maxRes
and the resolution along all other axes are set so as to ensure
a uniform cell size. All resolutions are rounded up to an even number.
The points of the features are assumed to be given in local
coordinates. The application has the option of specifying the
position and orientation of the grid center using the optional
argument TCL
. If this is non-null, then the widths
are set large enough about the specified center to accomodate
the features. Otherwise, if TCL
is null, then
the grid is fit in an axis-aligned manner to the local coordinate
frame (i.e., InterpolatingGridBase.getOrientation(maspack.matrix.RotationMatrix3d)
returns the identity), with
a center and widths computed to best fit the features.
To ensure that the grid is not degenerate, all widths are adjusted if necessary to ensure that they are at least 0.05 of the maximum width. Widths are then further adjusted by multiplying by
(1 + 2*marginFrac)to provide a margin around the features. Finally, if
maxRes > 0
,
widths may be grown to ensure uniform
cell size.features
- features used to compute the distance fieldmarginFrac
- specifies the fractional amount that the
grid should be grown in each direction to better contain the featuresTCL
- optional - if non-null, specifies the pose of the grid centermaxRes
- if > 0
,
specifies the resolution along the longest
width, with resolutions along other widths set to ensure uniform
cell sizesigned
- if true
, indicates that the field should be
signed. At present, signed fields can only be computed if all features
are faces (i.e., Face
).public void computeFromMesh(PolygonalMesh mesh, double marginFrac, RigidTransform3d TCL, int maxRes, boolean signed)
public void fitToFeatures(java.util.List<? extends Feature> features, double marginFrac, RigidTransform3d TCL, int maxRes)
computeFromFeatures(java.util.List<? extends maspack.geometry.Feature>, double, maspack.matrix.RigidTransform3d, int, boolean)
;
however, no distance field is computed and the features are not stored.features
- features used to fit the gridmarginFrac
- specifies the fractional amount that the
grid should be grown in each direction to better contain the featuresTCL
- optional - if non-null, specifies the pose of the grid centermaxRes
- if > 0
,
specifies the resolution along the longest
width, with resolutions along other widths set to ensure uniform
cell sizepublic void fitToFeatures(java.util.Collection<java.util.List<? extends Feature>> featureSets, double marginFrac, RigidTransform3d TCL, int maxRes)
computeFromFeatures(java.util.List<? extends maspack.geometry.Feature>, double, maspack.matrix.RigidTransform3d, int, boolean)
; however, no distance field is computed and the
features are not stored.featureSets
- lists of features used to fit the gridmarginFrac
- specifies the fractional amount that the
grid should be grown in each direction to better contain the featuresTCL
- optional - if non-null, specifies the pose of the grid centermaxRes
- if > 0
,
specifies the resolution along the longest
width, with resolutions along other widths set to ensure uniform
cell sizepublic void computeFromFeaturesOBB(java.util.List<? extends Feature> features, double marginFrac, int maxRes, boolean signed)
getFeatures()
. The existing grid
resolution is unchanged, unless maxRes
exceeds 0, in
which case the resolution along the longest width is set to
maxRes
and the resolution along all other axes are set so as
to ensure a uniform cell size. All resolutions are rounded up to an even
number.
The points of the features are assumed to be given in local
coordinates. The grid is fit to an oriented bounding box (OBB) with
respect to this frame; InterpolatingGridBase.getCenter(maspack.matrix.Vector3d)
and InterpolatingGridBase.getOrientation(maspack.matrix.RotationMatrix3d)
will
return the center and orientation of this box. To ensure that the grid is
not degenerate, all widths are adjusted if necessary to ensure that they
are at least 0.05 of the maximum width. Widths are then further adjusted
by multiplying by
(1 + 2*marginFrac)to provide a margin around the features.
features
- features used to compute the distance fieldmarginFrac
- specifies the fractional amount that the
grid should be grown in each direction to better contain the featuresmaxRes
- if >
0, specfies the resolution along the longest
width, with resolutions along other widths set to ensure uniform
cell sizesigned
- if true
, indicates that the field should be
signed. At present, signed fields can only be computed if all features
are faces (i.e., Face
).public void fitToFeaturesOBB(java.util.List<? extends Feature> features, double marginFrac, int maxRes)
computeFromFeaturesOBB(java.util.List<? extends maspack.geometry.Feature>, double, int, boolean)
; however, no distance field is computed and the
features are not stored.features
- features used to fit the gridmarginFrac
- specifies the fractional amount that the
grid should be grown in each direction to better contain the featuresmaxRes
- if >
0, specfies the resolution along the longest
width, with resolutions along other widths set to ensure uniform
cell sizepublic void fitToFeaturesOBB(java.util.Collection<java.util.List<? extends Feature>> featureSets, double marginFrac, int maxRes)
computeFromFeaturesOBB(java.util.List<? extends maspack.geometry.Feature>, double, int, boolean)
; however, no distance field is computed and the
features are not stored.featureSets
- lists of features used to fit the gridmarginFrac
- specifies the fractional amount that the
grid should be grown in each direction to better contain the featuresmaxRes
- if >
0, specfies the resolution along the longest
width, with resolutions along other widths set to ensure uniform
cell sizepublic void computeUnion(PolygonalMesh mesh)
public void computeIntersection(PolygonalMesh mesh)
public void computeDifference01(PolygonalMesh mesh)
public void computeDifference10(PolygonalMesh mesh)
public double getLocalDistance(Point3d point)
ScalarGridBase.OUTSIDE_GRID
is returned.point
- point at which to calculate the normal and distance
(local coordinates).OUTSIDE_GRID
.public double getWorldDistance(Point3d point)
ScalarGridBase.OUTSIDE_GRID
is returned.point
- point at which to calculate the normal and distance
(world coordinates).OUTSIDE_GRID
.public double getLocalDistanceAndNormal(Vector3d norm, double x, double y, double z)
getLocalDistanceAndNormal(Vector3d,Point3d)
.
If the point lies outside
the grid volume, ScalarGridBase.OUTSIDE_GRID
is returned.norm
- returns the normal (local coordinates)x
- x coordinate of point (local coordinates).y
- y coordinate of point (local coordinates).z
- z coordinate of point (local coordinates).OUTSIDE_GRID
.public double getWorldDistanceAndNormal(Vector3d norm, double x, double y, double z)
getLocalDistanceAndNormal(Vector3d,Point3d)
.
If the point lies outside
the grid volume, ScalarGridBase.OUTSIDE_GRID
is returned.norm
- returns the normal (world coordinates)x
- x coordinate of point (world coordinates).y
- y coordinate of point (world coordinates).z
- z coordinate of point (world coordinates).OUTSIDE_GRID
.public double getLocalDistanceAndGradient(Vector3d grad, double x, double y, double z)
getLocalDistanceAndGradient(Vector3d,Point3d)
.
If the point lies outside
the grid volume, ScalarGridBase.OUTSIDE_GRID
is returned.grad
- returns the gradient (local coordinates)x
- x coordinate of point (local coordinates).y
- y coordinate of point (local coordinates).z
- z coordinate of point (local coordinates).OUTSIDE_GRID
.public double getWorldDistanceAndGradient(Vector3d grad, double x, double y, double z)
getLocalDistanceAndGradient(Vector3d,Point3d)
.
If the point lies outside
the grid volume, ScalarGridBase.OUTSIDE_GRID
is returned.grad
- returns the gradient (world coordinates)x
- x coordinate of point (world coordinates).y
- y coordinate of point (world coordinates).z
- z coordinate of point (world coordinates).OUTSIDE_GRID
.public double getLocalDistanceAndNormal(Vector3d norm, Point3d point)
ScalarGridBase.OUTSIDE_GRID
is
returned.norm
- returns the normal (local coordinates)point
- point at which to calculate the normal and distance
(local coordinates).OUTSIDE_GRID
.public double getWorldDistanceAndNormal(Vector3d norm, Point3d point)
getLocalDistanceAndNormal(Vector3d,Point3d)
.
If the point lies outside
the grid volume, ScalarGridBase.OUTSIDE_GRID
is returned.norm
- returns the normal (world coordinates)point
- point at which to calculate the normal and distance
(world coordinates).OUTSIDE_GRID
.public double getLocalDistanceAndNormal(Vector3d norm, Matrix3d Dnrm, Point3d point)
getLocalDistanceAndNormal(Vector3d,Point3d)
. If Dnrm
is
non-null
, the method also calculated the normal derivative.
If the point lies outside the grid volume, ScalarGridBase.OUTSIDE_GRID
is
returned.norm
- returns the normal (local coordinates)Dnrm
- if non-null, returns the normal derivative (local coordinates)point
- point at which to calculate the normal and distance
(local coordinates).OUTSIDE_GRID
.public double getLocalDistanceAndGradient(Vector3d grad, Point3d point)
ScalarGridBase.OUTSIDE_GRID
is
returned.grad
- returns the gradient direction (local coordinates)point
- point at which to calculate the gradient and distance
(local coordinates).OUTSIDE_GRID
.public double getWorldDistanceAndGradient(Vector3d grad, Point3d point)
ScalarGridBase.OUTSIDE_GRID
is returned.grad
- returns the gradient (world coordinates)point
- point at which to calculate the gradient and distance
(world coordinates).public Feature getNearestLocalFeature(Point3d nearest, Point3d point)
getFeatures()
returns null
), or if the point lies outside the grid,
null
is returned.nearest
- returns the nearest point on the feature (local coordinates)point
- point for which to find nearest feature (local coordinates)public Feature getNearestWorldFeature(Point3d nearest, Point3d point)
getFeatures()
returns null
), or if the point lies outside
the grid, null
is returned.nearest
- returns the nearest point on the feature (world coordinates)point
- point for which to find nearest feature (world coordinates)public boolean isSigned()
true
if this grid is signed.public Feature[] getFeatures()
computeFromFeatures(java.util.List<? extends maspack.geometry.Feature>, double, maspack.matrix.RigidTransform3d, int, boolean)
, computeFromFeaturesOBB(java.util.List<? extends maspack.geometry.Feature>, double, int, boolean)
, or one of the associated constructors. They
will also be associated with the field if they were explicitly specified
via a call to setDistancesAndFeatures(double[], java.util.List<? extends maspack.geometry.Feature>, int[], boolean)
.
If no features are associated with this grid,
null
is returned.null
if
there are none.public void clearFeatures()
public Feature getClosestFeature(int idx)
idx
.
This assumes that the distance field is associated with features, as
described for getFeatures()
. If this is not the case,
null
is returned.idx
- vertex indexnull
if there are no featurespublic void smooth(int n)
n
- number of times to apply the smoothing operationpublic void smooth(double lambda, double mu)
lambda
- >
0, fraction of gradient to shrinkmu
- <
0, fraction of gradient to expandpublic void smooth(double lambda, double mu, int iters)
lambda
- >
0, fraction of gradient to shrinkmu
- <
0, fraction of gradient to expanditers
- number of applicationspublic void smooth()
public PolygonalMesh createQuadDistanceSurface(double val, int res)
val
.val
- iso surface valueres
- multiplies the resolution of this grid to obtain the
resolution of the grid used to create the mesh.public Vector3d getRenderVector(int xi, int yj, int zk)
public double getRenderVectorScale()
InterpolatingGridBase
InterpolatingGridBase.getRenderVector(maspack.matrix.Vector3d, int, int, int)
to determine the length of the
actual rendered vector. Returning 0 disables vector rendering.getRenderVectorScale
in class InterpolatingGridBase
public double getQuadDistance(Point3d point)
getQuadDistanceAndGradient(maspack.matrix.Vector3d, maspack.matrix.Matrix3d, maspack.matrix.Point3d)
.
If the point lies outside the grid volume, ScalarGridBase.OUTSIDE_GRID
is returned.point
- point at which to calculate the normal and distance
(local coordinates).OUTSIDE_GRID
.public java.lang.String getQuadCellInfo(Point3d point)
public double getWorldQuadDistance(Point3d point)
getQuadDistanceAndGradient(maspack.matrix.Vector3d, maspack.matrix.Matrix3d, maspack.matrix.Point3d)
.
If the point lies outside the grid volume, ScalarGridBase.OUTSIDE_GRID
is returned.point
- point at which to calculate the normal and distance
(world coordinates).OUTSIDE_GRID
.public DistanceGrid.TetDesc getQuadTet(Point3d point)
public double getQuadDistanceAndGradient(Vector3d grad, Matrix3d dgrad, Point3d point)
ScalarGridBase.OUTSIDE_GRID
is returned.grad
- returns the gradient direction (local coordinates)point
- point at which to calculate the gradient and distance
(local coordinates).dgrad
- if non-null, returns the derivative of the gradient (local
coordinates)OUTSIDE_GRID
.public double getWorldQuadDistanceAndGradient(Vector3d grad, Matrix3d dgrad, Point3d point)
getLocalDistanceAndGradient(Vector3d,Point3d)
.
If the point lies outside
the grid volume, ScalarGridBase.OUTSIDE_GRID
is returned.grad
- returns the gradient (world coordinates)dgrad
- if non-null, returns the derivative of the gradient (world
coordinates)point
- point at which to calculate the gradient and distance
(world coordinates).OUTSIDE_GRID
.public void computeQuadCoefsStub(double[] a, DistanceGrid.TetDesc tdesc)
public void computeQuadCoefs(double[] a, DistanceGrid.TetDesc tdesc)
public boolean findQuadSurfaceTangent(Point3d pt, Point3d p0, Point3d pa, Vector3d nrm)
pt
such that the line segment
pa-pt
is tangent to the quadratic zero distance surface and
also lies in the plane defined by point p0
and normal
nrm
. The search starts by finding a point on the surface
that is close to p0
and proceeds from there.pt
- returns the tangent pointp0
- point on the planepa
- origin for the tangent segmentnrm
- normal defining the planepublic boolean findQuadSurfaceIntersection(Point3d pi, Point3d p0, Point3d pa, Vector3d nrm)
pi
between
the quadratic zero surface and a ray (p0, pa)
lying in the
plane defined by point p0
and normal nrm
.
pa
is projected onto the plane internally so the caller does
not need to ensure this.pi
- returns the intersection pointp0
- point on the plane and ray originpa
- projected onto the plane to form the outer point defining the raynrm
- normal defining the planepublic boolean epsilonEquals(DistanceGrid grid, double tol)
true
if this grid equals another within a prescribed
tolerance. The grids are equal if the resolution, widths, center,
orientation, distances and signed property are equal. Specifying tol = 0
requires exact equality.grid
- grid to compare againsttol
- floating point tolerance (absolute)public void scaleDistance(double s)
scaleDistance
in class InterpolatingGridBase
s
- scaling factorpublic void markOutsideQuadtets(double mindist)
> mindist
as
being "outside".public void checkGradient()
public java.util.ArrayList<DistanceGrid.GridEdge> findSurfaceIntersectingEdges()
public java.util.ArrayList<DistanceGrid.TetDesc> findSurfaceEdgeTets(java.util.ArrayList<DistanceGrid.GridEdge> edges)