Simbody  3.5
SimTK::ContactGeometry::Ellipsoid Class Reference

This ContactGeometry subclass represents an ellipsoid centered at the origin, with its principal axes pointing along the x, y, and z axes and half dimensions a,b, and c (all > 0) along those axes, respectively. More...

+ Inheritance diagram for SimTK::ContactGeometry::Ellipsoid:

Public Member Functions

 Ellipsoid (const Vec3 &radii)
 Construct an Ellipsoid given its three principal half-axis dimensions a,b,c (all positive) along the local x,y,z directions respectively. More...
 
const Vec3getRadii () const
 Obtain the three half-axis dimensions a,b,c used to define this ellipsoid. More...
 
void setRadii (const Vec3 &radii)
 Set the three half-axis dimensions a,b,c (all positive) used to define this ellipsoid, overriding the current radii and recalculating the principal curvatures at a cost of about 30 flops. More...
 
const Vec3getCurvatures () const
 For efficiency we precalculate the principal curvatures whenever the ellipsoid radii are set; this avoids having to repeatedly perform these three expensive divisions at runtime. More...
 
UnitVec3 findUnitNormalAtPoint (const Vec3 &P) const
 Given a point P =(x,y,z) on the ellipsoid surface, return the unique unit outward normal to the ellipsoid at that point. More...
 
Vec3 findPointWithThisUnitNormal (const UnitVec3 &n) const
 Given a unit direction n, find the unique point P on the ellipsoid surface at which the outward-facing normal is n. More...
 
Vec3 findPointInSameDirection (const Vec3 &Q) const
 Given a direction d defined by the vector Q-O for an arbitrary point in space Q=(x,y,z)!=O, find the unique point P on the ellipsoid surface that is in direction d from the ellipsoid origin O. More...
 
void findParaboloidAtPoint (const Vec3 &Q, Transform &X_EP, Vec2 &k) const
 Given a point Q on the surface of the ellipsoid, find the approximating paraboloid at Q in a frame P where OP=Q, Pz is the outward-facing unit normal to the ellipsoid at Q, Px is the direction of maximum curvature and Py is the direction of minimum curvature. More...
 
void findParaboloidAtPointWithNormal (const Vec3 &Q, const UnitVec3 &n, Transform &X_EP, Vec2 &k) const
 If you already have both a point and the unit normal at that point, this will save about 40 flops by trusting that you have provided the correct normal; be careful – no one is going to check that you got this right. More...
 
const Impl & getImpl () const
 Internal use only. More...
 
Impl & updImpl ()
 Internal use only. More...
 
- Public Member Functions inherited from SimTK::ContactGeometry
 ContactGeometry ()
 Base class default constructor creates an empty handle. More...
 
 ContactGeometry (const ContactGeometry &src)
 Copy constructor makes a deep copy. More...
 
ContactGeometryoperator= (const ContactGeometry &src)
 Copy assignment makes a deep copy. More...
 
 ~ContactGeometry ()
 Base class destructor deletes the implementation object. Note that this is not virtual; handles should consist of just a pointer to the implementation. More...
 
DecorativeGeometry createDecorativeGeometry () const
 Generate a DecorativeGeometry that matches the shape of this ContactGeometry. More...
 
Vec3 findNearestPoint (const Vec3 &position, bool &inside, UnitVec3 &normal) const
 Given a point, find the nearest point on the surface of this object. More...
 
Vec3 projectDownhillToNearestPoint (const Vec3 &pointQ) const
 Given a query point Q, find the nearest point P on the surface of this object, looking only down the local gradient. More...
 
bool trackSeparationFromLine (const Vec3 &pointOnLine, const UnitVec3 &directionOfLine, const Vec3 &startingGuessForClosestPoint, Vec3 &newClosestPointOnSurface, Vec3 &closestPointOnLine, Real &height) const
 Track the closest point between this implicit surface and a given line, or the point of deepest penetration if the line intersects the surface. More...
 
bool intersectsRay (const Vec3 &origin, const UnitVec3 &direction, Real &distance, UnitVec3 &normal) const
 Determine whether this object intersects a ray, and if so, find the intersection point. More...
 
