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