Notes on Mercator’s Projection

The following work is my math essay assignment (an essay on any topic related to the Math 30 curriculum). I found it interesting enough that I’ll post it here on my math blog. Yay.


How do you make a world map? We find world maps everywhere — but the making of world maps is more complicated than most people realize. The earth is round, and a map is flat. How should one represent the surface of the earth, a sphere, on a flat sheet of paper?

The problem of displaying the surface of a sphere on a map has concerned cartographers, or map makers, for centuries. Over the years, cartographers have come up with a multitude of map projections, or systematic methods of transferring the surface of the earth onto flat paper.

There are many types of map projections, but one of the most common types is the cylindrical projection. Imagine a semi-transparent, hollow globe with a light inside of it. Now wrap a piece of blank map paper around it, like a cylinder:

Let the earth be represented as a unit sphere. The longitude (east-west location) shall be denoted as \lambda, and the latitude (north-south location) as \theta. A is the center of the earth, and B a point on the equator.

Now a light source is placed at point A; a point D on the surface of the globe will be projected in a straight line to point C on the cylinder such that points A, D, and C are collinear. Using some trigonometry, we have \tan \theta = \frac{BC}{AB}; since AB=1, the length of BC is \tan \theta.

After repeating this procedure for every point on the sphere, the cylinder is unrolled into a flat sheet. The result is the following map:

Here, if a point on the earth is known (we know its latitude and longitude), we know its position on the map. The position coordinates of the map, (x,y) can be described in terms of \lambda and \theta in the following equations:

x = \lambda

y = \tan \theta

This projection has several major drawbacks: most notably, as the latitude \theta increases toward 90^\circ, the slope of \tan(\theta) becomes increasingly steep, and the y-coordinates of the point of the map increases towards infinity. At high latitudes, the map is very distorted.

The map is not equal-area — in an equal-area map, shapes of equal area on the sphere have equal area on the projection. Looking at the map we created, this is clearly not the case. Neither is the map conformal, or shape preserving. In a conformal map, the angle between any two lines must be the same as their counterparts on the map.

Mercator’s Projection

One of the most well known cylindrical map projections is Mercator’s projection. Invented by Belgian cartographer Gerardus Mercator in 1569, Mercator addresses some of the problems in older cylindrical projections. Unlike our projection, Mercator’s projection is conformal, and it has been used for navigation for centuries:

However, the equations for Mercator’s projection are rather complicated. The position on Mercator’s projection, when given the latitude and longitude, is:

x = \lambda

y = \ln ( \sec \theta + \tan \theta )

This is a strange and seemingly random formula. But when we use this formula, we get a conformal map. How does this work?

Deriving Mercator’s Projection

Before deriving the formula for Mercator’s Projection, we take a step back and ask ourselves — what was the projection for in the first place?

Unlike most maps today which are used to decorate walls, Mercator’s map was more than that — it was used for nautical navigation. Suppose you were to travel by ship from Vancouver to Honolulu, equipped with nothing more than a compass. Which direction should you travel?

Now if you have a copy of Mercator’s map, this becomes simple. You draw a straight line from Vancouver to Honolulu, and with a protractor you find that the direction is 40^\circ south to the parallel. With a compass, it would be easy for you to maintain this constant heading for the entire duration of the trip.

(Footnote — The straight line on Mercator’s map is not straight on the Earth — it is actually a curved line, called a rhumb line. That is, after taking an initial bearing, one proceeds along the same bearing (relative to the north pole) for the duration of the trip. In the example, you would have to adjust your heading intermittently so that you are always heading at the direction 40^\circ south to the parallel.)

This technique is possible because it relies on one crucial property of Mercator’s map: that it is conformal. It would fail if we tried the same technique on the cylindrical projection given in the introduction.

However, notice that on Mercator’s map, all parallels have the same length: a parallel is simply a line across the width of the entire map. But on Earth, a sphere, parallels are longer when they are closer to the equator and shorter when near the poles:

Let the Earth be a unit sphere with center A. D is a point on the surface of the Earth and B is its corresponding point on the equator. Let C be a point such that CD is parallel to AB and AC is perpendicular to AB. As AB = 1, the circumference of the Earth at the equator is 2 \pi. But the circumference at latitude \theta is only 2\pi \cos \theta.

It is easy to see why: as AB || CD in the diagram, it follows that angles \angle DAB and \angle ADC are equal. Next, \cos \angle ADC = \frac{CD}{AD} = \frac{CD}{1} = CD, therefore CD = \cos \theta. The result follows.