void getBoundingSphere (Vec3 &center, Real &radius) const
 Get a bounding sphere which completely encloses this object. More...
 
bool isSmooth () const
 Returns true if this is a smooth surface, meaning that it can provide meaningful curvature information and continuous derivatives with respect to its parameterization. More...
 
void calcCurvature (const Vec3 &point, Vec2 &curvature, Rotation &orientation) const
 Compute the principal curvatures and their directions, and the surface normal, at a given point on a smooth surface. More...
 
const FunctiongetImplicitFunction () const
 Our smooth surfaces define a function f(P)=0 that provides an implicit representation of the surface. More...
 
Real calcSurfaceValue (const Vec3 &point) const
 Calculate the value of the implicit surface function, at a given point. More...
 
UnitVec3 calcSurfaceUnitNormal (const Vec3 &point) const
 Calculate the implicit surface outward facing unit normal at the given point. More...
 
Vec3 calcSurfaceGradient (const Vec3 &point) const
 Calculate the gradient of the implicit surface function, at a given point. More...
 
Mat33 calcSurfaceHessian (const Vec3 &point) const
 Calculate the hessian of the implicit surface function, at a given point. More...
 
Real calcGaussianCurvature (const Vec3 &gradient, const Mat33 &Hessian) const
 For an implicit surface, return the Gaussian curvature at the point p whose implicit surface function gradient g(p) and Hessian H(p) are supplied. More...
 
Real calcGaussianCurvature (const Vec3 &point) const
 This signature is for convenience; use the other one to save time if you already have the gradient and Hessian available for this point. More...
 
Real calcSurfaceCurvatureInDirection (const Vec3 &point, const UnitVec3 &direction) const
 For an implicit surface, return the curvature k of the surface at a given point p in a given direction tp. More...
 
void calcSurfacePrincipalCurvatures (const Vec3 &point, Vec2 &curvature, Rotation &R_SP) const
 For an implicit surface at a given point p, return the principal curvatures and principal curvature directions, using only the implicit function and its derivatives. More...
 
bool isConvex () const
 Returns true if this surface is known to be convex. More...
 
Vec3 calcSupportPoint (UnitVec3 direction) const
 Given a direction expressed in the surface's frame S, return the point P on the surface that is the furthest in that direction (or one of those points if there is more than one). More...
 
ContactGeometryTypeId getTypeId () const
 ContactTrackerSubsystem uses this id for fast identification of specific surface shapes. More...
 
 ContactGeometry (ContactGeometryImpl *impl)
 Internal use only. More...
 
bool isOwnerHandle () const
 Internal use only. More...
 
bool isEmptyHandle () const
 Internal use only. More...
 
bool hasImpl () const
 Internal use only. More...
 
const ContactGeometryImpl & getImpl () const
 Internal use only. More...
 
ContactGeometryImpl & updImpl ()
 Internal use only. More...
 
void initGeodesic (const Vec3 &xP, const Vec3 &xQ, const Vec3 &xSP, const GeodesicOptions &options, Geodesic &geod) const
 Given two points, find a geodesic curve connecting them. More...
 
void continueGeodesic (const Vec3 &xP, const Vec3 &xQ, const Geodesic &prevGeod, const GeodesicOptions &options, Geodesic &geod) const
 Given the current positions of two points P and Q moving on this surface, and the previous geodesic curve G' connecting prior locations P' and Q' of those same two points, return the geodesic G between P and Q that is closest in length to the previous one. More...
 
void makeStraightLineGeodesic (const Vec3 &xP, const Vec3 &xQ, const UnitVec3 &defaultDirectionIfNeeded, const GeodesicOptions &options, Geodesic &geod) const
 Produce a straight-line approximation to the (presumably short) geodesic between two points on this implicit surface. More...
 
void shootGeodesicInDirectionUntilLengthReached (const Vec3 &xP, const UnitVec3 &tP, const Real &terminatingLength, const GeodesicOptions &options, Geodesic &geod) const
 Compute a geodesic curve starting at the given point, starting in the given direction, and terminating at the given length. More...
 
void calcGeodesicReverseSensitivity (Geodesic &geodesic, const Vec2 &initSensitivity=Vec2(0, 1)) const
 Given an already-calculated geodesic on this surface connecting points P and Q, fill in the sensitivity of point P with respect to a change of tangent direction at Q. More...
 
