B-spline

From Rosetta Code
B-spline is a draft programming task. It is not yet considered ready to be promoted as a complete task, for reasons that should be found in its talk page.
Task

Generate a B-spline curve with a list of 12 points and plot or save image.

Coordinates of control points:

   start=171,171 
   1 185,111, 
   2 202,109, 
   3 202,189 
   4 328,160 
   5 208,254 
   6 241,330 
   7 164,252 
   8 69,278 
   9 139,208 
   10 72,148 
   end=168,172


Rules!!!!

Do not use third party libraries or functions

See also



ALGOL 68[edit]

Translation of: Lua
Which is
Translation of: Wren
Suppresses unused parts of the plot.
BEGIN # construct a B-Spline                                                 #

    # mode to hold a B Spline                                                #
    MODE BSPLINE = STRUCT( FLEX[ 1 : 0, 1 : 2 ]INT control points
                         , INT n
                         , INT k
                         , FLEX[ 1 : 0 ]INT t
                         );
    PROC  uniform knot vector = ( INT lb, ub )[]INT:
          BEGIN
             [ lb : ub ]INT range; 
             FOR n FROM lb TO ub DO range[ n ] := n OD;
	     range
          END # uniform knot vector # ;

    PROC calculate bspline = ( REF BSPLINE bs, INT i, k, x )REAL:
	 IF   k = 1
         THEN ABS ( ( t OF bs )[ i ] <= x AND x < ( t OF bs )[ i + 1 ] )
	 ELSE
              PROC helper = ( REF BSPLINE bs, INT i, k, x )REAL:
                   IF ( t OF bs )[ i + k ] /= ( t OF bs )[ i ]
                   THEN ( x - ( t OF bs )[ i ] ) / ( ( t OF bs )[ i + k ] - ( t OF bs )[ i ] )
                   ELSE 0
                   FI # helper # ;
              (     helper( bs, i,     k - 1, x ) ) * calculate bspline( bs, i,     k - 1, x )
            + ( 1 - helper( bs, i + 1, k - 1, x ) ) * calculate bspline( bs, i + 1, k - 1, x )
         FI # calculate bspline # ;

    PROC round = ( REAL n )INT: ENTIER( n + 0.5 );

    PROC get bspline points = ( REF BSPLINE bs )[,]INT:
         BEGIN
            INT p from = ( t OF bs )[     k OF bs ];
            INT p to   = ( t OF bs )[ 1 + n OF bs ] - 1;
            [ p from : p to, 1 : 2 ]INT points;
            FOR x FROM p from TO p to DO
                REAL sum x := 0;
                REAL sum y := 0;
		FOR i TO n OF bs DO
                    REAL f = calculate bspline( bs, i, k OF bs, x );
                    sum x +:= f * ( control points OF bs )[ i, 1 ];
                    sum y +:= f * ( control points OF bs )[ i, 2 ]
                OD;
		points[ x, 1 ] := round( sum x );
                points[ x, 2 ] := round( sum y )
            OD;
            points
         END # get bspline points # ;

    PROC raytrace = ( INT x0, y0, x2, y2, REF[,]BOOL plot, PROC( REF[,]BOOL, INT, INT )VOID visit )VOID:
         BEGIN
            INT x    := x0;
            INT y    := y0;
	    INT dx   := ABS ( x2 - x );
            INT dy   := ABS ( y2 - y );
	    INT n     = 1 + dx + dy;
	    INT dir x = IF x2 > x THEN 1 ELSE -1 FI;
            INT dir y = IF y2 > y THEN 1 ELSE -1 FI;
	    INT err  := dx - dy;
            dx *:= 2;
            dy *:= 2;
            FOR i TO n DO
                visit( plot, x, y );
                IF   err > 0
                THEN x +:= dir x; err -:= dy
                ELSE y +:= dir y; err +:= dx
                FI
            OD
         END # raytrace # ;

    PROC plot line = ( REF[,]BOOL plot, INT x1, y1, x2, y2 )VOID:
	raytrace( x1, y1, x2, y2, plot
                , ( REF[,]BOOL plot, INT x, INT y )VOID: IF   x >= 0
                                                         AND  y >= 0
                                                         AND  x < 1 UPB plot
                                                         AND  y < 2 UPB plot
                                                         THEN plot[ x + 1, y + 1 ] := TRUE
		                                         FI
                );

    PROC plot bspline = ( REF BSPLINE bs, REF[,]BOOL plot, REAL scale x, scale y )VOID:
         IF k OF bs > n OF bs OR k OF bs < 1 THEN
            print( ( "k (= ", whole( k OF bs, 0 ), ") can't be more than ", whole( n OF bs, 0 ), " or less than 1." ) );
            stop
         ELSE
            [,]INT points = get bspline points( bs );
            # Plot the curve. #
            FOR i FROM 1 LWB points TO 1 UPB points - 1 DO
                INT p1x = points[ i,     1 ], p1y = points[ i,     2 ];
                INT p2x = points[ i + 1, 1 ], p2y = points[ i + 1, 2 ];
                plot line( plot
                         , round( p1x * scale x ), round( p1y * scale y )
                         , round( p2x * scale x ), round( p2y * scale y )
                         )
            OD
	 FI # plot bspline # ;

    # print the plot - outputs @ or blank depending on whether the point is plotted or not #
    PROC print plot = ( [,]BOOL plot )VOID:
         FOR row FROM 1 LWB plot
         TO BEGIN # find the highest used row #
                INT max row := 1 UPB plot;
                WHILE IF max row < 1 LWB plot
                      THEN FALSE
                      ELSE
                          BOOL empty row := TRUE;
                          FOR column FROM 2 LWB plot TO 2 UPB plot
                          WHILE empty row := NOT plot[ column, max row ]
                          DO
                              SKIP
                          OD;
                          empty row
                      FI
                DO
                    max row -:= 1
                OD;
                max row
            END
         DO
             INT max column := 2 UPB plot;
             WHILE IF max column < 2 LWB plot THEN FALSE ELSE NOT plot[ max column, row ] FI
             DO
                 max column -:= 1
             OD;
             FOR column FROM 2 LWB plot TO max column DO
                 print( ( IF plot[ column, row ] THEN "@" ELSE " " FI ) )
             OD;
             print( ( newline ) )
         OD # print plot # ;

    # task #
    [,]INT  control points
        =   ( ( 171, 171 ), ( 185, 111 ), ( 202, 109 ), ( 202, 189 ), ( 328, 160 ), ( 208, 254 )
            , ( 241, 330 ), ( 164, 252 ), (  69, 278 ), ( 139, 208 ), (  72, 148 ), ( 168, 172 )
            );
    INT     k = 4;  # Polynomial degree is one less than this i.e. cubic. #
    BSPLINE bs
        := BSPLINE( control points
                  , UPB control points
                  , k
                  , uniform knot vector( 1, UPB control points + k )
                  );
 
    REAL scale x = 0.4; # Since we print the plot to the console as text let's scale things appropriately. #
    REAL scale y = 0.2;
    [ 1 : 350, 1 : 350 ]BOOL plot;
    FOR r FROM 1 LWB plot TO 1 UPB plot DO FOR c FROM 2 LWB plot TO 2 UPB plot DO plot[ c, r ] := FALSE OD OD;
    plot bspline( bs, plot, scale x, scale y );
    print plot( plot )
