1+ package us .ihmc .robotics .kinematics ;
2+
3+ import us .ihmc .euclid .matrix .interfaces .RotationMatrixBasics ;
4+ import us .ihmc .euclid .rotationConversion .RotationMatrixConversion ;
5+ import us .ihmc .euclid .tools .EuclidCoreTools ;
6+
7+ public class TrigonometricApproximation
8+ {
9+ public static double approximateArgumentForAddedArcCosAngles (double cosTheta1 , double cosTheta2 )
10+ {
11+ return cosTheta1 * cosTheta2 - Math .sqrt ((1.0 - cosTheta1 * cosTheta1 ) * (1.0 - cosTheta2 * cosTheta2 ));
12+ }
13+
14+ public static double approximateArgumentForSubtractedArcCosAngles (double cosTheta1 , double cosTheta2 )
15+ {
16+ // Warning: If solution is expected to be negative, it'll be out of the range of arccos()
17+ if (cosTheta1 > cosTheta2 )
18+ return Double .NaN ;
19+ return cosTheta1 * cosTheta2 + Math .sqrt ((1.0 - cosTheta1 * cosTheta1 ) * (1.0 - cosTheta2 * cosTheta2 ));
20+ }
21+
22+ public static double aTan2ApproximationOfArcCos (double cosTheta )
23+ {
24+ return 2.0 * Math .atan2 (Math .sqrt (1 - cosTheta * cosTheta ), 1 + cosTheta );
25+ }
26+
27+ public static double taylorArcTan2ApproximationOfArcCos (double cosTheta , int seriesDepth )
28+ {
29+ return 2.0 * taylorApproximationOfArcTan (Math .sqrt (1 - cosTheta * cosTheta ) / (1 + cosTheta ), seriesDepth );
30+ }
31+
32+ public static double taylorApproximationOfArcCos (double x , int seriesDepth )
33+ {
34+ if (seriesDepth < 1 )
35+ throw new IllegalArgumentException ("Depth must be at least 1" );
36+ double value = Math .PI / 2.0 - x ;
37+ int topFactorial = 1 ;
38+ int bottomFactorial = 1 ;
39+ int bottomConstExponential = 1 ;
40+ double valueExponential = x ;
41+ double xSquared = x * x ;
42+ int rhsDenominator = 1 ;
43+ for (int n = 1 ; n <= seriesDepth ; n ++)
44+ {
45+ topFactorial *= (2 * n ) * (2 * n - 1 );
46+ bottomFactorial *= n ;
47+ bottomConstExponential *= 4 ;
48+ rhsDenominator += 2 ;
49+ valueExponential *= xSquared ;
50+ double factorial = topFactorial / ((double ) (bottomFactorial * bottomFactorial ));
51+ double fullDenominator = bottomConstExponential * rhsDenominator ;
52+ value -= factorial * valueExponential / fullDenominator ;
53+ }
54+ return value ;
55+ }
56+
57+ /**
58+ * https://proofwiki.org/wiki/Power_Series_Expansion_for_Real_Arctangent_Function
59+ */
60+ public static double taylorApproximationOfArcTan (double x , int seriesDepth )
61+ {
62+ if (seriesDepth < 1 )
63+ throw new IllegalArgumentException ("Depth must be at least 1" );
64+ double xSquared = x * x ;
65+ double numerator = x ;
66+ double ret = x ;
67+ int denominator = 1 ;
68+ for (int i = 1 ; i <= seriesDepth ; i ++)
69+ {
70+ numerator *= -xSquared ;
71+ denominator += 2 ;
72+ ret += numerator / denominator ;
73+ }
74+ return ret ;
75+ }
76+
77+ public static double taylorApproximationOfCos (double x , int seriesDepth )
78+ {
79+ double xSquared = x * x ;
80+ int denominator = 1 ;
81+ double numerator = 1.0 ;
82+ double ret = 1.0 ;
83+ int doubleI = 0 ;
84+ for (int i = 1 ; i <= seriesDepth ; i ++)
85+ {
86+ numerator *= -xSquared ;
87+ doubleI += 2 ;
88+ denominator *= doubleI * (doubleI - 1 );
89+ ret += numerator / denominator ;
90+ }
91+ return ret ;
92+ }
93+
94+ public static double taylorApproximationOfSin (double x , int seriesDepth )
95+ {
96+ double xSquared = x * x ;
97+ int denominator = 1 ;
98+ double numerator = x ;
99+ double ret = x ;
100+ int doubleI = 1 ;
101+ for (int i = 1 ; i <= seriesDepth ; i ++)
102+ {
103+ numerator *= -xSquared ;
104+ doubleI += 2 ;
105+ denominator *= doubleI * (doubleI - 1 );
106+ ret += numerator / denominator ;
107+ }
108+ return ret ;
109+ }
110+
111+ public static void computePitchMatrix (double pitch , RotationMatrixBasics matrixToPack , int depth )
112+ {
113+ if (EuclidCoreTools .isAngleZero (pitch , RotationMatrixConversion .EPS ))
114+ {
115+ matrixToPack .setToZero ();
116+ }
117+ else
118+ {
119+ double sinPitch = TrigonometricApproximation .taylorApproximationOfSin (pitch , depth );
120+ double cosPitch = TrigonometricApproximation .taylorApproximationOfCos (pitch , depth );
121+ matrixToPack .setUnsafe (cosPitch , 0.0 , sinPitch , 0.0 , 1.0 , 0.0 , -sinPitch , 0.0 , cosPitch );
122+ }
123+ }
124+
125+ public static void computeRollMatrix (double roll , RotationMatrixBasics matrixToPack , int depth )
126+ {
127+ if (EuclidCoreTools .isAngleZero (roll , RotationMatrixConversion .EPS ))
128+ {
129+ matrixToPack .setToZero ();
130+ }
131+ else
132+ {
133+ double sinRoll = TrigonometricApproximation .taylorApproximationOfSin (roll , depth );
134+ double cosRoll = TrigonometricApproximation .taylorApproximationOfCos (roll , depth );
135+ matrixToPack .setUnsafe (1.0 , 0.0 , 0.0 , 0.0 , cosRoll , -sinRoll , 0.0 , sinRoll , cosRoll );
136+ }
137+ }
138+ }
0 commit comments