// works for up to 10 decimals!
string ftos_decimals(float number, float decimals)
{
+ // inhibit stupid negative zero
+ if(number == 0)
+ number = 0;
// we have sprintf...
return sprintf("%.*f", decimals, number);
}
return valstr;
}
+float dotproduct(vector a, vector b)
+{
+ return a_x * b_x + a_y * b_y + a_z * b_z;
+}
+
vector cross(vector a, vector b)
{
return
{
float created_saved_value;
entity e;
+
+ created_saved_value = FALSE;
if not(tmp_cvar || tmp_value)
{
{
float eps;
eps = (max(a, -a) + max(c, -c)) * 0.001;
+ if(a > c)
+ eps = -eps;
return b == median(a - eps, b, c + eps);
}
return v;
}
+vector solve_shotdirection(vector myorg, vector myvel, vector eorg, vector evel, float spd, float newton_style)
+{
+ vector ret;
+
+ // make origin and speed relative
+ eorg -= myorg;
+ if(newton_style)
+ evel -= myvel;
+
+ // now solve for ret, ret normalized:
+ // eorg + t * evel == t * ret * spd
+ // or, rather, solve for t:
+ // |eorg + t * evel| == t * spd
+ // eorg^2 + t^2 * evel^2 + 2 * t * (eorg * evel) == t^2 * spd^2
+ // t^2 * (evel^2 - spd^2) + t * (2 * (eorg * evel)) + eorg^2 == 0
+ vector solution = solve_quadratic(evel * evel - spd * spd, 2 * (eorg * evel), eorg * eorg);
+ // p = 2 * (eorg * evel) / (evel * evel - spd * spd)
+ // q = (eorg * eorg) / (evel * evel - spd * spd)
+ if(!solution_z) // no real solution
+ {
+ // happens if D < 0
+ // (eorg * evel)^2 < (evel^2 - spd^2) * eorg^2
+ // (eorg * evel)^2 / eorg^2 < evel^2 - spd^2
+ // spd^2 < ((evel^2 * eorg^2) - (eorg * evel)^2) / eorg^2
+ // spd^2 < evel^2 * (1 - cos^2 angle(evel, eorg))
+ // spd^2 < evel^2 * sin^2 angle(evel, eorg)
+ // spd < |evel| * sin angle(evel, eorg)
+ return '0 0 0';
+ }
+ else if(solution_x > 0)
+ {
+ // both solutions > 0: take the smaller one
+ // happens if p < 0 and q > 0
+ ret = normalize(eorg + solution_x * evel);
+ }
+ else if(solution_y > 0)
+ {
+ // one solution > 0: take the larger one
+ // happens if q < 0 or q == 0 and p < 0
+ ret = normalize(eorg + solution_y * evel);
+ }
+ else
+ {
+ // no solution > 0: reject
+ // happens if p > 0 and q >= 0
+ // 2 * (eorg * evel) / (evel * evel - spd * spd) > 0
+ // (eorg * eorg) / (evel * evel - spd * spd) >= 0
+ //
+ // |evel| >= spd
+ // eorg * evel > 0
+ //
+ // "Enemy is moving away from me at more than spd"
+ return '0 0 0';
+ }
+
+ // NOTE: we always got a solution if spd > |evel|
+
+ if(newton_style == 2)
+ ret = normalize(ret * spd + myvel);
+
+ return ret;
+}
+
+vector get_shotvelocity(vector myvel, vector mydir, float spd, float newton_style, float mi, float ma)
+{
+ if(!newton_style)
+ return spd * mydir;
+
+ if(newton_style == 2)
+ {
+ // true Newtonian projectiles with automatic aim adjustment
+ //
+ // solve: |outspeed * mydir - myvel| = spd
+ // outspeed^2 - 2 * outspeed * (mydir * myvel) + myvel^2 - spd^2 = 0
+ // outspeed = (mydir * myvel) +- sqrt((mydir * myvel)^2 - myvel^2 + spd^2)
+ // PLUS SIGN!
+ // not defined?
+ // then...
+ // myvel^2 - (mydir * myvel)^2 > spd^2
+ // velocity without mydir component > spd
+ // fire at smallest possible spd that works?
+ // |(mydir * myvel) * myvel - myvel| = spd
+
+ vector solution = solve_quadratic(1, -2 * (mydir * myvel), myvel * myvel - spd * spd);
+
+ float outspeed;
+ if(solution_z)
+ outspeed = solution_y; // the larger one
+ else
+ {
+ //outspeed = 0; // slowest possible shot
+ outspeed = solution_x; // the real part (that is, the average!)
+ //dprint("impossible shot, adjusting\n");
+ }
+
+ outspeed = bound(spd * mi, outspeed, spd * ma);
+ return mydir * outspeed;
+ }
+
+ // real Newtonian
+ return myvel + spd * mydir;
+}
+
void check_unacceptable_compiler_bugs()
{
if(cvar("_allow_unacceptable_compiler_bugs"))
while(acl)
{
t = car(acl); acl = cdr(acl);
+
d = 1;
if(substring(t, 0, 1) == "-")
{
}
else if(substring(t, 0, 1) == "+")
t = substring(t, 1, strlen(t) - 1);
+
if(substring(t, -1, 1) == "*")
{
t = substring(t, 0, strlen(t) - 1);
- s = substring(s, 0, strlen(t));
+ s = substring(str, 0, strlen(t));
}
else
s = str;
}
to_execute_next_frame = strzone(s);
}
+
+float cubic_speedfunc(float startspeedfactor, float endspeedfactor, float x)
+{
+ return
+ ((( startspeedfactor + endspeedfactor - 2
+ ) * x - 2 * startspeedfactor - endspeedfactor + 3
+ ) * x + startspeedfactor
+ ) * x;
+}
+
+float cubic_speedfunc_is_sane(float startspeedfactor, float endspeedfactor)
+{
+ if(startspeedfactor < 0 || endspeedfactor < 0)
+ return FALSE;
+
+ /*
+ // if this is the case, the possible zeros of the first derivative are outside
+ // 0..1
+ We can calculate this condition as condition
+ if(se <= 3)
+ return TRUE;
+ */
+
+ // better, see below:
+ if(startspeedfactor <= 3 && endspeedfactor <= 3)
+ return TRUE;
+
+ // if this is the case, the first derivative has no zeros at all
+ float se = startspeedfactor + endspeedfactor;
+ float s_e = startspeedfactor - endspeedfactor;
+ if(3 * (se - 4) * (se - 4) + s_e * s_e <= 12) // an ellipse
+ return TRUE;
+
+ // Now let s <= 3, s <= 3, s+e >= 3 (triangle) then we get se <= 6 (top right corner).
+ // we also get s_e <= 6 - se
+ // 3 * (se - 4)^2 + (6 - se)^2
+ // is quadratic, has value 12 at 3 and 6, and value < 12 in between.
+ // Therefore, above "better" check works!
+
+ return FALSE;
+
+ // known good cases:
+ // (0, [0..3])
+ // (0.5, [0..3.8])
+ // (1, [0..4])
+ // (1.5, [0..3.9])
+ // (2, [0..3.7])
+ // (2.5, [0..3.4])
+ // (3, [0..3])
+ // (3.5, [0.2..2.3])
+ // (4, 1)
+}
+
+#ifndef MENUQC
+vector cliptoplane(vector v, vector p)
+{
+ return v - (v * p) * p;
+}
+
+vector solve_cubic_pq(float p, float q)
+{
+ float D, u, v, a;
+ D = q*q/4.0 + p*p*p/27.0;
+ if(D < 0)
+ {
+ // irreducibilis
+ a = 1.0/3.0 * acos(-q/2.0 * sqrt(-27.0/(p*p*p)));
+ u = sqrt(-4.0/3.0 * p);
+ // a in range 0..pi/3
+ // cos(a)
+ // cos(a + 2pi/3)
+ // cos(a + 4pi/3)
+ return
+ u *
+ (
+ '1 0 0' * cos(a + 2.0/3.0*M_PI)
+ +
+ '0 1 0' * cos(a + 4.0/3.0*M_PI)
+ +
+ '0 0 1' * cos(a)
+ );
+ }
+ else if(D == 0)
+ {
+ // simple
+ if(p == 0)
+ return '0 0 0';
+ u = 3*q/p;
+ v = -u/2;
+ if(u >= v)
+ return '1 1 0' * v + '0 0 1' * u;
+ else
+ return '0 1 1' * v + '1 0 0' * u;
+ }
+ else
+ {
+ // cardano
+ u = cbrt(-q/2.0 + sqrt(D));
+ v = cbrt(-q/2.0 - sqrt(D));
+ return '1 1 1' * (u + v);
+ }
+}
+vector solve_cubic_abcd(float a, float b, float c, float d)
+{
+ // y = 3*a*x + b
+ // x = (y - b) / 3a
+ float p, q;
+ vector v;
+ p = (9*a*c - 3*b*b);
+ q = (27*a*a*d - 9*a*b*c + 2*b*b*b);
+ v = solve_cubic_pq(p, q);
+ v = (v - b * '1 1 1') * (1.0 / (3.0 * a));
+ if(a < 0)
+ v += '1 0 -1' * (v_z - v_x); // swap x, z
+ return v;
+}
+
+vector findperpendicular(vector v)
+{
+ vector p;
+ p_x = v_z;
+ p_y = -v_x;
+ p_z = v_y;
+ return normalize(cliptoplane(p, v));
+}
+
+vector W_CalculateSpread(vector forward, float spread, float spreadfactor, float spreadstyle)
+{
+ float sigma;
+ vector v1, v2;
+ float dx, dy, r;
+ float sstyle;
+ spread *= spreadfactor; //g_weaponspreadfactor;
+ if(spread <= 0)
+ return forward;
+ sstyle = spreadstyle; //autocvar_g_projectiles_spread_style;
+
+ if(sstyle == 0)
+ {
+ // this is the baseline for the spread value!
+ // standard deviation: sqrt(2/5)
+ // density function: sqrt(1-r^2)
+ return forward + randomvec() * spread;
+ }
+ else if(sstyle == 1)
+ {
+ // same thing, basically
+ return normalize(forward + cliptoplane(randomvec() * spread, forward));
+ }
+ else if(sstyle == 2)
+ {
+ // circle spread... has at sigma=1 a standard deviation of sqrt(1/2)
+ sigma = spread * 0.89442719099991587855; // match baseline stddev
+ v1 = findperpendicular(forward);
+ v2 = cross(forward, v1);
+ // random point on unit circle
+ dx = random() * 2 * M_PI;
+ dy = sin(dx);
+ dx = cos(dx);
+ // radius in our dist function
+ r = random();
+ r = sqrt(r);
+ return normalize(forward + (v1 * dx + v2 * dy) * r * sigma);
+ }
+ else if(sstyle == 3) // gauss 3d
+ {
+ sigma = spread * 0.44721359549996; // match baseline stddev
+ // note: 2D gaussian has sqrt(2) times the stddev of 1D, so this factor is right
+ v1 = forward;
+ v1_x += gsl_ran_gaussian(sigma);
+ v1_y += gsl_ran_gaussian(sigma);
+ v1_z += gsl_ran_gaussian(sigma);
+ return v1;
+ }
+ else if(sstyle == 4) // gauss 2d
+ {
+ sigma = spread * 0.44721359549996; // match baseline stddev
+ // note: 2D gaussian has sqrt(2) times the stddev of 1D, so this factor is right
+ v1_x = gsl_ran_gaussian(sigma);
+ v1_y = gsl_ran_gaussian(sigma);
+ v1_z = gsl_ran_gaussian(sigma);
+ return normalize(forward + cliptoplane(v1, forward));
+ }
+ else if(sstyle == 5) // 1-r
+ {
+ sigma = spread * 1.154700538379252; // match baseline stddev
+ v1 = findperpendicular(forward);
+ v2 = cross(forward, v1);
+ // random point on unit circle
+ dx = random() * 2 * M_PI;
+ dy = sin(dx);
+ dx = cos(dx);
+ // radius in our dist function
+ r = random();
+ r = solve_cubic_abcd(-2, 3, 0, -r) * '0 1 0';
+ return normalize(forward + (v1 * dx + v2 * dy) * r * sigma);
+ }
+ else if(sstyle == 6) // 1-r^2
+ {
+ sigma = spread * 1.095445115010332; // match baseline stddev
+ v1 = findperpendicular(forward);
+ v2 = cross(forward, v1);
+ // random point on unit circle
+ dx = random() * 2 * M_PI;
+ dy = sin(dx);
+ dx = cos(dx);
+ // radius in our dist function
+ r = random();
+ r = sqrt(1 - r);
+ r = sqrt(1 - r);
+ return normalize(forward + (v1 * dx + v2 * dy) * r * sigma);
+ }
+ else if(sstyle == 7) // (1-r) (2-r)
+ {
+ sigma = spread * 1.224744871391589; // match baseline stddev
+ v1 = findperpendicular(forward);
+ v2 = cross(forward, v1);
+ // random point on unit circle
+ dx = random() * 2 * M_PI;
+ dy = sin(dx);
+ dx = cos(dx);
+ // radius in our dist function
+ r = random();
+ r = 1 - sqrt(r);
+ r = 1 - sqrt(r);
+ return normalize(forward + (v1 * dx + v2 * dy) * r * sigma);
+ }
+ else
+ error("g_projectiles_spread_style must be 0 (sphere), 1 (flattened sphere), 2 (circle), 3 (gauss 3D), 4 (gauss plane), 5 (linear falloff), 6 (quadratic falloff), 7 (stronger falloff)!");
+ return '0 0 0';
+ /*
+ * how to derive falloff functions:
+ * rho(r) := (2-r) * (1-r);
+ * a : 0;
+ * b : 1;
+ * rhor(r) := r * rho(r);
+ * cr(t) := integrate(rhor(r), r, a, t);
+ * scr(t) := integrate(rhor(r) * r^2, r, a, t);
+ * variance : scr(b) / cr(b);
+ * solve(cr(r) = rand * cr(b), r), programmmode:false;
+ * sqrt(0.4 / variance), numer;
+ */
+}
+#endif