void shootGeodesicInDirectionUntilPlaneHit (const Vec3 &xP, const UnitVec3 &tP, const Plane &terminatingPlane, const GeodesicOptions &options, Geodesic &geod) const
 Compute a geodesic curve starting at the given point, starting in the given direction, and terminating when it hits the given plane. More...
 
void calcGeodesic (const Vec3 &xP, const Vec3 &xQ, const Vec3 &tPhint, const Vec3 &tQhint, Geodesic &geod) const
 Utility method to find geodesic between P and Q using split geodesic method with initial shooting directions tPhint and -tQhint. More...
 
void calcGeodesicUsingOrthogonalMethod (const Vec3 &xP, const Vec3 &xQ, const Vec3 &tPhint, Real lengthHint, Geodesic &geod) const
 Utility method to find geodesic between P and Q using the orthogonal method, with initial direction tPhint and initial length lengthHint. More...
 
void calcGeodesicUsingOrthogonalMethod (const Vec3 &xP, const Vec3 &xQ, Geodesic &geod) const
 This signature makes a guess at the initial direction and length and then calls the other signature. More...
 
Vec2 calcSplitGeodError (const Vec3 &P, const Vec3 &Q, const UnitVec3 &tP, const UnitVec3 &tQ, Geodesic *geod=0) const
 Utility method to calculate the "geodesic error" between one geodesic shot from P in the direction tP and another geodesic shot from Q in the direction tQ. More...
 
void shootGeodesicInDirectionUntilLengthReachedAnalytical (const Vec3 &xP, const UnitVec3 &tP, const Real &terminatingLength, const GeodesicOptions &options, Geodesic &geod) const
 Analytically compute a geodesic curve starting at the given point, starting in the given direction, and terminating at the given length. More...
 
void shootGeodesicInDirectionUntilPlaneHitAnalytical (const Vec3 &xP, const UnitVec3 &tP, const Plane &terminatingPlane, const GeodesicOptions &options, Geodesic &geod) const
 Analytically compute a geodesic curve starting at the given point, starting in the given direction, and terminating when it hits the given plane. More...
 
void calcGeodesicAnalytical (const Vec3 &xP, const Vec3 &xQ, const Vec3 &tPhint, const Vec3 &tQhint, Geodesic &geod) const
 Utility method to analytically find geodesic between P and Q with initial shooting directions tPhint and tQhint. More...
 
Vec2 calcSplitGeodErrorAnalytical (const Vec3 &P, const Vec3 &Q, const UnitVec3 &tP, const UnitVec3 &tQ, Geodesic *geod=0) const
 Utility method to analytically calculate the "geodesic error" between one geodesic shot from P in the direction tP and another geodesic shot from Q in the direction tQ. More...
 
const PlanegetPlane () const
 Get the plane associated with the geodesic hit plane event handler. More...
 
void setPlane (const Plane &plane) const
 Set the plane associated with the geodesic hit plane event handler. More...
 
const GeodesicgetGeodP () const
 Get the geodesic for access by visualizer. More...
 
const GeodesicgetGeodQ () const
 Get the geodesic for access by visualizer. More...
 
const int getNumGeodesicsShot () const
 Get the plane associated with the geodesic hit plane event handler. More...
 
void addVizReporter (ScheduledEventReporter *reporter) const
 Get the plane associated with the geodesic hit plane event handler. More...
 

Static Public Member Functions

static bool isInstance (const ContactGeometry &geo)
 Return true if the supplied ContactGeometry object is an Ellipsoid. More...
 
static const EllipsoidgetAs (const ContactGeometry &geo)
 Cast the supplied ContactGeometry object to a const Ellipsoid. More...
 
static EllipsoidupdAs (ContactGeometry &geo)
 Cast the supplied ContactGeometry object to a writable Ellipsoid. More...
 
static ContactGeometryTypeId classTypeId ()
 Obtain the unique id for Ellipsoid contact geometry. More...
 