END
Output:

leading blank lines removed...

                                                                          @@@@
                                                                             @@@@
                                                                                @@
                                                                                 @@
                                                                                  @@
                                                                                   @@
                                                                                    @@
                                                                                     @@
                                                                                      @@
                                                                                       @@
                                                                                        @@@@@@@@
                                                                                               @@@@@@@@@@@@@@
                                                                                                            @@@@@@@@
                                                                                                                 @@
                                                                                                               @@@
                                                                                                              @@
                                                                                                            @@@
                                                                                                           @@
                                              @                                                          @@@
                                             @@                                                         @@
                                            @@                                                        @@@
                                            @                                                        @@
                                           @@                                                      @@@
                                          @@                                                      @@
                                         @@                                                     @@@
                                         @                                                     @@
                                        @@                                                    @@
                                       @@                                                    @@
                                       @@@@@@@                                               @
                                             @@@@@@@@@@@@@@                                 @@
                                                          @@@@@@@@@                         @
                                                                  @@@@                      @
                                                                     @@@@@                 @@
                                                                         @@@@@             @
                                                                             @@@@         @@
                                                                                @@@@@     @
                                                                                    @@@@ @@
                                                                                       @@@

Julia[edit]

Choose BSpline D of 2, ie degree 1.

using Graphics, Plots

