kcl_lib/std/
utils.rs

1use std::f64::consts::PI;
2
3use kittycad_modeling_cmds::{shared::Angle, units::UnitLength};
4
5use super::args::TyF64;
6use crate::execution::types::NumericType;
7
8pub(crate) fn untype_point(p: [TyF64; 2]) -> ([f64; 2], NumericType) {
9    let (x, y, ty) = NumericType::combine_eq_coerce(p[0].clone(), p[1].clone());
10    ([x, y], ty)
11}
12
13pub(crate) fn untype_array<const N: usize>(p: [TyF64; N]) -> ([f64; N], NumericType) {
14    let (vec, ty) = NumericType::combine_eq_array(&p);
15    (
16        vec.try_into()
17            .unwrap_or_else(|v: Vec<f64>| panic!("Expected a Vec of length {} but it was {}", N, v.len())),
18        ty,
19    )
20}
21
22pub(crate) fn point_to_mm(p: [TyF64; 2]) -> [f64; 2] {
23    [p[0].to_mm(), p[1].to_mm()]
24}
25
26pub(crate) fn untyped_point_to_mm(p: [f64; 2], units: UnitLength) -> [f64; 2] {
27    [
28        crate::execution::types::adjust_length(units, p[0], UnitLength::Millimeters).0,
29        crate::execution::types::adjust_length(units, p[1], UnitLength::Millimeters).0,
30    ]
31}
32
33pub(crate) fn point_to_len_unit(p: [TyF64; 2], len: UnitLength) -> [f64; 2] {
34    [p[0].to_length_units(len), p[1].to_length_units(len)]
35}
36
37/// Precondition, `p` must be in `len` units (this function does no conversion).
38pub(crate) fn point_to_typed(p: [f64; 2], len: UnitLength) -> [TyF64; 2] {
39    [TyF64::new(p[0], len.into()), TyF64::new(p[1], len.into())]
40}
41
42pub(crate) fn point_3d_to_mm(p: [TyF64; 3]) -> [f64; 3] {
43    [p[0].to_mm(), p[1].to_mm(), p[2].to_mm()]
44}
45
46/// Get the distance between two points.
47pub(crate) fn distance(a: Coords2d, b: Coords2d) -> f64 {
48    ((b[0] - a[0]).powi(2) + (b[1] - a[1]).powi(2)).sqrt()
49}
50
51/// Get the angle between these points
52pub(crate) fn between(a: Coords2d, b: Coords2d) -> Angle {
53    let x = b[0] - a[0];
54    let y = b[1] - a[1];
55    normalize(Angle::from_radians(libm::atan2(y, x)))
56}
57
58/// Normalize the angle
59pub(crate) fn normalize(angle: Angle) -> Angle {
60    let deg = angle.to_degrees();
61    let result = ((deg % 360.0) + 360.0) % 360.0;
62    Angle::from_degrees(if result > 180.0 { result - 360.0 } else { result })
63}
64
65/// Gives the ▲-angle between from and to angles (shortest path)
66///
67/// Sign of the returned angle denotes direction, positive means counterClockwise 🔄
68/// # Examples
69///
70/// ```
71/// use std::f64::consts::PI;
72///
73/// use kcl_lib::std::utils::Angle;
74///
75/// assert_eq!(
76///     Angle::delta(Angle::from_radians(PI / 8.0), Angle::from_radians(PI / 4.0)),
77///     Angle::from_radians(PI / 8.0)
78/// );
79/// ```
80pub(crate) fn delta(from_angle: Angle, to_angle: Angle) -> Angle {
81    let norm_from_angle = normalize_rad(from_angle.to_radians());
82    let norm_to_angle = normalize_rad(to_angle.to_radians());
83    let provisional = norm_to_angle - norm_from_angle;
84
85    if provisional > -PI && provisional <= PI {
86        return Angle::from_radians(provisional);
87    }
88    if provisional > PI {
89        return Angle::from_radians(provisional - 2.0 * PI);
90    }
91    if provisional < -PI {
92        return Angle::from_radians(provisional + 2.0 * PI);
93    }
94    Angle::default()
95}
96
97pub(crate) fn normalize_rad(angle: f64) -> f64 {
98    let draft = angle % (2.0 * PI);
99    if draft < 0.0 { draft + 2.0 * PI } else { draft }
100}
101
102fn calculate_intersection_of_two_lines(line1: &[Coords2d; 2], line2_angle: f64, line2_point: Coords2d) -> Coords2d {
103    let line2_point_b = [
104        line2_point[0] + libm::cos(line2_angle.to_radians()) * 10.0,
105        line2_point[1] + libm::sin(line2_angle.to_radians()) * 10.0,
106    ];
107    intersect(line1[0], line1[1], line2_point, line2_point_b)
108}
109
110fn intersect(p1: Coords2d, p2: Coords2d, p3: Coords2d, p4: Coords2d) -> Coords2d {
111    let slope = |p1: Coords2d, p2: Coords2d| (p1[1] - p2[1]) / (p1[0] - p2[0]);
112    let constant = |p1: Coords2d, p2: Coords2d| p1[1] - slope(p1, p2) * p1[0];
113    let get_y = |for_x: f64, p1: Coords2d, p2: Coords2d| slope(p1, p2) * for_x + constant(p1, p2);
114
115    if p1[0] == p2[0] {
116        return [p1[0], get_y(p1[0], p3, p4)];
117    }
118    if p3[0] == p4[0] {
119        return [p3[0], get_y(p3[0], p1, p2)];
120    }
121
122    let x = (constant(p3, p4) - constant(p1, p2)) / (slope(p1, p2) - slope(p3, p4));
123    let y = get_y(x, p1, p2);
124    [x, y]
125}
126
127pub(crate) fn intersection_with_parallel_line(
128    line1: &[Coords2d; 2],
129    line1_offset: f64,
130    line2_angle: f64,
131    line2_point: Coords2d,
132) -> Coords2d {
133    calculate_intersection_of_two_lines(&offset_line(line1_offset, line1[0], line1[1]), line2_angle, line2_point)
134}
135
136fn offset_line(offset: f64, p1: Coords2d, p2: Coords2d) -> [Coords2d; 2] {
137    if p1[0] == p2[0] {
138        let direction = (p1[1] - p2[1]).signum();
139        return [[p1[0] + offset * direction, p1[1]], [p2[0] + offset * direction, p2[1]]];
140    }
141    if p1[1] == p2[1] {
142        let direction = (p2[0] - p1[0]).signum();
143        return [[p1[0], p1[1] + offset * direction], [p2[0], p2[1] + offset * direction]];
144    }
145    let x_offset = offset / libm::sin(libm::atan2(p1[1] - p2[1], p1[0] - p2[0]));
146    [[p1[0] + x_offset, p1[1]], [p2[0] + x_offset, p2[1]]]
147}
148
149pub(crate) fn get_y_component(angle: Angle, x: f64) -> Coords2d {
150    let normalised_angle = ((angle.to_degrees() % 360.0) + 360.0) % 360.0; // between 0 and 360
151    let y = x * libm::tan(normalised_angle.to_radians());
152    let sign = if normalised_angle > 90.0 && normalised_angle <= 270.0 {
153        -1.0
154    } else {
155        1.0
156    };
157    [x * sign, y * sign]
158}
159
160pub(crate) fn get_x_component(angle: Angle, y: f64) -> Coords2d {
161    let normalised_angle = ((angle.to_degrees() % 360.0) + 360.0) % 360.0; // between 0 and 360
162    let x = y / libm::tan(normalised_angle.to_radians());
163    let sign = if normalised_angle > 180.0 && normalised_angle <= 360.0 {
164        -1.0
165    } else {
166        1.0
167    };
168    [x * sign, y * sign]
169}
170
171pub(crate) fn arc_center_and_end(
172    from: Coords2d,
173    start_angle: Angle,
174    end_angle: Angle,
175    radius: f64,
176) -> (Coords2d, Coords2d) {
177    let start_angle = start_angle.to_radians();
178    let end_angle = end_angle.to_radians();
179
180    let center = [
181        -(radius * libm::cos(start_angle) - from[0]),
182        -(radius * libm::sin(start_angle) - from[1]),
183    ];
184
185    let end = [
186        center[0] + radius * libm::cos(end_angle),
187        center[1] + radius * libm::sin(end_angle),
188    ];
189
190    (center, end)
191}
192
193// Calculate the center of 3 points using an algebraic method
194// Handles if 3 points lie on the same line (collinear) by returning the average of the points (could return None instead..)
195pub(crate) fn calculate_circle_center(p1: [f64; 2], p2: [f64; 2], p3: [f64; 2]) -> [f64; 2] {
196    let (x1, y1) = (p1[0], p1[1]);
197    let (x2, y2) = (p2[0], p2[1]);
198    let (x3, y3) = (p3[0], p3[1]);
199
200    // Compute the determinant d = 2 * (x1*(y2-y3) + x2*(y3-y1) + x3*(y1-y2))
201    // Visually d is twice the area of the triangle formed by the points,
202    // also the same as: cross(p2 - p1, p3 - p1)
203    let d = 2.0 * (x1 * (y2 - y3) + x2 * (y3 - y1) + x3 * (y1 - y2));
204
205    // If d is nearly zero, the points are collinear, and a unique circle cannot be defined.
206    if d.abs() < f64::EPSILON {
207        return [(x1 + x2 + x3) / 3.0, (y1 + y2 + y3) / 3.0];
208    }
209
210    // squared lengths
211    let p1_sq = x1 * x1 + y1 * y1;
212    let p2_sq = x2 * x2 + y2 * y2;
213    let p3_sq = x3 * x3 + y3 * y3;
214
215    // This formula is derived from the circle equations:
216    //   (x - cx)^2 + (y - cy)^2 = r^2
217    // All 3 points will satisfy this equation, so we have 3 equations. Radius can be eliminated
218    // by subtracting one of the equations from the other two and the remaining 2 equations can
219    // be solved for cx and cy.
220    [
221        (p1_sq * (y2 - y3) + p2_sq * (y3 - y1) + p3_sq * (y1 - y2)) / d,
222        (p1_sq * (x3 - x2) + p2_sq * (x1 - x3) + p3_sq * (x2 - x1)) / d,
223    ]
224}
225
226pub struct CircleParams {
227    pub center: Coords2d,
228    pub radius: f64,
229}
230
231pub fn calculate_circle_from_3_points(points: [Coords2d; 3]) -> CircleParams {
232    let center = calculate_circle_center(points[0], points[1], points[2]);
233    CircleParams {
234        center,
235        radius: distance(center, points[1]),
236    }
237}
238
239#[cfg(test)]
240mod tests {
241    // Here you can bring your functions into scope
242    use std::f64::consts::TAU;
243
244    use approx::assert_relative_eq;
245    use pretty_assertions::assert_eq;
246
247    use super::{Angle, calculate_circle_center, get_x_component, get_y_component};
248
249    static EACH_QUAD: [(i32, [i32; 2]); 12] = [
250        (-315, [1, 1]),
251        (-225, [-1, 1]),
252        (-135, [-1, -1]),
253        (-45, [1, -1]),
254        (45, [1, 1]),
255        (135, [-1, 1]),
256        (225, [-1, -1]),
257        (315, [1, -1]),
258        (405, [1, 1]),
259        (495, [-1, 1]),
260        (585, [-1, -1]),
261        (675, [1, -1]),
262    ];
263
264    #[test]
265    fn test_get_y_component() {
266        let mut expected = Vec::new();
267        let mut results = Vec::new();
268
269        for &(angle, expected_result) in EACH_QUAD.iter() {
270            let res = get_y_component(Angle::from_degrees(angle as f64), 1.0);
271            results.push([res[0].round() as i32, res[1].round() as i32]);
272            expected.push(expected_result);
273        }
274
275        assert_eq!(results, expected);
276
277        let result = get_y_component(Angle::zero(), 1.0);
278        assert_eq!(result[0] as i32, 1);
279        assert_eq!(result[1] as i32, 0);
280
281        let result = get_y_component(Angle::from_degrees(90.0), 1.0);
282        assert_eq!(result[0] as i32, 1);
283        assert!(result[1] > 100000.0);
284
285        let result = get_y_component(Angle::from_degrees(180.0), 1.0);
286        assert_eq!(result[0] as i32, -1);
287        assert!((result[1] - 0.0).abs() < f64::EPSILON);
288
289        let result = get_y_component(Angle::from_degrees(270.0), 1.0);
290        assert_eq!(result[0] as i32, -1);
291        assert!(result[1] < -100000.0);
292    }
293
294    #[test]
295    fn test_get_x_component() {
296        let mut expected = Vec::new();
297        let mut results = Vec::new();
298
299        for &(angle, expected_result) in EACH_QUAD.iter() {
300            let res = get_x_component(Angle::from_degrees(angle as f64), 1.0);
301            results.push([res[0].round() as i32, res[1].round() as i32]);
302            expected.push(expected_result);
303        }
304
305        assert_eq!(results, expected);
306
307        let result = get_x_component(Angle::zero(), 1.0);
308        assert!(result[0] > 100000.0);
309        assert_eq!(result[1] as i32, 1);
310
311        let result = get_x_component(Angle::from_degrees(90.0), 1.0);
312        assert!((result[0] - 0.0).abs() < f64::EPSILON);
313        assert_eq!(result[1] as i32, 1);
314
315        let result = get_x_component(Angle::from_degrees(180.0), 1.0);
316        assert!(result[0] < -100000.0);
317        assert_eq!(result[1] as i32, 1);
318
319        let result = get_x_component(Angle::from_degrees(270.0), 1.0);
320        assert!((result[0] - 0.0).abs() < f64::EPSILON);
321        assert_eq!(result[1] as i32, -1);
322    }
323
324    #[test]
325    fn test_arc_center_and_end() {
326        let (center, end) = super::arc_center_and_end([0.0, 0.0], Angle::zero(), Angle::from_degrees(90.0), 1.0);
327        assert_eq!(center[0].round(), -1.0);
328        assert_eq!(center[1], 0.0);
329        assert_eq!(end[0].round(), -1.0);
330        assert_eq!(end[1], 1.0);
331
332        let (center, end) = super::arc_center_and_end([0.0, 0.0], Angle::zero(), Angle::from_degrees(180.0), 1.0);
333        assert_eq!(center[0].round(), -1.0);
334        assert_eq!(center[1], 0.0);
335        assert_eq!(end[0].round(), -2.0);
336        assert_eq!(end[1].round(), 0.0);
337
338        let (center, end) = super::arc_center_and_end([0.0, 0.0], Angle::zero(), Angle::from_degrees(180.0), 10.0);
339        assert_eq!(center[0].round(), -10.0);
340        assert_eq!(center[1], 0.0);
341        assert_eq!(end[0].round(), -20.0);
342        assert_eq!(end[1].round(), 0.0);
343    }
344
345    #[test]
346    fn test_calculate_circle_center() {
347        const EPS: f64 = 1e-4;
348
349        // Test: circle center = (4.1, 1.9)
350        let p1 = [1.0, 2.0];
351        let p2 = [4.0, 5.0];
352        let p3 = [7.0, 3.0];
353        let center = calculate_circle_center(p1, p2, p3);
354        assert_relative_eq!(center[0], 4.1, epsilon = EPS);
355        assert_relative_eq!(center[1], 1.9, epsilon = EPS);
356
357        // Tests: Generate a few circles and test its points
358        let center = [3.2, 0.7];
359        let radius_array = [0.001, 0.01, 0.6, 1.0, 5.0, 60.0, 500.0, 2000.0, 400_000.0];
360        let points_array = [[0.0, 0.33, 0.66], [0.0, 0.1, 0.2], [0.0, -0.1, 0.1], [0.0, 0.5, 0.7]];
361
362        let get_point = |radius: f64, t: f64| {
363            let angle = t * TAU;
364            [
365                center[0] + radius * libm::cos(angle),
366                center[1] + radius * libm::sin(angle),
367            ]
368        };
369
370        for radius in radius_array {
371            for point in points_array {
372                let p1 = get_point(radius, point[0]);
373                let p2 = get_point(radius, point[1]);
374                let p3 = get_point(radius, point[2]);
375                let c = calculate_circle_center(p1, p2, p3);
376                assert_relative_eq!(c[0], center[0], epsilon = EPS);
377                assert_relative_eq!(c[1], center[1], epsilon = EPS);
378            }
379        }
380
381        // Test: Equilateral triangle
382        let p1 = [0.0, 0.0];
383        let p2 = [1.0, 0.0];
384        let p3 = [0.5, 3.0_f64.sqrt() / 2.0];
385        let center = calculate_circle_center(p1, p2, p3);
386        assert_relative_eq!(center[0], 0.5, epsilon = EPS);
387        assert_relative_eq!(center[1], 1.0 / (2.0 * 3.0_f64.sqrt()), epsilon = EPS);
388
389        // Test: Collinear points (should return the average of the points)
390        let p1 = [0.0, 0.0];
391        let p2 = [1.0, 0.0];
392        let p3 = [2.0, 0.0];
393        let center = calculate_circle_center(p1, p2, p3);
394        assert_relative_eq!(center[0], 1.0, epsilon = EPS);
395        assert_relative_eq!(center[1], 0.0, epsilon = EPS);
396
397        // Test: Points forming a circle with radius = 1
398        let p1 = [0.0, 0.0];
399        let p2 = [0.0, 2.0];
400        let p3 = [2.0, 0.0];
401        let center = calculate_circle_center(p1, p2, p3);
402        assert_relative_eq!(center[0], 1.0, epsilon = EPS);
403        assert_relative_eq!(center[1], 1.0, epsilon = EPS);
404
405        // Test: Integer coordinates
406        let p1 = [0.0, 0.0];
407        let p2 = [0.0, 6.0];
408        let p3 = [6.0, 0.0];
409        let center = calculate_circle_center(p1, p2, p3);
410        assert_relative_eq!(center[0], 3.0, epsilon = EPS);
411        assert_relative_eq!(center[1], 3.0, epsilon = EPS);
412        // Verify radius (should be 3 * sqrt(2))
413        let radius = ((center[0] - p1[0]).powi(2) + (center[1] - p1[1]).powi(2)).sqrt();
414        assert_relative_eq!(radius, 3.0 * 2.0_f64.sqrt(), epsilon = EPS);
415    }
416}
417
418pub(crate) type Coords2d = [f64; 2];
419
420pub fn is_points_ccw_wasm(points: &[f64]) -> i32 {
421    // CCW is positive as that the Math convention
422
423    let mut sum = 0.0;
424    for i in 0..(points.len() / 2) {
425        let point1 = [points[2 * i], points[2 * i + 1]];
426        let point2 = [points[(2 * i + 2) % points.len()], points[(2 * i + 3) % points.len()]];
427        sum += (point2[0] + point1[0]) * (point2[1] - point1[1]);
428    }
429    sum.signum() as i32
430}
431
432pub(crate) fn is_points_ccw(points: &[Coords2d]) -> i32 {
433    let flattened_points: Vec<f64> = points.iter().flat_map(|&p| vec![p[0], p[1]]).collect();
434    is_points_ccw_wasm(&flattened_points)
435}
436
437fn get_slope(start: Coords2d, end: Coords2d) -> (f64, f64) {
438    let slope = if start[0] - end[0] == 0.0 {
439        f64::INFINITY
440    } else {
441        (start[1] - end[1]) / (start[0] - end[0])
442    };
443
444    let perp_slope = if slope == f64::INFINITY { 0.0 } else { -1.0 / slope };
445
446    (slope, perp_slope)
447}
448
449fn get_angle(point1: Coords2d, point2: Coords2d) -> f64 {
450    let delta_x = point2[0] - point1[0];
451    let delta_y = point2[1] - point1[1];
452    let angle = libm::atan2(delta_y, delta_x);
453
454    let result = if angle < 0.0 { angle + 2.0 * PI } else { angle };
455    result * (180.0 / PI)
456}
457
458fn delta_angle(from_angle: f64, to_angle: f64) -> f64 {
459    let norm_from_angle = normalize_rad(from_angle);
460    let norm_to_angle = normalize_rad(to_angle);
461    let provisional = norm_to_angle - norm_from_angle;
462
463    if provisional > -PI && provisional <= PI {
464        provisional
465    } else if provisional > PI {
466        provisional - 2.0 * PI
467    } else if provisional < -PI {
468        provisional + 2.0 * PI
469    } else {
470        provisional
471    }
472}
473
474fn deg2rad(deg: f64) -> f64 {
475    deg * (PI / 180.0)
476}
477
478fn get_mid_point(
479    center: Coords2d,
480    arc_start_point: Coords2d,
481    arc_end_point: Coords2d,
482    tan_previous_point: Coords2d,
483    radius: f64,
484    obtuse: bool,
485) -> Coords2d {
486    let angle_from_center_to_arc_start = get_angle(center, arc_start_point);
487    let angle_from_center_to_arc_end = get_angle(center, arc_end_point);
488    let delta_ang = delta_angle(
489        deg2rad(angle_from_center_to_arc_start),
490        deg2rad(angle_from_center_to_arc_end),
491    );
492    let delta_ang = delta_ang / 2.0 + deg2rad(angle_from_center_to_arc_start);
493    let shortest_arc_mid_point: Coords2d = [
494        libm::cos(delta_ang) * radius + center[0],
495        libm::sin(delta_ang) * radius + center[1],
496    ];
497    let opposite_delta = delta_ang + PI;
498    let longest_arc_mid_point: Coords2d = [
499        libm::cos(opposite_delta) * radius + center[0],
500        libm::sin(opposite_delta) * radius + center[1],
501    ];
502
503    let rotation_direction_original_points = is_points_ccw(&[tan_previous_point, arc_start_point, arc_end_point]);
504    let rotation_direction_points_on_arc = is_points_ccw(&[arc_start_point, shortest_arc_mid_point, arc_end_point]);
505    if rotation_direction_original_points != rotation_direction_points_on_arc && obtuse {
506        longest_arc_mid_point
507    } else {
508        shortest_arc_mid_point
509    }
510}
511
512fn intersect_point_n_slope(point1: Coords2d, slope1: f64, point2: Coords2d, slope2: f64) -> Coords2d {
513    let x = if slope1.abs() == f64::INFINITY {
514        point1[0]
515    } else if slope2.abs() == f64::INFINITY {
516        point2[0]
517    } else {
518        (point2[1] - slope2 * point2[0] - point1[1] + slope1 * point1[0]) / (slope1 - slope2)
519    };
520    let y = if slope1.abs() != f64::INFINITY {
521        slope1 * x - slope1 * point1[0] + point1[1]
522    } else {
523        slope2 * x - slope2 * point2[0] + point2[1]
524    };
525    [x, y]
526}
527
528/// Structure to hold input data for calculating tangential arc information.
529pub struct TangentialArcInfoInput {
530    /// The starting point of the arc.
531    pub arc_start_point: Coords2d,
532    /// The ending point of the arc.
533    pub arc_end_point: Coords2d,
534    /// The point from which the tangent is drawn.
535    pub tan_previous_point: Coords2d,
536    /// Flag to determine if the arc is obtuse. Obtuse means it flows smoothly from the previous segment.
537    pub obtuse: bool,
538}
539
540/// Structure to hold the output data from calculating tangential arc information.
541pub struct TangentialArcInfoOutput {
542    /// The center point of the arc.
543    pub center: Coords2d,
544    /// The midpoint on the arc.
545    pub arc_mid_point: Coords2d,
546    /// The radius of the arc.
547    pub radius: f64,
548    /// Start angle of the arc in radians.
549    pub start_angle: f64,
550    /// End angle of the arc in radians.
551    pub end_angle: f64,
552    /// If the arc is counter-clockwise.
553    pub ccw: i32,
554    /// The length of the arc.
555    pub arc_length: f64,
556}
557
558// tanPreviousPoint and arcStartPoint make up a straight segment leading into the arc (of which the arc should be tangential). The arc should start at arcStartPoint and end at, arcEndPoint
559// With this information we should everything we need to calculate the arc's center and radius. However there is two tangential arcs possible, that just varies on their direction
560// One is obtuse where the arc smoothly flows from the straight segment, and the other would be acute that immediately cuts back in the other direction. The obtuse boolean is there to control for this.
561pub fn get_tangential_arc_to_info(input: TangentialArcInfoInput) -> TangentialArcInfoOutput {
562    let (_, perp_slope) = get_slope(input.tan_previous_point, input.arc_start_point);
563    let tangential_line_perp_slope = perp_slope;
564
565    // Calculate the midpoint of the line segment between arcStartPoint and arcEndPoint
566    let mid_point: Coords2d = [
567        (input.arc_start_point[0] + input.arc_end_point[0]) / 2.0,
568        (input.arc_start_point[1] + input.arc_end_point[1]) / 2.0,
569    ];
570
571    let slope_mid_point_line = get_slope(input.arc_start_point, mid_point);
572
573    let center: Coords2d;
574    let radius: f64;
575
576    if tangential_line_perp_slope == slope_mid_point_line.0 {
577        // can't find the intersection of the two lines if they have the same gradient
578        // but in this case the center is the midpoint anyway
579        center = mid_point;
580        radius =
581            ((input.arc_start_point[0] - center[0]).powi(2) + (input.arc_start_point[1] - center[1]).powi(2)).sqrt();
582    } else {
583        center = intersect_point_n_slope(
584            mid_point,
585            slope_mid_point_line.1,
586            input.arc_start_point,
587            tangential_line_perp_slope,
588        );
589        radius =
590            ((input.arc_start_point[0] - center[0]).powi(2) + (input.arc_start_point[1] - center[1]).powi(2)).sqrt();
591    }
592
593    let arc_mid_point = get_mid_point(
594        center,
595        input.arc_start_point,
596        input.arc_end_point,
597        input.tan_previous_point,
598        radius,
599        input.obtuse,
600    );
601
602    let start_angle = libm::atan2(
603        input.arc_start_point[1] - center[1],
604        input.arc_start_point[0] - center[0],
605    );
606    let end_angle = libm::atan2(input.arc_end_point[1] - center[1], input.arc_end_point[0] - center[0]);
607    let ccw = is_points_ccw(&[input.arc_start_point, arc_mid_point, input.arc_end_point]);
608
609    let arc_mid_angle = libm::atan2(arc_mid_point[1] - center[1], arc_mid_point[0] - center[0]);
610    let start_to_mid_arc_length = radius
611        * delta(Angle::from_radians(start_angle), Angle::from_radians(arc_mid_angle))
612            .to_radians()
613            .abs();
614    let mid_to_end_arc_length = radius
615        * delta(Angle::from_radians(arc_mid_angle), Angle::from_radians(end_angle))
616            .to_radians()
617            .abs();
618    let arc_length = start_to_mid_arc_length + mid_to_end_arc_length;
619
620    TangentialArcInfoOutput {
621        center,
622        radius,
623        arc_mid_point,
624        start_angle,
625        end_angle,
626        ccw,
627        arc_length,
628    }
629}
630
631#[cfg(test)]
632mod get_tangential_arc_to_info_tests {
633    use approx::assert_relative_eq;
634
635    use super::*;
636
637    fn round_to_three_decimals(num: f64) -> f64 {
638        (num * 1000.0).round() / 1000.0
639    }
640
641    #[test]
642    fn test_basic_case() {
643        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
644            tan_previous_point: [0.0, -5.0],
645            arc_start_point: [0.0, 0.0],
646            arc_end_point: [4.0, 0.0],
647            obtuse: true,
648        });
649        assert_relative_eq!(result.center[0], 2.0);
650        assert_relative_eq!(result.center[1], 0.0);
651        assert_relative_eq!(result.arc_mid_point[0], 2.0);
652        assert_relative_eq!(result.arc_mid_point[1], 2.0);
653        assert_relative_eq!(result.radius, 2.0);
654        assert_relative_eq!(result.start_angle, PI);
655        assert_relative_eq!(result.end_angle, 0.0);
656        assert_eq!(result.ccw, -1);
657    }
658
659    #[test]
660    fn basic_case_with_arc_centered_at_0_0_and_the_tangential_line_being_45_degrees() {
661        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
662            tan_previous_point: [0.0, -4.0],
663            arc_start_point: [2.0, -2.0],
664            arc_end_point: [-2.0, 2.0],
665            obtuse: true,
666        });
667        assert_relative_eq!(result.center[0], 0.0);
668        assert_relative_eq!(result.center[1], 0.0);
669        assert_relative_eq!(round_to_three_decimals(result.arc_mid_point[0]), 2.0);
670        assert_relative_eq!(round_to_three_decimals(result.arc_mid_point[1]), 2.0);
671        assert_relative_eq!(result.radius, (2.0f64 * 2.0 + 2.0 * 2.0).sqrt());
672        assert_relative_eq!(result.start_angle, -PI / 4.0);
673        assert_relative_eq!(result.end_angle, 3.0 * PI / 4.0);
674        assert_eq!(result.ccw, 1);
675    }
676
677    #[test]
678    fn test_get_tangential_arc_to_info_moving_arc_end_point() {
679        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
680            tan_previous_point: [0.0, -4.0],
681            arc_start_point: [2.0, -2.0],
682            arc_end_point: [2.0, 2.0],
683            obtuse: true,
684        });
685        let expected_radius = (2.0f64 * 2.0 + 2.0 * 2.0).sqrt();
686        assert_relative_eq!(round_to_three_decimals(result.center[0]), 0.0);
687        assert_relative_eq!(result.center[1], 0.0);
688        assert_relative_eq!(result.arc_mid_point[0], expected_radius);
689        assert_relative_eq!(round_to_three_decimals(result.arc_mid_point[1]), -0.0);
690        assert_relative_eq!(result.radius, expected_radius);
691        assert_relative_eq!(result.start_angle, -PI / 4.0);
692        assert_relative_eq!(result.end_angle, PI / 4.0);
693        assert_eq!(result.ccw, 1);
694    }
695
696    #[test]
697    fn test_get_tangential_arc_to_info_moving_arc_end_point_again() {
698        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
699            tan_previous_point: [0.0, -4.0],
700            arc_start_point: [2.0, -2.0],
701            arc_end_point: [-2.0, -2.0],
702            obtuse: true,
703        });
704        let expected_radius = (2.0f64 * 2.0 + 2.0 * 2.0).sqrt();
705        assert_relative_eq!(result.center[0], 0.0);
706        assert_relative_eq!(result.center[1], 0.0);
707        assert_relative_eq!(result.radius, expected_radius);
708        assert_relative_eq!(round_to_three_decimals(result.arc_mid_point[0]), 0.0);
709        assert_relative_eq!(result.arc_mid_point[1], expected_radius);
710        assert_relative_eq!(result.start_angle, -PI / 4.0);
711        assert_relative_eq!(result.end_angle, -3.0 * PI / 4.0);
712        assert_eq!(result.ccw, 1);
713    }
714
715    #[test]
716    fn test_get_tangential_arc_to_info_acute_moving_arc_end_point() {
717        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
718            tan_previous_point: [0.0, -4.0],
719            arc_start_point: [2.0, -2.0],
720            arc_end_point: [-2.0, -2.0],
721            obtuse: false,
722        });
723        let expected_radius = (2.0f64 * 2.0 + 2.0 * 2.0).sqrt();
724        assert_relative_eq!(result.center[0], 0.0);
725        assert_relative_eq!(result.center[1], 0.0);
726        assert_relative_eq!(result.radius, expected_radius);
727        assert_relative_eq!(round_to_three_decimals(result.arc_mid_point[0]), -0.0);
728        assert_relative_eq!(result.arc_mid_point[1], -expected_radius);
729        assert_relative_eq!(result.start_angle, -PI / 4.0);
730        assert_relative_eq!(result.end_angle, -3.0 * PI / 4.0);
731        // would be cw if it was obtuse
732        assert_eq!(result.ccw, -1);
733    }
734
735    #[test]
736    fn test_get_tangential_arc_to_info_obtuse_with_wrap_around() {
737        let arc_end = libm::cos(std::f64::consts::PI / 4.0) * 2.0;
738        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
739            tan_previous_point: [2.0, -4.0],
740            arc_start_point: [2.0, 0.0],
741            arc_end_point: [0.0, -2.0],
742            obtuse: true,
743        });
744        assert_relative_eq!(result.center[0], -0.0);
745        assert_relative_eq!(result.center[1], 0.0);
746        assert_relative_eq!(result.radius, 2.0);
747        assert_relative_eq!(result.arc_mid_point[0], -arc_end);
748        assert_relative_eq!(result.arc_mid_point[1], arc_end);
749        assert_relative_eq!(result.start_angle, 0.0);
750        assert_relative_eq!(result.end_angle, -PI / 2.0);
751        assert_eq!(result.ccw, 1);
752    }
753
754    #[test]
755    fn test_arc_length_obtuse_cw() {
756        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
757            tan_previous_point: [-1.0, -1.0],
758            arc_start_point: [-1.0, 0.0],
759            arc_end_point: [0.0, -1.0],
760            obtuse: true,
761        });
762        let circumference = 2.0 * PI * result.radius;
763        let expected_length = circumference * 3.0 / 4.0; // 3 quarters of a circle circle
764        assert_relative_eq!(result.arc_length, expected_length);
765    }
766
767    #[test]
768    fn test_arc_length_acute_cw() {
769        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
770            tan_previous_point: [-1.0, -1.0],
771            arc_start_point: [-1.0, 0.0],
772            arc_end_point: [0.0, 1.0],
773            obtuse: true,
774        });
775        let circumference = 2.0 * PI * result.radius;
776        let expected_length = circumference / 4.0; // 1 quarters of a circle circle
777        assert_relative_eq!(result.arc_length, expected_length);
778    }
779
780    #[test]
781    fn test_arc_length_obtuse_ccw() {
782        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
783            tan_previous_point: [1.0, -1.0],
784            arc_start_point: [1.0, 0.0],
785            arc_end_point: [0.0, -1.0],
786            obtuse: true,
787        });
788        let circumference = 2.0 * PI * result.radius;
789        let expected_length = circumference * 3.0 / 4.0; // 1 quarters of a circle circle
790        assert_relative_eq!(result.arc_length, expected_length);
791    }
792
793    #[test]
794    fn test_arc_length_acute_ccw() {
795        let result = get_tangential_arc_to_info(TangentialArcInfoInput {
796            tan_previous_point: [1.0, -1.0],
797            arc_start_point: [1.0, 0.0],
798            arc_end_point: [0.0, 1.0],
799            obtuse: true,
800        });
801        let circumference = 2.0 * PI * result.radius;
802        let expected_length = circumference / 4.0; // 1 quarters of a circle circle
803        assert_relative_eq!(result.arc_length, expected_length);
804    }
805}
806
807pub(crate) fn get_tangent_point_from_previous_arc(
808    last_arc_center: Coords2d,
809    last_arc_ccw: bool,
810    last_arc_end: Coords2d,
811) -> Coords2d {
812    let angle_from_old_center_to_arc_start = get_angle(last_arc_center, last_arc_end);
813    let tangential_angle = angle_from_old_center_to_arc_start + if last_arc_ccw { -90.0 } else { 90.0 };
814    // What is the 10.0 constant doing???
815    [
816        libm::cos(tangential_angle.to_radians()) * 10.0 + last_arc_end[0],
817        libm::sin(tangential_angle.to_radians()) * 10.0 + last_arc_end[1],
818    ]
819}