6

I have a set of Long and Lat co-ordinates making up a five point polygon which I have converted in excel using the following formula:

=((X1*Y2)-(X2*Y1)+(X2*Y3)-(X3*Y2)+(X3*Y4)-(X4*Y3)+(X4*Y5)-(X5*Y4)+(X5*Y1)-(X1*Y5))*0.5

*Replacing X,Y with cell references of course.

This gives me a negative number but how do I convert this number to square miles or kilometres?

Improve this question
4
  • Thank you a lot for this, I have spent some time getting the right calculations - I Hope! However would you mind checking my work? It's allways nice to have a greater mind then mine give it a once over. I used the following site to easily get my co-ordinates : itouchmap.com/latlong.html I then input into the above formula, then worked out the area seperately, then to find the square miles I multiplied by 3.861E-07 or 0.00000003861 to get the result. However my result is 1789415 square miles and I am positive my area is not that big, any idea where I am going wrong? Thank you.
    – user2660
    Commented Apr 13, 2011 at 15:51
  • I just tried the formula I posted to try and approximate the area of my home city, and it worked fine. You need to make sure your last point is the same as your first point (so a rectangle is actually a polygon with 5 points). I suspect that is your issue.
    – Sasa Ivetic
    Commented Apr 13, 2011 at 17:35
  • You're probably not closing the ring: In the formula, (x5, y5) must be the same as (x1, y1). BTW, if your polygon doesn't circle a pole, you can omit the "2+" terms in the formula: they all cancel out. (Including them actually degrades numerical precision near the equator.)
    – whuber
    Commented Apr 13, 2011 at 20:46
  • BTW, please use the Post answer button only for actual answers. You should modify your original question to add additional information.
    – whuber
    Commented Apr 13, 2011 at 20:47

3 Answers 3

Reset to default
11

Try this formula (assuming your source is WGS1984, if not then you'll need to adjust the ellipsoid used by the second line):

area = rad(x2 - x1) * (2 + sin(rad(y1)) + sin(rad(y2))) + rad(x3 - x2) * (2 + sin(rad(y2)) + sin(rad(y3))) + rad(x4 - x3) * (2 + sin(rad(y3)) + sin(rad(y4))) + rad(x5 - x4) * (2 + sin(rad(y4)) + sin(rad(y5)))

area = abs(area * 6378137.0 * 6378137.0 / 2.0)

rad() is a function that converts Degrees to Radians (i.e. Degrees * PI / 180)

Source: OpenLayers LinearRing

This will result in an area returned in square meters.

Improve this answer
10
  • @whuber: Thanks for noticing the missing parenthesis!
    – Sasa Ivetic
    Commented Apr 13, 2011 at 14:39
  • 3
    +1 Quite effective. This formula applies a cylindrical equal-area projection (for a sphere, not an ellipsoid) to the vertices and then computes the area of the polygon formed by the projected vertices. This makes it a good approximation for small polygons: it is approximately as accurate as the spherical datum. Note that it will fail for large polygons (whose edges span more than a few degrees apiece) or for ones that cross the longitudinal cut (usually at +-180 degrees). The latter problem is fixable if 'x2-x1' etc. are computed modulo 360 degrees.
    – whuber
    Commented Apr 13, 2011 at 14:40
  • @whuber: Quick and effective, but like you mentioned it is NOT an exact area, instead it is an approximation. Very important point to remember.
    – Sasa Ivetic
    Commented Apr 13, 2011 at 14:48
  • 1
    I'm honestly not sure if it does have a name. Today I quickly looked it up in OpenLayers code (though I've used it before), and they cite a paper by Chamberlain and Duquette which you can get here: trs-new.jpl.nasa.gov/dspace/bitstream/2014/40409/1/07-03.pdf
    – Sasa Ivetic
    Commented Apr 13, 2011 at 15:18
  • 1
    @whuber - are you suggesting this isn't rocket science? :)
    – Kirk Kuykendall
    Commented Apr 14, 2011 at 13:15
3

I think Sasa's answer above is pretty solid.

A quick Measure: A rough and ready sanity check to I use to double check complicated math is the rough width / height of the shape. (ie its bounding box) At most inhabited latitudes .01 degrees lat lng is approximately .5 to .7 miles or roughly 1km. So a shape .01 x .01 would be about 1km+/- or .25 to .5 sq miles. This math will go haywire at the poles and international date line, so its just a rough guide. Depending on the type of shape it should be some rough percentage of the overall bounding box.

Example Below is the same measure (approximate hand drawn) in Boston and Atlanta for comparison.

Atlanta GA enter image description here

Boston, MA

enter image description here

Improve this answer
2
  • 1
    Not a bad approach, but it could use refinement: .01 degree latitude is 1.1 km while .01 degree of longitude is the same amount multiplied by the cosine of the latitude. It's easy to memorize a few cosines, especially for the latitudes where one tends to work.
    – whuber
    Commented Apr 13, 2011 at 20:01
  • Good idea. I like the simplicity the cosine conversion.
    – Glenn
    Commented Apr 14, 2011 at 10:22
1

(Somewhat elderly topic, I realize, but an interesting and useful one.)

If the earth were a sphere, you could get an analytically exact answer by splitting the region into triangles and computing the area of each triangle. You'd do that by computing the angle at each vertex, summing it up, and seeing by how much it exceeds 180 degrees = pi radians. If the "spherical excess" is theta degrees, the area is theta * r^2, where r=earth radius in your desired units.

This has two flaws. If the area is small, it'll be hard to avoid roundoff error (theta+pi, the total of the three vertex angles, will be very close to pi; subtract to get theta, and you're subtracting two nearly-equal quantities.) Also, the earth is not a sphere.

Best I can see for small areas would be to compute the sides of the triangles using, say, the method of Vicenty. (This allows one to compute quite precise distances and angles on an ellipsoid.) Then compute the area of a planar triangle with those sides. The error will be vanishingly low for small triangles.

For larger triangles... well, the best I can come up with is a recursive method that splits them into lots of small triangles :

area of a triangle :

If all three sides are below some threshhold in length, compute the area of a planar triangle with those sides, and you're done.

Otherwise, find the midpoint of the longest side. (The Vicenty algorithm makes that possible.) Say that point M lies on AB, opposite vertex C. Create two triangles, AMC and CMB, each with nearly half the area of the starting triangle ABC. Return the sum of the areas of AMC and CMB.

Improve this answer

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.