On the map, the length of the equator and the parallel at latitude \theta are equal. But on Earth, the parallel at latitude \theta is smaller than the equator by a factor of \cos \theta; in the projection any line at latitude \theta ends up being stretched horizontally by a factor of \frac{1}{\cos \theta} or \sec \theta.

Finally, in order to satisfy the property of conformality, any part of the map stretched horizontally by some factor, say k, must be stretched vertically by the same factor k. Doing so preserves the angles of that part of the map:

Imagine a very small piece of land at latitude \theta. When it is represented on the map, it is stretched horizontally by \sec \theta. Hence to preserve conformity, it must also be stretched vertically by \sec \theta.

But this can only work if the piece of land has no area. But as all pieces of land, no matter how small, has some area, the latitude of the piece of land is not \theta, but goes from \theta to \theta + \Delta \theta for some small number \Delta \theta. Similarly, the space this piece of land occupies on the map is not just y, but goes from y + \Delta y for some small number \Delta y. Now as \Delta \theta approaches 0, the ratio of \Delta y to \Delta \theta approaches \sec \theta:

\lim_{\Delta \theta \to 0} \frac{\Delta y}{\Delta \theta} = \sec \theta

Or equivalently,

\frac{dy}{d \theta} = \sec \theta

Solving for dy and integrating both sides gives:

\begin{array}{rl}y & = \int \sec \theta \, d \theta \\ & = \int \left( \sec \theta \cdot \frac{\sec \theta + \tan \theta}{\sec \theta + \tan \theta} \right) d \theta \\ & = \int \frac{\sec^2 \theta + \sec \theta \tan \theta}{\sec \theta + \tan \theta} \, d \theta\end{array}

If we let u be the denominator, u = \sec \theta + \tan \theta, it turns out that the numerator is the derivative of u, or du = \sec^2 \theta + \sec \theta \tan \theta. Thus we can substitute u into the integral:

\begin{array}{rl} y & = \int \frac{du}{u} \\ & = \ln |u| + C \\ & = \ln | \sec \theta + \tan \theta | + C \end{array}

We defined the map so that the equator lies on the x-axis, therefore when \theta = 0, y=0. Substituting and solving for C, we find that C=0. Also, since \sec \theta + \tan \theta > 0 for -\frac{\pi}{2} < \theta < \frac{\pi}{2}, we can remove the absolute value brackets. Therefore,

y = \ln ( \sec \theta + \tan \theta )

This is what is desired.

An Alternative Map Projection

Mercator’s Projection is merely a conformal map, not a perfect map. It has several drawbacks — the scale varies greatly from place to place: as the latitude increases, the map gets stretched more and more. This severely distorts larger figures, especially areas near the poles. At a latitude greater than about 70^\circ, the map is practically unusable.

The Gall-Peters Projection is an alternative cylindrical projection:

When a map is used in navigation, it is important that angles are preserved in the map. But when a map is used to present statistical data, it’s often more important to preserve areas in the map: otherwise one may be deceived into thinking one country is bigger than another when they are actually the same size.

The Gall-Peters projection does this, trading conformity for equal-area — two countries of the same area will have the same area on the map, no matter where they are on the map. This map projection is produced by the following equations:

x = \lambda \cos 45^\circ = \frac{\lambda}{\sqrt{2}}

y = \frac{\sin \theta}{\cos 45^\circ} = \sqrt{2} \sin \theta

A discussion on why this set of equations produces an equal-area map is beyond the scope of this essay.


It is a tricky matter to accurately depict the surface of the earth on a flat sheet of paper. The method of cylindrical projection produces several considerably different maps.

The first map we considered uses only simple trigonometry, but the map produced is not very useful. The second map we investigate, Mercator’s Projection, preserves angles on the map. The last map we consider, the Gall-Peters projection, preserves areas on the map.

Each of these map projections are both similar and very different. They are similar in that they are all produced in some way by wrapping a cylinder around a spherical Earth. But the mathematics needed to produce them differ.

Throwing a rock off a cliff: Calculating the optimal angle

You want to throw a rock far, as far as you can possibly throw. Using all your energy, you throw the rock into the sky, and it lands some distance away a few seconds later. But at what angle should you throw it?

It is easy to see that the optimal angle to throw a rock when you are on flat ground is 45^\circ. That is, to make the rock travel as far as possible, one should throw the rock halfway between vertical and horizontal. Any other angle and you do worse than throwing it at a 45^\circ angle.