Point(t::Tuple) = Vec2(Float64(t[1]), Float64(t[2]))
const controlpoints = Point.([(171, 171), (185, 111), (202, 109), (202, 189), (328, 160),
    (208, 254), (241, 330), (164,252), (69, 278), (139, 208), (72, 148), (168, 172)])
plt = plot(map(a -> a.x, controlpoints), map(a -> a.y, controlpoints))
savefig(plt, "BSplineplot.png")

Lua[edit]

Translation of: Wren
local function Range(from, to)
	local range = {}
	for n = from, to do  table.insert(range, n)  end
	return range
end

local function Bspline(controlPoints, k)
	return {
		controlPoints = controlPoints,
		n = #controlPoints,
		k = k,
		t = Range(1, #controlPoints+k), -- Use a uniform knot vector, delta=1.
	}
end

local function helper(bspline, i, k, x)
	return (bspline.t[i+k] ~= bspline.t[i])
	   and (x - bspline.t[i]) / (bspline.t[i+k] - bspline.t[i])
	   or  0
end

local function calculateBspline(bspline, i, k, x)
	if k == 1 then
		return (bspline.t[i] <= x and x < bspline.t[i+1]) and 1 or 0
	end
	return (  helper(bspline, i  , k-1, x)) * calculateBspline(bspline, i  , k-1, x)
	     + (1-helper(bspline, i+1, k-1, x)) * calculateBspline(bspline, i+1, k-1, x)
end

local function round(n)
	return math.floor(n+.5)
end

local function getBsplinePoints(bspline)
	local points = {}

	for x = bspline.t[bspline.k], bspline.t[bspline.n+1]-1 do
		local sumX = 0
		local sumY = 0

		for i = 1, bspline.n do
			local f = calculateBspline(bspline, i, bspline.k, x)
			sumX    = sumX + f * bspline.controlPoints[i].x
			sumY    = sumY + f * bspline.controlPoints[i].y
		end
		table.insert(points, {x=round(sumX), y=round(sumY)})
	end

	return points
end

local function Plot(unscaledWidth,unscaledHeight, scaleX,scaleY)
	local plot = {
		width  = round(unscaledWidth  * scaleX),
		height = round(unscaledHeight * scaleY),
		scaleX = scaleX,
		scaleY = scaleY,
	}
	for row = 1, plot.height do
		plot[row] = {}
	end
	return plot
end

local function raytrace(x,y, x2,y2, visit)
	local dx   = math.abs(x2 - x)
	local dy   = math.abs(y2 - y)
	local n    = 1 + dx + dy
	local dirX = (x2 > x) and 1 or -1
	local dirY = (y2 > y) and 1 or -1
	local err  = dx - dy
	dx, dy     = 2*dx, 2*dy

	for n = 1, n do
		visit(x, y)
		if err > 0 then  x, err = x+dirX, err-dy
		else             y, err = y+dirY, err+dx  end
	end
end

local function plotLine(plot, x1,y1, x2,y2)
	raytrace(x1,y1, x2,y2, function(x, y)
		if x >= 0 and y >= 0 and x < plot.width and y < plot.height then
			plot[y+1][x+1] = true
		end
	end)
end

local function plotBspline(bspline, plot)
	if bspline.k > bspline.n or bspline.k < 1 then
		error("k (= "..bspline.k..") can't be more than "..bspline.n.." or less than 1.")
	end
	local points = getBsplinePoints(bspline)

	-- Plot the curve.
	for i = 1, #points-1 do
		local p1 = points[i]
		local p2 = points[i+1]
		plotLine(plot,
			round(p1.x*plot.scaleX), round(p1.y*plot.scaleY),
			round(p2.x*plot.scaleX), round(p2.y*plot.scaleY)
		)
	end
end

local function printPlot(plot)
	for row = 1, plot.height do
		for column = 1, plot.width do
			io.write(plot[row][column] and "@" or " ")
		end
		io.write("\n")
	end
end

local controlPoints = {
	{x=171, y=171}, {x=185, y=111}, {x=202, y=109}, {x=202, y=189}, {x=328, y=160}, {x=208, y=254},
	{x=241, y=330}, {x=164, y=252}, {x= 69, y=278}, {x=139, y=208}, {x= 72, y=148}, {x=168, y=172},
}
local k       = 4 -- Polynomial degree is one less than this i.e. cubic.
local bspline = Bspline(controlPoints, k)

local scaleX = .4 -- Since we print the plot to the console as text let's scale things appropriately.
local scaleY = .2
local plot   = Plot(350,350, scaleX,scaleY)
plotBspline(bspline, plot)

printPlot(plot)
Output:
                                   @@@@
                                      @@@@
                                         @@
                                          @@
                                           @@
                                            @@
                                             @@
                                              @@
                                               @@
                                                @@
                                                 @@@@@@@@
                                                        @@@@@@@@@@@@@@
                                                                     @@@@@@@@
                                                                          @@
                                                                        @@@
                                                                       @@
                                                                     @@@
                                                                    @@
       @                                                          @@@
      @@                                                         @@
     @@                                                        @@@
     @                                                        @@
    @@                                                      @@@
   @@                                                      @@
  @@                                                     @@@
  @                                                     @@
 @@                                                    @@
@@                                                    @@
@@@@@@@                                               @
      @@@@@@@@@@@@@@                                 @@
                   @@@@@@@@@                         @
                           @@@@                      @
                              @@@@@                 @@
                                  @@@@@             @
                                      @@@@         @@
                                         @@@@@     @
                                             @@@@ @@
                                                @@@

Mathematica/Wolfram Language[edit]

Graphics[
 BSplineCurve[{{171, 171}, {185, 111}, {202, 109}, {202, 189}, {328, 
    160}, {208, 254}, {241, 330}, {164, 252}, {69, 278}, {139, 
    208}, {72, 148}, {168, 172}}, SplineClosed -> True, 
  SplineDegree -> 2]]
Output:

Outputs a graphical representation of a B-spline.

Perl[edit]

Translation of: Raku
use strict;
use warnings;
use Class::Struct;
use Cairo;

{ package Line;
    struct( A => '@', B => '@');
}

my ($WIDTH,  $HEIGHT,  $W_LINE,  $CURVE_F,  $DETACHED,            $OUTPUT  ) =
   (   400,      400,        2,      0.25,          0,  'run/b-spline.png' );

my @pt = (
   [171, 171], [185, 111], [202, 109], [202, 189], [328, 160], [208, 254],
   [241, 330], [164, 252], [ 69, 278], [139, 208], [ 72, 148], [168, 172]
);
my $cnt = @pt;

sub angle {
    my($g) = @_;
    atan2 $g->B->[1] - $g->A->[1], $g->B->[0] - $g->A->[0]
}

sub control_points {
    my($g, $l) = @_;

    my $h = Line->new;
    my $lgt = sqrt( ($g->B->[0] - $l->A->[0])**2 + ($g->B->[1] - $l->A->[1])**2 );

    @{$h->B} = @{$l->A};
    @{$h->A} = ($g->B->[0] - $lgt * cos(angle $g) , $g->B->[1] - $lgt * sin(angle $g));
    my $a    = angle $h;
    my @p1   = ($g->B->[0] + $lgt * cos($a) * $CURVE_F,  $g->B->[1] + $lgt * sin($a) * $CURVE_F);

    @{$h->A} = @{$g->B};
    @{$h->B} = ($l->A->[0] + $lgt * cos(angle $l) , $l->A->[1] + $lgt * sin(angle $l));
    $a       = angle $h;
    my @p2   = ($l->A->[0] - $lgt * cos($a) * $CURVE_F,  $l->A->[1] - $lgt * sin($a) * $CURVE_F);

    \@p1, \@p2
}

my $surf = Cairo::ImageSurface->create ('argb32', $WIDTH, $HEIGHT);
my $cr   = Cairo::Context->create ($surf);
$cr->set_line_width($W_LINE);
$cr->move_to($pt[$DETACHED - 1 + $cnt][0], $pt[$DETACHED - 1 + $cnt][1]);

my Line ($g,$l);
for my $j ($DETACHED..$cnt-1) {
    $g = Line->new( A=>$pt[($j + $cnt - 2) % $cnt], B=>$pt[($j + $cnt - 1) % $cnt]);
    $l = Line->new( A=>$pt[($j + $cnt + 0) % $cnt], B=>$pt[($j + $cnt + 1) % $cnt]);
    my($p1,$p2) = control_points($g, $l);
    $cr->curve_to($$p1[0], $$p1[1], $$p2[0], $$p2[1], $pt[$j][0], $pt[$j][1]);
}
$cr->stroke;
$surf->write_to_png($OUTPUT);

Output: b-spline.png (offsite image)

Phix[edit]

Translation of: Wren
Library: Phix/pGUI
Library: Phix/online

You can run this online here.

--
-- demo\rosetta\B-spline.exw
-- =========================
--
--  Use +/- to change the order between k = 1 and k = 4.
--
with javascript_semantics
include pGUI.e
include IupGraph.e

constant ctrl_points = {{171, 171}, {185, 111}, {202, 109}, {202, 189}, {328, 160}, {208, 254},
                        {241, 330}, {164, 252}, { 69, 278}, {139, 208}, { 72, 148}, {168, 172}}
integer k = 2, n
sequence t

function w(integer i, k, x)  // B-spline helper function
    return iff(t[i+k]!=t[i] ? (x-t[i])/(t[i+k]-t[i]) : 0 )
end function

function b(integer i, k, x)  // B-spline function
    if k==1 then return iff(t[i]<=x and x<t[i+1] ? 1 : 0) end if
    return w(i,k-1,x)*b(i,k-1,x) + (1-w(i+1,k-1,x))*b(i+1,k-1,x)
end function

function b_spline(Ihandle graph)
    n = length(ctrl_points)
    t = tagset(n+1+k) // use a uniform knot vector, delta = 1
    assert(k<=n+1 and k>=1,"k (= %d) cannot be more than %d or less than 1.",{k,n+1})

    sequence px = {}, py = {}
    for x=t[k] to t[n+1] do
        atom sumX = 0,
             sumY = 0
        for i=1 to n do
            atom f = b(i,k,x)
            sumX += f*ctrl_points[i][1]
            sumY += f*ctrl_points[i][2]
        end for
        px &= round(sumX)
        py &= round(sumY)
    end for

    integer xtick = 40,
            ytick = 40,
            xmin = trunc(min(px)/xtick)*xtick,
            xmax = ceil(max(px)/xtick)*xtick,
            ymin = trunc(min(py)/ytick)*ytick,
            ymax = ceil(max(py)/ytick)*ytick

    IupSetInt(graph,"XTICK",xtick)
    IupSetInt(graph,"XMIN",xmin)
    IupSetInt(graph,"XMAX",xmax)
    IupSetInt(graph,"YTICK",ytick)
    IupSetInt(graph,"YMIN",ymin)
    IupSetInt(graph,"YMAX",ymax)
    sequence graphdata = {{px,py,CD_BLUE}}
    return graphdata
end function

procedure set_title(Ihandle dlg)
    IupSetStrAttribute(dlg, "TITLE", "B-spline curve (order k = %d)",{k})
end procedure

function key_cb(Ihandle dlg, atom c)
    if c=K_ESC then return IUP_CLOSE end if
    if c='+' then k = min(k+1,4) end if
    if c='-' then k = max(k-1,1) end if
    set_title(dlg)
    IupRedraw(dlg)
    return IUP_IGNORE
end function

procedure main()
    IupOpen()
    Ihandle graph = IupGraph(b_spline,`RASTERSIZE=600x600`)
    Ihandle dlg = IupDialog(graph)
    IupSetCallback(dlg, "KEY_CB", Icallback("key_cb"))
    set_title(dlg)
    IupShow(dlg)
    if platform()!=JS then
        IupMainLoop()
        IupClose()
    end if
end procedure

main()

Processing[edit]

//Aamrun, 26th June 2022

int abscissae[]={171,185,202,202,328,208,241,164,69,139,72,168};
int ordinates[] = {171,111,109,189 ,160,254 ,330 ,252,278 ,208 ,148 ,172};

void setup() {
  size(450, 450);
  background(255);
  smooth();

  noFill();
  stroke(0);
  
  beginShape();
  for(int i=0;i<abscissae.length;i++){
    curveVertex(abscissae[i], ordinates[i]);
  }
  endShape();

  fill(255, 0, 0);
  noStroke();
  for (int i = 0; i < abscissae.length; i ++) {
    ellipse(abscissae[i], ordinates[i], 3, 3);
  }
}

Raku[edit]

A minimal translation of this C program, by Bernhard R. Fischer.

# 20211112 Raku programming solution

use Cairo;

# class point_t { has Num ($.x,$.y) is rw } # get by with two element lists
class line_t  { has ($.A,$.B) is rw }

my (\WIDTH,  \HEIGHT,  \W_LINE,  \CURVE_F,  \DETACHED,          \OUTPUT  ) = 
       400,      400,        2,      0.25,          0,  './b-spline.png' ;

my \cnt = #`(Number of points) ( my \pt = [
   [171, 171], [185, 111], [202, 109], [202, 189], [328, 160], [208, 254],
   [241, 330], [164, 252], [ 69, 278], [139, 208], [ 72, 148], [168, 172], ]
).elems;

