Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

linestring getLength not returning expected value #3533

Closed
alexandre-melard opened this issue Apr 9, 2015 · 48 comments
Closed

linestring getLength not returning expected value #3533

alexandre-melard opened this issue Apr 9, 2015 · 48 comments
Labels

Comments

@alexandre-melard
Copy link

I am using getLength to retrieve the linestring length.

For the same segment:
1- when using google map measure tool, I get 228m
2- when using IGN geoportail measure tool, I get 228m
3- when I use e.feature.getGeometry().getLength() I get 330m

Here are the flat coordinates:
e.feature.getGeometry().getFlatCoordinates() :
[571382.4214041593, 5723486.068714521, 571593.8175605105, 5723741.65502785]

in 4326:
[5.132815622245775, 45.644023326845485, 5.134714626228319, 45.64562844964627]

When I check the coordinates position on either ol3 or google map, I get the same points. The difference must come from the calcul...

Did I miss something and I should not use the getLength method? Please give me some direction if you think this is not an issue.

@alexandre-melard
Copy link
Author

I've done many tests and I get a 1.44 factor all the time between the distances. Maybe it's a unit problem. But I though the default unit in ol3 was meter...

@alexandre-melard
Copy link
Author

well for the moment, I'll just divide my result by 1.44

@bartvde
Copy link
Member

bartvde commented Apr 10, 2015

I am not sure but maybe the other frameworks use geodesic measurements? Check the measure example in ol3 which has a geodesic option.

@bartvde
Copy link
Member

bartvde commented Apr 10, 2015

@bartvde bartvde closed this as completed Apr 10, 2015
@bill-chadwick
Copy link
Contributor

I think Your coordinate spacing is too small for geodesic to make a difference. Do you get a better result with a latitude near the equator?

@alexandre-melard
Copy link
Author

@bartvde Crossposting is a good way to rush people, I don't expect everybody to read the ol3 issues.

I've changed my length calcul to use geodesic, I get an exact answer now.

    var wgs84Sphere = new ol.Sphere(6378137);
    var coordinates = e.feature.getGeometry().getCoordinates();
    var length = 0;
    var sourceProj = map.getView().getProjection();
    for (var i = 0, ii = coordinates.length - 1; i < ii; ++i) {
      var c1 = ol.proj.transform(coordinates[i], sourceProj, 'EPSG:4326');
      var c2 = ol.proj.transform(coordinates[i + 1], sourceProj, 'EPSG:4326');
      length += wgs84Sphere.haversineDistance(c1, c2);
    }

Now, I don't understand why the getLength doesn't return the same value.
And I am also a bit worried by the fact that all ol.Sphere functions (even the constructor) are experimental.
(see http://openlayers.org/en/v3.4.0/apidoc/ol.Sphere.html?unstable=true)

@bartvde
Copy link
Member

bartvde commented Apr 10, 2015

No cross-posting is actually quite rude, and increases your chances of getting ignored.
Most developers use github notifications and are subscribed to the openlayers-3 tag on stack overflow, they will get notified twice.

Please only open up a github issue once you've confirmed and are sure it's a bug in the library.

@alexandre-melard
Copy link
Author

Ok, fine, I was not aware of your organisation, so is it better to post on stackexchange first, and then post an issue here if the bug is confirmed? I didn't mean to be rude but just to help by alerting on a suspect behavior

@marcjansen
Copy link
Member

@mylen Yes, questions go to stackoverflow, and if a bug / misbehaviour is confirmed there, open an issue here. If you are quite sure you have discovered a bug, you may skip the stackoverflow question, but of course it wouldn't hurt to ask there first.

@bill-chadwick
Copy link
Contributor

Mylen, you get the factor of 1.4 on EPSG:3857 (Web Mercator) because 1 / cos (latitude of 45.6 N) is about that and your points are around 45.6 degrees North.

With Web Mercator, scale changes with latitude (see the Wikipedia article on the Mercator Projection for more details). The method ol.geom.LineString.prototype.getLength will only work (give accurate answers in metres) for short lines near the equator.

Currently, the comment / API doc for the method is this

/**

  • @return {number} Length (on projected plane).
  • @api stable
    */

I think that is not telling the whole story - IMHO a stronger health warning is needed. It should perhaps say that the returned units are not always metres. There are other object methods in ol.geom that will have similar issues e.g. Polygon.getArea

Perhaps this issue should be reopened until the docs are sorted.

Personally I am not convinced that the methods should be 'api stable' as they are so dangerous for use with Web Merctor. Even relative comparisons of lengths / areas are not possible with geometries at significantly different latitudes.

@marcjansen
Copy link
Member

Pull requests are always welcome.

@bill-chadwick
Copy link
Contributor

I could do that. Its a fair fit with my pending graticule work.
On 10 Apr 2015 16:42, "Marc Jansen" notifications@github.com wrote:

Pull requests are always welcome.


Reply to this email directly or view it on GitHub
#3533 (comment).

@probins
Copy link
Contributor

probins commented Apr 10, 2015

I am not convinced that the methods should be 'api stable' as they are so dangerous for use with Web Merctor

I think the assumption is (and IMO should be) that those using the API know the limitations of using global projections like Mercator. The method is 'stable' in that it's not going to change. Whether it's useful or not is another matter :-) IMO OL2 was more helpful in this respect in having getGeodesicLength/Area methods (Vincenty rather than Haversine) as well as simple planar ones. Haversine is exported currently. Vincenty is available in the source (under ellipsoid/), but this is not currently exported in the standard build.

