Game Development Community

Need some help with ballistic math

by David Erenger · in Technical Issues · 10/21/2005 (6:07 am) · 55 replies

Hi!

I'm working on an rts and have a problem with the math..

I have bullet ballistics without drag and I need to calculate the angle to fire for the projectile to hit a point..

If I understand the math right there is two angles that will result in a hit like this.

img489.imageshack.us/img489/7773/ballistic6co.jpg
I need both angles, the long for my artillery and the second for directhit units like tanks and rifles.

I have made an aproximation that works pretty well for direct hit bullets at short ranges but it doasnt do the trick for me anymore :/...

So I thougt I could calculate the real ekvation for a.
I can use this ekvation:

y = x tan(a) - g / (2Vo^2 cos(a)^2) * x^2

The x, y, Vo and g are all given, the only thing I need is a.

I have tried a couple of hours to get out the a but without 100% succes :(

So if there is anyone with some greater math skills then me please help me with this...


Thanks all
#21
10/26/2005 (8:56 am)
It's worth noting here that Mauchly and Eckert invented the modern computer, the ENIAC, during WW2 to solve this problem because the computations needed to fill out the gunnery tables for all the new weapons being built were to time consuming for even 100 people working full time. --Google tidbit
#22
10/26/2005 (12:28 pm)
Ok, after checking a little more I'm convinced the problem has to be solved iteratively, if you don't want to specify the maximum height or apogee of the projectile. Can't be done in closed form. Even then it's kind of hard to setup and write the solver. But here's how for the simpler case, where the landing elevation is the same as the launch elevation.

The range as a function of angle is

R(a) = 2(V0^2)sin(a)cos(a)/g;

Assume you want to solve for the angle, limited to a1 and a2 degrees, at a range of 2200 meters. Define

f(a) = R(a) - 2200;

then do this

f0 = f(a1);
if( f(a2) * f0 ) > 0 {
   error( 'no root in interval');
}

da = (a2-a1)/100;
a = a1 + da;
while ( a < a2 ) {
   if ( f(a)*f0 < 0 ) {
      aleft = a-da;
      aright = a;
      break;
   }
   a = a + da;
}

a will be within the interval [aleft,aright] that contains the solution. Change the 100 to something smaller to speed it up as long as the projectile lands close enough for your game.
#23
10/26/2005 (1:18 pm)
If you ask me the problem is this line right here:

t = x/(Vo * cos(a))

here t was replaced with an expression in terms of angle a but neither variable was ever specified. later on, even though an equation in 1 variable was reached, this ambiguity was still somehow propogated into the equation creating something which cannot be solved for a. I think a travel time, Vx or Vy value needs to be specified.
#24
10/26/2005 (1:26 pm)
Quote:I think a travel time, Vx or Vy value needs to be specified

Well, I did that when I proposed the alternative solution of specifying Vx as some useful value like 200 m/sec. That tie the travel time down. But absent that, there doesn't seem to be a direct solution to the problem. I realized that when I found that the gunnery people back in WW2 had teams of 100's of mathematicians, before the ENIAC, doing the iterative tables to 100 msec resolution.
#25
10/26/2005 (5:38 pm)
I am attempting to do exactly the same thing, David.

I will try to utilize some of these equations if I have some time.. I will share a solution if I find one.
#26
10/28/2005 (7:16 am)
Just for fun, I ran this problem thru Mathematica

height = (V^2*Sin[A]^2)/(2*G) ;

B = 5;
R = 2200;
V = 200;
G = 9.82;
H = height;

Solve[V*Sin[A]/G + (2(H - B)/G)^.5 -  R/(V*Cos[A]) == 0, A]
which produced the following:
Solve::"ifun": "Inverse functions are being used by Solve, so some 
solutions may not be found."

{{A -> -2.35848 - 0.701304 \[ImaginaryI]}, {A -> -2.35848 + 
        0.701304 \[ImaginaryI]}, {A -> 0.407581}, {A -> 1.16779}}

So it has real roots at angles of 23.33 degrees (.407 radians) and 67 degrees (1.167 radians).

So Mathematica can do it, but we can't!

:o)

Edit:
Solving this symbolically, w/o defining B,R, and V, produces about 10 pages of output. A lot of that has to do with the precision at which I was printing but it is daunting nonetheless. However, I can probably reduce this to a table of constants and a useful algorithm which reduces the precision to something appropriate for a game. I'll take a look at it tonight, before my Mathematica 4-day license expires, and grab the output for posterity and the GG community which seems to do a lot of this cannon firing sort of thing. :0
#27
10/28/2005 (8:48 am)
Tim a table of constants sounds like a good idea. you could even describe the problem in detail and divulge your results as a resource. from what Thomas said, others have obviously run into this problem as well. the only problem would be coming up with a reasonable domain. the simplest thing to do would be to plot the graph as:

y=ASin(a) - B/Cos(a)

but how to come up with a reasonable domain of values for A and B...

any kind of resource pertaining to this stuff would be very cool.
#28
10/28/2005 (9:06 am)
I'll see if I can produce a resource that describes something like

%angle = getFiringAngle(V,R,B);

Need to do it quick though while I have Mathematica to play with. My machine config has changed since the last time I had the Student Version installed and now I'm stuck with only a 4 day password because I'm not a student anymore. Rats! Gotta go enroll for a class everytime I upgrade my computer.
#29
10/29/2005 (9:18 am)
Good god, please tell me Mathmatica generated that for you! ;)
#30
10/29/2005 (12:40 pm)
Oh man that is the craziest equation i've ever seen. ( I don't solve many equations) hehe.

