9105269ff32ebbc00d9c54e832557943412cedd4
[xonotic/xonotic-data.pk3dir.git] / qcsrc / lib / warpzone / mathlib.qc
1 #include "mathlib.qh"
2 #if defined(CSQC)
3 #elif defined(MENUQC)
4 #elif defined(SVQC)
5 #endif
6
7 int fpclassify(float e)
8 {
9         if(isnan(e))
10                 return FP_NAN;
11         if(isinf(e))
12                 return FP_INFINITE;
13         if(e == 0)
14                 return FP_ZERO;
15         return FP_NORMAL;
16 }
17 bool isfinite(float e)
18 {
19         return !(isnan(e) || isinf(e));
20 }
21 bool isinf(float e)
22 {
23         return (e != 0) && (e + e == e);
24 }
25 bool isnan(float e)
26 {
27         // the sane way to detect NaN is broken because of a compiler bug
28         // (works with constants but breaks when assigned to variables)
29         // use conversion to string instead
30
31         //float f = e;
32         //return (e != f);
33         return ftos(e) == "-nan";
34 }
35 bool isnormal(float e)
36 {
37         return isfinite(e);
38 }
39 bool signbit(float e)
40 {
41         return (e < 0);
42 }
43
44 float acosh(float e)
45 {
46         return log(e + sqrt(e*e - 1));
47 }
48 float asinh(float e)
49 {
50         return log(e + sqrt(e*e + 1));
51 }
52 float atanh(float e)
53 {
54         return 0.5 * log((1+e) / (1-e));
55 }
56 float cosh(float e)
57 {
58         return 0.5 * (exp(e) + exp(-e));
59 }
60 float sinh(float e)
61 {
62         return 0.5 * (exp(e) - exp(-e));
63 }
64 float tanh(float e)
65 {
66         return sinh(e) / cosh(e);
67 }
68
69 float exp(float e)
70 {
71         return pow(M_E, e);
72 }
73 float exp2(float e)
74 {
75         return pow(2, e);
76 }
77 float expm1(float e)
78 {
79         return exp(e) - 1;
80 }
81
82 vector frexp(float e)
83 {
84         vector v;
85         v.z = 0;
86         v.y = ilogb(e) + 1;
87         v.x = e / pow(2, v.y);
88         return v;
89 }
90 int ilogb(float e)
91 {
92         return floor(log2(fabs(e)));
93 }
94 float ldexp(float x, int e)
95 {
96         return x * pow(2, e);
97 }
98 float logn(float e, float base)
99 {
100         return log(e) / log(base);
101 }
102 float log10(float e)
103 {
104         return log(e) * M_LOG10E;
105 }
106 float log1p(float e)
107 {
108         return log(e + 1);
109 }
110 float log2(float e)
111 {
112         return log(e) * M_LOG2E;
113 }
114 float logb(float e)
115 {
116         return floor(log2(fabs(e)));
117 }
118 vector modf(float f)
119 {
120         return '1 0 0' * (f - trunc(f)) + '0 1 0' * trunc(f);
121 }
122
123 float scalbn(float e, int n)
124 {
125         return e * pow(2, n);
126 }
127
128 float cbrt(float e)
129 {
130         return copysign(pow(fabs(e), (1.0/3.0)), e);
131 }
132 float hypot(float e, float f)
133 {
134         return sqrt(e*e + f*f);
135 }
136
137 float erf(float e)
138 {
139         // approximation taken from wikipedia
140         float f;
141         f = e*e;
142         return copysign(sqrt(1 - exp(-f * (1.273239544735163 + 0.14001228868667 * f) / (1 + 0.14001228868667 * f))), e);
143 }
144 float erfc(float e)
145 {
146         return 1.0 - erf(e);
147 }
148 vector lgamma(float e)
149 {
150         // TODO improve accuracy
151         if(!isfinite(e))
152                 return fabs(e) * '1 0 0' + copysign(1, e) * '0 1 0';
153         if(e < 1 && e == floor(e))
154                 return nan("gamma") * '1 1 1';
155         if(e < 0.1)
156         {
157                 vector v;
158                 v = lgamma(1.0 - e);
159                 // reflection formula:
160                 // gamma(1-z) * gamma(z) = pi / sin(pi*z)
161                 // lgamma(1-z) + lgamma(z) = log(pi) - log(sin(pi*z))
162                 // sign of gamma(1-z) = sign of gamma(z) * sign of sin(pi*z)
163                 v.z = sin(M_PI * e);
164                 v.x = log(M_PI) - log(fabs(v.z)) - v.x;
165                 if(v.z < 0)
166                         v.y = -v.y;
167                 v.z = 0;
168                 return v;
169         }
170         if(e < 1.1)
171                 return lgamma(e + 1) - log(e) * '1 0 0';
172         e -= 1;
173         return (0.5 * log(2 * M_PI * e) + e * (log(e) - 1)) * '1 0 0' + '0 1 0';
174 }
175 float tgamma(float e)
176 {
177         vector v = lgamma(e);
178         return exp(v.x) * v.y;
179 }
180
181 /**
182  * Pythonic mod:
183  * TODO: %% operator?
184  *
185  *  1 %  2 ==  1
186  * -1 %  2 ==  1
187  *  1 % -2 == -1
188  * -1 % -2 == -1
189  */
190 float pymod(float e, float f)
191 {
192         return e - f * floor(e / f);
193 }
194
195 float nearbyint(float e)
196 {
197         return rint(e);
198 }
199 float trunc(float e)
200 {
201         return (e>=0) ? floor(e) : ceil(e);
202 }
203
204 float fmod(float e, float f)
205 {
206         return e - f * trunc(e / f);
207 }
208 float remainder(float e, float f)
209 {
210         return e - f * rint(e / f);
211 }
212 vector remquo(float e, float f)
213 {
214         vector v;
215         v.z = 0;
216         v.y = rint(e / f);
217         v.x = e - f * v.y;
218         return v;
219 }
220
221 float copysign(float e, float f)
222 {
223         return fabs(e) * ((f>0) ? 1 : -1);
224 }
225 /// Always use `isnan` function to compare because `float x = nan(); x == x;` gives true
226 float nan(string tag)
227 {
228         return sqrt(-1);
229 }
230 float nextafter(float e, float f)
231 {
232         // TODO very crude
233         if(e == f)
234                 return nan("nextafter");
235         if(e > f)
236                 return -nextafter(-e, -f);
237         // now we know that e < f
238         // so we need the next number > e
239         float d, a, b;
240         d = max(fabs(e), 0.00000000000000000000001);
241         a = e + d;
242         do
243         {
244                 d *= 0.5;
245                 b = a;
246                 a = e + d;
247         }
248         while(a != e);
249         return b;
250 }
251 float nexttoward(float e, float f)
252 {
253         return nextafter(e, f);
254 }
255
256 float fdim(float e, float f)
257 {
258         return max(e-f, 0);
259 }
260 float fmax(float e, float f)
261 {
262         return max(e, f);
263 }
264 float fmin(float e, float f)
265 {
266         return min(e, f);
267 }
268 float fma(float e, float f, float g)
269 {
270         return e * f + g;
271 }
272
273 int isgreater(float e, float f)
274 {
275         return e > f;
276 }
277 int isgreaterequal(float e, float f)
278 {
279         return e >= f;
280 }
281 int isless(float e, float f)
282 {
283         return e < f;
284 }
285 int islessequal(float e, float f)
286 {
287         return e <= f;
288 }
289 int islessgreater(float e, float f)
290 {
291         return e < f || e > f;
292 }
293 int isunordered(float e, float f)
294 {
295         return !(e < f || e == f || e > f);
296 }