@alexandre-melard
Copy link
Author

@probins, when you make an assumption :

the assumption is (and IMO should be) that those using the API know the limitations of using global projections like Mercator

You mean that the usage of openlayers is targeted for expert? I'm not sure that only geographers are using it .
If there is a problem with length calculation in the Mercator projection, why is it used as default on the map, could you give me an example of a use of getLength in a mercator projection?

@bill-chadwick
Copy link
Contributor

I agree with Mylen. More of a worry to me is that there is no good/easy solution.

The Geometry classes do not record the projection of their points, so that length and area methods returning answers in metres and square metres (or other projection units) can not be implemented by them without significant changes. Even if they did know their own projection, we would not want per projection calculations in the geometry classes. IMHO the length / area maths belongs in the projection classes.

Either the length and area methods in ol.geom.Geometry classes would have to take an ol.proj.ProjectionLike parameter so that appropriate maths can be done, or else the area and length methods are moved elsewhere - e.g. to the ol.proj.Projection objects.

The latter make most sense to me. An ol.proj. Projection could be required to provide functions to calculate a Geometry length / perimeter, in projection units (often metres), and its area, in square projection units (often square metres). Having to pass a ProjectionLike into ol.geom.Geometry class length and area methods when the developer already knows (or should) the projection of the coordinates in the Geometry, would seem odd.

For 'local' projections like UTM and the UK National Grid (which have very near uniform scale across their extent) , simple maths as already exists in ol.geom.flat.lengthflatgeom and ol.geom.flat.areaflatgeom, will suffice. Perhaps some of the old OL2 code could be reused for the 'global' Web Mercator (3857) and Lat/Lon (4326) length and area calculations.

@probins
Copy link
Contributor

probins commented Apr 11, 2015

but won't the methods in ol.Sphere/ol.Ellipsoid, as used in the measure example, suffice? This is independent of projection. Do you actually need per projection calculations? I'd agree that it's potentially confusing for people who aren't aware of these issues (99+% of the population) that getLength isn't appropriate for 3857, but I don't think it's the job of API docs to discuss these issues. They could usefully point people to the measure example when using 3857 and similar projections. And, as I say, IMO the API would be better if there were (stable) geodesic methods on the geometries as in OL2.

@ahocevar
Copy link
Member

Re-opening because this is turning into a good discussion.

@ahocevar ahocevar reopened this Apr 13, 2015
@bill-chadwick
Copy link
Contributor

Thanks for re-opening. There are several solution options to think about.

@probins
Copy link
Contributor

probins commented Apr 13, 2015

the text in the measure example could be improved too. "Earth curvature is not taken into account" isn't really correct. The whole point of projections is that they do take curvature into account. The issue is that they have different priorities. Mercator isn't intended for measuring distance or area, so getLength/Area isn't appropriate here. On the other hand, I would be pretty sure that getLength, for example, on the OSGB projection that @bill-chadwick mentions will be more accurate than Haversine.

