Converting Longitude and Latitude Coordinates to Square Miles?

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.

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

Boston, MA

(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.