Good luck and great job so far!
#31
10/29/2005 (2:01 pm)
It just looks weird because it contains the Taylor series expansions (or the like) for functions (like sine and cosine) that we never see. That's just the way that software like Mathematica and Maple do it so that equations can be fully manipulated and reduced. All we care about is the end effect.
#32
10/30/2005 (7:36 am)
Looks like I will have to give up trying to find the closed form because I just can't get past the "implicit root' derivations that MM spits out in the intermediate solution form. They use some internal iterative solver for that which I am unable to easily reproduce and the likelyhood that I will is slim. But I dug out my old copy of "Numerical Recipes in C" to have a look. Wow I haven't touched that in a long time but it used to be on my work desk all the time.

Anyway there are some solvers there that would do it, like Newton's Method, etc, but this has gotten way out of scope. The simple iterative solving scheme proposed in psuedo code above is more appropriate for this.

On the positive side, I learned a lot about how Mathematica generates solutions that can be directly implemented in C. The bad news is that this one had implicit roots buried in the each term. I tried a few other problems and they didn't generate implicit roots. But then I know the original problem is hard.

More bad news is that my MM license is expiring now, (first time I have to reboot), so I can't rely on following up any further work with cross checks using it. Rats! But I'll keep thinking about it in the long term mostly for the usefulness of having a it as a code(like) generator.
#33
10/30/2005 (1:07 pm)
Oh well. You tried.. I will also be keeping this on my back burner brain segment. I need to set aside some time to figure it out...
#34
10/31/2005 (12:28 am)
Hi again, a friend of mine came up with a good teori for how to solve this.

a = (( sin^-1 * ( ((g*cos(aTan(y/x))*sqrt(x^2+y^2)) / Vo^2 ) ) / 2) + aTan(y/x)

dont have the time to explain it nor test it now, I'm comming back later tonight...

fell free to test it, might work :)
#35
10/31/2005 (4:26 am)
Hi David. the equation is very encouraging but it doesn't look like it is written down correctly. Is the first sin^-1*(.... actually asin(.....?
#36
10/31/2005 (7:40 pm)
David, I came back tonight to have a look at the new equation and noticed that it seems to reduce to something more compact but not necessarily useful. The part, cos(atan(y/x)) * sqrt(x^2+y^2) is just 'x', in fact. Why?
Cos(atan(y/x) is defined to be x/R, where R is the distance to x,y. And, of course, R = sqrt(x^2+y^2). so that part of the equation is x/R * R, or just plain 'x'.

So a = asin(g*x/(2*V^2)) + atan(y/x), which doesn't seem to be what we are looking for. Besides, we know the correct equation has two solutions, (two real roots) for each muzzle velocity, one fired at low trajectory (rifle shot) and one fired in a high trajectory (mortar mode). So the equation needs to be quadratic.

For my part, I'll be looking for a good C quadratic solver we can plug into this problem.
#37
11/01/2005 (9:33 am)
Two things spring to mind:

1. From a realism point of view, few weapons ever fire at angles > 45 degrees when shooting at ground targets. Mortars are the exception, and they generall never fire at < 45 degrees.

2. If the issue is that your GUI and/or AI requires a working gunsight or sense of the elevation to use, you should do what people who make gunsights actually do: fire your simulated gun at N different angles of elevation and observe the range attained. Take those range,elevation pairs and build a spline datamember that can relate a given range (within the envelope of the system) to the elevation that will reach it. Math for this should not be fancy AT ALL.

Note in particular that when you discard the idea of solving the formulas listed above and instead use splines, you can put drag and any other effects into play again and it will all "just work". If variable factors (say, wind or air density, or a target at dissimilar elevation) cause a firing scenario to differ from that the sight is designed for, you could have the AI employ a binary search to bracket and hit the target.

tone
#38
11/01/2005 (10:01 am)
@Tim, about the ekvation I wrote down yesterday it was just a teory I was given that I tried to whrite down very fast without checking it at all, sorry about that. but.....

great news everyone, I have got the answers from the math professor I asked, here it is

a = φ/2 + (1/2)arcsin((y + g(x/v0)^2)/√(x^2 + y^2)) + n*Pi
a = Pi/2 + φ/2 − (1/2)arcsin((y + g(x/v0)^2)/√(x^2 + y^2)) + n*Pi

the answers we need are when n = 0. becouse φ is in [0,Pi/2) so, φ = arctan(y/x).

need to try it out first but a guess he would know :)

@Anthony, I did some thinking about using some kind databank, but it would be a great job building it, say I have 5 different mortars with a range from 100- 4000 meters, and the altitude differs with a couple of hundered meters, that would be a lot of data, I plus that I have 20+ units with different ammo...
#39
11/01/2005 (10:24 am)
I can now confirm that it works perfekt, here it is again easier to read maby.

fastest way:
a =aTan(y/x)/2 + (1/2)aSin((y + g(x/v0)^2)/sqrt(x^2 + y^2))

High (?) way:
a = Pi/2 + aTan(y/x)/2 - (1/2)aSin((y + g(x/v0)^2)/sqrt(x^2 + y^2))
#40
11/01/2005 (11:35 am)
David, this seems really good. But when I run the same data as used above with Mathematica for which I got

{A -> 0.407581}, {A -> 1.16779}},

through Matlab, I get

B = 5;
R = 2200;
V = 200;
G = 9.82;
A =atan(B/R)/2 + (1/2)*asin((B + G*(R/V)^2)/sqrt(R^2 + B^2));
A1 = A
A2 = 3.1415/2 - A

A1 = 0.2878
A2 = 1.2830

So I wonder which one is really correct?

OTH, yours is "best" because we can actually use it. The proof is in where the shells land in the game, I suppose. But substituting the original expressions back into your new solution should produce an identity if it is correct.