@bill-chadwick
Copy link
Contributor

I am coming round to the idea that the Geometry classes should know the projection of their coordinates. We could add an optional projectionLike parameter to the Geometry constructors (defaulting to EPSG:3857) and store the projectionLike as a class member. This new member would be updated if the Geometry is transformed using ol.geom.Geometry.prototype.transform but would need a setter method call after using ol.geom.Geometry.prototype.applyTransform (or the applyTransform method should be removed).

If the Geometry classes know their projection, then the getLength and getArea methods could use ol.Sphere / ol.Ellipsoid methods for projections where ol.proj.Projection.prototype.isGlobal returns true and the existing ol.geom.flat.lengthflatgeom and ol.geom.flat.areaflatgeom if isGlobal returns false. The returned length units would normally be metres and the area units metres x metres. However for non-global projections with other distance units, like feet, then the returns would be in e.g. feet or feet x feet.

@ahocevar
Copy link
Member

We have discussed to let geometries know their projection in the past, but decided against it for several reasons.

@probins
Copy link
Contributor

probins commented Apr 13, 2015

alternatively, you could move that decision into userland, and pass in a parameter: say, projected plane as default, or 'h' for Haversine, 'v' for Vincenty: getLength('h')

Another oddity of getLength at the moment is that if your view is in 4326 you get a 'length' in decimal degrees - I struggle to think how you could use that :-)

@bill-chadwick
Copy link
Contributor

Vincenty and Haversine work on 4326 coords, how would we know the projection to transform from ?

@probins
Copy link
Contributor

probins commented Apr 13, 2015

true. 'Pass in 2 parameters ...'

@ahocevar
Copy link
Member

true. 'Pass in 2 parameters ...'

Makes sense to me, and would be consistent with other functions.

@bill-chadwick
Copy link
Contributor

Should the params should be optional and the defaults be Haversine and 3857 ?

@bill-chadwick
Copy link
Contributor

or perhaps default to Haversine for 3857 and 4326 and Projected Plane for non 'global' projections.

@ahocevar
Copy link
Member

Haversine as default yes, but no default for the projection. So projection mandatory, and method optional.

@bill-chadwick
Copy link
Contributor

I did a bit on this last night and plan to have a PR ready in a few days. I will add some new length and area method tests.

I propose calculation methods of 'CARTESIAN', 'HAVERSINE' and 'VINCENTY' though 'CARTESIAN', 'SPHERICAL' and 'ELLIPSOIDAL' might do instead (any opinions as to which are best?). The functions do not need a projection for the 'CARTESIAN' method though they could return NaN if projection units are degrees (4326) and the method is 'CARTESIAN'.

Calculation method will default to 'HAVERSINE' for global projections and 'CARTESIAN' for non global projections.

@probins
Copy link
Contributor

probins commented Apr 14, 2015

another issue here is that of units. A really friendly function would give you the ability to define the units 'give me the length in miles'. That could be kept in userland though. Units can also be 'pixels' in which case only cartesian.

IMO 'haversine' and 'vincenty' is better than spherical/ellipsoid, as it's more specific. And I would prefer an abbreviated form like 'h' and 'v' - because I don't like typing :-)

@probins
Copy link
Contributor

probins commented Apr 14, 2015

no default for the projection. So projection mandatory

this will always be the view projection though, won't it? Perhaps add a note in the docs that this will normally be view.getProjection()

@bill-chadwick
Copy link
Contributor

I guess we would want an enum of geometry distance units and a default units parameter value of metres. We can't really reuse the ol.proj.Units enum as it includes degrees. We would have acres, hectares etc for geometry area units.

@bill-chadwick
Copy link
Contributor

To avoid bloat, I think perhaps the length and area methods should always return metres and metres x metres, with any unit conversions being done by the user as required. A unit parameter could always be added later if deemed necessary.

@probins
Copy link
Contributor

probins commented Apr 14, 2015

projection mandatory

that also makes this a breaking change, so there'll have to be release notes for this too.

@probins
Copy link
Contributor

probins commented Apr 14, 2015

the text in the measure example could be improved

I was going to submit a PR for this, but the example should be changed to reflect this new method anyway (it will be simpler), so I'll hold off on this.

@bill-chadwick
Copy link
Contributor

AARGHH!

For convenience, I would like to build ol3 on Windows but the build system is so fragile it drives me mad.

There is some issue with generate-info.js and ('child_process').spawn that is very unreliable. Sometimes (occasionally) the build works other (most) times it fails without generating info.json.
A fix would be MOST welcome.

Incidentally we will need to look at / update ol.interaction.Draw as it calls LineString.getLength

@ahocevar
Copy link
Member

Improvements to the build system are an ongoing effort and should land in master soon. As a result, building on Windows should be possible again soon.

@bill-chadwick
Copy link
Contributor

Looks like ol.interaction.Draw is going to unravel my plans a bit.

ol.interaction.Draw uses geometry.setRadius(sketchLineGeom.getLength()); to set the radius of a circle. We can use event.map.getView().getProjection() to find the projection to use for an updated lineString.getLength method in ol.interaction.Draw.prototype.modifyDrawing_ but ...

For the draw interaction to work on a 4326 View I think we would have to allow the use of the Cartesian length calculation method for the 4326 projection (and other projections with units of degrees). To draw circles with non degrees projections requires that lineString.getLength must always return its result in projection units.

@bill-chadwick
Copy link
Contributor

http://www.andre-simon.de/doku/ansifilter/en/ansifilter.php is useful for looking at test output on Windows.

@bill-chadwick
Copy link
Contributor

I have committed the getLength changes to a branch here

bill-chadwick@49ca111

the getArea change is still TODO.

I have added several tests for the updated getLength method

Feedback welcome.

To think about some more ...

4326 geometry simplification (Doouglas Peuker) at extreme latitudes (where x coords are smaller than y coords)

ol.geom.Circle.setRadius - 'The radius is in the units of the projection.' is a misleading comment for 3857

@bill-chadwick
Copy link
Contributor

For 3857, we can calculate LineString length as a summation of line segment rhumbline distances (on the 3857 sphere). I believe this should be the default method for 3857 as it is the length of the lines as plotted on the map, rather than a summation of great circle segment length.

see

http://en.wikipedia.org/wiki/Mercator_projection - Formulae for distance section.

This should have better performance than using Haversine or Vincenty.

I think the maths for rhumb line length should go in epsg3857projection.js

@probins
Copy link
Contributor

probins commented Apr 17, 2015

I'm not sure this is worth it: I think these days build size and latency are more the issue than performing calculations. And an advantage of having Haversine as default is it should produce the same result as measuring the same line on e.g. GMaps.

@bill-chadwick
Copy link
Contributor

I will document that Haversine and Vincenty are great circle methods. Calculation of Rhumb line (loxodrome) length from lat/lon and EPSG:3857 coords uses very similar maths. So for completeness, I do plan to add another Geometry Length Calculation method for Rhumb Lines on the sphere, though will leave out Ellipsoidal Rhumb line maths for now.

see http://www.movable-type.co.uk/scripts/latlong.html

which I plan to use for data for the tests.

The existing measure example at

http://openlayers.org/en/v3.4.0/examples/measure.html?q=measure

is really rather wrong in several respects.

@bill-chadwick
Copy link
Contributor

I plan to add a rhumb line length method to ol.Sphere where there is currently the comment

// FIXME add rhumb lines

@bill-chadwick
Copy link
Contributor

I have added spherical Rhumb Line as a length calculation method. I now need to add rhumb line tests for sphere.js where the new rhumb length maths is.

@bill-chadwick
Copy link
Contributor

PR here
#3609

@bill-chadwick
Copy link
Contributor

Updated PR here
#3609

@stale
Copy link

stale bot commented May 22, 2019

This issue has been automatically marked as stale because it has not had recent activity. It will be closed if no further activity occurs. Thank you for your contributions.

@stale stale bot added the stale label May 22, 2019
@stale stale bot closed this as completed May 29, 2019
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
Labels
Projects
None yet
Development

No branches or pull requests

6 participants