Bitmap/Bézier curves/Quadratic

Bitmap/Bézier curves/Quadratic
You are encouraged to solve this task according to the task description, using any language you may know.

Using the data storage type defined on this page for raster images, and the draw_line function defined in this one, draw a quadratic bezier curves (definition on Wikipedia).

Ada

<lang ada> procedure Quadratic_Bezier

         (  Picture    : in out Image;
            P1, P2, P3 : Point;
            Color      : Pixel;
            N          : Positive := 20
         )  is
  Points : array (0..N) of Point;

begin

  for I in Points'Range loop
     declare
        T : constant Float := Float (I) / Float (N);
        A : constant Float := (1.0 - T)**2;
        B : constant Float := 2.0 * T * (1.0 - T);
        C : constant Float := T**2;
     begin
        Points (I).X := Positive (A * Float (P1.X) + B * Float (P2.X) + C * Float (P3.X));
        Points (I).Y := Positive (A * Float (P1.Y) + B * Float (P2.Y) + C * Float (P3.Y));
     end;
  end loop;
  for I in Points'First..Points'Last - 1 loop
     Line (Picture, Points (I), Points (I + 1), Color);
  end loop;

end Quadratic_Bezier; </lang> The following test <lang ada>

  X : Image (1..16, 1..16);

begin

  Fill (X, White);
  Quadratic_Bezier (X, (8, 2), (13, 8), (2, 15), Black);
  Print (X);

</lang> should produce;

              H
             H
             H
            H
           H
          HH
 HH      H
  HH  HHH
    HH

C

<lang c>#include <math.h>

/* number of segments for the curve */

1. define N_SEG 20

1. define plot(x, y) put_pixel_clip(img, x, y, r, g, b)
2. define line(x0,y0,x1,y1) draw_line(img, x0,y0,x1,y1, r,g,b)

void quad_bezier(

       image img,
       unsigned int x1, unsigned int y1,
       unsigned int x2, unsigned int y2,
       unsigned int x3, unsigned int y3,
       color_component r,
       color_component g,
       color_component b )

{

   unsigned int i;
   double pts[N_SEG+1][2];
   for (i=0; i <= N_SEG; ++i)
   {
       double t = (double)i / (double)N_SEG;
       double a = pow((1.0 - t), 2.0);
       double b = 2.0 * t * (1.0 - t);
       double c = pow(t, 2.0);
       double x = a * x1 + b * x2 + c * x3;
       double y = a * y1 + b * y2 + c * y3;
       pts[i][0] = x;
       pts[i][1] = y;
   }

1. if 0

   /* draw only points */
   for (i=0; i <= N_SEG; ++i)
   {
       plot( pts[i][0],
             pts[i][1] );
   }

1. else

   /* draw segments */
   for (i=0; i < N_SEG; ++i)
   {
       int j = i + 1;
       line( pts[i][0], pts[i][1],
             pts[j][0], pts[j][1] );
   }

1. endif

}

1. undef plot
2. undef line</lang>

OCaml

<lang ocaml>let quad_bezier ~img ~color

       ~p1:(_x1, _y1)
       ~p2:(_x2, _y2)
       ~p3:(_x3, _y3) =
 let (x1, y1, x2, y2, x3, y3) =
   (float _x1, float _y1, float _x2, float _y2, float _x3, float _y3)
 in
 let bz t =
   let a = (1.0 -. t) ** 2.0
   and b = 2.0 *. t *. (1.0 -. t)
   and c = t ** 2.0
   in
   let x = a *. x1 +. b *. x2 +. c *. x3
   and y = a *. y1 +. b *. y2 +. c *. y3
   in
   (int_of_float x, int_of_float y)
 in
 let rec loop _t acc =
   if _t > 20 then acc else
   begin
     let t = (float _t) /. 20.0 in
     let x, y = bz t in
     loop (succ _t) ((x,y)::acc)
   end
 in
 let pts = loop 0 [] in

 (*
 (* draw only points *)
 List.iter (fun (x, y) -> put_pixel img color x y) pts;
 *)

 (* draw segments *)
 let line = draw_line ~img ~color in
 let by_pair li f =
   let rec aux prev = function
     | [] -> ()
     | x::xs ->
         f prev x;
         aux x xs
   in
   aux (List.hd li) (List.tl li)
 in
 by_pair pts (fun p0 p1 -> line ~p0 ~p1);

</lang>

Vedit macro language

This implementation uses de Casteljau's algorithm to recursively split the Bezier curve into two smaller segments until the segment is short enough to be approximated with a straight line. The advantage of this method is that only integer calculations are needed, and the most complex operations are addition and shift right. (I have used multiplication and division here for clarity.)

Constant recursion depth is used here. Recursion depth of 5 seems to give accurate enough result in most situations. In real world implementations, some adaptive method is often used to decide when to stop recursion.

// Daw a Cubic bezier curve
//  #20, #30 = Start point
//  #21, #31 = Control point 1
//  #22, #32 = Control point 2
//  #23, #33 = end point
//  #40 = depth of recursion

:CUBIC_BEZIER:
if (#40 > 0) {
    #24 = (#20+#21)/2;              #34 = (#30+#31)/2
    #26 = (#22+#23)/2;              #36 = (#32+#33)/2
    #27 = (#20+#21*2+#22)/4;        #37 = (#30+#31*2+#32)/4
    #28 = (#21+#22*2+#23)/4;        #38 = (#31+#32*2+#33)/4
    #29 = (#20+#21*3+#22*3+#23)/8;  #39 = (#30+#31*3+#32*3+#33)/8
    Num_Push(20,40)
        #21 = #24; #31 = #34    // control 1
        #22 = #27; #32 = #37    // control 2
        #23 = #29; #33 = #39    // end point
        #40--
        Call("CUBIC_BEZIER")    // Draw "left" part
    Num_Pop(20,40)
    Num_Push(20,40)
        #20 = #29; #30 = #39    // start point
        #21 = #28; #31 = #38    // control 1
        #22 = #26; #32 = #36    // control 2
        #40--
        Call("CUBIC_BEZIER")    // Draw "right" part
    Num_Pop(20,40)
} else {
    #1=#20; #2=#30; #3=#23; #4=#33
    Call("DRAW_LINE")
}
return