Simbody  3.8
ContactGeometry.h
Go to the documentation of this file.
1 #ifndef SimTK_SIMMATH_CONTACT_GEOMETRY_H_
2 #define SimTK_SIMMATH_CONTACT_GEOMETRY_H_
3 
4 /* -------------------------------------------------------------------------- *
5  * Simbody(tm): SimTKmath *
6  * -------------------------------------------------------------------------- *
7  * This is part of the SimTK biosimulation toolkit originating from *
8  * Simbios, the NIH National Center for Physics-Based Simulation of *
9  * Biological Structures at Stanford, funded under the NIH Roadmap for *
10  * Medical Research, grant U54 GM072970. See https://simtk.org/home/simbody. *
11  * *
12  * Portions copyright (c) 2008-12 Stanford University and the Authors. *
13  * Authors: Peter Eastman, Michael Sherman *
14  * Contributors: Ian Stavness, Andreas Scholz *
15  * *
16  * Licensed under the Apache License, Version 2.0 (the "License"); you may *
17  * not use this file except in compliance with the License. You may obtain a *
18  * copy of the License at http://www.apache.org/licenses/LICENSE-2.0. *
19  * *
20  * Unless required by applicable law or agreed to in writing, software *
21  * distributed under the License is distributed on an "AS IS" BASIS, *
22  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. *
23  * See the License for the specific language governing permissions and *
24  * limitations under the License. *
25  * -------------------------------------------------------------------------- */
26 
31 #include "SimTKcommon.h"
36 
37 #include <cassert>
38 #include <functional>
39 
40 namespace SimTK {
41 
46 
47 class ContactGeometryImpl;
48 class OBBTreeNodeImpl;
49 class OBBTree;
50 class Plane;
51 
52 
53 
54 //==============================================================================
55 // CONTACT GEOMETRY
56 //==============================================================================
112 public:
113 class HalfSpace;
114 class Sphere;
115 class Ellipsoid;
116 class Torus;
117 class SmoothHeightMap;
118 class Cylinder;
119 class Brick;
120 class TriangleMesh;
121 
122 // TODO
123 class Cone;
124 
128  Real arcLength = NaN;
129 
131  Vec3 point {NaN};
133  UnitVec3 tangent {NaN, NaN, NaN};
134 
136  Real jacobiRot = NaN;
138  Real jacobiTrans = NaN;
139 
141  Real jacobiRotDot = NaN;
143  Real jacobiTransDot = NaN;
144 };
145 
147 ContactGeometry() : impl(0) {}
156 
159 
164 bool isSurfaceDefined(const Vec3& point) const;
165 
176 Vec3 findNearestPoint(const Vec3& position, bool& inside, UnitVec3& normal) const;
177 
193 
258 bool trackSeparationFromLine(const Vec3& pointOnLine,
259  const UnitVec3& directionOfLine,
260  const Vec3& startingGuessForClosestPoint,
261  Vec3& newClosestPointOnSurface,
262  Vec3& closestPointOnLine,
263  Real& height) const;
264 
267 {
268  Converged,
269  MaxIterationsExceeded,
270  PointFallsBelowSurface,
271 };
272 
306  const Vec3& pointA,
307  const Vec3& pointB,
308  int maxIterations,
309  Real tolerance,
310  Vec3& nearestPointOnLine) const;
311 
323 bool intersectsRay(const Vec3& origin, const UnitVec3& direction,
324  Real& distance, UnitVec3& normal) const;
325 
331 void getBoundingSphere(Vec3& center, Real& radius) const;
332 
336 bool isSmooth() const;
337 
354 void calcCurvature(const Vec3& point, Vec2& curvature,
355  Rotation& orientation) const;
356 
368 
374 Real calcSurfaceValue(const Vec3& point) const;
375 
386 UnitVec3 calcSurfaceUnitNormal(const Vec3& point) const;
387 
393 Vec3 calcSurfaceGradient(const Vec3& point) const;
394 
400 Mat33 calcSurfaceHessian(const Vec3& point) const;
401 
431  const Mat33& Hessian) const;
432 
436 Real calcGaussianCurvature(const Vec3& point) const {
437  return calcGaussianCurvature(calcSurfaceGradient(point),
438  calcSurfaceHessian(point));
439 }
440 
449 Real calcSurfaceCurvatureInDirection(const Vec3& point, const UnitVec3& direction) const;
450 
459 Real calcSurfaceTorsionInDirection(const Vec3& point, const UnitVec3& direction) const;
460 
469 void calcSurfacePrincipalCurvatures(const Vec3& point, Vec2& curvature,
470  Rotation& R_SP) const;
471 
474 bool isConvex() const;
475 
481 Vec3 calcSupportPoint(UnitVec3 direction) const;
482 
486 
537 static Vec2 evalParametricCurvature(const Vec3& P, const UnitVec3& nn,
538  const Vec3& dPdu, const Vec3& dPdv,
539  const Vec3& d2Pdu2, const Vec3& d2Pdv2,
540  const Vec3& d2Pdudv,
541  Transform& X_EP);
542 
616 static void combineParaboloids(const Rotation& R_SP1, const Vec2& k1,
617  const UnitVec3& x2, const Vec2& k2,
618  Rotation& R_SP, Vec2& k);
619 
624 static void combineParaboloids(const Rotation& R_SP1, const Vec2& k1,
625  const UnitVec3& x2, const Vec2& k2,
626  Vec2& k);
627 
628 
642 void initGeodesic(const Vec3& xP, const Vec3& xQ, const Vec3& xSP,
643  const GeodesicOptions& options, Geodesic& geod) const;
644 
648 
675  const Vec3& initialPointApprox,
676  const Vec3& initialTangentApprox,
677  Real finalArcLength,
678  int numberOfKnotPoints,
679  const std::function<void(const ContactGeometry::GeodesicKnotPoint&)>&
680  geodesicKnotPointsSink) const;
681 
705  const Vec3& initialPointApprox,
706  const Vec3& initialTangentApprox,
707  Real finalArcLength,
708  Real initialIntegratorStepSize,
709  Real integratorAccuracy,
710  Real constraintTolerance,
711  int maxIterations,
712  const std::function<void(const ContactGeometry::GeodesicKnotPoint&)>&
713  geodesicKnotPointsSink) const;
714 
757 // XXX if xP and xQ are the exact end-points of prevGeod; then geod = prevGeod;
758 void continueGeodesic(const Vec3& xP, const Vec3& xQ, const Geodesic& prevGeod,
759  const GeodesicOptions& options, Geodesic& geod) const;
760 
787 void makeStraightLineGeodesic(const Vec3& xP, const Vec3& xQ,
788  const UnitVec3& defaultDirectionIfNeeded,
789  const GeodesicOptions& options, Geodesic& geod) const;
790 
791 
801 // XXX what to do if tP is not in the tangent plane at P -- project it?
803  (const Vec3& xP, const UnitVec3& tP, const Real& terminatingLength,
804  const GeodesicOptions& options, Geodesic& geod) const;
805 
820  (Geodesic& geodesic,
821  const Vec2& initSensitivity = Vec2(0,1)) const; // j, jdot at end point
822 
823 
833 // XXX what to do if tP is not in the tangent plane at P -- project it?
834 // XXX what to do if we don't hit the plane
836  const Plane& terminatingPlane, const GeodesicOptions& options,
837  Geodesic& geod) const;
838 
839 
842 void calcGeodesic(const Vec3& xP, const Vec3& xQ,
843  const Vec3& tPhint, const Vec3& tQhint, Geodesic& geod) const;
844 
847 void calcGeodesicUsingOrthogonalMethod(const Vec3& xP, const Vec3& xQ,
848  const Vec3& tPhint, Real lengthHint, Geodesic& geod) const;
849 
852 void calcGeodesicUsingOrthogonalMethod(const Vec3& xP, const Vec3& xQ,
853  Geodesic& geod) const
854 {
855  const Vec3 r_PQ = xQ - xP;
856  const Real lengthHint = r_PQ.norm();
857  const UnitVec3 n = calcSurfaceUnitNormal(xP);
858  // Project r_PQ into the tangent plane.
859  const Vec3 t_PQ = r_PQ - (~r_PQ*n)*n;
860  const Real tLength = t_PQ.norm();
861  const UnitVec3 tPhint =
862  tLength != 0 ? UnitVec3(t_PQ/tLength, true)
863  : n.perp(); // some arbitrary perpendicular to n
864  calcGeodesicUsingOrthogonalMethod(xP, xQ, Vec3(tPhint), lengthHint, geod);
865 }
866 
867 
875 Vec2 calcSplitGeodError(const Vec3& P, const Vec3& Q,
876  const UnitVec3& tP, const UnitVec3& tQ,
877  Geodesic* geod=0) const;
878 
879 
880 
891 // XXX what to do if tP is not in the tangent plane at P -- project it?
893  (const Vec3& xP, const UnitVec3& tP, const Real& terminatingLength,
894  const GeodesicOptions& options, Geodesic& geod) const;
895 
896 
907 // XXX what to do if tP is not in the tangent plane at P -- project it?
908 // XXX what to do if we don't hit the plane
910  const Plane& terminatingPlane, const GeodesicOptions& options,
911  Geodesic& geod) const;
912 
913 
918 void calcGeodesicAnalytical(const Vec3& xP, const Vec3& xQ,
919  const Vec3& tPhint, const Vec3& tQhint, Geodesic& geod) const;
920 
930  const UnitVec3& tP, const UnitVec3& tQ,
931  Geodesic* geod=0) const;
932 
943 const Plane& getPlane() const;
946 void setPlane(const Plane& plane) const;
948 const Geodesic& getGeodP() const;
950 const Geodesic& getGeodQ() const;
951 const int getNumGeodesicsShot() const;
957 explicit ContactGeometry(ContactGeometryImpl* impl);
958 bool isOwnerHandle() const;
959 bool isEmptyHandle() const;
960 bool hasImpl() const {return impl != 0;}
962 const ContactGeometryImpl& getImpl() const {assert(impl); return *impl;}
964 ContactGeometryImpl& updImpl() {assert(impl); return *impl; }
965 
966 protected:
967 ContactGeometryImpl* impl;
968 };
969 
970 
971 
972 //==============================================================================
973 // HALF SPACE
974 //==============================================================================
979 public:
986 
989 
991 static bool isInstance(const ContactGeometry& geo)
992 { return geo.getTypeId()==classTypeId(); }
994 static const HalfSpace& getAs(const ContactGeometry& geo)
995 { assert(isInstance(geo)); return static_cast<const HalfSpace&>(geo); }
998 { assert(isInstance(geo)); return static_cast<HalfSpace&>(geo); }
999 
1002 
1003 class Impl;
1004 const Impl& getImpl() const;
1005 Impl& updImpl();
1006 };
1007 
1008 
1009 
1010 //==============================================================================
1011 // CYLINDER
1012 //==============================================================================
1017 public:
1018 explicit Cylinder(Real radius);
1019 Real getRadius() const;
1020 void setRadius(Real radius);
1021 
1023 static bool isInstance(const ContactGeometry& geo)
1024 { return geo.getTypeId()==classTypeId(); }
1026 static const Cylinder& getAs(const ContactGeometry& geo)
1027 { assert(isInstance(geo)); return static_cast<const Cylinder&>(geo); }
1030 { assert(isInstance(geo)); return static_cast<Cylinder&>(geo); }
1031 
1034 
1035 class Impl;
1036 const Impl& getImpl() const;
1037 Impl& updImpl();
1038 };
1039 
1040 
1041 
1042 //==============================================================================
1043 // SPHERE
1044 //==============================================================================
1048 public:
1049 explicit Sphere(Real radius);
1050 Real getRadius() const;
1051 void setRadius(Real radius);
1052 
1054 static bool isInstance(const ContactGeometry& geo)
1055 { return geo.getTypeId()==classTypeId(); }
1057 static const Sphere& getAs(const ContactGeometry& geo)
1058 { assert(isInstance(geo)); return static_cast<const Sphere&>(geo); }
1061 { assert(isInstance(geo)); return static_cast<Sphere&>(geo); }
1062 
1065 
1066 class Impl;
1067 const Impl& getImpl() const;
1068 Impl& updImpl();
1069 };
1070 
1071 
1072 
1073 //==============================================================================
1074 // ELLIPSOID
1075 //==============================================================================
1098 public:
1102 explicit Ellipsoid(const Vec3& radii);
1105 const Vec3& getRadii() const;
1111 void setRadii(const Vec3& radii);
1112 
1118 const Vec3& getCurvatures() const;
1119 
1133 
1142 
1152 
1175 void findParaboloidAtPoint(const Vec3& Q, Transform& X_EP, Vec2& k) const;
1176 
1183  Transform& X_EP, Vec2& k) const;
1184 
1186 static bool isInstance(const ContactGeometry& geo)
1187 { return geo.getTypeId()==classTypeId(); }
1189 static const Ellipsoid& getAs(const ContactGeometry& geo)
1190 { assert(isInstance(geo)); return static_cast<const Ellipsoid&>(geo); }
1193 { assert(isInstance(geo)); return static_cast<Ellipsoid&>(geo); }
1194 
1197 
1198 class Impl;
1199 const Impl& getImpl() const;
1200 Impl& updImpl();
1201 };
1202 
1203 
1204 
1205 //==============================================================================
1206 // SMOOTH HEIGHT MAP
1207 //==============================================================================
1220 public:
1224 explicit SmoothHeightMap(const BicubicSurface& surface);
1225 
1229 
1232 const OBBTree& getOBBTree() const;
1233 
1235 static bool isInstance(const ContactGeometry& geo)
1236 { return geo.getTypeId()==classTypeId(); }
1238 static const SmoothHeightMap& getAs(const ContactGeometry& geo)
1239 { assert(isInstance(geo)); return static_cast<const SmoothHeightMap&>(geo); }
1242 { assert(isInstance(geo)); return static_cast<SmoothHeightMap&>(geo); }
1243 
1246 
1247 class Impl;
1248 const Impl& getImpl() const;
1249 Impl& updImpl();
1250 };
1251 
1252 
1253 //==============================================================================
1254 // BRICK
1255 //==============================================================================
1259 public:
1262 explicit Brick(const Vec3& halfLengths);
1266 const Vec3& getHalfLengths() const;
1268 void setHalfLengths(const Vec3& halfLengths);
1269 
1271 const Geo::Box& getGeoBox() const;
1272 
1274 static bool isInstance(const ContactGeometry& geo)
1275 { return geo.getTypeId()==classTypeId(); }
1277 static const Brick& getAs(const ContactGeometry& geo)
1278 { assert(isInstance(geo)); return static_cast<const Brick&>(geo); }
1281 { assert(isInstance(geo)); return static_cast<Brick&>(geo); }
1282 
1285 
1286 class Impl;
1287 const Impl& getImpl() const;
1288 Impl& updImpl();
1289 };
1290 
1291 
1292 //==============================================================================
1293 // TRIANGLE MESH
1294 //==============================================================================
1316 : public ContactGeometry {
1317 public:
1318 class OBBTreeNode;
1329 TriangleMesh(const ArrayViewConst_<Vec3>& vertices, const ArrayViewConst_<int>& faceIndices, bool smooth=false);
1338 explicit TriangleMesh(const PolygonalMesh& mesh, bool smooth=false);
1340 int getNumEdges() const;
1342 int getNumFaces() const;
1344 int getNumVertices() const;
1348 const Vec3& getVertexPosition(int index) const;
1354 int getFaceEdge(int face, int edge) const;
1359 int getFaceVertex(int face, int vertex) const;
1364 int getEdgeFace(int edge, int face) const;
1369 int getEdgeVertex(int edge, int vertex) const;
1374 void findVertexEdges(int vertex, Array_<int>& edges) const;
1377 const UnitVec3& getFaceNormal(int face) const;
1380 Real getFaceArea(int face) const;
1386 Vec3 findPoint(int face, const Vec2& uv) const;
1391 Vec3 findCentroid(int face) const;
1396 UnitVec3 findNormalAtPoint(int face, const Vec2& uv) const;
1407 Vec3 findNearestPoint(const Vec3& position, bool& inside, UnitVec3& normal) const;
1420 Vec3 findNearestPoint(const Vec3& position, bool& inside, int& face, Vec2& uv) const;
1421 
1431 Vec3 findNearestPointToFace(const Vec3& position, int face, Vec2& uv) const;
1432 
1433 
1445 bool intersectsRay(const Vec3& origin, const UnitVec3& direction, Real& distance, UnitVec3& normal) const;
1459 bool intersectsRay(const Vec3& origin, const UnitVec3& direction, Real& distance, int& face, Vec2& uv) const;
1463 
1467 
1469 static bool isInstance(const ContactGeometry& geo)
1470 { return geo.getTypeId()==classTypeId(); }
1472 static const TriangleMesh& getAs(const ContactGeometry& geo)
1473 { assert(isInstance(geo)); return static_cast<const TriangleMesh&>(geo); }
1476 { assert(isInstance(geo)); return static_cast<TriangleMesh&>(geo); }
1477 
1480 
1481 class Impl;
1482 const Impl& getImpl() const;
1483 Impl& updImpl();
1484 };
1485 
1486 
1487 
1488 //==============================================================================
1489 // TRIANGLE MESH :: OBB TREE NODE
1490 //==============================================================================
1496 public:
1497 OBBTreeNode(const OBBTreeNodeImpl& impl);
1502 bool isLeafNode() const;
1511 const Array_<int>& getTriangles() const;
1515 int getNumTriangles() const;
1516 
1517 private:
1518 const OBBTreeNodeImpl* impl;
1519 };
1520 
1521 //==============================================================================
1522 // TORUS
1523 //==============================================================================
1530 public:
1531 Torus(Real torusRadius, Real tubeRadius);
1533 void setTorusRadius(Real radius);
1535 void setTubeRadius(Real radius);
1536 
1538 static bool isInstance(const ContactGeometry& geo)
1539 { return geo.getTypeId()==classTypeId(); }
1541 static const Torus& getAs(const ContactGeometry& geo)
1542 { assert(isInstance(geo)); return static_cast<const Torus&>(geo); }
1545 { assert(isInstance(geo)); return static_cast<Torus&>(geo); }
1546 
1549 
1550 class Impl;
1551 const Impl& getImpl() const;
1552 Impl& updImpl();
1553 };
1554 
1555 
1556 
1557 
1558 //==============================================================================
1559 // GEODESIC EVALUATOR helper classes
1560 //==============================================================================
1561 
1562 
1566 class Plane {
1567 public:
1568  Plane() : m_normal(1,0,0), m_offset(0) { }
1569  Plane(const Vec3& normal, const Real& offset)
1570  : m_normal(normal), m_offset(offset) { }
1571 
1572  Real getDistance(const Vec3& pt) const {
1573  return ~m_normal*pt - m_offset;
1574  }
1575 
1576  Vec3 getNormal() const {
1577  return m_normal;
1578  }
1579 
1580  Real getOffset() const {
1581  return m_offset;
1582  }
1583 
1584 private:
1585  Vec3 m_normal;
1586  Real m_offset;
1587 }; // class Plane
1588 
1589 
1595 public:
1597  : TriggeredEventHandler(Stage::Position) { }
1598 
1599  explicit GeodHitPlaneEvent(const Plane& aplane)
1600  : TriggeredEventHandler(Stage::Position) {
1601  plane = aplane;
1602  }
1603 
1604  // event is triggered if distance of geodesic endpoint to plane is zero
1605  Real getValue(const State& state) const override {
1606  if (!enabled) {
1607  return 1;
1608  }
1609  Vec3 endpt(&state.getQ()[0]);
1610  Real dist = plane.getDistance(endpt);
1611 // std::cout << "dist = " << dist << std::endl;
1612  return dist;
1613  }
1614 
1615  // This method is called whenever this event occurs.
1616  void handleEvent(State& state, Real accuracy, bool& shouldTerminate) const override {
1617  if (!enabled) {
1618  return;
1619  }
1620 
1621  // This should be triggered when geodesic endpoint to plane is zero.
1622  Vec3 endpt;
1623  const Vector& q = state.getQ();
1624  endpt[0] = q[0]; endpt[1] = q[1]; endpt[2] = q[2];
1625  shouldTerminate = true;
1626  }
1627 
1628  void setPlane(const Plane& aplane) const {
1629  plane = aplane;
1630  }
1631 
1632  const Plane& getPlane() const {
1633  return plane;
1634  }
1635 
1636  const void setEnabled(bool enabledFlag) {
1637  enabled = enabledFlag;
1638  }
1639 
1640  const bool isEnabled() {
1641  return enabled;
1642  }
1643 
1644 private:
1645  mutable Plane plane;
1646  bool enabled;
1647 
1648 }; // class GeodHitPlaneEvent
1649 
1654 public:
1655  PathDecorator(const Vector& x, const Vec3& O, const Vec3& I, const Vec3& color) :
1656  m_x(x), m_O(O), m_I(I), m_color(color) { }
1657 
1658  virtual void generateDecorations(const State& state,
1659  Array_<DecorativeGeometry>& geometry) override {
1660 // m_system.realize(state, Stage::Position);
1661 
1662  Vec3 P, Q;
1663  P[0] = m_x[0]; P[1] = m_x[1]; P[2] = m_x[2];
1664  Q[0] = m_x[3]; Q[1] = m_x[4]; Q[2] = m_x[5];
1665 
1666  geometry.push_back(DecorativeSphere(Real(.05)).setColor(Black).setTransform(m_O));
1667  geometry.push_back(DecorativeSphere(Real(.05)).setColor(Black).setTransform(P));
1668  geometry.push_back(DecorativeSphere(Real(.05)).setColor(Black).setTransform(Q));
1669  geometry.push_back(DecorativeSphere(Real(.05)).setColor(Black).setTransform(m_I));
1670 
1671  geometry.push_back(DecorativeLine(m_O,P)
1672  .setColor(m_color)
1673  .setLineThickness(2));
1674  geometry.push_back(DecorativeLine(Q,m_I)
1675  .setColor(m_color)
1676  .setLineThickness(2));
1677 
1678  }
1679 
1680 private:
1681  const Vector& m_x; // x = ~[P Q]
1682  const Vec3& m_O;
1683  const Vec3& m_I;
1684  const Vec3& m_color;
1685  Rotation R_plane;
1686  Vec3 offset;
1687 }; // class DecorationGenerator
1688 
1689 
1694 public:
1695  PlaneDecorator(const Plane& plane, const Vec3& color) :
1696  m_plane(plane), m_color(color) { }
1697 
1698  virtual void generateDecorations(const State& state,
1699  Array_<DecorativeGeometry>& geometry) override {
1700 // m_system.realize(state, Stage::Position);
1701 
1702  // draw plane
1703  R_plane.setRotationFromOneAxis(UnitVec3(m_plane.getNormal()),
1705  offset = 0;
1706  offset[0] = m_plane.getOffset();
1707  geometry.push_back(
1708  DecorativeBrick(Vec3(Real(.01),1,1))
1709  .setTransform(Transform(R_plane, R_plane*offset))
1710  .setColor(m_color)
1711  .setOpacity(Real(.2)));
1712  }
1713 
1714 private:
1715  const Plane& m_plane;
1716  const Vec3& m_color;
1717  Rotation R_plane;
1718  Vec3 offset;
1719 }; // class DecorationGenerator
1720 
1721 
1722 } // namespace SimTK
1723 
1724 #endif // SimTK_SIMMATH_CONTACT_GEOMETRY_H_
This file defines the BicubicSurface class, and the BicubicFunction class that uses it to create a tw...
This file defines the Geodesic class.
Includes internal headers providing declarations for the basic SimTK Core classes,...
This is the header file that every Simmath compilation unit should include first.
#define SimTK_SIMMATH_EXPORT
Definition: SimTKmath/include/simmath/internal/common.h:64
This Array_ helper class is the base class for ArrayView_ which is the base class for Array_; here we...
Definition: Array.h:324
void push_back(const T &value)
This method increases the size of the Array by one element at the end and initializes that element by...
Definition: Array.h:2399
This class will create a smooth surface that approximates a two-argument function F(X,...
Definition: BicubicSurface.h:158
This is a unique integer type for quickly identifying specific types of contact geometry for fast loo...
This ContactGeometry subclass represents a rectangular solid centered at the origin.
Definition: ContactGeometry.h:1258
const Vec3 & getHalfLengths() const
Get the half-dimensions of this Brick, expressed in its own frame.
static Brick & updAs(ContactGeometry &geo)
Cast the supplied ContactGeometry object to a writable brick.
Definition: ContactGeometry.h:1280
void setHalfLengths(const Vec3 &halfLengths)
Change the shape or size of this brick by setting its half-dimensions.
Brick(const Vec3 &halfLengths)
Create a brick-shaped contact shape of the given half-dimensions, expressed in the brick's local fram...
static ContactGeometryTypeId classTypeId()
Obtain the unique id for Brick contact geometry.
Impl & updImpl()
Internal use only.
const Impl & getImpl() const
Internal use only.
static bool isInstance(const ContactGeometry &geo)
Return true if the supplied ContactGeometry object is a brick.
Definition: ContactGeometry.h:1274
static const Brick & getAs(const ContactGeometry &geo)
Cast the supplied ContactGeometry object to a const brick.
Definition: ContactGeometry.h:1277
const Geo::Box & getGeoBox() const
Get the Geo::Box object used to represent this Brick.
This ContactGeometry subclass represents a cylinder centered at the origin, with radius r in the x-y ...
Definition: ContactGeometry.h:1016
static bool isInstance(const ContactGeometry &geo)
Return true if the supplied ContactGeometry object is a cylinder.
Definition: ContactGeometry.h:1023
static Cylinder & updAs(ContactGeometry &geo)
Cast the supplied ContactGeometry object to a writable cylinder.
Definition: ContactGeometry.h:1029
static const Cylinder & getAs(const ContactGeometry &geo)
Cast the supplied ContactGeometry object to a const cylinder.
Definition: ContactGeometry.h:1026
static ContactGeometryTypeId classTypeId()
Obtain the unique id for Cylinder contact geometry.
Impl & updImpl()
Internal use only.
const Impl & getImpl() const
Internal use only.
This ContactGeometry subclass represents an ellipsoid centered at the origin, with its principal axes...
Definition: ContactGeometry.h:1097
Ellipsoid(const Vec3 &radii)
Construct an Ellipsoid given its three principal half-axis dimensions a,b,c (all positive) along the ...
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 ...
static const Ellipsoid & getAs(const ContactGeometry &geo)
Cast the supplied ContactGeometry object to a const Ellipsoid.
Definition: ContactGeometry.h:1189
static bool isInstance(const ContactGeometry &geo)
Return true if the supplied ContactGeometry object is an Ellipsoid.
Definition: ContactGeometry.h:1186
const Vec3 & getCurvatures() const
For efficiency we precalculate the principal curvatures whenever the ellipsoid radii are set; this av...
Impl & updImpl()
Internal use only.
static ContactGeometryTypeId classTypeId()
Obtain the unique id for Ellipsoid contact geometry.
const Impl & getImpl() const
Internal use only.
void setRadii(const Vec3 &radii)
Set the three half-axis dimensions a,b,c (all positive) used to define this ellipsoid,...
Vec3 findPointInSameDirection(const Vec3 &Q) const
Given a direction d defined by the vector Q-O for an arbitrary point in space Q=(x,...
static Ellipsoid & updAs(ContactGeometry &geo)
Cast the supplied ContactGeometry object to a writable Ellipsoid.
Definition: ContactGeometry.h:1192
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 ...
Vec3 findPointWithThisUnitNormal(const UnitVec3 &n) const
Given a unit direction n, find the unique point P on the ellipsoid surface at which the outward-facin...
const Vec3 & getRadii() const
Obtain the three half-axis dimensions a,b,c used to define this ellipsoid.
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 ellip...
This ContactGeometry subclass represents an object that occupies the entire half-space x>0.
Definition: ContactGeometry.h:978
HalfSpace()
Create a HalfSpace for contact, with surface passing through the origin and outward normal -XAxis (in...
Impl & updImpl()
Internal use only.
static ContactGeometryTypeId classTypeId()
Obtain the unique id for HalfSpace contact geometry.
const Impl & getImpl() const
Internal use only.
static const HalfSpace & getAs(const ContactGeometry &geo)
Cast the supplied ContactGeometry object to a const halfspace.
Definition: ContactGeometry.h:994
static bool isInstance(const ContactGeometry &geo)
Return true if the supplied ContactGeometry object is a halfspace.
Definition: ContactGeometry.h:991
UnitVec3 getNormal() const
Return the HalfSpace outward normal in its own frame as a unit vector.
static HalfSpace & updAs(ContactGeometry &geo)
Cast the supplied ContactGeometry object to a writable halfspace.
Definition: ContactGeometry.h:997
This ContactGeometry subclass represents a smooth surface fit through a set of sampled points using b...
Definition: ContactGeometry.h:1219
const Impl & getImpl() const
Internal use only.
const OBBTree & getOBBTree() const
(Advanced) Return a reference to the oriented bounding box tree for this surface.
static ContactGeometryTypeId classTypeId()
Obtain the unique id for SmoothHeightMap contact geometry.
const BicubicSurface & getBicubicSurface() const
Return a reference to the BicubicSurface object being used by this SmoothHeightMap.
Impl & updImpl()
Internal use only.
static bool isInstance(const ContactGeometry &geo)
Return true if the supplied ContactGeometry object is a SmoothHeightMap.
Definition: ContactGeometry.h:1235
static SmoothHeightMap & updAs(ContactGeometry &geo)
Cast the supplied ContactGeometry object to a writable SmoothHeightMap.
Definition: ContactGeometry.h:1241
static const SmoothHeightMap & getAs(const ContactGeometry &geo)
Cast the supplied ContactGeometry object to a const SmoothHeightMap.
Definition: ContactGeometry.h:1238
SmoothHeightMap(const BicubicSurface &surface)
Create a SmoothHeightMap surface using an already-existing BicubicSurface.
This ContactGeometry subclass represents a sphere centered at the origin.
Definition: ContactGeometry.h:1047
Impl & updImpl()
Internal use only.
static Sphere & updAs(ContactGeometry &geo)
Cast the supplied ContactGeometry object to a writable sphere.
Definition: ContactGeometry.h:1060
static bool isInstance(const ContactGeometry &geo)
Return true if the supplied ContactGeometry object is a sphere.
Definition: ContactGeometry.h:1054
const Impl & getImpl() const
Internal use only.
static const Sphere & getAs(const ContactGeometry &geo)
Cast the supplied ContactGeometry object to a const sphere.
Definition: ContactGeometry.h:1057
static ContactGeometryTypeId classTypeId()
Obtain the unique id for Sphere contact geometry.
This ContactGeometry subclass represents a torus centered at the origin with the axial direction alig...
Definition: ContactGeometry.h:1529
void setTubeRadius(Real radius)
const Impl & getImpl() const
Internal use only.
Torus(Real torusRadius, Real tubeRadius)
Impl & updImpl()
Internal use only.
static bool isInstance(const ContactGeometry &geo)
Return true if the supplied ContactGeometry object is a torus.
Definition: ContactGeometry.h:1538
static const Torus & getAs(const ContactGeometry &geo)
Cast the supplied ContactGeometry object to a const torus.
Definition: ContactGeometry.h:1541
static ContactGeometryTypeId classTypeId()
Obtain the unique id for Torus contact geometry.
void setTorusRadius(Real radius)
static Torus & updAs(ContactGeometry &geo)
Cast the supplied ContactGeometry object to a writable torus.
Definition: ContactGeometry.h:1544
This class represents a node in the Oriented Bounding Box Tree for a TriangleMesh.
Definition: ContactGeometry.h:1495
const OBBTreeNode getFirstChildNode() const
Get the first child node.
bool isLeafNode() const
Get whether this is a leaf node.
const OrientedBoundingBox & getBounds() const
Get the OrientedBoundingBox which encloses all triangles in this node or its children.
const Array_< int > & getTriangles() const
Get the indices of all triangles contained in this node.
const OBBTreeNode getSecondChildNode() const
Get the second child node.
int getNumTriangles() const
Get the number of triangles inside this node.
This ContactGeometry subclass represents an arbitrary shape described by a mesh of triangular faces.
Definition: ContactGeometry.h:1316
Vec3 findNearestPoint(const Vec3 &position, bool &inside, int &face, Vec2 &uv) const
Given a point, find the nearest point on the surface of this object.
int getFaceEdge(int face, int edge) const
Get the index of one of the edges of a face.
int getEdgeVertex(int edge, int vertex) const
Get the index of one of the vertices shared by an edge.
Impl & updImpl()
Internal use only.
Real getFaceArea(int face) const
Get the area of a face.
static ContactGeometryTypeId classTypeId()
Obtain the unique id for TriangleMesh contact geometry.
const UnitVec3 & getFaceNormal(int face) const
Get the normal vector for a face.
const Impl & getImpl() const
Internal use only.
static TriangleMesh & updAs(ContactGeometry &geo)
Cast the supplied ContactGeometry object to a writable triangle mesh.
Definition: ContactGeometry.h:1475
const Vec3 & getVertexPosition(int index) const
Get the position of a vertex in the mesh.
int getNumEdges() const
Get the number of edges in the mesh.
OBBTreeNode getOBBTreeNode() const
Get the OBBTreeNode which forms the root of this mesh's Oriented Bounding Box Tree.
Vec3 findCentroid(int face) const
Calculate the location of a face's centroid, that is, the point uv=(1/3,1/3) which is the average of ...
Vec3 findNearestPointToFace(const Vec3 &position, int face, Vec2 &uv) const
Given a point and a face of this object, find the point of the face that is nearest the given point.
UnitVec3 findNormalAtPoint(int face, const Vec2 &uv) const
Calculate the normal vector at a point on the surface.
Vec3 findPoint(int face, const Vec2 &uv) const
Calculate the location of a point on the surface, in the local frame of the TriangleMesh.
static bool isInstance(const ContactGeometry &geo)
Return true if the supplied ContactGeometry object is a triangle mesh.
Definition: ContactGeometry.h:1469
bool intersectsRay(const Vec3 &origin, const UnitVec3 &direction, Real &distance, UnitVec3 &normal) const
Determine whether this mesh intersects a ray, and if so, find the intersection point.
TriangleMesh(const ArrayViewConst_< Vec3 > &vertices, const ArrayViewConst_< int > &faceIndices, bool smooth=false)
Create a TriangleMesh.
Vec3 findNearestPoint(const Vec3 &position, bool &inside, UnitVec3 &normal) const
Given a point, find the nearest point on the surface of this object.
int getEdgeFace(int edge, int face) const
Get the index of one of the faces shared by an edge.
TriangleMesh(const PolygonalMesh &mesh, bool smooth=false)
Create a TriangleMesh based on a PolygonalMesh object.
int getNumVertices() const
Get the number of vertices in the mesh.
int getNumFaces() const
Get the number of faces in the mesh.
int getFaceVertex(int face, int vertex) const
Get the index of one of the vertices of a face.
static const TriangleMesh & getAs(const ContactGeometry &geo)
Cast the supplied ContactGeometry object to a const triangle mesh.
Definition: ContactGeometry.h:1472
bool intersectsRay(const Vec3 &origin, const UnitVec3 &direction, Real &distance, int &face, Vec2 &uv) const
Determine whether this mesh intersects a ray, and if so, find what face it hit.
PolygonalMesh createPolygonalMesh() const
Generate a PolygonalMesh from this TriangleMesh; useful mostly for debugging because you can create a...
void findVertexEdges(int vertex, Array_< int > &edges) const
Find all edges that intersect a vertex.
A ContactGeometry object describes the shape of all or part of the boundary of a solid object,...
Definition: ContactGeometry.h:111
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.
Definition: ContactGeometry.h:852
ContactGeometryImpl * impl
Internal use only.
Definition: ContactGeometry.h:967
bool isConvex() const
Returns true if this surface is known to be convex.
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 t...
UnitVec3 calcSurfaceUnitNormal(const Vec3 &point) const
Calculate the implicit surface outward facing unit normal at the given point.
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...
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 di...
const int getNumGeodesicsShot() const
Get the plane associated with the geodesic hit plane event handler
void getBoundingSphere(Vec3 &center, Real &radius) const
Get a bounding sphere which completely encloses this object.
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 ...
DecorativeGeometry createDecorativeGeometry() const
Generate a DecorativeGeometry that matches the shape of this ContactGeometry.
void setPlane(const Plane &plane) const
Set the plane associated with the geodesic hit plane event handler
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 sensitivi...
const Plane & getPlane() const
Get the plane associated with the geodesic hit plane event handler
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...
bool isSurfaceDefined(const Vec3 &point) const
This function checks to see if a point is within the defined bounds of this particular ContactGeometr...
const Function & getImplicitFunction() const
Our smooth surfaces define a function f(P)=0 that provides an implicit representation of the surface.
bool hasImpl() const
Internal use only.
Definition: ContactGeometry.h:960
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 ...
Vec3 calcSurfaceGradient(const Vec3 &point) const
Calculate the gradient of the implicit surface function, at a given point.
Real calcSurfaceTorsionInDirection(const Vec3 &point, const UnitVec3 &direction) const
For an implicit surface, return the geodesic torsion tau of the surface at a given point p in a given...
~ContactGeometry()
Base class destructor deletes the implementation object. Note that this is not virtual; handles shoul...
Vec3 findNearestPoint(const Vec3 &position, bool &inside, UnitVec3 &normal) const
Given a point, find the nearest point on the surface of this object.
bool isSmooth() const
Returns true if this is a smooth surface, meaning that it can provide meaningful curvature informatio...
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 penet...
NearestPointOnLineResult calcNearestPointOnLineImplicitly(const Vec3 &pointA, const Vec3 &pointB, int maxIterations, Real tolerance, Vec3 &nearestPointOnLine) const
Compute the nearest point on the line spanned by points pA-pB to the surface.
bool isAnalyticFormAvailable() const
Use this function to determine if you can compute a geodesic analytically using shootGeodesicInDirect...
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.
Mat33 calcSurfaceHessian(const Vec3 &point) const
Calculate the hessian of the implicit surface function, at a given point.
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,...
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 f...
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 dir...
NearestPointOnLineResult
Exit flag returned by calcNearestPointOnLineImplicitly.
Definition: ContactGeometry.h:267
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...
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...
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 c...
ContactGeometry(const ContactGeometry &src)
Copy constructor makes a deep copy.
ContactGeometry()
Base class default constructor creates an empty handle.
Definition: ContactGeometry.h:147
Real calcGaussianCurvature(const Vec3 &point) const
This signature is for convenience; use the other one to save time if you already have the gradient an...
Definition: ContactGeometry.h:436
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 i...
bool isEmptyHandle() const
Internal use only.
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 ...
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,...
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,...
const ContactGeometryImpl & getImpl() const
Internal use only.
Definition: ContactGeometry.h:962
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,...
const Geodesic & getGeodP() const
Get the geodesic for access by visualizer.
bool isOwnerHandle() const
Internal use only.
void shootGeodesicInDirectionImplicitly(const Vec3 &initialPointApprox, const Vec3 &initialTangentApprox, Real finalArcLength, Real initialIntegratorStepSize, Real integratorAccuracy, Real constraintTolerance, int maxIterations, const std::function< void(const ContactGeometry::GeodesicKnotPoint &)> &geodesicKnotPointsSink) const
Given an approximate initial point and an approximate initial tangent direction compute the geodesic ...
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 d...
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,...
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 directio...
ContactGeometry(ContactGeometryImpl *impl)
Internal use only.
Real calcSurfaceValue(const Vec3 &point) const
Calculate the value of the implicit surface function, at a given point.
void addVizReporter(ScheduledEventReporter *reporter) const
Get the plane associated with the geodesic hit plane event handler
ContactGeometryTypeId getTypeId() const
ContactTrackerSubsystem uses this id for fast identification of specific surface shapes.
const Geodesic & getGeodQ() const
Get the geodesic for access by visualizer.
ContactGeometryImpl & updImpl()
Internal use only.
Definition: ContactGeometry.h:964
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.
void shootGeodesicInDirectionAnalytically(const Vec3 &initialPointApprox, const Vec3 &initialTangentApprox, Real finalArcLength, int numberOfKnotPoints, const std::function< void(const ContactGeometry::GeodesicKnotPoint &)> &geodesicKnotPointsSink) const
Given an approximate initial point and an approximate initial tangent direction compute the geodesic ...
ContactGeometry & operator=(const ContactGeometry &src)
Copy assignment makes a deep copy.
Definition: CoordinateAxis.h:194
A DecorationGenerator is used to define geometry that may change over the course of a simulation.
Definition: DecorationGenerator.h:45
This defines a rectangular solid centered at the origin and aligned with the local frame axes.
Definition: DecorativeGeometry.h:426
This is the client-side interface to an implementation-independent representation of "Decorations" su...
Definition: DecorativeGeometry.h:86
A line between two points.
Definition: DecorativeGeometry.h:306
This defines a sphere centered at the origin.
Definition: DecorativeGeometry.h:367
A 3d rectangular box aligned with an unspecified frame F and centered at that frame's origin.
Definition: Geo_Box.h:48
A event handler to terminate integration when geodesic hits the plane.
Definition: ContactGeometry.h:1594
GeodHitPlaneEvent()
Definition: ContactGeometry.h:1596
const bool isEnabled()
Definition: ContactGeometry.h:1640
void setPlane(const Plane &aplane) const
Definition: ContactGeometry.h:1628
const Plane & getPlane() const
Definition: ContactGeometry.h:1632
GeodHitPlaneEvent(const Plane &aplane)
Definition: ContactGeometry.h:1599
const void setEnabled(bool enabledFlag)
Definition: ContactGeometry.h:1636
Real getValue(const State &state) const override
Get the value of the event trigger function for a State.
Definition: ContactGeometry.h:1605
void handleEvent(State &state, Real accuracy, bool &shouldTerminate) const override
This method is invoked to handle the event.
Definition: ContactGeometry.h:1616
This class stores options for calculating geodesics.
Definition: Geodesic.h:311
This class stores a geodesic curve after it has been determined.
Definition: Geodesic.h:51
This class represents a small matrix whose size is known at compile time, containing elements of any ...
Definition: Mat.h:97
TODO.
Definition: OBBTree.h:100
This class represents a rectangular box with arbitrary position and orientation.
Definition: OrientedBoundingBox.h:42
This class generates decoration for contact points and straight line path segments.
Definition: ContactGeometry.h:1653
virtual void generateDecorations(const State &state, Array_< DecorativeGeometry > &geometry) override
This will be called every time a new State is about to be visualized.
Definition: ContactGeometry.h:1658
PathDecorator(const Vector &x, const Vec3 &O, const Vec3 &I, const Vec3 &color)
Definition: ContactGeometry.h:1655
This class generates decoration for a plane.
Definition: ContactGeometry.h:1693
virtual void generateDecorations(const State &state, Array_< DecorativeGeometry > &geometry) override
This will be called every time a new State is about to be visualized.
Definition: ContactGeometry.h:1698
PlaneDecorator(const Plane &plane, const Vec3 &color)
Definition: ContactGeometry.h:1695
A simple plane class.
Definition: ContactGeometry.h:1566
Real getOffset() const
Definition: ContactGeometry.h:1580
Vec3 getNormal() const
Definition: ContactGeometry.h:1576
Real getDistance(const Vec3 &pt) const
Definition: ContactGeometry.h:1572
Plane(const Vec3 &normal, const Real &offset)
Definition: ContactGeometry.h:1569
Plane()
Definition: ContactGeometry.h:1568
This class provides a description of a mesh made of polygonal faces (not limited to triangles).
Definition: PolygonalMesh.h:79
Rotation_ & setRotationFromOneAxis(const UnitVec3P &uvec, CoordinateAxis axis)
Calculate R_AB by knowing one of B's unit vectors expressed in A.
ScheduledEventReporter is a subclass of EventReporter for events that occur at a particular time that...
Definition: EventReporter.h:72
This class is basically a glorified enumerated type, type-safe and range checked but permitting conve...
Definition: Stage.h:66
This object is intended to contain all state information for a SimTK::System, except topological info...
Definition: State.h:280
const Vector & getQ(SubsystemIndex) const
Per-subsystem access to the global shared variables.
TriggeredEventHandler is a subclass of EventHandler for events that occur when some condition is sati...
Definition: EventHandler.h:109
UnitVec< P, 1 > perp() const
Return a new unit vector perpendicular to this one but otherwise arbitrary.
Definition: UnitVec.h:181
CNT< ScalarNormSq >::TSqrt norm() const
Definition: Vec.h:610
Vec< 2 > Vec2
This is the most common 2D vector type: a column of 2 Real values stored consecutively in memory (pac...
Definition: SmallMatrix.h:126
Vec< 3 > Vec3
This is the most common 3D vector type: a column of 3 Real values stored consecutively in memory (pac...
Definition: SmallMatrix.h:129
const Real NaN
This is the IEEE "not a number" constant for this implementation of the default-precision Real type; ...
const Complex I
We only need one complex constant, i = sqrt(-1). For the rest just multiply the real constant by i,...
This is the top-level SimTK namespace into which all SimTK names are placed to avoid collision with o...
Definition: Assembler.h:37
const Vec3 Black
RGB=( 0, 0, 0)
UnitVec< Real, 1 > UnitVec3
Definition: UnitVec.h:41
SimTK_Real Real
This is the default compiled-in floating point type for SimTK, either float or double.
Definition: SimTKcommon/include/SimTKcommon/internal/common.h:607
Transform_< Real > Transform
Definition: Transform.h:44
SimTK_DEFINE_UNIQUE_INDEX_TYPE(AssemblyConditionIndex)
The state of a geodesic at a (knot) point along the geodesic.
Definition: ContactGeometry.h:126