Skip to content

Commit d6861a9

Browse files
fontanffontanf
authored
Improve arc arc intersection
1 parent f24f899 commit d6861a9

2 files changed

Lines changed: 86 additions & 117 deletions

File tree

‎src/elements_intersections.cpp‎

Lines changed: 85 additions & 116 deletions
Original file line numberDiff line numberDiff line change
@@ -490,136 +490,105 @@ ShapeElementIntersectionsOutput compute_arc_arc_intersections(
490490
}
491491

492492
std::vector<Point> points;
493-
ElementPos p_size = 0;
494493
LengthDbl radius_1 = distance(arc.start, arc.center);
495494
LengthDbl radius_2 = distance(arc_2.start, arc_2.center);
495+
496+
LengthDbl xm = arc.center.x;
497+
LengthDbl ym = arc.center.y;
498+
LengthDbl xm2 = arc_2.center.x;
499+
LengthDbl ym2 = arc_2.center.y;
500+
LengthDbl a = 2 * (xm2 - xm);
501+
LengthDbl b = 2 * (ym2 - ym);
502+
LengthDbl c = rsq - (xm * xm) - (ym * ym)
503+
- r2sq + (xm2 * xm2) + (ym2 * ym2);
504+
//std::cout << "a " << a << " b " << b << " c " << c << std::endl;
505+
506+
LengthDbl c_prime = c - a * xm - b * ym;
507+
508+
// No intersection.
509+
if (strictly_lesser(rsq * (a * a + b * b), c_prime * c_prime))
510+
return {};
511+
512+
std::vector<Point> intersections;
513+
LengthDbl discriminant = rsq * (a * a + b * b) - c_prime * c_prime;
514+
//std::cout << "discriminant " << discriminant << std::endl;
515+
if (discriminant < 0)
516+
discriminant = 0;
517+
LengthDbl denom = a * a + b * b;
518+
LengthDbl eta_1 = (a * c_prime + b * std::sqrt(discriminant)) / denom;
519+
LengthDbl eta_2 = (a * c_prime - b * std::sqrt(discriminant)) / denom;
520+
LengthDbl teta_1 = (b * c_prime - a * std::sqrt(discriminant)) / denom;
521+
LengthDbl teta_2 = (b * c_prime + a * std::sqrt(discriminant)) / denom;
522+
points = {{xm + eta_1, ym + teta_1}, {xm + eta_2, ym + teta_2}};
523+
//std::cout << "p0 " << points[0].to_string() << std::endl;
524+
//std::cout << "p1 " << points[1].to_string() << std::endl;
525+
526+
Point middle = 0.5 * (points[0] + points[1]);
527+
if (equal(distance(middle, arc.center), radius_1)
528+
|| equal(distance(middle, arc.center), radius_2)) {
529+
points = {middle};
530+
}
531+
496532
bool circle_1_contains_arc_2_start = equal(distance(arc_2.start, arc.center), radius_1);
533+
if (circle_1_contains_arc_2_start) {
534+
if (points.size() == 1
535+
|| squared_distance(arc_2.start, points[0]) < squared_distance(arc_2.start, points[1])) {
536+
points[0] = arc_2.start;
537+
} else {
538+
points[1] = arc_2.start;
539+
}
540+
}
497541
bool circle_1_contains_arc_2_end = equal(distance(arc_2.end, arc.center), radius_1);
542+
if (circle_1_contains_arc_2_end) {
543+
if (points.size() == 1
544+
|| squared_distance(arc_2.end, points[0]) < squared_distance(arc_2.end, points[1])) {
545+
points[0] = arc_2.end;
546+
} else {
547+
points[1] = arc_2.end;
548+
}
549+
}
498550
bool circle_2_contains_arc_1_start = equal(distance(arc.start, arc_2.center), radius_2);
551+
if (circle_2_contains_arc_1_start) {
552+
if (points.size() == 1
553+
|| squared_distance(arc.start, points[0]) < squared_distance(arc.start, points[1])) {
554+
points[0] = arc.start;
555+
} else {
556+
points[1] = arc.start;
557+
}
558+
}
499559
bool circle_2_contains_arc_1_end = equal(distance(arc.end, arc_2.center), radius_2);
500-
#ifdef ELEMENTS_INTERSECTIONS_ENABLE_DEBUG
501-
std::cout << "circle_1_contains_arc_2_start " << circle_1_contains_arc_2_start << std::endl;
502-
std::cout << "circle_1_contains_arc_2_end " << circle_1_contains_arc_2_end << std::endl;
503-
std::cout << "circle_2_contains_arc_1_start " << circle_2_contains_arc_1_start << std::endl;
504-
std::cout << "circle_2_contains_arc_1_end " << circle_2_contains_arc_1_end << std::endl;
505-
#endif
506-
if (points.size() < 2 && circle_1_contains_arc_2_start)
507-
if (points.empty() || !(equal(arc_2.start, points.back())))
508-
points.push_back(arc_2.start);
509-
if (points.size() < 2 && circle_1_contains_arc_2_end)
510-
if (points.empty() || !(equal(arc_2.end, points.back())))
511-
points.push_back(arc_2.end);
512-
if (points.size() < 2 && circle_2_contains_arc_1_start)
513-
if (points.empty() || !(equal(arc.start, points.back())))
514-
points.push_back(arc.start);
515-
if (points.size() < 2 && circle_2_contains_arc_1_end)
516-
if (points.empty() || !(equal(arc.end, points.back())))
517-
points.push_back(arc.end);
518-
519-
#ifdef ELEMENTS_INTERSECTIONS_ENABLE_DEBUG
520-
std::cout << "points.size() " << points.size() << std::endl;
521-
#endif
560+
if (circle_2_contains_arc_1_end) {
561+
if (points.size() == 1
562+
|| squared_distance(arc.end, points[0]) < squared_distance(arc.end, points[1])) {
563+
points[0] = arc.end;
564+
} else {
565+
points[1] = arc.end;
566+
}
567+
}
522568

