
Geometry : 
Subject: Re: [geometry] Help regarding line intersection on sphere issue
From: Adam Wulkiewicz (adam.wulkiewicz_at_[hidden])
Date: 20150526 10:27:11
Hi Stephen,
Stephen Hardy wrote:
> On 17 May 2015, at 10:50 pm, Stephen Hardy <stephen.hardy_at_[hidden]> wrote:
>
>> Hi,
>>
>> I have been working on a code deriving plate tectonic motions on the earth’s surface based on geophysical evidence.
>>
>> I recently tried replacing my spherical geometry code with boost::geometry, using the spherical_equatorial coordinate system. The code became much smaller  this is great.
>>
>> However, I had a few regression issues, and I’ve traced them down to some peculiarities around the poles.
>>
>> The following is as simple as I can get it  intersect two lines, both crossing the pole, intersecting at the pole.
>>
>> ===
>> https://gist.github.com/sjhardy/cf3d61c0bbdd334eaf4f
>> ===
>>
>> #include <iostream>
>> #include <boost/geometry.hpp>
>>
>> int main(int argc, const char *argv[])
>> {
>>
>> typedef boost::geometry::model::point<double, 2, boost::geometry::cs::spherical_equatorial<boost::geometry::degree> > spoint;
>> typedef boost::geometry::model::linestring<spoint> sline;
>>
>> sline sline1, sline2;
>>
>> boost::geometry::append(sline1, spoint(90.0, 80.0));
>> boost::geometry::append(sline1, spoint(90.0, 80.0));
>> boost::geometry::append(sline2, spoint(180.0, 80.0));
>> boost::geometry::append(sline2, spoint(0.0, 80.0));
>>
>> typedef std::vector<spoint> pointvec;
>> pointvec res3;
>> boost::geometry::intersection(sline1, sline2, res3);
>> for(pointvec::iterator p = res3.begin(); p != res3.end(); p++)
>> std::cout << p>get<0>() << " " << p>get<1>() << "\n";
>>
>> return 1;
>> }
>>
>> ===
>>
>> This returns
>>
>> $ ./boostgeomtestcase
>> 90 80
>> 0 80
>>
>> rather than (0 90) or equivalent, as expected. I’m using boost 1.58.0.
>>
>> * Is there something that I am doing wrong above that is giving me the incorrect result?
>> * Is this a known issue (I could not find a bug for it)?
>>
>> Many thanks,
>> Stephen
>>
>>
>
> Hi,
>
> Further to the above, I’ve now found Ticket #9162 ( https://svn.boost.org/trac/boost/ticket/9162 ) from Sep ’13  "intersects returns incorrect result in spherical equatorial". This relates to point in polygon tests failing.
>
> I can confirm that this test still fails on boost 1.58.0.
Thanks for the info about the failure!
I think an issue mentioned in this ticket is a different one than the
one you're experiencing. The ticket is related to the spherical side
and/or winding strategies. But in your case I think that the
segment/segment intersection strategy should be blamed. In particular,
in this case as in all other cases cart_intersect strategy is used in
get_turns(). So the cartesian intersection of the segments of
linestrings is found. AFAIR noncartesian intersection of segments isn't
properly supported right now. It's even probable that it's wrong that
this code compiles. Anyway, could you please create a separate ticket?
Or maybe would you like to contribute to the library and add the
spherical CS support for intersection of 2 segments?
Regards,
Adam
Geometry list run by mateusz at loskot.net