 # Geometry :

Subject: [geometry] Help regarding line intersection on sphere issue
From: Stephen Hardy (Stephen.Hardy_at_[hidden])
Date: 2015-05-17 08:50:16

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.

#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

________________________________

The information in this e-mail may be confidential and subject to legal professional privilege and/or copyright. National ICT Australia Limited accepts no liability for any damage caused by this email or its attachments.

Geometry list run by mateusz at loskot.net