Arc Projectiles: Computing Angles

Given a target $\langle \Delta x, \Delta y \rangle$ away and an initial speed $v$, can you find the launch angle $\theta$ of a projectile that will hit the target?

In this post, I'll go over the step-by-step process in attempting to solve this problem. I'll focus mainly on math, but since I initially encountered this problem while working on a video game, I will include some code at the end!

Some basic experience with physics and calculus is recommended for this post.

Derivation

Getting to a solution was much harder than I expected, and involved lots of failed attempts. In the following sections, I will highlight the main parts of the solution, including some of the things I tried that didn't work. If you just want the solution, skip ahead to Final Steps.

Kinematics

My first instinct was to use kinematics to solve for the launch angle.

I'm going to assume air resistance is negligible, so once the projectile is launched, it will have an acceleration of $\vec{a} = \langle 0, -g \rangle$. Using this, we have two kinematic equations:

$\begin{aligned} \Delta y &= v\sin(\theta)t - \frac{1}{2}gt^2 \\ \Delta x &= v\cos(\theta)t \end{aligned}$

We have two unknowns and two equations, so we should be good to go to solve it! I used the second equation to express $t$ in terms of $\theta$:

$\begin{aligned} \Delta x &= v\cos(\theta)t \\ t &= \frac{\Delta x}{v\cos(\theta)} \end{aligned}$

Plugging $t$ into the first equation:

$\Delta y = \Delta x\tan(\theta) - \frac{g\Delta x^2}{2v^2\cos^2(\theta)}$

Hmm. The equation we land on doesn't seem to have an algebraic solution for $\theta$. At this point, my gut told me that this was okay: I'd just use solve for it using some approximation formula. If I couldn't get a closed-form solution, I'd settle for a computed answer.

Approximation

In order to perform numeric approximation, we're first going to make a single-variable function:

$f(\theta) = -\Delta y + \Delta x\tan(\theta) - \frac{g\Delta x^2}{2v^2\cos^2(\theta)}$

Thus, when $f(\theta) = 0$, we have a solution to the kinematic equation. First, I tried to use Newton's method to find a sufficient answer. If you're unfamiliar, the main statement is this:

$\theta_{n+1} = \theta_n - \frac{f(\theta_n)}{f'(\theta_n)} \quad \text{where} \,\, f(\theta) = 0$

As $n \to \infty$, the approximation gets better. Newton's method is great (and efficient) for solving equations using code! But notice that we need the derivative of $f$. Naively differentiating f resulted in this:

$f'(\theta) = \Delta x\sec^2(\theta)\left(1 - \frac{g\Delta x\tan(\theta)}{v^2}\right)$

With this, I wrote up a quick Python script to approximate $\theta$ using $f$, $f'$, and Newton's Method. But I kept getting incorrect answers!¹ I turned to Desmos to try to figure out why, and here's what I found:

The red line is the graph of my $f'$, and the green one is the graph of the derivative of $f$ that is generated by Desmos. They should be the exact same, but they're not! I'm still not sure why (my knowledge of calculus isn't deep enough, but I'm sure there's a reason I can't naively differentiate like I did before).

Secant Method

So, I attempted a different method of approximation, called the secant method. It's very similar to Newton's method, but does not rely on knowing the derivative of the function.

After typing up an implementation in Python, I started to get correct angles!

Issues

But there's a huge problem with approximation that I ignored at the beginning: since we're working with a trigonometric function, there are 2 distinct roots of $f$, and an infinite number of translations of those roots. To illustrate what I mean, here's Desmos graph of $f$ for some hard-coded constants $\Delta x$, $\Delta y$, and $v$:

All numeric approximation methods take an initial guess. Depending on the guess, I'd get one of the solutions on the graph, and have absolutely no information about where the other solutions are.²

Because of this, numeric approximation wasn't going to cut it.

Substitution

Let's take a look at our original equation again:

$\Delta y = \Delta x\tan(\theta) - \frac{g\Delta x^2}{2v^2\cos^2(\theta)}$

The reason this is algebraically unsolvable stems from the $\tan$ and $\cos$. I remembered that a neat trick with these types of equations is substitution. Ideally, I could swap out $\tan$ or $\cos$ with some variable, which we'll call $z$, that might make the equation solvable, and then do inverse trigonometry on $z$ to get $\theta$ back out.

It's possible to express $\cos^2(x)$ in terms of $\tan$:

$\cos^2(x) = \frac{1}{1 + \tan^2(x)}$

Using this, we have:

$\Delta y = \Delta x\tan(\theta) - \frac{g\Delta x^2(1 + \tan^2(\theta))}{2v^2}$

Let $z = \tan(\theta)$.

$\Delta y = \Delta xz - \frac{g\Delta x^2(1 + z^2)}{2v^2}$

And suddenly, we're left with a dead-simple quadratic equation for $z$! The wonders of trigonometric identities and algebra have led us to something we can actually solve.

Final Steps

For good measure, here's the final expression for $z$ using the quadratic formula:

$\begin{aligned} a &= -\frac{g\Delta x^2}{2v^2} \\ z &= \frac{-\Delta x \pm \sqrt{\Delta x^2 - 4a(a - \Delta y)}}{2a} \end{aligned}$

$z = \tan(\theta) \Rightarrow \theta = \tan^{-1}(z)$! But remember that we get two answers because of the $\pm$. So which one do we use? Well, it depends on what direction we're shooting in, which is indicated by the sign of $\Delta x$.

We also have a handy physical meaning for the discriminant: if there are no solutions (i.e. the discriminant < 0), then our target is too far for our initial speed.

Wrap-Up

In this post, I walked through my trial and error process for deriving a solution to a classic physics problem. Although it took me a while, it was fun to figure it out on my own.

I hope you learned something from this post! Below, I put two code implementations: one with Python, and one with GDScript.³

Code

Python

import math


def firing_angle(v: float, dx: float, dy: float) -> float | None:
    GRAVITY = 9.81

    a = -((GRAVITY * dx**2) / (2 * v**2))
    b = dx
    c = dy + a
    discriminant = b**2 - (4 * a * c)
    if discriminant < 0:
        return None
    root = math.sqrt(discriminant)
    angle1 = math.atan((-b - root) / (2 * a))
    angle2 = math.atan((-b + root) / (2 * a))

    if dx > 0:
        return angle1
    else:
        return angle2

GDScript

Godot's coordinate system is a bit strange: the y-axis points downwards, which should change our answers. The fix isn't too bad, we just have to mess around with the quadrants a little!

static func firing_angle(v: float, dx: float, dy: float) -> float:
	# Units are in pixels, so use 980
	const GRAVITY = 980.0

	var a := -((GRAVITY * pow(dx, 2)) / (2 * pow(v, 2)))
	var b := dx
	var c := dy + a
	var discriminant := pow(b, 2) - (4 * a * c)
	# This is an assertion, but you can change it to whatever you want
	assert(discriminant > 0)
	var root := sqrt(discriminant)
	var angle1 := atan((-b - root) / (2 * a))
	var angle2 := atan((-b + root) / (2 * a))

	if dx > 0:
		# angle1 should be in the 1st quadrant
		assert(angle1 > 0)
		# Negate to put into Godot's coordinate system (-y points upwards)
		return -angle1
	else:
		# angle2 should be in the 4th quadrant
		assert(angle2 < 0)
		# Translate to the 2nd quadrant (for Godot's coordinate system)
		return -PI - angle2