1 #include "calculations.qh"
3 // =============================
4 // Explosion Force Calculation
5 // =============================
7 float explosion_calcpush_getmultiplier(vector explosion_v, vector target_v)
10 a = explosion_v * (explosion_v - target_v);
13 // target is too fast to be hittable by this
16 a /= (explosion_v * explosion_v);
17 // we know we can divide by this, or above a would be == 0
23 vector explosion_calcpush(vector explosion_v, float explosion_m, vector target_v, float target_m, float elasticity)
25 // solution of the equations:
26 // v' = v + a vp // central hit
27 // m*v' + mp*vp' = m*v + mp*vp // conservation of momentum
28 // m*v'^2 + mp*vp'^2 = m*v^2 + mp*vp^2 // conservation of energy (ELASTIC hit)
29 // -> a = 0 // case 1: did not hit
30 // -> a = 2*mp*(vp^2 - vp.v) / ((m+mp) * vp^2) // case 2: did hit
31 // // non-elastic hits are somewhere between these two
33 // this would be physically correct, but we don't do that
34 return explosion_v * explosion_calcpush_getmultiplier(explosion_v, target_v) * (
38 target_m + explosion_m
40 ); // note: this factor is at least 0, at most 2
44 // simplified formula, tuned so that if the target has velocity 0, we get exactly the original force
45 vector damage_explosion_calcpush(vector explosion_f, vector target_v, float speedfactor)
47 // if below 1, the formulas make no sense (and would cause superjumps)
54 // speedfactor * (1 + e) * m / (1 + m) == 1
55 m = 1 / ((1 + 0) * speedfactor - 1);
57 v = explosion_calcpush(explosion_f * speedfactor, m, target_v, 1, 0);
58 // the factor we then get is:
60 LOG_INFOF("MASS: %f\nv: %v -> %v\nENERGY BEFORE == %f + %f = %f\nENERGY AFTER >= %f\n",
62 target_v, target_v + v,
63 target_v * target_v, m * explosion_f * speedfactor * explosion_f * speedfactor, target_v * target_v + m * explosion_f * speedfactor * explosion_f * speedfactor,
64 (target_v + v) * (target_v + v));
67 return explosion_f * explosion_calcpush_getmultiplier(explosion_f * speedfactor, target_v);
71 // =========================
72 // Shot Spread Calculation
73 // =========================
75 vector cliptoplane(vector v, vector p)
77 return v - (v * p) * p;
80 vector solve_cubic_pq(float p, float q)
83 D = q*q/4.0 + p*p*p/27.0;
87 a = 1.0/3.0 * acos(-q/2.0 * sqrt(-27.0/(p*p*p)));
88 u = sqrt(-4.0/3.0 * p);
96 '1 0 0' * cos(a + 2.0/3.0*M_PI)
98 '0 1 0' * cos(a + 4.0/3.0*M_PI)
111 return '1 1 0' * v + '0 0 1' * u;
113 return '0 1 1' * v + '1 0 0' * u;
118 u = cbrt(-q/2.0 + sqrt(D));
119 v = cbrt(-q/2.0 - sqrt(D));
120 return '1 1 1' * (u + v);
123 vector solve_cubic_abcd(float a, float b, float c, float d)
130 q = (27*a*a*d - 9*a*b*c + 2*b*b*b);
131 v = solve_cubic_pq(p, q);
132 v = (v - b * '1 1 1') * (1.0 / (3.0 * a));
134 v += '1 0 -1' * (v.z - v.x); // swap x, z
138 vector findperpendicular(vector v)
144 return normalize(cliptoplane(p, v));
148 int W_GunAlign(entity this, int preferred_align)
150 entity own = this.owner;
151 // using wasfreed, as we don't actually clear .gunaligns yet
152 if(!own.gunaligns[preferred_align] || wasfreed(own.gunaligns[preferred_align]) || own.gunaligns[preferred_align] == this)
154 own.gunaligns[preferred_align] = this;
155 return preferred_align; // fall back if the good one is already choosable
158 for(int j = 4; j > 0; --j) // start from left and try the others
160 if(!own.gunaligns[j] || wasfreed(own.gunaligns[j]) || own.gunaligns[j] == this)
162 own.gunaligns[j] = this;
167 own.gunaligns[preferred_align] = this;
168 return preferred_align; // no other choice
171 int W_GunAlign(entity this, int preferred_align)
173 return this.m_gunalign > 0 ? this.m_gunalign : preferred_align;
178 int W_GetGunAlignment(entity player)
180 int gunalign = STAT(GUNALIGN, player);
181 if(gunalign < 1 || gunalign > 4)
182 gunalign = 3; // default value
189 vector W_CalculateSpread(vector forward, float spread, float spreadfactor, float spreadstyle)
192 vector v1 = '0 0 0', v2;
194 spread *= spreadfactor; //g_weaponspreadfactor;
202 // this is the baseline for the spread value!
203 // standard deviation: sqrt(2/5)
204 // density function: sqrt(1-r^2)
205 return forward + randomvec() * spread;
209 // same thing, basically
210 return normalize(forward + cliptoplane(randomvec() * spread, forward));
214 // circle spread... has at sigma=1 a standard deviation of sqrt(1/2)
215 sigma = spread * 0.89442719099991587855; // match baseline stddev
216 v1 = findperpendicular(forward);
217 v2 = cross(forward, v1);
218 // random point on unit circle
219 dx = random() * 2 * M_PI;
222 // radius in our dist function
225 return normalize(forward + (v1 * dx + v2 * dy) * r * sigma);
229 sigma = spread * 0.44721359549996; // match baseline stddev
230 // note: 2D gaussian has sqrt(2) times the stddev of 1D, so this factor is right
232 v1_x += gsl_ran_gaussian(sigma);
233 v1_y += gsl_ran_gaussian(sigma);
234 v1_z += gsl_ran_gaussian(sigma);
239 sigma = spread * 0.44721359549996; // match baseline stddev
240 // note: 2D gaussian has sqrt(2) times the stddev of 1D, so this factor is right
241 v1_x = gsl_ran_gaussian(sigma);
242 v1_y = gsl_ran_gaussian(sigma);
243 v1_z = gsl_ran_gaussian(sigma);
244 return normalize(forward + cliptoplane(v1, forward));
248 sigma = spread * 1.154700538379252; // match baseline stddev
249 v1 = findperpendicular(forward);
250 v2 = cross(forward, v1);
251 // random point on unit circle
252 dx = random() * 2 * M_PI;
255 // radius in our dist function
257 r = solve_cubic_abcd(-2, 3, 0, -r) * '0 1 0';
258 return normalize(forward + (v1 * dx + v2 * dy) * r * sigma);
262 sigma = spread * 1.095445115010332; // match baseline stddev
263 v1 = findperpendicular(forward);
264 v2 = cross(forward, v1);
265 // random point on unit circle
266 dx = random() * 2 * M_PI;
269 // radius in our dist function
273 return normalize(forward + (v1 * dx + v2 * dy) * r * sigma);
275 case 7: // (1-r) (2-r)
277 sigma = spread * 1.224744871391589; // match baseline stddev
278 v1 = findperpendicular(forward);
279 v2 = cross(forward, v1);
280 // random point on unit circle
281 dx = random() * 2 * M_PI;
284 // radius in our dist function
288 return normalize(forward + (v1 * dx + v2 * dy) * r * sigma);
291 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)!");
296 * how to derive falloff functions:
297 * rho(r) := (2-r) * (1-r);
300 * rhor(r) := r * rho(r);
301 * cr(t) := integrate(rhor(r), r, a, t);
302 * scr(t) := integrate(rhor(r) * r^2, r, a, t);
303 * variance : scr(b) / cr(b);
304 * solve(cr(r) = rand * cr(b), r), programmmode:false;
305 * sqrt(0.4 / variance), numer;