But what about throwing a rock off a cliff into, say, the ocean? Is 45^\circ still the optimal angle, or is it different? Perhaps it varies?

Suppose you are standing on the edge of a cliff overlooking the sea. It is h meters above the sea, and you can throw a rock with a velocity v, at an angle \theta between 0^\circ and 90^\circ. You are trying to maximize d, the distance thrown.

We assume that there is no friction or air resistance, and we let g be the gravitational constant, 9.81 m/s^2.

When throwing a rock at angle \theta, the velocity of the rock is a vector with a horizontal component of v \cos \theta and a vertical component of v \sin \theta.

Then the distance traveled can be expressed partially in terms of \theta:

d = vt \cos \theta .

The next equation we can get by substituting our values into one of the motion equations:

-h = (v \sin \theta)t - \frac{1}{2}gt^2

The unknown in the two equations is t, or the total time traveled. From the first equation we have

t = \frac{d}{v \cos \theta} .

Substituting into the second equation, we get

-h = (v \sin \theta) (\frac{d}{v \cos \theta}) - \frac{1}{2}g(\frac{d}{v \cos \theta})^2

-h = d \tan \theta - \frac{gd^2 \sec^2 \theta}{2v^2}

From here we want to find the maximum for d, or when \frac{dd}{d \theta}=0. To solve for \frac{dd}{d \theta}, we first differentiate both sides of the equation with respect to \theta:

0 = \frac{d}{d \theta} (d \tan \theta) - \frac{g}{2v^2} \frac{d}{d \theta} (d^2 \sec^2 \theta)

0 = d \sec^2 \theta + \tan \theta \frac{dd}{d \theta} - \frac{g}{2v^2}(\frac{dd^2}{d \theta} \sec^2 \theta + 2d^2 \tan \theta \sec^2 \theta)

As \frac{dd^2}{d \theta} = 2d \frac{dd}{d \theta}, we can isolate \frac{dd}{d \theta}:

\frac{dd}{d \theta} = \frac{\frac{gd^2 \tan \theta \sec^2 \theta}{v^2}-d \sec^2 \theta}{\tan \theta - \frac{gd \sec^2 \theta}{v^2}}

Since \frac{dd}{d \theta} = 0, we eliminate the denominator and are left with

0 = \frac{gd^2 \tan \theta \sec^2 \theta}{v^2} - d \sec^2 \theta

0 = d(\frac{gd \tan \theta \sec^2 \theta}{v^2} - \sec^2 \theta)

Ignoring the trivial case where d=0, we get

gd \tan \theta \sec^2 \theta = v^2 \sec^2 \theta

d = \frac{v^2}{g \tan \theta}

We have the optimal distance; we can substitute this back into a previous equation:

-h = (\frac{v^2}{g \tan \theta}) \tan \theta - \frac{g (\frac{v^2}{g \tan \theta})^2 \sec^2 \theta}{2v^2}

-h = \frac{v^2}{g} - \frac{v^2 \sec^2 \theta}{2g \tan^2 \theta}

-h = \frac{v^2}{g}(1 - \frac{1}{2\sin^2 \theta})

Finally we solve for \theta:

\frac{1}{2 \sin^2 \theta} = \frac{hg}{v^2} + 1

\sin^2 \theta = \frac{v^2}{2hg+2v^2}

\theta = \arcsin(\frac{v}{\sqrt{2hg+2v^2}})

This is the final formula for \theta in terms of v, h, and g.

Investigations on star polygons with a star polygon / stellar number generator

Note 9/3/2010: Redone all of the diagrams

If you are just looking for a download link to the program, here it is. To run, save as Stellar.jar (Java is required)

As school has been over for more than a week now, I’ve had time to (among other things) study geometry. I’m just spending an hour or two every day working through various exercises in Geometry Revisited.

This blog post is long overdue; it had been sitting in my drafts section with a few sentences written since the middle of April. Maybe it’s still useful now, maybe not. Anyways.

I myself am not in IB, because of my hatred of being forced to do copious and tedious amounts of homework. Many of my friends, however, take IB math. One famous project in IB math (not only in my school, but also in other schools) is a large project known as the math portfolio. From what I’ve heard, this consists of a dozen or so pages of writeup on some extended problem.