sub angle(\g) { atan2(g.B.[1] - g.A.[1], g.B.[0] - g.A.[0]) }

sub control_points(\g, \l, @p1, @p2){

#`[ This function calculates the control points. It takes two lines g and l as
 * arguments but it takes three lines into account for calculation. This is
 * line g (P0/P1), line h (P1/P2), and line l (P2/P3). The control points being
 * calculated are actually those for the middle line h, this is from P1 to P2.
 * Line g is the predecessor and line l the successor of line h.
 * @param g Pointer to first line (P0 to P1)
 * @param l Pointer to third line (P2 to P3)
 * @param p1 Pointer to memory of first control point. 
 * @param p2 Pointer to memory of second control point. ]
 
   my \h = $ = line_t.new;

   my \lgt = sqrt([+]([ g.B.[0]-l.A.[0], g.B.[1]-l.A.[1] ]>>²));#length of P1 to P2

   h.B = l.A.clone;  # end point of 1st tangent
   # start point of tangent at same distance as end point along 'g'
   h.A = g.B.[0] - lgt * cos(angle g) , g.B.[1] - lgt * sin(angle g);

   my $a = angle h ; # angle of tangent
   # 1st control point on tangent at distance 'lgt * CURVE_F'
   @p1 = g.B.[0] + lgt * cos($a) * CURVE_F,  g.B.[1] + lgt * sin($a) * CURVE_F;

   h.A = g.B.clone; # start point of 2nd tangent
   # end point of tangent at same distance as start point along 'l'
   h.B = l.A.[0] + lgt * cos(angle l) , l.A.[1] + lgt * sin(angle l);

   $a = angle h; # angle of tangent
   # 2nd control point on tangent at distance 'lgt * CURVE_F'
   @p2 = l.A.[0] - lgt * cos($a) * CURVE_F,  l.A.[1] - lgt * sin($a) * CURVE_F;
}


