This post is the third in a series of posts on spherical trigonometry. We first looked at the analog of the Pythagorean theorem on a sphere, then the analog of the law of cosines on a sphere. Now we look at the analog of the law of sines on a sphere.

As before we denote the sides of a triangle with a, b, and c and denote the angle opposite a side by the capital version of that letter. So A is the angle opposite side a and so forth.

The law of sines on a sphere says

sin a / sin A = sin b / sin B = sin c / sin C.

Small triangles

Suppose for a moment that we’re looking at a relatively small triangle, one for which the length of the sides is small relative to the radius of the sphere. Then we can apply the small angle approximation

sin x ≈ x

to each side [1] and find

a / sin A ≈ b / sin B ≈ c / sin C.

In other words, the law of sines from plane trigonometry holds approximately in spherical trigonometry for small triangles. This is not surprising since we know we can ignore the curvature of the earth when working over relatively small areas, but it is reassuring.

Example: LAX to IAH

In the previous post we found the air distance between LAX (Los Angeles) and IAH (Houston) airports by making a spherical triangle with vertices at LAX, IAH, and the north pole. We labeled the side from LAX to the north pole a and the side from IAH to the north pole b. We found the length of side c using the law of cosines.

Now suppose we’re at LAX and about to fly to IAH. What direction should we fly? Since Los Angeles and Houston are approximately at the same latitude, and they’re not terribly far apart relative to the size of the earth, we know that we’ll roughly fly due east. Let’s find the flight direction more precisely.

The angle B in our triangle gives us the angle between our heading and a heading of due north. We know b from the latitude of IAH, we know C from the difference in longitude of the two airports, and we found c in the previous post. Therefore we can find B from the law of sines:

sin b / sin B = sin c / sin C

and we know everything but B. So

sin B = sin b sin C / sin c = sin 60.02° sin 23.07° / sin 19.92°.

This tells us sin B = 0.9960. So what is B? If we ask software for the arcsine of 0.9960 we will get 84.88°. This would say we head 5.12° north of due east to get from Los Angeles to Houston. But Houston is at a lower latitude than Los Angeles and so that doesn’t seem right.

In the first post in this series I mentioned that a flight between Houston and Wuhan would fly over Alaska even though the two cities are at essentially the same latitude. Wuhan is at a slightly higher latitude than Houston, but a flight leaving Wuhan for Houston would head somewhat north. So it is possible that a flight from a higher latitude to a lower latitude would head north [2].

But Los Angeles is much closer to Houston than Wuhan is, and the difference in latitude is greater. A flight from LA to Houston would head mostly east and a little south. So where did we go wrong?

The arcsine function in any software package returns a solution to the equation

sin(x) = y,

but it might not return the solution we need. The angle 84.88° is a valid solution to the equation

sin(B) = 0.9960

but it’s not the appropriate solution in our context. We need the solution 95.12°. Our flight heada 5.12° south of east. This assumes our flight takes a great circle path; airlines may have reasons not to take the shortest possible path between two airports.

***

[1] The approximation assumes x is measured in radians. If we measure angles in degrees or any other units, the right side of the approximation will need to be multiplied by a constant. But the same constant will appear in all the numerators and denominators of the law of sines and so the approximation holds in any units.

[2] Here’s another illustration of the same point. Draw a circle around the north pole, say a circle a mile in diameter, and suppose you want to walk from longitude 0 to longitude 180°. Following a path of constant latitude is walking around the circle. Heading north instead is walking across the diameter of the circle. If your destination is a little bit outside the circle on the other side, it’s still shorter to walk across the circle, i.e. head north, even though your destination is at a lower latitude, i.e. south.

You might think that if zip codes are close, then the regions they represent are close. Or that if zip codes are consecutive, then their regions touch. Neither of these are true. I explore how far they are from being true in the next post.

But these statements could have been true [1]. It’s possible to lay a grid on a region and number the squares in the grid in such a way that consecutive numbers correspond to contiguous squares, and nearby numbers correspond to nearby squares. There are geocoding systems that accomplish this using Hilbert curves. For example, Google’s S2 system is based on Hilbert curves.

A Hilbert curve is a curve that winds through a square, coming arbitrarily close to every point in the square. A Hilbert curve is a fractal, defined as the limit of an iterative process. We aren’t concerned with the limit because we only want to carry out a finite number of steps in the construction. But the fact that the limit exists tells us that with enough steps we can get any resolution we want.

To illustrate Hilbert curves and how they could be used to label grids, we will use a Hilbert curve to tour a chessboard. We will number the squares of a chess board so that consecutive numbers correspond to neighboring squares.

Here’s Python code to produce our numbering.

    from hilbertcurve.hilbertcurve import HilbertCurve

    p = 3 # 2^3 squares on a side
    n = 2 # two dimensional cube, i.e. square

    hc = HilbertCurve(p, n)

    for row in range(2**p):
        for col in range(2**p):
            d = hc.distance_from_point([row, col])
            print(format(d,"2d"), " ", end="")
        print()

This produces the following output.

     0   1  14  15  16  19  20  21  
     3   2  13  12  17  18  23  22  
     4   7   8  11  30  29  24  25  
     5   6   9  10  31  28  27  26  
    58  57  54  53  32  35  36  37  
    59  56  55  52  33  34  39  38  
    60  61  50  51  46  45  40  41  
    63  62  49  48  47  44  43  42  

Here’s Python code to draw lines rather than print numbers.

    N = (2**p)**n
    points = hc.points_from_distances(range(N))

    for i in range(1, N):
        x0, y0 = points[i-1]
        x1, y1 = points[i]
        plt.plot([x0, x1], [y0, y1])
    plt.gca().set_aspect('equal')    
    plt.show()

This produces the following image.

There’s a small problem. Following the square numberings above produces a Hilbert curve, and the plotting code produces and Hilbert curve, but they’re not the same Hilbert curve. The implicit axes created by our sequence of print statements doesn’t match the geometry in the plotting code. To get the curves to match we need to do a sort of change of coordinates. We can do this by changing the call to plt.plot to

    plt.plot([y0, y1], [8-x0, 8-x1])

This produces the following plot, matching the numbering above.

Related posts

• Maidenhead geocode system
• Decoding a grid square
• Mathematics and radio
• Online Maidenhead calculator

[1] It could have been true in the sense that it was theoretically possible. It might not have been practically possible. Postal codes have different design constraints than geocodes.