I first heard about the existence of such a project when a few of my friends asked me for help on it. The project was related to ‘stellar numbers’; a picture of the assignment can be found here. The problem itself seems rather trivial, thus I don’t feel it necessary to discuss it in this blog post.

My investigations here, although inspired by this assignment, is only marginally related to the assignment itself. The end product to be handed in for grades was supposed to include diagrams. However, it isn’t so easy to draw diagrams of stars (see above image) either by hand or on the computer. I saw many of my friends hastily put together their diagrams of stars using Powerpoint or MS Word, using the line and circle tools and a heck lot of copy-and-pasting. The result, unsurprisingly, is messy and inexact. One friend of mine wrote a computer program to output Latex code to draw the stars; the details of how this is done is something I consider more difficult and interesting than the project itself.

Keeping along with my theme of geometry this week, it is useful to figure out the geometry of star polygons before writing programs to generate them.

Some geometric constructions

Consider the five sider star (fig.1):

We can see that the five corners of the star, ABCDE, are the same distance from the center K of the star. This means the five points lie on a circle; furthermore the distances between the corners are the same.

A similar phenomenon occurs with the five points between the corners, FGHIJ. They, too, lie on a circle centered around K.

Let us call the radius of the larger circle r, and the radius of the smaller circle i. Any star can be built around the framework of two circles. For example, our five pointed star (fig.2):

There is, however, one more characteristic present in star polygons that we neglected in the above diagram. It is that we chose the values of r and i pretty much arbitrarily.

Notice that in (fig.1), the sides EF and GB are parallel. It can be shown geometrically that if parallel, EF and GB are also concurrent. In (fig.2), EF and GB are not parallel nor concurrent.

Although any values of i and r for i<r are enough to produce a star, it is necessary to further adjust the ratio between i and r in order to produce a correct star polygon.

Let us denote n to mean the number of points of a star, and d to denote the density of a star.

What do we mean by density? In order to define density, we will take a step back and take a look at alternative ways that star polygons can be constructed.

At present we are looking at star polygons as points on two distinct circles joined together with alternating line segments. Here is a completely different way of presenting it (fig.3):

Points ABCDE are spread equally on a circle. For the five points, every two points are connected by a line. We define density to be this number: the density for this star is 2.

It is possible, perhaps easier, to implement a drawing routine this way. However, by doing this we would have to erase the lines FGHIJ, which seems straightforward but is troublesome in implementation. We are going to stick to the two circles paradigm.

If d=1 , then our ‘star’ would be degenerate; it would merely be a n-sided polygon. On the other hand, we assume that the inner circle has positive radius, thus,

d < \frac{n}{2}

For any integer n with n \geq 3 , the maximum value of d is given by:

d \leq \lfloor \frac{n-1}{2} \rfloor

Next is finding a way to express the ratio between the radii in terms of n and d.

Here we combined the previous two diagrams, and added some new lines. Since our star contains so much symmetry, what follows can be extended to all corners of the star, and also to stars with different densities and number of corners.

Let us call the angle \angle EKA to be x, which is \frac{1}{n} of the circle, or \frac{2 \pi}{n} radians.

Let y be the angle such that

x:x+y = 1:d

.. whence here d=2 , therefore x=y=72^\circ .

As \triangle EKB is isosceles, \angle BEK = \angle EBK , and \angle BEK + \angle EBK + x + y = 180^\circ ; therefore

\angle BEK = \frac{1}{2}\left[ 180-(x+y) \right] = 90-\frac{1}{2}(x+y)

.. whence here \angle BEK = 18^\circ .

Next, \angle EKL = \frac{1}{2}x , in our case 36^\circ , thus the measure of \angle EFK can be found:

\angle EFK = 180 - \left[ 90-\frac{1}{2}(x+y) \right] - \frac{1}{2}x = 90 + \frac{1}{2}y

In our case \angle EFK = 126^\circ.

By the sine law in \triangle EKF,

\frac{r}{\sin (90 + \frac{1}{2}y)} = \frac{i}{\sin \left[90-\frac{1}{2}(x+y) \right]}


i = \frac{r \sin \left[90-\frac{1}{2}(x+y) \right]}{\sin (90 + \frac{1}{2}y)}

Simplifying, we get

i = \frac{r \cos \left[ \frac{1}{2}(x+y) \right]}{\cos (\frac{1}{2}y)}

We defined that x = \frac{2 \pi}{n} and that x+y = \frac{2 \pi d}{n} , thus we write y:

y=\frac{2 \pi (d-1)}{n}