523569
ShapeElementIntersectionsOutput output;
524-
if (points.size() == 2) {
525-
for (const Point& p: points)
526-
if (arc.in_circular_arc_cone(p) && arc_2.in_circular_arc_cone(p))
570+
for (const Point& p: points) {
571+
// Check if any intersection coincides with an arc_2 endpoint
572+
if (equal(p, arc.start)) {
573+
if (arc_2.contains(p))
574+
output.improper_intersections.push_back(arc.start);
575+
} else if (equal(p, arc.end)) {
576+
if (arc_2.contains(p))
577+
output.improper_intersections.push_back(arc.end);
578+
} else if (equal(p, arc_2.start)) {
579+
if (arc.contains(p))
580+
output.improper_intersections.push_back(arc_2.start);
581+
} else if (equal(p, arc_2.end)) {
582+
if (arc.contains(p))
583+
output.improper_intersections.push_back(arc_2.end);
584+
} else if (arc.contains(p) && arc_2.contains(p)) {
585+
if (points.size() == 1) {
527586
output.improper_intersections.push_back(p);
528-
529-
} else if (points.size() == 1) {
530-
Point u = arc_2.center - arc.center;
531-
LengthDbl l2 = dot_product(u, u);
532-
533-
LengthDbl alpha = (l2 + radius_1 * radius_1 - radius_2 * radius_2) / (2 * l2);
534-
Point m = arc.center + alpha * u;
535-
Point p = 2 * m - points[0];
536-
537-
#ifdef ELEMENTS_INTERSECTIONS_ENABLE_DEBUG
538-
std::cout << "p " << p.to_string() << std::endl;
539-
#endif
540-
541-
if (arc.in_circular_arc_cone(points[0]) && arc_2.in_circular_arc_cone(points[0]))
542-
output.improper_intersections.push_back(points[0]);
543-
if (equal(p, points[0]))
544-
return output;
545-
if (arc.in_circular_arc_cone(p) && arc_2.in_circular_arc_cone(p))
546-
output.proper_intersections.push_back(p);
547-
548-
} else {
549-
LengthDbl xm = arc.center.x;
550-
LengthDbl ym = arc.center.y;
551-
LengthDbl xm2 = arc_2.center.x;
552-
LengthDbl ym2 = arc_2.center.y;
553-
LengthDbl a = 2 * (xm2 - xm);
554-
LengthDbl b = 2 * (ym2 - ym);
555-
LengthDbl c = rsq - (xm * xm) - (ym * ym)
556-
- r2sq + (xm2 * xm2) + (ym2 * ym2);
557-
//std::cout << "a " << a << " b " << b << " c " << c << std::endl;
558-
559-
LengthDbl c_prime = c - a * xm - b * ym;
560-
561-
// No intersection.
562-
if (strictly_lesser(rsq * (a * a + b * b), c_prime * c_prime))
563-
return {};
564-
565-
std::vector<Point> intersections;
566-
LengthDbl discriminant = rsq * (a * a + b * b) - c_prime * c_prime;
567-
//std::cout << "discriminant " << discriminant << std::endl;
568-
if (discriminant < 0)
569-
discriminant = 0;
570-
LengthDbl denom = a * a + b * b;
571-
LengthDbl eta_1 = (a * c_prime + b * std::sqrt(discriminant)) / denom;
572-
LengthDbl eta_2 = (a * c_prime - b * std::sqrt(discriminant)) / denom;
573-
LengthDbl teta_1 = (b * c_prime - a * std::sqrt(discriminant)) / denom;
574-
LengthDbl teta_2 = (b * c_prime + a * std::sqrt(discriminant)) / denom;
575-
points = {{xm + eta_1, ym + teta_1}, {xm + eta_2, ym + teta_2}};
576-
//std::cout << "p0 " << points[0].to_string() << std::endl;
577-
//std::cout << "p1 " << points[1].to_string() << std::endl;
578-
579-
// Single intersection point.
580-
if (equal(points[0], points[1])) {
581-
Point p;
582-
p.x = (points[0].x + points[1].x) / 2.0;
583-
p.y = (points[0].y + points[1].y) / 2.0;
584-
//std::cout << "p " << p.to_string() << std::endl;
585-
586-
if (equal(p, arc.start)) {
587-
p = arc.start;
588-
} else if (equal(p, arc.end)) {
589-
p = arc.end;
590-
} else if (equal(p, arc_2.start)) {
591-
p = arc_2.start;
592-
} else if (equal(p, arc_2.end)) {
593-
p = arc_2.end;
594-
}
595-
if (arc.contains(p) && arc_2.contains(p)) {
596-
return {{}, {p}, {}};
597587
} else {
598-
return {};
599-
}
600-
}
601-
602-
for (const Point& p: points) {
603-
//std::cout << "p " << p.to_string() << std::endl;
604-
// Check if any intersection coincides with an arc endpoint
605-
if (equal(p, arc.start)) {
606-
if (arc_2.contains(p))
607-
output.improper_intersections.push_back(arc.start);
608-
} else if (equal(p, arc.end)) {
609-
if (arc_2.contains(p))
610-
output.improper_intersections.push_back(arc.end);
611-
} else if (equal(p, arc_2.start)) {
612-
if (arc.contains(p))
613-
output.improper_intersections.push_back(arc_2.start);
614-
} else if (equal(p, arc_2.end)) {
615-
if (arc.contains(p))
616-
output.improper_intersections.push_back(arc_2.end);
617-
} else if (arc.contains(p) && arc_2.contains(p)) {
618588
output.proper_intersections.push_back(p);
619589
}
620590
}
621591
}
622-
623592
return output;
624593
}
625594

‎test/elements_intersections_test.cpp‎

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -349,7 +349,7 @@ INSTANTIATE_TEST_SUITE_P(
349349
}, {
350350
build_circular_arc({50.637459075944477, 829.58819029134031}, {50.605297582203768, 829.80611801559166}, {50.476672281381227, 829.67579862074012}, ShapeElementOrientation::Anticlockwise),
351351
build_circular_arc({50.480314241423272, 829.49272922597515}, {50.595723240021599, 829.81568159347239}, {50.4749160608335, 829.67675553731203}, ShapeElementOrientation::Anticlockwise),
352-
{{}, {}, {{50.63249994923487, 829.77197737420897}}},
352+
{{}, {}, {{50.632484210633152, 829.7719768687432}}},
353353
}, {
354354
build_line_segment({601.58208879290987, 795.1347195779681}, {601.5802365263462, 795.13583615892458}),
355355
build_circular_arc({601.58200217728086, 795.13771832729913}, {601.58206032271914, 795.13571917270087}, {601.58203125, 795.13671875}, ShapeElementOrientation::Anticlockwise),

0 commit comments

Comments
 (0)