2

I am trying to measure different angles on 2D map using lat/lon coordinates taken from Google Maps, and Java library JTS.

I have noticed a weird thing: my measurements are always distorted, except angles that are perpendicular to meridians and parallels. It gets as high as 25 degrees when the angle sides are at an angle of 45 degrees to the parallalels, example here:

enter image description here

The result is 68.97 and the library I use uses the law of cosines to measure the angle. I also double checked the result using the method provided here, and the difference was still there.

Here's the Java code:

double point2X = 2.35956;
double point2Y = 48.86054;

double point1X = 2.35922;
double point1Y = 48.86215;

double fixedX = 2.3582;
double fixedY = 48.86124;


double angle1 = Math.atan2(point1Y - fixedY, point1X - fixedX);
double angle2 = Math.atan2(point2Y - fixedY, point2X - fixedX);

return Math.toDegrees(angle1 - angle2);

Or as done using JTS:

Coordinate point1 = new Coordinate(2.35922, 48.86215);
Coordinate point2 = new Coordinate(2.35956, 48.86054);
Coordinate fixed = new Coordinate(2.3582, 48.86124);

return Math.toDegrees(Angle.angleBetween(point1, fixed, point2));

And here is a JS fiddle showing the calculations:

https://jsfiddle.net/hrg04v3n/

I found out that Web Mercator is non-conformal, thus it doesn't preserve angles, but does it get as bad for even such small areas as presented above?

  • 1
    You've only provided a portion of the necessary code in the question. Where is the JTS code? – Vince Oct 7 '18 at 13:23
  • This doesn't have anything to do with Web Mercator, since your code is all in degrees (GCS). Your Law of Cosines code is spherical while the JTS code is likely spheroidal, which will make a slight difference. I don't have a code library immediately available, but the fact that cosine (48.86) = 0.6579 indicates that your X component is ~2/3 of Y in degrees (and Y is ~3/2 exaggerated in Mercator), which makes the calculations likely correct -- If the map looks square in Mercator at that latitude, it's actually rectangular (which is what non-conformal means). – Vince Oct 7 '18 at 16:49
  • @Vince do you mean, by saying that "it's actually rectangular" that I shouldn't expect 90 degrees angle there? But even checking on Google Earth - the streets there are exactly perpendicular. – Broccoli Oct 7 '18 at 19:11
2

Geographic (lat/lon) coordinates are already angular measures, so you can't compute a simple arctangent with the coordinates based on atan2(delta_y,delta_x) as if they were Cartesian values.

Using a port of the US Geodetic Survey FORTRAN code to solve the inverse (aka reverse) geodetic problem, and using the coordinate values in your Java code (assuming WGS84 datum), I got the following:

p1                     |p2                      |gcbearing
POINT (2.3582 48.86124)|POINT (2.35922 48.86215)|36.4845
POINT (2.3582 48.86124)|POINT (2.35956 48.86054)|127.9567

The difference between the bearings is 91.4722 degrees.

Reprojecting the points to WGS84 UTM Zone 31N and using Cartesian bearing, I got:

p1                            |p2                            |bearing
POINT (452926.994 5412229.237)|POINT (453002.658 5412329.763)|36.9681
POINT (452926.994 5412229.237)|POINT (453026.088 5412150.583)|128.4401

The difference here is 91.4720 degrees (since Paris is on the left side of the UTM cylinder the bearings are slightly different, as well)

Finally, reprojecting into a custom WGS84 datum flavor of the Europe Albers Equal Area Conic projection (central meridian 10.0, standard parallels 43.0,62.0), I got:

p1                            |p2                             |bearing
POINT ( -553711.962 2098386.8)|POINT (-553627.667 2098480.823)|41.8774
POINT ( -553711.962 2098386.8)|POINT (-553621.972 2098298.257)|134.5356

with a difference of 92.6582 degrees.

You cannot rely on a Web Mercator canvas for whether the intersection "looks like" a right angle because of the significant Y-axis distortion in Web Mercator, even at only 49N.

I don't have any experience in JTS, but it apparently only has Cartesian (not geodetic) functions, so you'd need to reproject to compute bearing changes (based on the above I'd recommend using UTM for relatively small areas).

  • 1
    JTS is purely cartesian, so you would need to use a local projection – Ian Turton Oct 9 '18 at 7:38

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.