And that brings us to the final formula:

i = \frac{r \cos (\frac{\pi d}{n})}{\cos \left[ \frac{\pi (d-1)}{n} \right]}

Evaluating it for our star, we find that if r=1, then d \approx 0.38197.

Implementation in Java

Now that we have the math done, it is a fairly easy programming exercise to implement this as a java application. Here’s what I came up with after ~2 hours of hacking in Netbeans (animated gif):

This may be useful to future people doing the same assignment (if it ever gets assigned again), or with modifications it may be useful for other purposes too.

To use it, adjust the settings, take a screenshot, crop it, and it may be necessary to further edit it in photoshop.

Without further ado, here’s the application as a jar file:

(broken link removed)

Or the source as a Netbeans project: (37kb)

Do whatever you want with them, and have fun!

Random Math Problems (2)

Okay, these problems are not really random, they’re from the same section of a textbook:

This book is called Geometry Revisited (by Coxeter and Greitzer). It’s often recommended as a olympiad level geometry textbook. I would recommend this book too (although I haven’t gone through the entire book yet).

Section 1.1 of this book proves a theorem referred to as the Extended Law of Sines:

\frac{a}{\sin A} = \frac{b}{\sin B} = \frac{c}{\sin C} = 2R

where R is the circumradius of triangle ABC.

The extended Law of Sines is often given in a truncated form, without mention to the circumradius.

After the proof the book gives four exercises, which I shall rephrase and solve here.

Problem 1

Show that in any triangle \triangle ABC,

a = b \cos C + c \cos B

From this, deduce the ‘addition formula’ for sines:

\sin (B+C) = \sin B \cos C + \sin C \cos B


From the Law of Cosines:

\begin{array}{rcl} c^2 &=& a^2 + b^2 - 2ab \cos C \\ b^2 &=& a^2 + c^2 - 2ac \cos B \end{array}

Adding the two:

\begin{array}{rcl} c^2 + b^2 &=& 2a^2 + c^2 + b^2 - 2ab \cos C - 2ac \cos B \\ 2a^2 &=& 2ab \cos C + 2ac \cos B \\ a &=& b \cos C + c \cos B \end{array}

as desired.

We now prove the addition formula. From the law of sines:

\begin{array}{rcl} \frac{a}{\sin A} &=& 2R \\ a &=& 2R \sin A \end{array}


\begin{array}{rcl} b &=& 2R \sin B \\ c &=& 2R \sin C \end{array}

Substituting into our previous equation:

\begin{array}{rcl} a &=& b \cos C + c \cos B \\ 2R \sin A &=& 2R \sin B \cos C + 2R \sin C \cos B \\ \sin A &=& \sin B \cos C + \sin C \cos B \end{array}

Since \triangle ABC is a triangle, A+B+C = 180. Thus, A = 180 - (B+C).

Note that \sin \theta = \sin (180 - \theta). Thus,

\sin (B+C) = \sin A = \sin B \cos C + \sin C \cos B

as desired.

Problem 2

Show that in any triangle \triangle ABC,

a (\sin B - \sin C) + b (\sin C - \sin A) + c(\sin A -\sin B) = 0


Using the law of sines again:

\begin{array}{rcl} \frac{a}{\sin A} &=& \frac{b}{\sin B} \\ a \sin B &=& b \sin A \\ a \sin B - b \sin A &=& 0 \end{array}


\begin{array}{rcl} b \sin C - c \sin B &=& 0 \\ c \sin A - a \sin C &=& 0 \end{array}

Adding the three equations:

\begin{array}{rcl} a \sin B - b \sin A + b \sin C - c \sin B + c \sin A - a \sin C &=& 0 \\ a (\sin B - \sin C) + b(\sin C - \sin A) + c(\sin A - \sin B) &=& 0 \end{array}

as desired.

Problem 3

Show that in any triangle \triangle ABC,

[ABC] = \frac{abc}{4R}

(although the book uses (ABC), I’ve accustomed to writing [ABC] on my blog to denote area of \triangle ABC)


Let h be the height of \triangle ABC (from A perpendicular to BC), so that [ABC] = \frac{1}{2} ah:

Using trigonometry,

h = b \sin C

so we can write the area as

[ABC] = \frac{1}{2} ab \sin C

From the extended Law of Sines,

\begin{array}{rcl} \frac{c}{\sin C} &=& 2R \\ \sin C &=& \frac{c}{2R} \end{array}