given Cairo::Image.create(Cairo::FORMAT_ARGB32, WIDTH, HEIGHT) {
   given Cairo::Context.new($_) {
      my line_t ($g,$l);
      my (@p1,@p2);

      .line_width = W_LINE;
      .move_to(pt[DETACHED - 1 + cnt].[0], pt[DETACHED - 1 + cnt].[1]);

      for DETACHED..^cnt -> \j { 
         $g = line_t.new: A=>pt[(j + cnt - 2) % cnt], B=>pt[(j + cnt - 1) % cnt];
         $l = line_t.new: A=>pt[(j + cnt + 0) % cnt], B=>pt[(j + cnt + 1) % cnt];
        
         # Calculate controls points for points pt[j-1] and pt[j].
         control_points($g, $l, @p1, @p2);
         
         .curve_to(@p1[0], @p1[1], @p2[0], @p2[1], pt[j].[0], pt[j].[1]);
      }
      .stroke;
   };
   .write_png(OUTPUT) and die # C return
}

Output: (Offsite image file)

Wren[edit]

Library: DOME

In the absence of any clarification on what to use (see Talk page), the following uses a degree of 3 (i.e order k = 4) and a uniform knot vector from 1 to 16 (as there are 12 control points) with a delta of 1.

If one uses a value for k of 1, then the script will simply plot the control points as in the Julia example.

