DistRectangle3AlignedBox3.h 6.1 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141
  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/DCPQuery.h>
  9. #include <Mathematics/LCPSolver.h>
  10. #include <Mathematics/AlignedBox.h>
  11. #include <Mathematics/Rectangle.h>
  12. #include <Mathematics/Vector3.h>
  13. // Compute the distance between a rectangle and an aligned box in 3D. The
  14. // algorithm is based on using an LCP solver for the convex quadratic
  15. // programming problem. For details, see
  16. // https://www.geometrictools.com/Documentation/ConvexQuadraticProgramming.pdf
  17. namespace WwiseGTE
  18. {
  19. template <typename Real>
  20. class DCPQuery<Real, Rectangle3<Real>, AlignedBox3<Real>>
  21. {
  22. public:
  23. struct Result
  24. {
  25. bool queryIsSuccessful;
  26. // These members are valid only when queryIsSuccessful is true;
  27. // otherwise, they are all set to zero.
  28. Real distance, sqrDistance;
  29. std::array<Real, 2> rectangleParameter;
  30. std::array<Real, 3> boxParameter;
  31. Vector3<Real> closestPoint[2];
  32. // The number of iterations used by LCPSolver regardless of
  33. // whether the query is successful.
  34. int numLCPIterations;
  35. };
  36. // The default maximum iterations is 81 (n = 9, maxIterations = n*n).
  37. // If the solver fails to converge, try increasing the maximum number
  38. // of iterations.
  39. void SetMaxLCPIterations(int maxLCPIterations)
  40. {
  41. mLCP.SetMaxIterations(maxLCPIterations);
  42. }
  43. Result operator()(Rectangle3<Real> const& rectangle, AlignedBox3<Real> const& box)
  44. {
  45. Result result;
  46. // Translate the rectangle and aligned box so that the aligned
  47. // box becomes a canonical box.
  48. Vector3<Real> K = box.max - box.min;
  49. Vector3<Real> V = rectangle.center - box.min;
  50. // Convert the oriented rectangle to a regular one (origin at a
  51. // corner).
  52. Vector3<Real> scaledE0 = rectangle.axis[0] * rectangle.extent[0];
  53. Vector3<Real> scaledE1 = rectangle.axis[1] * rectangle.extent[1];
  54. Vector3<Real> E0 = scaledE0 * (Real)2;
  55. Vector3<Real> E1 = scaledE1 * (Real)2;
  56. V -= scaledE0 + scaledE1;
  57. // Compute quantities to initialize q and M in the LCP.
  58. Real dotVE0 = Dot(V, E0);
  59. Real dotVE1 = Dot(V, E1);
  60. Real dotE0E0 = Dot(E0, E0);
  61. Real dotE1E1 = Dot(E1, E1);
  62. // The LCP has 5 variables and 5 (nontrivial) inequality
  63. // constraints.
  64. std::array<Real, 10> q =
  65. {
  66. -V[0], -V[1], -V[2], dotVE0, dotVE1, K[0], K[1], K[2], (Real)1, (Real)1
  67. };
  68. std::array<std::array<Real, 10>, 10> M;
  69. M[0] = { (Real)1, (Real)0, (Real)0, -E0[0], -E1[0], (Real)1, (Real)0, (Real)0, (Real)0, (Real)0 };
  70. M[1] = { (Real)0, (Real)1, (Real)0, -E0[1], -E1[1], (Real)0, (Real)1, (Real)0, (Real)0, (Real)0 };
  71. M[2] = { (Real)0, (Real)0, (Real)1, -E0[2], -E1[2], (Real)0, (Real)0, (Real)1, (Real)0 , (Real)0 };
  72. M[3] = { -E0[0], -E0[1], -E0[2], dotE0E0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)1, (Real)0 };
  73. M[4] = { -E1[0], -E1[1], -E1[2], (Real)0, dotE1E1, (Real)0, (Real)0, (Real)0, (Real)0, (Real)1 };
  74. M[5] = { (Real)-1, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0 };
  75. M[6] = { (Real)0, (Real)-1, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0 };
  76. M[7] = { (Real)0, (Real)0, (Real)-1, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0 };
  77. M[8] = { (Real)0, (Real)0, (Real)0, (Real)-1, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0 };
  78. M[9] = { (Real)0, (Real)0, (Real)0, (Real)0, (Real)-1, (Real)0, (Real)0, (Real)0, (Real)0, (Real)0 };
  79. std::array<Real, 10> w, z;
  80. if (mLCP.Solve(q, M, w, z))
  81. {
  82. result.queryIsSuccessful = true;
  83. Real t0 = (z[3] * (Real)2 - (Real)1) * rectangle.extent[0];
  84. Real t1 = (z[4] * (Real)2 - (Real)1) * rectangle.extent[1];
  85. result.rectangleParameter[0] = t0;
  86. result.rectangleParameter[1] = t1;
  87. result.closestPoint[0] = rectangle.center + t0 * rectangle.axis[0] + t1 * rectangle.axis[1];
  88. for (int i = 0; i < 3; ++i)
  89. {
  90. result.boxParameter[i] = z[i] + box.min[i];
  91. result.closestPoint[1][i] = result.boxParameter[i];
  92. }
  93. Vector3<Real> diff = result.closestPoint[1] - result.closestPoint[0];
  94. result.sqrDistance = Dot(diff, diff);
  95. result.distance = std::sqrt(result.sqrDistance);
  96. }
  97. else
  98. {
  99. // If you reach this case, the maximum number of iterations
  100. // was not specified to be large enough or there is a problem
  101. // due to floating-point rounding errors. If you believe the
  102. // latter is true, file a bug report.
  103. result.queryIsSuccessful = false;
  104. for (int i = 0; i < 2; ++i)
  105. {
  106. result.rectangleParameter[i] = (Real)0;
  107. }
  108. for (int i = 0; i < 3; ++i)
  109. {
  110. result.boxParameter[i] = (Real)0;
  111. result.closestPoint[0][i] = (Real)0;
  112. result.closestPoint[1][i] = (Real)0;
  113. }
  114. result.distance = (Real)0;
  115. result.sqrDistance = (Real)0;
  116. }
  117. result.numLCPIterations = mLCP.GetNumIterations();
  118. return result;
  119. }
  120. private:
  121. LCPSolver<Real, 10> mLCP;
  122. };
  123. }