Now we can rewrite the area again:

\begin{array}{rcl} [ABC] &=& \frac{1}{2}ab \frac{c}{2R} \\ &=& \frac{abc}{4R} \end{array}

as desired.

Problem 4

For a triangle \triangle ABC, let p be the radius of a circle going through A and tangent to BC at B. Let q be the radius of a circle going through A and tangent to BC at C.

Show that:

pq = R^2


Let M be the center of the circle with radius p, and N be the center of the circle with radius q. Also K and L are two points on circles M and N respectively:

Using the law of sines on circle N and M, we have:

\begin{array}{rcl} \frac{b}{\sin L} &=& 2q \\ \frac{c}{\sin K} &=& 2p \end{array}

Because BC is tangent to the two circles, \angle L = \angle C and \angle K = \angle B.


\begin{array}{rcl} \frac{b}{\sin C} = 2q  & \Rightarrow & b = 2q \sin C \\ \frac{c}{\sin B} = 2p & \Rightarrow & c = 2p \sin B \end{array}

Now applying the Law of Sines again in \triangle ABC:

\begin{array}{rcl} \frac{b}{\sin B} = 2R & \Rightarrow & \sin B = \frac{b}{2R} \\ \frac{c}{\sin C} = 2R & \Rightarrow & \sin C = \frac{c}{2R} \end{array}

Substituting these values back into the previous equation:

\begin{array}{rcl} b = 2q \frac{c}{2R} = \frac{qc}{R} & \Rightarrow & q = \frac{bR}{c} \\ c = 2p \frac{b}{2R} = \frac{pb}{R} & \Rightarrow & p = \frac{cR}{b}\end{array}

By multiplication,

pq = \frac{bR}{c} \frac{cR}{b} = R^2

as desired.

A Geometry Exercise

This post was inspired by a problem I came across when studying for a math contest, namely problem 23 in the 2007 Cayley Contest.

I think this is an interesting problem because the solution is not so obvious, and there are several incorrect approaches that I’ve tried (and given up on). I’ve modified the problem slightly to make it more interesting.

Here we have a rectangle, ABCD. It’s rotated to the right by \theta degrees, using A as the pivot. This forms AEFG as the new rectangle. We are given x and y, we are asked to find the area of the shaded region.

Can you come up with a general formula for the area of the shaded region, using the given variables x, y, and \theta?


I will use the notation (ABC) to denote the area of triangle \triangle ABC.

The solution is to draw a line parallel to BC through E, the blue line. This splits the lower shaded area into two triangles: \triangle AXE and \triangle EYH, and a rectangle: BCYX.

Since ABCD is the same as AEFG, the top shaded area is the same as the bottom shaded area.

Let’s find the area of triangle \triangle AXE.

The angle \angle BAE is equal to \theta. Of course, AE is equal to y. Using some basic trigonometry, AX = y \cos \theta, XE = y \sin \theta, and BX = y - y \cos \theta. Now we know the area of triangle \triangle AXE:

\begin{array}{rcl} (AXE) &=& \frac{1}{2} (y \cos \theta) (y \sin \theta) \\ &=& \frac{1}{4} y^2 \sin(2 \theta)\end{array}

Next we find the area of \triangle EYH. EY = x - y \sin \theta, and using trigonometry again, HY = \tan \theta (x - y \sin \theta).

\begin{array}{rcl} (EYH) &=& \frac{1}{2} (x - y \sin \theta) [\tan \theta (x - y \sin \theta)] \\ &=& \frac{1}{2} (x - y\sin \theta)^2 \tan \theta \end{array}

Finally, for the rectangle BCYX:

\begin{array}{rcl} (BCYX) &=& x (y - y \cos \theta) \end{array}

Now we add all this up and simplify (T is total):

\begin{array}{rcl} T &=& 2[(AXE) + (EYH) + (BCYX)] \\ &=& 2[\frac{1}{4} y^2 \sin(2 \theta) + \frac{1}{2}(x - y \sin \theta)^2 \tan \theta + x (y - y \cos \theta)] \\ &=& \tan \theta (x^2+y^2) - 2xy (\sec \theta - 1)\end{array}

It’s rather tedious to get from the second to the last step, but Wolfram Alpha could do it for me.

Notice that we assume that EF intersects DC. If \theta is large enough, then EF would intersect AD, and this formula would no longer work.

Of course it would also fail if y is greater than x, so my formula works for only a rather limited range of values.