CosEstimate.h 4.6 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139
  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/Math.h>
  9. // Minimax polynomial approximations to cos(x). The polynomial p(x) of
  10. // degree D has only even-power terms, is required to have constant term 1,
  11. // and p(pi/2) = cos(pi/2) = 0. It minimizes the quantity
  12. // maximum{|cos(x) - p(x)| : x in [-pi/2,pi/2]} over all polynomials of
  13. // degree D subject to the constraints mentioned.
  14. namespace WwiseGTE
  15. {
  16. template <typename Real>
  17. class CosEstimate
  18. {
  19. public:
  20. // The input constraint is x in [-pi/2,pi/2]. For example,
  21. // float x; // in [-pi/2,pi/2]
  22. // float result = CosEstimate<float>::Degree<4>(x);
  23. template <int D>
  24. inline static Real Degree(Real x)
  25. {
  26. return Evaluate(degree<D>(), x);
  27. }
  28. // The input x can be any real number. Range reduction is used to
  29. // generate a value y in [-pi/2,pi/2] and a sign s for which
  30. // cos(y) = s*cos(x). For example,
  31. // float x; // x any real number
  32. // float result = CosEstimate<float>::DegreeRR<3>(x);
  33. template <int D>
  34. inline static Real DegreeRR(Real x)
  35. {
  36. Real y, sign;
  37. Reduce(x, y, sign);
  38. Real poly = sign * Degree<D>(y);
  39. return poly;
  40. }
  41. private:
  42. // Metaprogramming and private implementation to allow specialization
  43. // of a template member function.
  44. template <int D> struct degree {};
  45. inline static Real Evaluate(degree<2>, Real x)
  46. {
  47. Real xsqr = x * x;
  48. Real poly;
  49. poly = (Real)GTE_C_COS_DEG2_C1;
  50. poly = (Real)GTE_C_COS_DEG2_C0 + poly * xsqr;
  51. return poly;
  52. }
  53. inline static Real Evaluate(degree<4>, Real x)
  54. {
  55. Real xsqr = x * x;
  56. Real poly;
  57. poly = (Real)GTE_C_COS_DEG4_C2;
  58. poly = (Real)GTE_C_COS_DEG4_C1 + poly * xsqr;
  59. poly = (Real)GTE_C_COS_DEG4_C0 + poly * xsqr;
  60. return poly;
  61. }
  62. inline static Real Evaluate(degree<6>, Real x)
  63. {
  64. Real xsqr = x * x;
  65. Real poly;
  66. poly = (Real)GTE_C_COS_DEG6_C3;
  67. poly = (Real)GTE_C_COS_DEG6_C2 + poly * xsqr;
  68. poly = (Real)GTE_C_COS_DEG6_C1 + poly * xsqr;
  69. poly = (Real)GTE_C_COS_DEG6_C0 + poly * xsqr;
  70. return poly;
  71. }
  72. inline static Real Evaluate(degree<8>, Real x)
  73. {
  74. Real xsqr = x * x;
  75. Real poly;
  76. poly = (Real)GTE_C_COS_DEG8_C4;
  77. poly = (Real)GTE_C_COS_DEG8_C3 + poly * xsqr;
  78. poly = (Real)GTE_C_COS_DEG8_C2 + poly * xsqr;
  79. poly = (Real)GTE_C_COS_DEG8_C1 + poly * xsqr;
  80. poly = (Real)GTE_C_COS_DEG8_C0 + poly * xsqr;
  81. return poly;
  82. }
  83. inline static Real Evaluate(degree<10>, Real x)
  84. {
  85. Real xsqr = x * x;
  86. Real poly;
  87. poly = (Real)GTE_C_COS_DEG10_C5;
  88. poly = (Real)GTE_C_COS_DEG10_C4 + poly * xsqr;
  89. poly = (Real)GTE_C_COS_DEG10_C3 + poly * xsqr;
  90. poly = (Real)GTE_C_COS_DEG10_C2 + poly * xsqr;
  91. poly = (Real)GTE_C_COS_DEG10_C1 + poly * xsqr;
  92. poly = (Real)GTE_C_COS_DEG10_C0 + poly * xsqr;
  93. return poly;
  94. }
  95. // Support for range reduction.
  96. inline static void Reduce(Real x, Real& y, Real& sign)
  97. {
  98. // Map x to y in [-pi,pi], x = 2*pi*quotient + remainder.
  99. Real quotient = (Real)GTE_C_INV_TWO_PI * x;
  100. if (x >= (Real)0)
  101. {
  102. quotient = (Real)((int)(quotient + (Real)0.5));
  103. }
  104. else
  105. {
  106. quotient = (Real)((int)(quotient - (Real)0.5));
  107. }
  108. y = x - (Real)GTE_C_TWO_PI * quotient;
  109. // Map y to [-pi/2,pi/2] with cos(y) = sign*cos(x).
  110. if (y > (Real)GTE_C_HALF_PI)
  111. {
  112. y = (Real)GTE_C_PI - y;
  113. sign = (Real)-1;
  114. }
  115. else if (y < (Real)-GTE_C_HALF_PI)
  116. {
  117. y = (Real)-GTE_C_PI - y;
  118. sign = (Real)-1;
  119. }
  120. else
  121. {
  122. sign = (Real)1;
  123. }
  124. }
  125. };
  126. }