SVG Elliptical Arc to Chord Approximation

arcs.cljs

;; NB there is a bug when switching between circles/ellipses.  Anyone spot it?

(ns perturb.arcs
  (:require [perturb.util
             :refer [cos
                     sin
                     PI
                     TAU
                     avg
                     rad
                     sq
                     sqrt
                     pairwise
                     atan
                     +'
                     -'
                     not-neg?
                     polar->cart]]))


(def end        1.000001)
(def angle-zero (rad 0.000001))

;; this assumes we're dealing with circles, but just adding an ry param ought to be fine

(defn- arc-points* [[x y] [x' y'] rx & [{alpha :alpha sweep?* :sweep? large? :large?}]]
  (let [ry      rx
        alpha   (or alpha 0)
        angle   (- PI alpha)
        sweep?  (cond-> sweep?* large? not)
        c       (cos angle)
        s       (sin angle)
        e       (/ rx ry)
        ax      (- (* x c) (* y s))
        ay      (* (+ (* y c) (* x s)) e)
        bx      (- (* x' c) (* y' s))
        by      (* (+ (* y' c) (* x' s)) e)
        sx      (avg ax bx)
        sy      (avg ay by)
        vx      (- ay by)
        vy      (- bx ax)
        l       (sqrt (max 0 (- (/ (sq rx) (+ (sq vx) (sq vy))) 0.25)))
        vx      (* vx l)
        vy      (* vy l)
        [sx sy] (pairwise (if sweep? + -) sx sy vx vy)
        a0      (atan (- ay sy) (- ax sx))
        a1      (atan (- by sy) (- bx sx))
        da      (- a1 a0)
        [a0 a1]
        (cond
          (<= (abs (- (abs da) PI)) angle-zero)
          (let [db     (- (avg a0 a1) (atan (- by ay) (- bx ax)))
                [db]   (drop-while #(< %1 (- PI)) (iterate (+' TAU) db))
                [db]   (drop-while #(< PI     %1) (iterate (-' TAU) db))
                sweep? (or (neg? db) sweep?*)]
            [(cond-> a0 (and sweep?     (neg? da)) (- TAU))
             (cond-> a1 (and sweep? (not-neg? da)) (- TAU))])
          large?
          [(cond-> a0 (and (neg?     da) (< (- PI) da)) (- TAU))
           (cond-> a1 (and (not-neg? da) (< da     PI)) (- TAU))]
          :else
          [(cond-> a0 (< da (- PI)) (- TAU))
           (cond-> a1 (< PI da)     (- TAU))])
        da (- a1 a0)]
    (fn point-calc [incr]
      (let [t (+ a0 (* da incr))
            x (+ sx (* rx (cos t)))
            y (+ (/ sy e) (* ry (sin t)))
            c (cos (- angle))
            s (sin (- angle))]
        [(- (* x c) (* y s))
         (+ (* x s) (* y c))]))))

(defn arc-points
  ([p1 p2 r & [opts]]
   (let [f (arc-points* p1 p2 r opts)]
     (into [] (map f (range 0 end (:incr opts 0.05)))))))

(defn approx-arc [[px py] {radius :r :as circle} t & [opts]]
  (let [[end-x end-y] (polar->cart circle t)]
    (arc-points [px py] [end-x end-y] radius opts)))

(defn approx-arcs [prev & params]
  (second
   (reduce
    (fn [[prev acc] [circle angle sweep? large?]]
      (let [pts  (approx-arc prev circle angle {:large? large? :sweep? sweep?})
            prev (last pts)]
        [prev (into acc pts)]))
    [prev []] params)))