@@ -676,43 +676,31 @@ mod tests {
676676}
677677
678678pub mod inscribe_circles_algorithms {
679- use core:: ops:: Range ;
680- use kurbo:: { ParamCurve , ParamCurveDeriv , ParamCurveExtrema } ;
681-
682- const ROUND_ACCURACY : f64 = 1e-5 ;
679+ use kurbo:: { ParamCurve , ParamCurveDeriv , common:: solve_itp} ;
683680
684681 #[ derive( Clone , Copy , Debug , PartialEq ) ]
685682 pub struct CircleInscription {
686- pub time_1 : f64 ,
687- pub time_2 : f64 ,
683+ pub t_first : f64 ,
684+ pub t_second : f64 ,
688685 pub theta : f64 ,
689- pub circle_centre1 : glam :: DVec2 ,
690- pub circle_centre2 : glam:: DVec2 ,
686+ pub radius_to_centre : f64 ,
687+ pub circle_centre : glam:: DVec2 ,
691688 }
692689
693690 /// Find the normalised tangent at a particular time. Avoid using for t=0 or t=1 due to errors.
694- fn tangent ( segment : kurbo:: PathSeg , t : f64 ) -> kurbo :: Vec2 {
691+ fn tangent ( segment : kurbo:: PathSeg , t : f64 ) -> glam :: DVec2 {
695692 let tangent = match segment {
696693 kurbo:: PathSeg :: Line ( line) => line. deriv ( ) . eval ( t) ,
697694 kurbo:: PathSeg :: Quad ( quad_bez) => quad_bez. deriv ( ) . eval ( t) ,
698695 kurbo:: PathSeg :: Cubic ( cubic_bez) => cubic_bez. deriv ( ) . eval ( t) ,
699696 }
700697 . to_vec2 ( )
701698 . normalize ( ) ;
699+ let tangent = glam:: DVec2 :: new ( tangent. x , tangent. y ) ;
702700 debug_assert ! ( tangent. is_finite( ) , "cannot round corner with NaN tangent" ) ;
703701 tangent
704702 }
705703
706- /// Rotate 90 degrees in one direction
707- fn offset_1 ( value : kurbo:: Vec2 , radius : f64 ) -> kurbo:: Vec2 {
708- kurbo:: Vec2 :: new ( -value. y , value. x ) * radius
709- }
710-
711- /// Rotate 90 degrees in one direction
712- fn offset_2 ( value : kurbo:: Vec2 , radius : f64 ) -> kurbo:: Vec2 {
713- kurbo:: Vec2 :: new ( value. y , -value. x ) * radius
714- }
715-
716704 /// Compute the tangent at t=0 for the path segment
717705 pub fn tangent_at_start ( segment : kurbo:: PathSeg ) -> kurbo:: Vec2 {
718706 let tangent = match segment {
@@ -735,81 +723,6 @@ pub mod inscribe_circles_algorithms {
735723 tangent
736724 }
737725
738- /// Resolve the bounding boxes offset by radius in either direciton.
739- fn offset_bounding_boxes ( segment : kurbo:: PathSeg , radius : f64 ) -> [ kurbo:: Rect ; 2 ] {
740- let [ start_tangent, end_tangent] = [ tangent_at_start ( segment) , -tangent_at_start ( segment. reverse ( ) ) ] ;
741-
742- let mut bbox1 = kurbo:: Rect :: from_points ( segment. start ( ) + offset_1 ( start_tangent, radius) , segment. end ( ) + offset_1 ( end_tangent, radius) ) ;
743- let mut bbox2 = kurbo:: Rect :: from_points ( segment. start ( ) + offset_2 ( start_tangent, radius) , segment. end ( ) + offset_2 ( end_tangent, radius) ) ;
744- // The extrema for the original curve should be the same as for the offset curve
745- for extremum in segment. extrema ( ) {
746- let value = segment. eval ( extremum) ;
747- let derivative = tangent ( segment, extremum) ;
748- bbox1 = bbox1. union_pt ( value + offset_1 ( derivative, radius) ) ;
749- bbox2 = bbox2. union_pt ( value + offset_2 ( derivative, radius) ) ;
750- }
751- debug_assert ! ( bbox1. is_finite( ) && bbox2. is_finite( ) , "a wild NaN appeared :(" ) ;
752- [ bbox1, bbox2]
753- }
754-
755- /// If the width and height both smaller than accuracy then we can end the recursion
756- fn rect_within_accuracy ( rect : kurbo:: Rect , accuracy : f64 ) -> bool {
757- rect. width ( ) . abs ( ) < accuracy && rect. height ( ) . abs ( ) < accuracy
758- }
759-
760- /// Resursively find position to inscribe circles
761- fn inscribe_internal ( segment1 : kurbo:: PathSeg , t1 : Range < f64 > , segment2 : kurbo:: PathSeg , t2 : Range < f64 > , radius : f64 ) -> Option < CircleInscription > {
762- let bbox1 = offset_bounding_boxes ( segment1. subsegment ( t1. clone ( ) ) , radius) ;
763- let bbox2 = offset_bounding_boxes ( segment2. subsegment ( t2. clone ( ) ) , radius) ;
764- let mid_t1 = ( t1. start + t1. end ) / 2. ;
765- let mid_t2 = ( t2. start + t2. end ) / 2. ;
766-
767- // Check if the bounding boxes overlap
768- let mut any_overlap = false ;
769- for i in 0 ..4usize {
770- let [ index_1, index_2] = [ i >> 1 , i & 1 ] ;
771- let [ first, second] = [ bbox1[ index_1] , bbox2[ index_2] ] ;
772-
773- // Ignore non overlapping
774- if !first. overlaps ( second) {
775- continue ;
776- }
777- // If the rects are small enough then complete the recursion
778- if rect_within_accuracy ( first, ROUND_ACCURACY ) && rect_within_accuracy ( second, ROUND_ACCURACY ) {
779- let tangents = [ ( segment1, mid_t1) , ( segment2, mid_t2) ] . map ( |( segment, t) | tangent ( segment, t) ) ;
780- let normal_1 = [ offset_1, offset_2] [ index_1] ( tangents[ 0 ] , 1. ) ;
781- let normal_2 = [ offset_1, offset_2] [ index_2] ( tangents[ 1 ] , 1. ) ;
782- let circle_centre_1 = segment1. eval ( mid_t1) + normal_1 * radius;
783- let circle_centre_2 = segment2. eval ( mid_t2) + normal_2 * radius;
784- return Some ( CircleInscription {
785- time_1 : mid_t1,
786- time_2 : mid_t2,
787- theta : normal_1. dot ( normal_2) . clamp ( -1. , 1. ) . acos ( ) ,
788- circle_centre1 : glam:: DVec2 :: new ( circle_centre_1. x , circle_centre_1. y ) ,
789- circle_centre2 : glam:: DVec2 :: new ( circle_centre_2. x , circle_centre_2. y ) ,
790- } ) ;
791- }
792- any_overlap = true ;
793- }
794- if !any_overlap {
795- return None ;
796- }
797-
798- let [ start_t1, end_t1] = [ t1. start , t1. end ] ;
799- let [ start_t2, end_t2] = [ t2. start , t2. end ] ;
800-
801- // Repeat checking the intersection with the combinations of the two halves of each curve
802- if let Some ( result) = None
803- . or_else ( || inscribe_internal ( segment1, start_t1..mid_t1, segment2, start_t2..mid_t2, radius) )
804- . or_else ( || inscribe_internal ( segment1, start_t1..mid_t1, segment2, mid_t2..end_t2, radius) )
805- . or_else ( || inscribe_internal ( segment1, mid_t1..end_t1, segment2, start_t2..mid_t2, radius) )
806- . or_else ( || inscribe_internal ( segment1, mid_t1..end_t1, segment2, mid_t2..end_t2, radius) )
807- {
808- return Some ( result) ;
809- }
810- None
811- }
812-
813726 /// Convert [`crate::subpath::Bezier`] to [`kurbo::PathSeg`]
814727 pub fn bezier_to_path_seg ( bezier : crate :: subpath:: Bezier ) -> kurbo:: PathSeg {
815728 let [ start, end] = [ ( bezier. start ( ) . x , bezier. start ( ) . y ) , ( bezier. end ( ) . x , bezier. end ( ) . y ) ] ;
@@ -834,20 +747,62 @@ pub mod inscribe_circles_algorithms {
834747 }
835748 }
836749
750+ /// Find the t value that is distance `radius` from the start
751+ fn distance_from_start ( seg : kurbo:: PathSeg , radius : f64 ) -> Option < f64 > {
752+ let r_squared = radius * radius;
753+ let final_distance = ( seg. end ( ) - seg. start ( ) ) . length_squared ( ) ;
754+ if final_distance < radius {
755+ return None ;
756+ }
757+ let evaluate = |t| ( seg. eval ( t) - seg. start ( ) ) . length_squared ( ) - r_squared;
758+ Some ( solve_itp ( evaluate, 0. , 1. , 1e-9 , 1 , 0.2 , evaluate ( 0. ) , evaluate ( 1. ) ) )
759+ }
760+
837761 /// Attemt to inscribe circle into the start of the [`kurbo::PathSeg`]s
838762 pub fn inscribe ( first : kurbo:: PathSeg , second : kurbo:: PathSeg , radius : f64 ) -> Option < CircleInscription > {
839- inscribe_internal ( first, 0.0 ..1. , second, 0.0 ..1. , radius)
763+ let [ t_first, t_second] = [ distance_from_start ( first, radius) ?, distance_from_start ( second, radius) ?] ;
764+
765+ let tangents = [ ( first, t_first) , ( second, t_second) ] . map ( |( segment, t) | tangent ( segment, t) ) ;
766+ let points = [ ( first, t_first) , ( second, t_second) ] . map ( |( segment, t) | segment. eval ( t) ) . map ( |x| glam:: DVec2 :: new ( x. x , x. y ) ) ;
767+
768+ let mut normals = tangents. map ( glam:: DVec2 :: perp) ;
769+ // Make sure the normals are pointing in the right direction
770+ normals[ 0 ] *= normals[ 0 ] . dot ( tangents[ 1 ] ) . signum ( ) ;
771+ normals[ 1 ] *= normals[ 1 ] . dot ( tangents[ 0 ] ) . signum ( ) ;
772+
773+ let mid = ( points[ 0 ] + points[ 1 ] ) / 2. ;
774+
775+ if normals[ 0 ] . abs_diff_eq ( glam:: DVec2 :: ZERO , 1e-6 ) || normals[ 1 ] . abs_diff_eq ( glam:: DVec2 :: ZERO , 1e-6 ) || mid. abs_diff_eq ( points[ 0 ] , 1e-6 ) {
776+ return None ;
777+ }
778+
779+ let radius_to_centre = ( mid - points[ 0 ] ) . length_squared ( ) / ( normals[ 0 ] . dot ( mid - points[ 0 ] ) ) ;
780+ let circle_centre = points[ 0 ] + normals[ 0 ] * radius_to_centre;
781+
782+ if radius_to_centre > radius * 10. {
783+ return None ; // Don't inscribe if it is a long way from the centre
784+ }
785+
786+ info ! ( "Points {points:?}\n tangents {tangents:?}\n normals {normals:?}\n centres {circle_centre}" ) ;
787+ return Some ( CircleInscription {
788+ t_first,
789+ t_second,
790+ theta : normals[ 0 ] . dot ( normals[ 1 ] ) . clamp ( -1. , 1. ) . acos ( ) ,
791+ radius_to_centre,
792+ circle_centre,
793+ } ) ;
840794 }
841795
842796 #[ cfg( test) ]
843797 mod inscribe_tests {
798+ const ROUND_ACCURACY : f64 = 1e-6 ;
844799 #[ test]
845800 fn test_perpendicular_lines ( ) {
846801 let l1 = kurbo:: PathSeg :: Line ( kurbo:: Line :: new ( ( 0. , 0. ) , ( 100. , 0. ) ) ) ;
847802 let l2 = kurbo:: PathSeg :: Line ( kurbo:: Line :: new ( ( 0. , 0. ) , ( 0. , 100. ) ) ) ;
848803
849804 let result = super :: inscribe ( l1, l2, 5. ) ;
850- assert ! ( result. unwrap( ) . circle_centre1 . abs_diff_eq( glam:: DVec2 :: new( 5. , 5. ) , super :: ROUND_ACCURACY * 10. ) , "{result:?}" ) ;
805+ assert ! ( result. unwrap( ) . circle_centre . abs_diff_eq( glam:: DVec2 :: new( 5. , 5. ) , ROUND_ACCURACY ) , "{result:?}" ) ;
851806 assert_eq ! ( result. unwrap( ) . theta, std:: f64 :: consts:: FRAC_PI_2 , "unexpected {result:?}" ) ;
852807 }
853808
@@ -857,8 +812,8 @@ pub mod inscribe_circles_algorithms {
857812 let l2 = kurbo:: PathSeg :: Line ( kurbo:: Line :: new ( ( 0. , 0. ) , ( 0. , 100. ) ) ) ;
858813
859814 let result = super :: inscribe ( l1, l2, 5. ) ;
860- let expected_centre = glam:: DVec2 :: new ( 5. , 5. + 5. * std :: f64:: consts:: SQRT_2 ) ;
861- assert ! ( result. unwrap( ) . circle_centre1 . abs_diff_eq( expected_centre, super :: ROUND_ACCURACY * 10. ) , "unexpected {result:?}" ) ;
815+ let expected_centre = glam:: DVec2 :: new ( 10. / core :: f64:: consts:: SQRT_2 - 5. , 5. ) ;
816+ assert ! ( result. unwrap( ) . circle_centre . abs_diff_eq( expected_centre, ROUND_ACCURACY ) , "unexpected {result:?}" ) ;
862817 assert_eq ! ( result. unwrap( ) . theta, std:: f64 :: consts:: FRAC_PI_4 * 3. , "unexpected {result:?}" ) ;
863818 }
864819
@@ -868,8 +823,8 @@ pub mod inscribe_circles_algorithms {
868823 let l2 = kurbo:: PathSeg :: Line ( kurbo:: Line :: new ( ( 0. , 0. ) , ( 40. , 30. ) ) ) ;
869824
870825 let result = super :: inscribe ( l1, l2, 5. ) ;
871- let expected_centre = glam:: DVec2 :: new ( 25. , 25. ) ;
872- assert ! ( result. unwrap( ) . circle_centre1 . abs_diff_eq( expected_centre, super :: ROUND_ACCURACY * 10. ) , "{result:?}" ) ;
826+ let expected_centre = glam:: DVec2 :: new ( 25. / 7. , 25. / 7 .) ;
827+ assert ! ( result. unwrap( ) . circle_centre . abs_diff_eq( expected_centre, ROUND_ACCURACY ) , "{result:?}" ) ;
873828 assert_eq ! ( result. unwrap( ) . theta, ( -24f64 / 25. ) . acos( ) , "{result:?}" ) ;
874829 }
875830
@@ -879,7 +834,7 @@ pub mod inscribe_circles_algorithms {
879834 let l2 = kurbo:: PathSeg :: Cubic ( kurbo:: CubicBez :: new ( ( 0. , 0. ) , ( 0. , 33. ) , ( 0. , 67. ) , ( 0. , 100. ) ) ) ;
880835
881836 let result = super :: inscribe ( l1, l2, 5. ) ;
882- assert ! ( result. unwrap( ) . circle_centre1 . abs_diff_eq( glam:: DVec2 :: new( 5. , 5. ) , super :: ROUND_ACCURACY * 10. ) , "{result:?}" ) ;
837+ assert ! ( result. unwrap( ) . circle_centre . abs_diff_eq( glam:: DVec2 :: new( 5. , 5. ) , ROUND_ACCURACY ) , "{result:?}" ) ;
883838 assert_eq ! ( result. unwrap( ) . theta, std:: f64 :: consts:: FRAC_PI_2 , "unexpected {result:?}" ) ;
884839 }
885840 }
0 commit comments