import "dome" for Window, Process
import "graphics" for Canvas, Color

class BSpline {
    construct new(width, height, cpoints, k) {
        Window.resize(width, height)
        Canvas.resize(width, height)
        Window.title = "B-spline curve"
        _p = cpoints
        _n = cpoints.count - 1
        _k = k
        _t = (1.._n + 1 + k).toList // use a uniform knot vector, delta = 1
    }

    // B-spline helper function
    w(i, k, x) { (_t[i+k] != _t[i]) ? (x - _t[i]) / (_t[i+k] - _t[i]) : 0 }
        
    // B-spline function
    b(i, k, x) {
        if (k == 1) return (_t[i] <= x  &&  x < _t[i + 1]) ? 1 : 0
        return w(i, k-1, x) * b(i, k-1, x) + (1 - w(i+1, k-1, x)) * b(i+1, k-1, x)
    }

    // B-spline points
    p() {
        var bpoints = []
        for (x in _t[_k-1]..._t[_n + 1]) {
            var sumX = 0
            var sumY = 0
            for (i in 0.._n) {
                var f = b(i, _k, x)
                sumX = sumX + f * _p[i][0]
                sumY = sumY + f * _p[i][1]
            }
            bpoints.add([sumX.round, sumY.round])
         }
         return bpoints
    }

    init() {
        if (_k > _n + 1 || _k < 1) {
            System.print("k (= %(_k)) can't be more than %(_n+1) or less than 1.")
            Process.exit()
        }
        var bpoints = p()
        // plot the curve
        for (i in 1...bpoints.count) {
            Canvas.line(bpoints[i-1][0], bpoints[i-1][1], bpoints[i][0], bpoints[i][1], Color.white)
        }
    }

    update() {}

    draw(alpha) {}
}

var cpoints = [
    [171, 171], [185, 111], [202, 109], [202, 189], [328, 160], [208, 254],
    [241, 330], [164, 252], [ 69, 278], [139, 208], [ 72, 148], [168, 172]
]
var k = 4 // polynomial degree is one less than this i.e. cubic
var Game = BSpline.new(400, 400, cpoints, k)