- Static Public Member Functions inherited from SimTK::ContactGeometry
static Vec2 evalParametricCurvature (const Vec3 &P, const UnitVec3 &nn, const Vec3 &dPdu, const Vec3 &dPdv, const Vec3 &d2Pdu2, const Vec3 &d2Pdv2, const Vec3 &d2Pdudv, Transform &X_EP)
 Calculate surface curvature at a point using differential geometry as suggested by Harris 2006, "Curvature of ellipsoids and other surfaces" Ophthal. More...
 
static void combineParaboloids (const Rotation &R_SP1, const Vec2 &k1, const UnitVec3 &x2, const Vec2 &k2, Rotation &R_SP, Vec2 &k)
 This utility method is useful for characterizing the relative geometry of two locally-smooth surfaces in contact, in a way that is useful for later application of Hertz compliant contact theory for generating forces. More...
 
static void combineParaboloids (const Rotation &R_SP1, const Vec2 &k1, const UnitVec3 &x2, const Vec2 &k2, Vec2 &k)
 This is a much faster version of combineParaboloids() for when you just need the curvatures of the difference paraboloid, but not the directions of those curvatures. More...
 

Additional Inherited Members

- Protected Attributes inherited from SimTK::ContactGeometry
ContactGeometryImpl * impl
 Internal use only. More...
 

Detailed Description

This ContactGeometry subclass represents an ellipsoid centered at the origin, with its principal axes pointing along the x, y, and z axes and half dimensions a,b, and c (all > 0) along those axes, respectively.

The implicit equation f(x,y,z)=0 of the ellipsoid surface is

    f(x,y,z) = Ax^2+By^2+Cz^2 - 1
    where A=1/a^2, B=1/b^2, C=1/c^2     

A,B, and C are the squares of the principal curvatures ka=1/a, kb=1/b, and kc=1/c.

The interior of the ellipsoid consists of all points such that f(x,y,z)<0 and points exterior satisfy f(x,y,z)>0. The region around any point (x,y,z) on an ellipsoid surface is locally an elliptic paraboloid with equation

    -2 z' = kmax x'^2 + kmin y'^2   

where z' is measured along the the outward unit normal n at (x,y,z), x' is measured along the the unit direction u of maximum curvature, and y' is measured along the unit direction v of minimum curvature. kmax,kmin are the curvatures with kmax >= kmin > 0. The signs of the mutually perpendicular vectors u and v are chosen so that (u,v,n) forms a right-handed coordinate system for the paraboloid.

Constructor & Destructor Documentation

Construct an Ellipsoid given its three principal half-axis dimensions a,b,c (all positive) along the local x,y,z directions respectively.

The curvatures (reciprocals of radii) are precalculated here at a cost of about 30 flops.

Member Function Documentation

const Vec3& SimTK::ContactGeometry::Ellipsoid::getRadii ( ) const

Obtain the three half-axis dimensions a,b,c used to define this ellipsoid.

void SimTK::ContactGeometry::Ellipsoid::setRadii ( const Vec3 radii)

Set the three half-axis dimensions a,b,c (all positive) used to define this ellipsoid, overriding the current radii and recalculating the principal curvatures at a cost of about 30 flops.

Parameters
[in]radiiThe three half-dimensions of the ellipsoid, in the ellipsoid's local x, y, and z directions respectively.
const Vec3& SimTK::ContactGeometry::Ellipsoid::getCurvatures ( ) const

For efficiency we precalculate the principal curvatures whenever the ellipsoid radii are set; this avoids having to repeatedly perform these three expensive divisions at runtime.

The curvatures are ka=1/a, kb=1/b, and kc=1/c so that the ellipsoid's implicit equation can be written Ax^2+By^2+Cz^2=1, with A=ka^2, etc.

UnitVec3 SimTK::ContactGeometry::Ellipsoid::findUnitNormalAtPoint ( const Vec3 P) const

Given a point P =(x,y,z) on the ellipsoid surface, return the unique unit outward normal to the ellipsoid at that point.

If P is not on the surface, the result is the same as for the point obtained by scaling the vector P - O until it just touches the surface. That is, we compute P'=findPointInThisDirection(P) and then return the normal at P'. Cost is about 40 flops regardless of whether P was initially on the surface.

