| 594 | } |
| 595 | |
| 596 | func intersectionLineEllipse(zs Intersections, l0, l1, center, radius Point, phi, theta0, theta1 float64) Intersections { |
| 597 | if Equal(radius.X, radius.Y) { |
| 598 | return intersectionLineCircle(zs, l0, l1, center, radius.X, theta0, theta1) |
| 599 | } else if l0.Equals(l1) { |
| 600 | return zs // zero-length Close |
| 601 | } |
| 602 | |
| 603 | // TODO: needs more testing |
| 604 | // TODO: intersection inconsistency due to numerical stability in finding tangent collisions for subsequent path segments (line -> ellipse), or due to the endpoint of a line not touching with another arc, but the subsequent segment does touch with its starting point |
| 605 | dira := l1.Sub(l0).Angle() |
| 606 | |
| 607 | // we take the ellipse center as the origin and counter-rotate by phi |
| 608 | l0 = l0.Sub(center).Rot(-phi, Origin) |
| 609 | l1 = l1.Sub(center).Rot(-phi, Origin) |
| 610 | |
| 611 | // line: cx + dy + e = 0 |
| 612 | c := l0.Y - l1.Y |
| 613 | d := l1.X - l0.X |
| 614 | e := l0.PerpDot(l1) |
| 615 | |
| 616 | // follow different code paths when line is mostly horizontal or vertical |
| 617 | horizontal := math.Abs(c) <= math.Abs(d) |
| 618 | |
| 619 | // ellipse: x^2/a + y^2/b = 1 |
| 620 | a := radius.X * radius.X |
| 621 | b := radius.Y * radius.Y |
| 622 | |
| 623 | // rewrite as a polynomial by substituting x or y to obtain: |
| 624 | // At^2 + Bt + C = 0, with t either x (horizontal) or y (!horizontal) |
| 625 | var A, B, C float64 |
| 626 | A = a*c*c + b*d*d |
| 627 | if horizontal { |
| 628 | B = 2.0 * a * c * e |
| 629 | C = a*e*e - a*b*d*d |
| 630 | } else { |
| 631 | B = 2.0 * b * d * e |
| 632 | C = b*e*e - a*b*c*c |
| 633 | } |
| 634 | |
| 635 | // find solutions |
| 636 | roots := []float64{} |
| 637 | r0, r1 := solveQuadraticFormula(A, B, C) |
| 638 | if !math.IsNaN(r0) { |
| 639 | roots = append(roots, r0) |
| 640 | if !math.IsNaN(r1) && !Equal(r0, r1) { |
| 641 | roots = append(roots, r1) |
| 642 | } |
| 643 | } |
| 644 | |
| 645 | for _, root := range roots { |
| 646 | // get intersection position with center as origin |
| 647 | var x, y, t0, t1 float64 |
| 648 | if horizontal { |
| 649 | x = root |
| 650 | y = -e/d - c*root/d |
| 651 | t0 = l0.X |
| 652 | t1 = l1.X |
| 653 | } else { |