
Geometry : 
Subject: Re: [geometry] Support for geographic coordinate system
From: Adam Wulkiewicz (adam.wulkiewicz_at_[hidden])
Date: 20141126 15:25:28
Hi Barend,
Barend Gehrels wrote:
> Adam Wulkiewicz wrote On 26112014 0:16:
>> Barend Gehrels wrote:
>>> Thanks! But...
>>> What is the conclusion?
>>>  SSF can only be used for spherical and not for geographic (all
>>> nonspheres)?
>>
>> It can be used but will give wrong results for spheroids.
>
> If I ask if it can be used, I mean of course with correct results. So
> your answer is no.
>
Well, yes and no. It's more about choosing a default for geographic CS.
Similar but not the same  distance Andoyer vs Vincenty vs haversine.
But in the case of distance often a distance between far Points is
measured. On the other hand, in the case of side calculation, I'm
guessing that it isn't very common to define Polygons containing
150kmlong segments.
>
>> SSF gives the same results as Vincenty for sphere/flattening=0.
>>
>>>  the method comparing azimuths (you mentioned is probably not good
>>> enough) is not sufficient?
>>>
>>
>> It's because the further the Point the calculation becomes less
>> accurate. A "real" crosstrack distance should probably be calculated
>> and compared with EPS (if needed), similar to side_by_cross_track (or
>> side_by_azimuth). However even there the distance nor radius isn't
>> taken into account and without it we can't calculate a value of XTD.
>> But maybe it's sufficient to do it this way... At least for doubles,
>> I'm guessing that for floats it'd be different.
>
> So for points up to about 0.6 meter left or right of a segment of 153
> km, there is a difference in ssf and Vincenty. And you state that
> Vincenty is correct and ssf is wrong, and you mention the Vincenty
> calculation is less accurate (but maybe not in this case). If doubles,
> for floats it is unknown. Besides that, Vincenty is very good, but it
> is an approximation.
>
> Just a question, I'm just curious  how do you know for sure then that
> the Vincenty approximation is correct here, and the ssf results are
> wrong? It is not that easy to verify...
>
Of course in all formulas some approximation of a globe is used. My
assessment is based on a fact that SSF assume thet a globe is a sphere
and Vincenty that it's a spheroid. So the approximation is closer to
reality. And the greater the flattening, the greater the difference
between SSF and Vincenty. I'm not considering precision, numerical
errors, stability etc.
As for the Vincenty formula here: http://www.icsm.gov.au/gda/gdav2.3.pdf
states:
"Vincenty's formulae (Vincenty, 1975) may be used for lines ranging from
a few cm to nearly 20,000 km, with millimetre accuracy."
Calculation done with SSF is probably also precise but it uses a sphere
model, that's all.
And as mentioned above, I feel that the above is more important for
calculating distances. Often the distance between far Points is
measured. I don't know if we need such great precision for the sides
calculation. Who defines a Polygon containing 153kmlong segments? And
is an error of 0.6m a big deal in this case?
Well, maybe on Saturn which has flattening=0.1 and is bigger ;)
>
>> Do you mean, when a point is on a segment? The epsilon is so small
>> that all test points are missing the segment (at least for double).
>>
>> Btw, in spherical_side_formula EPS isn't used in the final result
>> calculation:
>>
>> returndist>zero?1
>> : dist < zero ? 1
>> : 0;
>>
>> AFAIU dist should be compared with 0 using math::equals(). Do you
>> agree? Or is there a reason why it's implemented like this?
>
> Yes there is a reason: in the past we used it using FPinput a lot and
> had to fight with robustness issues. I remember that it has been using
> math::equals() but that is gone.
> However, all these robustness problems are gone for intersections so
> we might add the math::equals again. Thanks for pointing this out.
>
>>
>>>
>>>> Here is the code and the results:
>>>> https://github.com/awulkiew/testgeoside
>>>
>>> Can you give a coordinatepair where the deviation is large
>>> (probably easy to read from the graph)?
>>>
>>
>> For the Earth the width of the erroreous part around lon=11.3,
>> lat=17.3 is more or less 1.7 km but this is for this very long
>> segment (51 51, 51 51). For a 1 deg segment (11 17, 12 18) (153km
>> long) the width of the erroreous fragment is around 0.6m. So the
>> difference seems to be relativaly small. I didn't check it for
>> coordinates closer to poles.
>
> I assume lon=11,lat=17 and similarly. I tried that using a third point
> somewhere in the middle.
Yes, (lon1 lat1, lon2 lat2) in degrees. I'm guessing that it'd be hard
to find a point of difference becuase a segment is quite "curvy" when
raw coordinates are considered, as shown here:
https://github.com/awulkiew/testgeoside.
> However, vincenty_inverse is not in the library  you must use some
> future improvement... Can you point where I can find that?
It's in my fork awulkiew/geometry branch feature/geographic:
https://github.com/awulkiew/geometry/tree/feature/geographic
also in a PR:
https://github.com/boostorg/geometry/pull/183
So e.g. to run https://github.com/awulkiew/testgeoside besides GLUT
installed in the system (e.g. freeglutdev package AFAIR) you'd have to
check out this branch from my fork, e.g.:
In BG local clone>
git remote add awulkiew git_at_[hidden]:awulkiew/geometry
git fetch awulkiew
git checkout feature/geographic
git pull
Regards,
Adam
Geometry list run by mateusz at loskot.net