Parameters
[in]PA point on the ellipsoid surface, measured and expressed in the ellipsoid's local frame. See text for what happens if P is not actually on the ellipsoid surface.
Returns
The outward-facing unit normal at point P (or at the surface point pointed to by P).
See also
findPointInSameDirection()
Vec3 SimTK::ContactGeometry::Ellipsoid::findPointWithThisUnitNormal ( const UnitVec3 n) const

Given a unit direction n, find the unique point P on the ellipsoid surface at which the outward-facing normal is n.

Cost is about 40 flops.

Parameters
[in]nThe unit vector for which we want to find a match on the ellipsoid surface, expressed in the ellipsoid's local frame.
Returns
The point on the ellipsoid's surface at which the outward-facing normal is the same as n. The point is measured and expressed in the ellipsoid's local frame.
Vec3 SimTK::ContactGeometry::Ellipsoid::findPointInSameDirection ( const Vec3 Q) const

Given a direction d defined by the vector Q-O for an arbitrary point in space Q=(x,y,z)!=O, find the unique point P on the ellipsoid surface that is in direction d from the ellipsoid origin O.

That is, P=s*d for some scalar s > 0 such that f(P)=0. Cost is about 40 flops.

Parameters
[in]QA point in space measured from the ellipsoid origin but not the origin.
Returns
P, the intersection of the ray in the direction Q-O with the ellipsoid surface
void SimTK::ContactGeometry::Ellipsoid::findParaboloidAtPoint ( const Vec3 Q,
Transform X_EP,
Vec2 k 
) const

Given a point Q on the surface of the ellipsoid, find the approximating paraboloid at Q in a frame P where OP=Q, Pz is the outward-facing unit normal to the ellipsoid at Q, Px is the direction of maximum curvature and Py is the direction of minimum curvature.

k=(kmax,kmin) are the returned curvatures with kmax >= kmin > 0. The equation of the resulting paraboloid in the P frame is -2z = kmax*x^2 + kmin*y^2. Cost is about 260 flops; you can save a little time if you already know the normal at Q by using the other overloaded signature for this method.

Warning
It is up to you to make sure that Q is actually on the ellipsoid surface. If it is not you will quietly get a meaningless result.
Parameters
[in]QA point on the surface of this ellipsoid, measured and expressed in the ellipsoid's local frame.
[out]X_EPThe frame of the paraboloid P, measured and expressed in the ellipsoid local frame E. X_EP.p() is Q, X_EP.x() is the calculated direction of maximum curvature kmax; y() is the direction of minimum curvature kmin; z is the outward facing normal at Q.
[out]kThe maximum (k[0]) and minimum (k[1]) curvatures of the ellipsoid (and paraboloid P) at point Q.
See also
findParaboloidAtPointWithNormal()
void SimTK::ContactGeometry::Ellipsoid::findParaboloidAtPointWithNormal ( const Vec3 Q,
const UnitVec3 n,
Transform X_EP,
Vec2 k 
) const

If you already have both a point and the unit normal at that point, this will save about 40 flops by trusting that you have provided the correct normal; be careful – no one is going to check that you got this right.

The results are meaningless if the point and normal are not consistent. Cost is about 220 flops.

See also
findParaboloidAtPoint() for details
static bool SimTK::ContactGeometry::Ellipsoid::isInstance ( const ContactGeometry geo)
inlinestatic

Return true if the supplied ContactGeometry object is an Ellipsoid.

static const Ellipsoid& SimTK::ContactGeometry::Ellipsoid::getAs ( const ContactGeometry geo)
inlinestatic

Cast the supplied ContactGeometry object to a const Ellipsoid.

static Ellipsoid& SimTK::ContactGeometry::Ellipsoid::updAs ( ContactGeometry geo)
inlinestatic

Cast the supplied ContactGeometry object to a writable Ellipsoid.

static ContactGeometryTypeId SimTK::ContactGeometry::Ellipsoid::classTypeId ( )
static

Obtain the unique id for Ellipsoid contact geometry.

const Impl& SimTK::ContactGeometry::Ellipsoid::getImpl ( ) const

Internal use only.

Internal use only.

Impl& SimTK::ContactGeometry::Ellipsoid::updImpl ( )

Internal use only.


The documentation for this class was generated from the following file: