Torus3.h 7.8 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230
  1. // David Eberly, Geometric Tools, Redmond WA 98052
  2. // Copyright (c) 1998-2020
  3. // Distributed under the Boost Software License, Version 1.0.
  4. // https://www.boost.org/LICENSE_1_0.txt
  5. // https://www.geometrictools.com/License/Boost/LICENSE_1_0.txt
  6. // Version: 4.0.2019.08.13
  7. #pragma once
  8. #include <Mathematics/Vector3.h>
  9. // A torus with origin (0,0,0), outer radius r0 and inner radius r1 (with
  10. // (r0 >= r1) is defined implicitly as follows. The point P0 = (x,y,z) is on
  11. // the torus. Its projection onto the xy-plane is P1 = (x,y,0). The circular
  12. // cross section of the torus that contains the projection has radius r0 and
  13. // center P2 = r0*(x,y,0)/sqrt(x^2+y^2). The points triangle <P0,P1,P2> is a
  14. // right triangle with right angle at P1. The hypotenuse <P0,P2> has length
  15. // r1, leg <P1,P2> has length z and leg <P0,P1> has length
  16. // |r0 - sqrt(x^2+y^2)|. The Pythagorean theorem says
  17. // z^2 + |r0 - sqrt(x^2+y^2)|^2 = r1^2. This can be algebraically
  18. // manipulated to
  19. // (x^2 + y^2 + z^2 + r0^2 - r1^2)^2 - 4 * r0^2 * (x^2 + y^2) = 0
  20. //
  21. // A parametric form is
  22. // x = (r0 + r1 * cos(v)) * cos(u)
  23. // y = (r0 + r1 * cos(v)) * sin(u)
  24. // z = r1 * sin(v)
  25. // for u in [0,2*pi) and v in [0,2*pi).
  26. //
  27. // Generally, let the torus center be C with plane of symmetry containing C
  28. // and having directions D0 and D1. The axis of symmetry is the line
  29. // containing C and having direction N (the plane normal). The radius from
  30. // the center of the torus is r0 and the radius of the tube of the torus is
  31. // r1. A point P may be written as P = C + x*D0 + y*D1 + z*N, where matrix
  32. // [D0 D1 N] is orthonormal and has determinant 1. Thus, x = Dot(D0,P-C),
  33. // y = Dot(D1,P-C) and z = Dot(N,P-C). The implicit form is
  34. // [|P-C|^2 + r0^2 - r1^2]^2 - 4*r0^2*[|P-C|^2 - (Dot(N,P-C))^2] = 0
  35. // Observe that D0 and D1 are not present in the equation, which is to be
  36. // expected by the symmetry. The parametric form is
  37. // P(u,v) = C + (r0 + r1*cos(v))*(cos(u)*D0 + sin(u)*D1) + r1*sin(v)*N
  38. // for u in [0,2*pi) and v in [0,2*pi).
  39. //
  40. // In the class Torus3, the members are 'center' C, 'direction0' D0,
  41. // 'direction1' D1, 'normal' N, 'radius0' r0 and 'radius1' r1.
  42. namespace WwiseGTE
  43. {
  44. template <typename Real>
  45. class Torus3
  46. {
  47. public:
  48. // Construction and destruction. The default constructor sets center
  49. // to (0,0,0), direction0 to (1,0,0), direction1 to (0,1,0), normal
  50. // to (0,0,1), radius0 to 2 and radius1 to 1.
  51. Torus3()
  52. :
  53. center(Vector3<Real>::Zero()),
  54. direction0(Vector3<Real>::Unit(0)),
  55. direction1(Vector3<Real>::Unit(1)),
  56. normal(Vector3<Real>::Unit(2)),
  57. radius0((Real)2),
  58. radius1((Real)1)
  59. {
  60. }
  61. Torus3(Vector3<Real> const& inCenter, Vector3<Real> const& inDirection0,
  62. Vector3<Real> const& inDirection1, Vector3<Real> const& inNormal,
  63. Real inRadius0, Real inRadius1)
  64. :
  65. center(inCenter),
  66. direction0(inDirection0),
  67. direction1(inDirection1),
  68. normal(inNormal),
  69. radius0(inRadius0),
  70. radius1(inRadius1)
  71. {
  72. }
  73. // Evaluation of the surface. The function supports derivative
  74. // calculation through order 2; that is, maxOrder <= 2 is required.
  75. // If you want only the position, pass in maxOrder of 0. If you want
  76. // the position and first-order derivatives, pass in maxOrder of 1,
  77. // and so on. The output 'values' are ordered as: position X;
  78. // first-order derivatives dX/du, dX/dv; second-order derivatives
  79. // d2X/du2, d2X/dudv, d2X/dv2. The input array 'jet' must have enough
  80. // storage for the specified order.
  81. void Evaluate(Real u, Real v, unsigned int maxOrder, Vector3<Real>* jet) const
  82. {
  83. // Compute position.
  84. Real csu = std::cos(u);
  85. Real snu = std::sin(u);
  86. Real csv = std::cos(v);
  87. Real snv = std::sin(v);
  88. Real r1csv = radius1 * csv;
  89. Real r1snv = radius1 * snv;
  90. Real r0pr1csv = radius0 + r1csv;
  91. Vector3<Real> combo0 = csu * direction0 + snu * direction1;
  92. Vector3<Real> r0pr1csvcombo0 = r0pr1csv * combo0;
  93. Vector3<Real> r1snvnormal = r1snv * normal;
  94. jet[0] = center + r0pr1csvcombo0 + r1snvnormal;
  95. if (maxOrder >= 1)
  96. {
  97. // Compute first-order derivatives.
  98. Vector3<Real> combo1 = -snu * direction0 + csu * direction1;
  99. jet[1] = r0pr1csv * combo1;
  100. jet[2] = -r1snv * combo0 + r1csv * normal;
  101. if (maxOrder == 2)
  102. {
  103. // Compute second-order derivatives.
  104. jet[3] = -r0pr1csvcombo0;
  105. jet[4] = -r1snv * combo1;
  106. jet[5] = -r1csv * combo0 - r1snvnormal;
  107. }
  108. }
  109. }
  110. // Reverse lookup of parameters from position.
  111. void GetParameters(Vector3<Real> const& X, Real& u, Real& v) const
  112. {
  113. Vector3<Real> delta = X - center;
  114. // (r0 + r1*cos(v))*cos(u)
  115. Real dot0 = Dot(direction0, delta);
  116. // (r0 + r1*cos(v))*sin(u)
  117. Real dot1 = Dot(direction1, delta);
  118. // r1*sin(v)
  119. Real dot2 = Dot(normal, delta);
  120. // r1*cos(v)
  121. Real r1csv = std::sqrt(dot0 * dot0 + dot1 * dot1) - radius0;
  122. u = std::atan2(dot1, dot0);
  123. v = std::atan2(dot2, r1csv);
  124. }
  125. Vector3<Real> center, direction0, direction1, normal;
  126. Real radius0, radius1;
  127. public:
  128. // Comparisons to support sorted containers.
  129. bool operator==(Torus3 const& torus) const
  130. {
  131. return center == torus.center
  132. && direction0 == torus.direction0
  133. && direction1 == torus.direction1
  134. && normal == torus.normal
  135. && radius0 == torus.radius0
  136. && radius1 == torus.radius1;
  137. }
  138. bool operator!=(Torus3 const& torus) const
  139. {
  140. return !operator==(torus);
  141. }
  142. bool operator< (Torus3 const& torus) const
  143. {
  144. if (center < torus.center)
  145. {
  146. return true;
  147. }
  148. if (center > torus.center)
  149. {
  150. return false;
  151. }
  152. if (direction0 < torus.direction0)
  153. {
  154. return true;
  155. }
  156. if (direction0 > torus.direction0)
  157. {
  158. return false;
  159. }
  160. if (direction1 < torus.direction1)
  161. {
  162. return true;
  163. }
  164. if (direction1 > torus.direction1)
  165. {
  166. return false;
  167. }
  168. if (normal < torus.normal)
  169. {
  170. return true;
  171. }
  172. if (normal > torus.normal)
  173. {
  174. return false;
  175. }
  176. if (radius0 < torus.radius0)
  177. {
  178. return true;
  179. }
  180. if (radius0 > torus.radius0)
  181. {
  182. return false;
  183. }
  184. return radius1 < torus.radius1;
  185. }
  186. bool operator<=(Torus3 const& torus) const
  187. {
  188. return !torus.operator<(*this);
  189. }
  190. bool operator> (Torus3 const& torus) const
  191. {
  192. return torus.operator<(*this);
  193. }
  194. bool operator>=(Torus3 const& torus) const
  195. {
  196. return !operator<(torus);
  197. }
  198. };
  199. }