What kind of spline does `draw-spline` draw?

Alternately: how can I convert points that I liked with draw-spline to points consistent with curve-to?

I'm working on drawing some fire with racket/draw + pict (using the dc constructor for picts). Once I get the right outline, I can fill/scale/stack to get a layered effect.

I started with dc<%>'s draw-spline to draw (what I thought were) Bézier curves with a single control point (aka quadratic Bézier curves). When I realized I couldn't fill the resulting shape easily this way, I switched to dc-path%'s curve-to, which uses cubic Bézier curves (i.e., 2 control points). Fortunately, any quadratic Bézier can be converted to a cubic (algorithm - Convert a quadratic bezier to a cubic one - Stack Overflow).

But when I run this transformation on my initial code, I get slightly different curves, so I suspect draw-spline was never a quadratic Bézier!

You can see the discrepancies.

  • First row: draw-spline then curve-to shapes
  • Second row: curve-to on top of draw-spline, then vice-versa.

Fire shape attempts

(The following code was a live experiment; please excuse the abuse of list where structs might be better.)

#lang racket

(require racket/draw
         pict
         )

;; list of point and control points defining a path of quadratic bézier curves
;; (list pt cp)
;; (list p2 __)
;; ------------
;; Draw from point pt through control-point cp to point p2.
;; The last pair's control-point is ignored.
(define q-points
  (list
   (list (list 25 45)
         (list 18 42))
   (list (list 15 30)
         (list 16 23))
   (list (list 18 19)
         (list 22 12))
   (list (list 20 7)
         (list 29 11))
   (list (list 35 27)
         (list 37 28))
   (list (list 40 25)
         (list 41 21))
   (list (list 39 18)
         (list 46 27))
   (list (list 43 35)
         (list 39 43))
   (list (list 25 45) ; must be same as (first (first (first points)))
         #f)))

(define (fire-shape*)
  (dc (lambda (dc dx dy)
        (for ([pc1 q-points]
              [pc2 (cdr q-points)])
          (match-define (list (list p1x p1y) (list cx cy)) pc1)
          (match-define (list (list p2x p2y) _) pc2)
          (send dc draw-spline
                p1x p1y
                cx cy
                p2x p2y)))
      50 50))

;; convert a quadratic bézier curve to a cubic
;; https://stackoverflow.com/q/3162645/4400820
(define (quad->cubic p1 cp p2)
  (define (ct pt cp)
    (+ pt (* 2/3 (- cp pt))))
  (values
   p1
   (map ct p1 cp)
   (map ct p2 cp)
   p2))

;; cubic bézier points
(define c-points
  (append
   (for/list ([pc1 (in-list q-points)]
              [pc2 (in-list (cdr q-points))])
     (match-define (list p1 cp) pc1)
     (match-define (list p2 _) pc2)
     (define-values (p1* cp1 cp2 p2*) (quad->cubic p1 cp p2))
     (list p1* cp1 cp2))
   (list (list (car (last q-points)) #f #f))))

(define path
  (let ([p (new dc-path%)])
    (match-define
      (cons (cons (list x0 y0) _) _)
      c-points)
    (send p move-to x0 y0)
    (for ([pcc1 (in-list c-points)]
          [pcc2 (in-list (cdr c-points))])
      (match-define (list _ (list c1x c1y) (list c2x c2y)) pcc1)
      (match-define (list (list p2x p2y) _ _) pcc2)
      (send p curve-to
            c1x c1y
            c2x c2y
            p2x p2y))
    p))

(define (fire-shape border-color fill-color)
  (dc (lambda (dc dx dy)
        (define old-brush (send dc get-brush))
        (define old-pen (send dc get-pen))
        (send* dc
          (set-brush (new brush% [color fill-color]))
          (set-pen (new pen% [color border-color]))
          (draw-path path dx dy)
          (set-brush old-brush)
          (set-pen old-pen)))
      50 50))

(hc-append
 (fire-shape*)
 (fire-shape "black" "red"))

(hc-append
 (cc-superimpose
  (fire-shape*)
  (fire-shape "black" "red"))
 (cc-superimpose
  (fire-shape "black" "red")
  (fire-shape*)))
1 Like

draw-spline is implemented here: draw/dc.rkt at master · racket/draw · GitHub

The curve is drawn with cairo_curve_to, which says it draws a cubic spline: Paths: Cairo: A Vector Graphics Library

I'll have to study that a bit more; it looks like that code computes a bunch of midpoints and uses those for the cubic spline.

Is there a reason for that particular formulation? It doesn't appear to be the same as the quadratic-to-cubic formula I referenced, so the input points cannot be taken as if for a quadratic spline.

For these type of images MetaPict can be used to design and compute the Bezier curves.
I am sure @benknoble is already done at this point, but I'd like to show how
the curve construct in MetaPict can help to draw these type of curves.

Step 1
Get some points.
It's simple to draw the points on grid paper.
In this case, where ben has a bitmap, I used WebPlotDigizer to find the coordinates
for a points. Among the points I made sure to pick points with horizontal and vertical tangents
as well as points with cusps.

Step 2
Draw and label the points with metapict.

Step 3
Use curve to conect the points with ...
Using (curve A .. B .. C .. D) draws a "nice" curve from A to D through B and C.

Step 4
Insert directions of in going and outgoing tangents in points, where
I am dissatisfied with the automatic curve.

Step 5
Instead of drawing the curves, simply return them. This gives a list of Bezier
curves that can be used anywhere (with or without metapict).

#lang racket
(require metapict)

;; These points were read from the bitmap using WebPlotDigitizer
(define A (pt 45.5 191.3)) ; point with horizontal tangent (lowest point)
(define B (pt 30.0 180.0))
(define C (pt 23.8 160.2)) ; point with vertical tangent   (at the left)
(define D (pt 37.6 123.4)) ; point with vertical tangent   (at the right)
(define E (pt 32.9 113.7)) ; highest point
(define F (pt 54.1 132.3))
(define G (pt 63.3 155.8))
(define H (pt 66.7 156.6)) ; point with hor. tangent
(define I (pt 76.1 146.9))
(define J (pt 71.9 135.7))
(define K (pt 84.2 157.9)) ; point with ver. tangent
(define L (pt 74.8 180.6))

(define pts (list A B C D E F G H I J K L))

(define xmin (apply min (map pt-x pts)))
(define xmax (apply max (map pt-x pts)))
(define ymin (apply min (map pt-y pts)))
(define ymax (apply max (map pt-y pts)))
(define Δx   (- xmax xmin))
(define Δy   (- ymax ymin))

;; Set the smallest x-value to 0.
(set!-values (A B C D E F G H I J K L)
             (apply values
                    (map (λ (p) (pt (- (pt-x p) xmin) (- (pt-y p) ymin)))
                         pts)))

;; Flip the y-axis 
(set!-values (A B C D E F G H I J K L)
             (apply values
                    (map (λ (p) (pt (pt-x p) (+ 10 (- ymax (pt-y p)))))
                         pts)))

(displayln (list xmin xmax ymin ymax Δx Δy)) ; (23.8 84.2 113.7 191.3 60.4 77.7)

;; Size of the resulting pict
(set-curve-pict-size 500 500)

(with-window (window 0 110 0 110) ; the logical coordinates
  ; First make a plot so we can see each point
  (draw (dot-label "A" A)
        (dot-label "B" B)
        (dot-label "C" C (rt))
        (dot-label "D" D)
        (dot-label "E" E)
        (dot-label "F" F)
        (dot-label "G" G)
        (dot-label "H" H)
        (dot-label "I" I (rt))
        (dot-label "J" J)
        (dot-label "K" K)))

(with-window (window 0 110 0 110)
  ; Second, draw three curves through the points
  (draw (dot-label "A" A)
        (dot-label "B" B)
        (dot-label "C" C (rt))
        (dot-label "D" D)
        (dot-label "E" E)
        (dot-label "F" F)
        (dot-label "G" G)
        (dot-label "H" H)
        (dot-label "I" I (rt))
        (dot-label "J" J)
        (dot-label "K" K)

        ; Using `curve` with .. means: draw a nice curve through the points.
        (curve A .. B .. C .. D .. E)
        (curve A .. L .. K .. J)
        (curve E .. F .. G .. H .. I .. J)))

(with-window (window 0 110 0 110)
  (draw (dot-label "A" A)
        (dot-label "B" B)
        (dot-label "C" C (rt))
        (dot-label "D" D)
        (dot-label "E" E)
        (dot-label "F" F)
        (dot-label "G" G)
        (dot-label "H" H)
        (dot-label "I" I (rt))
        (dot-label "J" J)
        (dot-label "K" K)

        ; When two points are connected with .. as in P .. Q,
        ; we can specify the direction the curve leaves P and enters Q
        ; using  P out-dir .. in-dir Q.

        ; The directions up, down, left, right are predefined.
        ; Any vector can be used as a direction.
        ; The call (dir d) gives a vector where the angle is d (degrees) between the vector and the x-axis.        
        (curve A left .. B .. up C up .. up D up .. E)           ; C and D have vertical tangents 
        (curve E (dir -20) .. F .. G .. right H right .. I .. J) ; H has a horizontal tangent
        (curve A right  .. L .. up K up .. J)
        ))

;; The `curve` form returns the computed Bezier curves.

(with-window (window 0 110 0 110)
  (list
   (curve A left .. B .. up C up .. up D up .. E)
   (curve E (dir -20) .. F .. G .. right H right .. I .. J)
   (curve A right  .. L .. up K up .. J)))

The computed Bezier curves were:

(list
 (curve
  #f
  (list
   (bez (pt 45.5 10.0) (pt 38.70407852054976 10.0) (pt 33.57856711100168 15.382250645057347) (pt 30.0 21.30000000000001))
   (bez (pt 30.0 21.30000000000001) (pt 26.35435291213044 27.32867713064906) (pt 23.8 34.05090329181637) (pt 23.8 41.10000000000002))
   (bez (pt 23.8 41.10000000000002) (pt 23.8 54.63158945760202) (pt 37.6 64.36841054239801) (pt 37.6 77.9))
   (bez (pt 37.6 77.9) (pt 37.6 81.68214432046649) (pt 35.86828988758938 85.25609923331395) (pt 32.9 87.60000000000001))))
 (curve
  #f
  (list
   (bez (pt 32.9 87.60000000000001) (pt 42.25706459644486 84.194307006788) (pt 51.0592526208653 78.40240373110163) (pt 54.1 69.0))
   (bez (pt 54.1 69.0) (pt 56.83659739877486 60.538069301868184) (pt 55.03952966189786 49.45395459521617) (pt 63.3 45.5))
   (bez (pt 63.3 45.5) (pt 64.36327481408537 44.991053152537994) (pt 65.52095410162818 44.70000000000002) (pt 66.7 44.70000000000002))
   (bez (pt 66.7 44.70000000000002) (pt 71.77047551332339 44.70000000000002) (pt 75.54030298803202 49.16267560892743) (pt 76.1 54.400000000000006))
   (bez (pt 76.1 54.400000000000006) (pt 76.54746678411853 58.587138134646565) (pt 74.99031922725416 62.73953161961822) (pt 71.9 65.60000000000002))))
 (curve
  #f
  (list
   (bez (pt 45.5 10.0) (pt 56.317102828452065 9.999999999999998) (pt 67.12798211527789 13.10009552015201) (pt 74.8 20.700000000000017))
   (bez (pt 74.8 20.700000000000017) (pt 80.85103285446482 26.694155956014235) (pt 84.2 34.881680874212535) (pt 84.2 43.400000000000006))
   (bez (pt 84.2 43.400000000000006) (pt 84.2 52.42527449297748) (pt 79.5520002268251 60.81434725431754) (pt 71.9 65.60000000000002)))))

2 Likes

I did initially look at doing some of this with metapict, but I couldn't figure out how to size the images the way I wanted (let alone mess with colors and strokes, though that was an after thought).

For example, the whole thing needs to fit inside a (disk 50) when I'm done. Easier to just draw it at that size than to make guesses scaling, but I very much appreciate you providing the resulting points (especially with how nice the result is…). Still, at a glance the points seem to span about 20–80 in the x coordinate and 10–90 in the y coordinate. That's going to be too big :frowning: I may steal the idea of computing the points, however, if I can get it to the right size. Otherwise I'll keep my cruddy points and figure out how to translate them through draw-spline odd formula.

I think that code goes back to the original WxWindows code; see eg old-plt/src/wxmac/src/base/xfspline.cc at master · racket/old-plt · GitHub and old-plt/src/wxmac/src/base/xfspline.cc at 3ff80c07109a59067def52b677a901e00732eb72 · racket/old-plt · GitHub. I don't know why they did it that way, though.

The difficult part is getting the right shape.

The function crop/inked crops a pict to the bounding box of the inked part of the pict.
That is, all the surrounding whitespace is removed.
Since the flame is higher than the wide, we can get an approximate scale by using 50 / inked_height.
Drawing the resulting curve on top of a disk makes it easy to make adjustments.

(require metapict/crop)

(define T ((shifted 5 0) (scaled (* 0.85 (/ 50. 377.))))) 

(define p
  (with-window (window 0 110 0 110)
    (penwidth 2
      (draw (T (curve A left .. B .. up C up .. up D up .. E))
            (T (curve E (dir -20) .. F .. G .. right H right .. I .. J))
            (T (curve A right  .. L .. up K up .. J))))))

(define cropped (crop/inked p))
; (define h (pict-height cropped))
; (displayln h) ; without scaling, the height is 377

(set-curve-pict-size 50 50)
(with-window (window 0 50 0 50)
  (draw (color "red" (disk 50))
        cropped))

(list (T (curve A left .. B .. up C up .. up D up .. E))
      (T (curve E (dir -20) .. F .. G .. right H right .. I .. J))
      (T (curve A right  .. L .. up K up .. J)))

image

The curves are:

(list
 (curve
  #f
  (list
   (bez (pt 10.129310344827587 1.1273209549071619) (pt 9.363191875658792 1.1273209549071619) (pt 8.785382233998863 1.7340733485807356) (pt 8.381962864721485 2.401193633952256))
   (bez (pt 8.381962864721485 2.401193633952256) (pt 7.970981429086324 3.0808190399272815) (pt 7.683023872679046 3.838629681438185) (pt 7.683023872679046 4.633289124668438))
   (bez (pt 7.683023872679046 4.633289124668438) (pt 7.683023872679046 6.158733559543995) (pt 9.23872679045093 7.256385803851235) (pt 9.23872679045093 8.781830238726792))
   (bez (pt 9.23872679045093 8.781830238726792) (pt 9.23872679045093 9.20819929342129) (pt 9.043507480696416 9.611098719935924) (pt 8.708885941644564 9.87533156498674))))
 (curve
  #f
  (list
   (bez (pt 8.708885941644564 9.87533156498674) (pt 9.763727441243784 9.491400657263899) (pt 10.756016542139989 8.838467264116233) (pt 11.098806366047747 7.778514588859418))
   (bez (pt 11.098806366047747 7.778514588859418) (pt 11.407308725326079 6.824583409361799) (pt 11.204721513609176 5.575047931821452) (pt 12.135941644562335 5.129310344827587))
   (bez (pt 12.135941644562335 5.129310344827587) (pt 12.25580684243668 5.071935700219801) (pt 12.38631445442758 5.039124668435016) (pt 12.51923076923077 5.039124668435016))
   (bez (pt 12.51923076923077 5.039124668435016) (pt 13.09083609898208 5.039124668435016) (pt 13.51581664984446 5.542211441324711) (pt 13.578912466843502 6.132625994694962))
   (bez (pt 13.578912466843502 6.132625994694962) (pt 13.629356335079676 6.604650850722757) (pt 13.453815828006107 7.072758869585609) (pt 13.105437665782494 7.395225464190985))))
 (curve
  #f
  (list
   (bez (pt 10.129310344827587 1.1273209549071619) (pt 11.348745013817542 1.1273209549071617) (pt 12.567478089918595 1.4768012191152797) (pt 13.432360742705571 2.333554376657827))
   (bez (pt 13.432360742705571 2.333554376657827) (pt 14.11450635627256 3.009288138277467) (pt 14.492042440318304 3.9322849791884162) (pt 14.492042440318304 4.892572944297084))
   (bez (pt 14.492042440318304 4.892572944297084) (pt 14.492042440318304 5.910011050269345) (pt 13.968063686047923 6.855728801879299) (pt 13.105437665782494 7.395225464190985)))))

So the actual code

  1. Draws a line from point1 to (midpoint point1 control-point);
  2. Draws a cubic from (midpoint (midpoint point1 control-point) control-point) to (midpoint control-point point2) with control-point (midpoint control-point (midpoint control-point point2)); and
  3. Draws a line from (midpoint control-point point2) to point2

where (send dc draw-spline x1 y1 x2 y2 x3 y3) has point1 = (x1,y1), control-point = (x2,y2), and point2 = (x3,y3). I think.

At this point I'm just going to steal Soegaard's points (with suitable re-ordering) and call it a day. Alas, it comes out tiny and upside down. I suppose I will continue to experiment, then.

Here's a dc-path%-based conversion from the original draw-spline points:

(define (avg x . xs)
  (/ (apply + x xs)
     (add1 (length xs))))

(define better-path
  (let ([p (new dc-path%)])
    (match-define
      (cons (cons (list x0 y0) _) _)
      q-points)
    (send p move-to x0 y0)
    (for ([pc1 (in-list q-points)]
          [pc2 (in-list (cdr q-points))])
      (match-define (list (list p1x p1y) (list cx cy)) pc1)
      (match-define (list (list p2x p2y) _) pc2)

      (define x21 (avg p1x cx))
      (define y21 (avg p1y cy))
      (send p line-to x21 y21)

      (define x22 (avg cx p2x))
      (define y22 (avg cy p2y))
      (define xm1 (avg x21 cx))
      (define ym1 (avg y21 cy))
      (define xm2 (avg cx x22))
      (define ym2 (avg cy y22))
      (send p curve-to
            xm1 ym1
            xm2 ym2
            x22 y22)

      (send p line-to p2x p2y))
    p))

It matches exactly, which is no surprise.