74 static double computeEdgeDistance(
const geom::CoordinateXY& p,
const geom::CoordinateXY& p0,
const geom::CoordinateXY& p1);
81 precisionModel(initialPrecisionModel),
98 if(isInteriorIntersection(0)) {
101 if(isInteriorIntersection(1)) {
116 for(std::size_t i = 0; i < result; ++i) {
117 if(!(intPt[i].equals2D(*inputLines[inputLineIndex][0])
118 || intPt[i].equals2D(*inputLines[inputLineIndex][1]))) {
134 precisionModel = newPM;
143 void computeIntersection(
const geom::CoordinateXY& p,
const geom::CoordinateXY& p1,
const geom::CoordinateXY& p2);
146 static bool hasIntersection(
const geom::CoordinateXY& p,
const geom::CoordinateXY& p1,
const geom::CoordinateXY& p2);
153 POINT_INTERSECTION = 1,
156 COLLINEAR_INTERSECTION = 2
160 template<
typename C1,
typename C2>
162 const C2& p3,
const C2& p4)
164 inputLines[0][0] = &p1;
165 inputLines[0][1] = &p2;
166 inputLines[1][0] = &p3;
167 inputLines[1][1] = &p4;
168 result = computeIntersect(p1, p2, p3, p4);
175 std::string toString()
const;
185 return result != NO_INTERSECTION;
196 const geom::CoordinateXY*
199 return inputLines[segmentIndex][ptIndex];
219 const geom::CoordinateXYZM&
222 return intPt[intIndex];
243 for(std::size_t i = 0; i < result; ++i) {
244 if(intPt[i].equals2D(pt)) {
268 return hasIntersection() && isProperVar;
315 const geom::CoordinateXY* inputLines[2][2];
321 geom::CoordinateXYZM intPt[2];
327 std::size_t intLineIndex[2][2];
336 return result == COLLINEAR_INTERSECTION;
339 template<
typename C1,
typename C2>
340 uint8_t computeIntersect(
const C1& p1,
const C1& p2,
const C2& q1,
const C2& q2)
345 if(!geom::Envelope::intersects(p1, p2, q1, q2)) {
346 return NO_INTERSECTION;
352 int Pq1 = Orientation::index(p1, p2, q1);
353 int Pq2 = Orientation::index(p1, p2, q2);
355 if((Pq1 > 0 && Pq2 > 0) || (Pq1 < 0 && Pq2 < 0)) {
356 return NO_INTERSECTION;
359 int Qp1 = Orientation::index(q1, q2, p1);
360 int Qp2 = Orientation::index(q1, q2, p2);
362 if((Qp1 > 0 && Qp2 > 0) || (Qp1 < 0 && Qp2 < 0)) {
363 return NO_INTERSECTION;
369 bool collinear = Pq1 == 0 && Pq2 == 0 && Qp1 == 0 && Qp2 == 0;
371 return computeCollinearIntersection(p1, p2, q1, q2);
389 geom::CoordinateXYZM p;
390 double z = DoubleNotANumber;
391 double m = DoubleNotANumber;
393 if(Pq1 == 0 || Pq2 == 0 || Qp1 == 0 || Qp2 == 0) {
414 if (p1.equals2D(q1)) {
416 z = Interpolate::zGet(p1, q1);
417 m = Interpolate::mGet(p1, q1);
419 else if (p1.equals2D(q2)) {
421 z = Interpolate::zGet(p1, q2);
422 m = Interpolate::mGet(p1, q2);
424 else if (p2.equals2D(q1)) {
426 z = Interpolate::zGet(p2, q1);
427 m = Interpolate::mGet(p2, q1);
429 else if (p2.equals2D(q2)) {
431 z = Interpolate::zGet(p2, q2);
432 m = Interpolate::mGet(p2, q2);
439 z = Interpolate::zGetOrInterpolate(q1, p1, p2);
440 m = Interpolate::mGetOrInterpolate(q1, p1, p2);
444 z = Interpolate::zGetOrInterpolate(q2, p1, p2);
445 m = Interpolate::mGetOrInterpolate(q2, p1, p2);
449 z = Interpolate::zGetOrInterpolate(p1, q1, q2);
450 m = Interpolate::mGetOrInterpolate(p1, q1, q2);
454 z = Interpolate::zGetOrInterpolate(p2, q1, q2);
455 m = Interpolate::mGetOrInterpolate(p2, q1, q2);
459 p = intersection(p1, p2, q1, q2);
460 z = Interpolate::zInterpolate(p, p1, p2, q1, q2);
461 m = Interpolate::mInterpolate(p, p1, p2, q1, q2);
463 intPt[0] = geom::CoordinateXYZM(p.x, p.y, z, m);
465 std::cerr <<
" POINT_INTERSECTION; intPt[0]:" << intPt[0].toString() << std::endl;
467 return POINT_INTERSECTION;
473 return hasIntersection() && !isProperVar;
476 void computeIntLineIndex();
478 void computeIntLineIndex(std::size_t segmentIndex);
480 template<
typename C1,
typename C2>
481 uint8_t computeCollinearIntersection(
const C1& p1,
const C1& p2,
const C2& q1,
const C2& q2)
483 bool q1inP = geom::Envelope::intersects(p1, p2, q1);
484 bool q2inP = geom::Envelope::intersects(p1, p2, q2);
485 bool p1inQ = geom::Envelope::intersects(q1, q2, p1);
486 bool p2inQ = geom::Envelope::intersects(q1, q2, p2);
489 intPt[0] = zmGetOrInterpolateCopy(q1, p1, p2);
490 intPt[1] = zmGetOrInterpolateCopy(q2, p1, p2);
491 return COLLINEAR_INTERSECTION;
494 intPt[0] = zmGetOrInterpolateCopy(p1, q1, q2);
495 intPt[1] = zmGetOrInterpolateCopy(p2, q1, q2);
496 return COLLINEAR_INTERSECTION;
500 intPt[0] = zmGetOrInterpolateCopy(q1, p1, p2);
501 intPt[1] = zmGetOrInterpolateCopy(p1, q1, q2);
503 return (q1 == p1) && !q2inP && !p2inQ ? POINT_INTERSECTION : COLLINEAR_INTERSECTION;
507 intPt[0] = zmGetOrInterpolateCopy(q1, p1, p2);
508 intPt[1] = zmGetOrInterpolateCopy(p2, q1, q2);
510 return (q1 == p2) && !q2inP && !p1inQ ? POINT_INTERSECTION : COLLINEAR_INTERSECTION;
514 intPt[0] = zmGetOrInterpolateCopy(q2, p1, p2);
515 intPt[1] = zmGetOrInterpolateCopy(p1, q1, q2);
517 return (q2 == p1) && !q1inP && !p2inQ ? POINT_INTERSECTION : COLLINEAR_INTERSECTION;
521 intPt[0] = zmGetOrInterpolateCopy(q2, p1, p2);
522 intPt[1] = zmGetOrInterpolateCopy(p2, q1, q2);
523 return (q2 == p2) && !q1inP && !p1inQ ? POINT_INTERSECTION : COLLINEAR_INTERSECTION;
525 return NO_INTERSECTION;
537 template<
typename C1,
typename C2>
538 geom::CoordinateXYZM intersection (
const C1& p1,
const C1& p2,
const C2& q1,
const C2& q2)
const {
539 auto intPtOut = intersectionSafe(p1, p2, q1, q2);
556 if(! isInSegmentEnvelopes(intPtOut)) {
558 intPtOut = nearestEndpoint(p1, p2, q1, q2);
561 if(precisionModel !=
nullptr) {
578 bool isInSegmentEnvelopes(
const geom::CoordinateXY& pt)
const
580 geom::Envelope env0(*inputLines[0][0], *inputLines[0][1]);
581 geom::Envelope env1(*inputLines[1][0], *inputLines[1][1]);
582 return env0.contains(pt) && env1.contains(pt);
597 template<
typename C1,
typename C2>
598 geom::CoordinateXYZM intersectionSafe(
const C1& p1,
const C1& p2,
599 const C2& q1,
const C2& q2)
const
601 geom::CoordinateXYZM ptInt(Intersection::intersection(p1, p2, q1, q2));
602 if (ptInt.isNull()) {
604 ptInt =
static_cast<const C1&
>(nearestEndpoint(p1, p2, q1, q2));
628 static const geom::CoordinateXY& nearestEndpoint(
const geom::CoordinateXY& p1,
629 const geom::CoordinateXY& p2,
630 const geom::CoordinateXY& q1,
631 const geom::CoordinateXY& q2);
634 template<
typename C1,
typename C2>
635 static geom::CoordinateXYZM zmGetOrInterpolateCopy(
640 geom::CoordinateXYZM pCopy(p);
641 pCopy.z = Interpolate::zGetOrInterpolate(p, p1, p2);
642 pCopy.m = Interpolate::mGetOrInterpolate(p, p1, p2);