Some Bezier curve drawing

[Tourist Trap]

Hi all,

simply adapted quick and dirty from one Python (I guess) sample from that excellent page: https://pomax.github.io/bezierinfo/

Code: Select all

type XY
    as single   _x
    as single   _y
end type

'the last point of the array will be totally ignored
dim as XY points(7)
points(0)._x = 10   : points(0)._y = 10
points(1)._x = 80   : points(1)._y = 40
points(2)._x = 40   : points(2)._y = 190
points(3)._x = 190  : points(3)._y = 50
points(4)._x = 100  : points(4)._y = 170
points(5)._x = 180  : points(5)._y = 170
points(6)._x = 100  : points(6)._y = 0
points(7)._x = 0  : points(7)._y = 0

function drawCurve(points() as XY, byval t as single) as integer
  dim as XY newpoints(lBound(points) to uBound(points))
  for index as integer = lBound(points) to uBound(points)
    newpoints(index) = points(index)
  next index
  '
  if (uBound(points) - lBound(points) + 1)=1 then
    circle (points(0)._x, points(0)._y), 4
  else
    redim preserve newpoints(lBound(points) to uBound(points) - 1)
    for index as integer = lBound(newpoints) to uBound(newpoints) - 1
      var x = (1-t) * points(index)._x + t * points(index + 1)._x
      var y = (1-t) * points(index)._y + t * points(index + 1)._y
      newpoints(index) = type(x,y)
    next
    drawCurve(newpoints(), t)
  end if
end function

screenRes 400,400,32
line (points(0)._x, points(0)._y)-(points(1)._x, points(1)._y), rgb(255,0,0)
draw string (points(0)._x, points(0)._y), "0"
line (points(1)._x, points(1)._y)-(points(2)._x, points(2)._y), rgb(255,0,0)
draw string (points(1)._x, points(1)._y), "1"
line (points(2)._x, points(2)._y)-(points(3)._x, points(3)._y), rgb(255,0,0)
draw string (points(2)._x, points(2)._y), "2"
line (points(3)._x, points(3)._y)-(points(4)._x, points(4)._y), rgb(255,0,0)
draw string (points(3)._x, points(3)._y), "3"
line (points(4)._x, points(4)._y)-(points(5)._x, points(5)._y), rgb(255,0,0)
draw string (points(4)._x, points(4)._y), "4"
line (points(5)._x, points(5)._y)-(points(6)._x, points(6)._y), rgb(255,0,0)
draw string (points(5)._x, points(5)._y), "5"
line (points(6)._x, points(6)._y)-(points(7)._x, points(7)._y), rgb(255,0,0)
draw string (points(6)._x, points(6)._y), "6"

for t as single = 0 to 1 step .01
    drawCurve(points(), t)
next t

getKey()

'(eof)

[badidea]

It the code in that page python? I did not know python allows stuff like:
for(k=0; k<n; k++):

Probably JavaScript: http://paperjs.org/examples/smoothing/

[Tourist Trap]

It the code in that page python? I did not know python allows stuff like:
for(k=0; k<n; k++):

Maybe javascript then?

[counting_pine]

I don't know of any languages that feature the "factorial" postfix operator.

Code: Select all

function Bezier(n,t):
  sum = 0
  for(k=0; k<n; k++):
    sum += n!/(k!*(n-k)!) * (1-t)^(n-k) * t^(k)
  return sum

I think the above is either pseudocode or some mathematically inclined language.

[jj2007]

Why reinvent the wheel?

Code: Select all

GpStatus WINGDIPAPI GdipAddPathBeziers(
   GpPath *path, 
   GDIPCONST GpPointF *points, 
   INT count
);

[badidea]

Why reinvent the wheel?

I can think of 2 reasons:
A) It is nice to make your own wheel
B) To understand how a wheel works

[dodicat]

Thanks TT

Compare with catmull rom in 3D

Code: Select all

Type Point 
    As long x,y,z
End Type
#define Intrange(f,l) int(Rnd*(((l)+1)-(f))+(f))
#define map(a,b,x,c,d)  ((d)-(c))*((x)-(a))/((b)-(a))+(c)

function bspline(p() as point, byval t as single) as point '3D bezier
    dim as single w=1-t
    return type<point> _
    (p(2).x*w*w*w + 3*p(1).x*(w*w)*t + 3*p(4).x*(w)*(t*t) + p(3).x*(t*t*t), _
     p(2).y*w*w*w + 3*p(1).y*(w*w)*t + 3*p(4).y*(w)*(t*t) + p(3).y*(t*t*t), _
     p(2).z*w*w*w + 3*p(1).z*(w*w)*t + 3*p(4).z*(w)*(t*t) + p(3).z*(t*t*t))
end function

function cspline(p() As Point,t As Single) As Point '3D catmull
    #macro set(n)
    0.5 *(     (2 * P(2).n) +_
   (-1*P(1).n + P(3).n) * t +_
    (2*P(1).n - 5*P(2).n + 4*P(3).n - P(4).n) * t*t +_
    (-1*P(1).n + 3*P(2).n- 3*P(3).n + P(4).n) * t*t*t)
    #endmacro
    return type<point>(set(x),set(y),set(z))
end function

Sub GetSpline(v() As Point,outarray() As Point,arraysize As Integer=2000, _
    fn as function(() as point,as single) as point)
    Dim As Point p(1 To 4)
    Redim outarray(0)
    Dim As Single stepsize=(Ubound(v)-1)/(arraysize)
    If stepsize>1 Then stepsize=1
    For n As Integer=2 To Ubound(v)-2
        p(1)=v(n-1):p(2)=v(n):p(3)=v(n+1):p(4)=v(n+2)
        For t As Single=0 To 1 Step stepsize
            Redim Preserve outarray(1 To Ubound(outarray)+1)
            outarray(Ubound(outarray))=fn(p(),t)
        Next t
    Next n
End Sub  

Sub DrawCurve(a() As Point,col as ulong)
     dim as ubyte rd=Cptr(Ubyte Ptr,@col)[2] ,gr=Cptr(Ubyte Ptr,@col)[1],bl=Cptr(Ubyte Ptr,@col)[0], _
     alp=Cptr(Ubyte Ptr,@col)[3] 
    For z As Integer=Lbound(a) To Ubound(a)
        var rad=map(0,500,a(z).z,1,8)
        var f=map(lbound(a),ubound(a),z,.25,1)
        var c=rgba(rd*f,gr*f,bl*f,alp)
        circle(a(z).x,a(z).y),rad,c,,,,f
    Next z
End Sub

sub drawpoints(p() as point,col as ulong)
    for n as long=lbound(p) to ubound(p)
         var rad=map(0,500,p(n).z,1,12)
        circle (p(n).x,p(n).y),rad/2,col,,,,f
        draw string(p(n).x-4,p(n).y),str(n),rgb(0,200,0)
    next n
end sub

Sub SetUpPoints(p1() As Point,ypos As long,range As long)
    dim as integer xres,yres
    screeninfo xres,yres
    For n As long=1 To Ubound(p1) 
        Var xpos=map(lbound(p1),Ubound(p1),n,(-.2*xres),(xres+.2*xres))
        p1(n)=type<point>(xpos,Intrange((ypos-range),(ypos+range)),xpos)
    Next n
End Sub

screen 19,32,,64
color ,rgb(255,255,255)
do
    cls
redim as point p(1 to intrange(5,15))
redim as point catmul(),bezier()
setuppoints(p(),300,intrange(50,300))
getspline(p(),catmul(),,@cspline)
getspline(p(),bezier(),,@bspline)

drawcurve(bezier(),rgba(0,100,255,15))
drawcurve(catmul(),rgba(255,100,0,10))
drawpoints(p(),rgb(0,0,0))
sleep
loop until inkey=chr(27)
sleep


  

[Tourist Trap]

Why reinvent the wheel?

Code: Select all

GpStatus WINGDIPAPI GdipAddPathBeziers(
   GpPath *path, 
   GDIPCONST GpPointF *points, 
   INT count
);

Unfortunately this will be windows only. But please don't hesitate to share any additional idea. It will serve if not today then some other day for sure.

Thanks TT

Compare with catmull rom in 3D

It's really impressive. I already knew of your implementations of catmul, even if I may have missed the 3d version. If I can make any addition to my unfinishable game project, I will probably borrow this one to you, if you don't see any inconvenience :)

[jj2007]

@badidea & Tourist Trap: I like reinventing the wheel, too ;-)

Right now I am a little bit confused. The Windows function PolyBezier expects triples of points. I tried to use your data with this function, but the result looks different, see below (->source).

What is the math behind your version? Quadratic Bézier curves? Dodicat's 3D lines (beautiful indeed) use 4 points each, so I guess that is cubic Beziers, right?

Here is an article on cubic Bézier lines in Excel.

[dodicat]

Using windows and TT's points.

Code: Select all

 

#Include "windows.bi"

dim as point pts(7)={(10,10),(80,40),(40,190),(190,50),(100,170),(180,170),(100,0),(0,0)}
var k=3.5  'blow up a bit
for n as long=0 to 7
    pts(n).x=k*pts(n).x
    pts(n).y=k*pts(n).y
    next n
var Main_Win=CreateWindowEx(0,"#32770","TT's bezier points",WS_OVERLAPPEDWINDOW Or WS_VISIBLE,20,20,800,600,0,0,0,0)
                        
Dim As Any Ptr MyScreen = GetDC(main_win)

dim as msg msg
While GetMessage(@msg,Main_Win,0,0)
    TranslateMessage(@msg)
     DispatchMessage(@msg)
        Select Case msg.message
        Case 273,WM_QUIT
            End
        case else
        var p=polybezier(myscreen,@pts(0),7)
    End Select
  
Wend
  

[Tourist Trap]

What is the math behind your version? Quadratic Bézier curves? Dodicat's 3D lines (beautiful indeed) use 4 points each, so I guess that is cubic Beziers, right?

I don't know if it helps but , I took this from a paragraph named: de Casteljau's algorithm.

[jj2007]

Using windows and TT's points

Yes, that's exactly what I get. But apparently the math is a different one.

@Tourist Trap: https://en.wikipedia.org/wiki/De_Castel ... _algorithm (but math has always been my weakest point...)

[dodicat]

A quadratic curve would use a polynomial degree 2 to approximate
a*x^2 +b*x +c
The cubic would use degree 3
k*x^3 +a*x^2 +b*x +c

The splines in my code are cubic (t instead of x)
Here is an example of a quadratic curve going through three points.

Code: Select all

 

Sub GaussJordan(matrix() As Double,rhs() As Double,ans() As Double)
    Dim As Long n=Ubound(matrix,1)
    Redim ans(0):Redim ans(1 To n)
    Dim As Double b(1 To n,1 To n),r(1 To n)
    For c As Long=1 To n 'take copies
        r(c)=rhs(c)
        For d As Long=1 To n
            b(c,d)=matrix(c,d)
        Next d
    Next c
    #macro pivot(num)
    For p1 As Long  = num To n - 1
        For p2 As Long  = p1 + 1 To n  
            If Abs(b(p1,num))<Abs(b(p2,num)) Then
                Swap r(p1),r(p2)
                For g As Long=1 To n
                    Swap b(p1,g),b(p2,g)
                Next g
            End If
        Next p2
    Next p1
    #endmacro
    For k As Long=1 To n-1
        pivot(k)              'full pivoting 
        For row As Long =k To n-1
            If b(row+1,k)=0 Then Exit For
            Var f=b(k,k)/b(row+1,k)
            r(row+1)=r(row+1)*f-r(k)
            For g As Long=1 To n
                b((row+1),g)=b((row+1),g)*f-b(k,g)
            Next g
        Next row
    Next k
    'back substitute 
    For z As Long=n To 1 Step -1
        ans(z)=r(z)/b(z,z)
        For j As Long = n To z+1 Step -1
            ans(z)=ans(z)-(b(z,j)*ans(j)/b(z,z))
        Next j
    Next    z
End Sub

'Interpolate through points.
Sub Interpolate(x_values() As Double,y_values() As Double,p() As Double)
    dim as long U=Ubound(x_values),L=lbound(x_values),ctrA=0,ctrB=0
    Dim As Double matrix(1 To (U-L+1),1 To (U-L+1)),rhs(1 To (U-L+1))
    For a As Long=L To U
        ctrA+=1
        rhs(ctrA)=y_values(a)
        ctrB=0
        For b As Long=L To U
            ctrB+=1
            matrix(ctrA,ctrB)=x_values(a)^(ctrB-1)
        Next b
    Next a
    'Solve the linear equations
    GaussJordan(matrix(),rhs(),p())
End Sub

'Evaluate a polynomial at x
Function polyeval(Coefficients() As Double,Byval x As Double) As Double
    Dim As Double acc
    For i As Long=Ubound(Coefficients) To Lbound(Coefficients) Step -1
        acc=acc*x+Coefficients(i)
    Next i
    Return acc
End Function

'======================== SET UP POINTS =============== 

#macro setup()
Dim As Double x(1 To ...)={50,300, 600}
Dim As Double y(1 To ...)={250,400,100}
#endmacro
'====================================================

    setup()
    
    Redim As Double Poly(0)
    'Get the polynomial Poly()
    Interpolate(x(),y(),Poly())
    
    'print coefficients to console
    Print "Polynomial Coefficients:"
    Print
    For z As Long=1 To Ubound(Poly)
        If z=1 Then
            Print "constant term  ";Tab(20);Poly(z)
        Else
            Print Tab(8); "x^";z-1;"   ";Tab(20);Poly(z)
        End If
    Next z
    Screen 19
    
    windowtitle "esc to quit"
    'plot the interpolating polynomial
    Pset(0,polyeval(poly(),0))
    For x As Double=0 To 800 Step 800/10000
        Line -(x,polyeval(poly(),x))
    Next x
    
    'Mark the x~y points
    For z As Long=Lbound(x) To Ubound(x)
        Circle(x(z),y(z)),5,4,,,,f
        Draw String(x(z),y(z)-20),Str(x(z)) &","&str(y(z))
    Next z
  


Sleep 


 

[UEZ]

Here a short demo using GDIPlus:

Code: Select all

'Coded by UEZ
#Include "fbgfx.bi"

#Ifdef __Fb_64bit__
    #Inclib "gdiplus"
    #Include Once "win/gdiplus-c.bi"
#Else
    #Include Once "win/gdiplus.bi"
    Using gdiplus
#Endif

Using FB

Dim Shared gdipToken As ULONG_PTR
Dim Shared GDIp As GdiplusStartupInput 
GDIp.GdiplusVersion = 1
If GdiplusStartup(@gdipToken, @GDIp, NULL) <> 0 Then End

Const As Ushort iW = 400, iH = iW

Screencontrol SET_DRIVER_NAME, "GDI"
Screenres iW, iH, 32, 1, GFX_HIGH_PRIORITY Or GFX_NO_SWITCH

Windowtitle "GDI+ Bezier Demo"

'center windows by adding the taskbar to the calculation
Dim As Integer iDW, iDH
Screencontrol GET_DESKTOP_SIZE, iDW, iDH
Dim tWorkingArea As RECT
SystemParametersInfo(SPI_GETWORKAREA, null, @tWorkingArea, null)
Screencontrol SET_WINDOW_POS, (iDW - iW) \ 2, _
                              ((tWorkingArea.Bottom - iH) - (iDH - tWorkingArea.Bottom)) \ 2

'init GDI / GDI+ canvas, pens, brushes, etc. for drawing
Dim As HWND hHWND
Screencontrol(FB.GET_WINDOW_HANDLE, Cast(Integer, hHWND))

Dim As Any Ptr hDC = GetDC(hHWND), _
					 hHBitmap = CreateCompatibleBitmap(hDC, iW, iH), _
					 hDC_backbuffer = CreateCompatibleDC(hDC), _
					 hCanvas, hPen, hBrush, hPath
Var hObjOld = SelectObject(hDC_backbuffer, hHBitmap) 

GdipCreateFromHDC(hDC_backbuffer, @hCanvas)
GdipSetSmoothingMode(hCanvas, SmoothingModeAntiAlias)
GdipSetPixelOffsetMode(hCanvas, PixelOffsetModeHalf)

GdipCreatePen1(&hFFFFFFFF, 2, 2, @hPen)
GdipCreateSolidFill(&hFF404040, @hBrush)

GdipCreatePath(0, @hPath)

Dim As Single x1 = 10, y1 = 10, x2 = 150, y2 = 390, x3 = 300, y3 = 100, x4 = 390, y4 = 390, r = 10, r2 = r / 2

GdipAddPathEllipse(hPath, x1 - r2, y1 - r2, r, r)
GdipAddPathEllipse(hPath, x2 - r2, y2 - r2, r, r)
GdipAddPathEllipse(hPath, x3 - r2, y3 - r2, r, r)
GdipAddPathEllipse(hPath, x4 - r2, y4 - r2, r, r)
GdipFillPath(hCanvas, hBrush, hPath)

GdipResetPath(hPath)
GdipAddPathBezier(hPath, x1, y1, x2, y2, x3, y3, x4, y4)
GdipDrawPath(hCanvas, hPen, hPath)


Do
	BitBlt(hDC, 0, 0, iW, iH, hDC_backbuffer, 0, 0, SRCCOPY)
	Sleep(30, 1)
Loop Until Inkey = Chr(27)

'release resources
GdipDeletePath(hPath)
GdipDeletePen(hPen)
GdipDeleteBrush(hBrush)
SelectObject(hDC_backbuffer, hObjOld)
DeleteDC(hDC_backbuffer)
ReleaseDC(hHWND, hDC)
DeleteObject(hHBitmap)
GdipDeleteGraphics(hCanvas)
GdiplusShutdown(gdipToken)

[dodicat]

Here are TT's points with Windows Bezier on a gfx screen.
At last I have found a place for screensync.

If you don't like TT's points you can drag them around.
As usual, crap on 64 bits, so I choose -gen gas.
Thank you UEZ.

Code: Select all

screen 19,32

color,rgb(255,255,255)
dim  as integer xres,yres
screeninfo xres,yres
windowtitle "Drag the points with the mouse"

type point
    as long x,y
end type

Extern "windows" 
Declare Function GetDC Alias "GetDC" (Byval As Any Ptr) As Any Ptr
declare function PolyBezier(byval as any ptr, byval as POINT ptr, byval as integer) as boolean
End Extern
Declare Function timeBeginPeriod       Alias "timeBeginPeriod"(As Ulong=1) As Long
Declare Function timeEndPeriod         Alias "timeEndPeriod"  (As Ulong=1) As Long


sub drawpoints(p() as point,col as ulong)
    for n as long=lbound(p) to ubound(p)
        circle (p(n).x,p(n).y),10,col,,,,f
        draw string(p(n).x-4,p(n).y+16),str(n),rgb(0,0,0)
    next n
end sub

sub show(m as any ptr,p() as point)
    screenlock
    cls
 polybezier(m,@p(0),7)
 drawpoints(p(),rgb(200,0,0))
 screenunlock
 screensync
 sleep 1
    end sub

Dim As any ptr win
Screencontrol(2,Cast(integer,win))
Dim As any Ptr myscreen = GetDC(win)

'TT's points
dim as point pts(7)={(10,10),(80,40),(40,190),(190,50),(100,170),(180,170),(100,0),(0,0)}
var k=3.0  'blow up a bit
for n as long=0 to 7
    pts(n).x=k*pts(n).x
    pts(n).y=k*pts(n).y
next n

#define onscreen (mx>10) and (mx<xres-10) and (my>10) and (my<yres-10)
#define incircle(cx,cy,radius,x,y) (cx-x)*(cx-x) +(cy-y)*(cy-y)<= radius*radius

#macro mouse(m)
Dim As Long x=mx,y=my,dx,dy
While mb = 1
    show(myscreen,pts()):sleep 1,1
    Getmouse mx,my,,mb
    If onscreen Then
        If mx<>x Or my<>y  Then
            dx = mx - x
            dy = my - y
            x = mx
            y = my
            pts(m).x=x+dx
            pts(m).y=y+dy
        End If
    End If
Wend
#endmacro

dim as integer mx,my,mb
timebeginperiod
do
    getmouse mx,my,,mb
    show(myscreen,pts())
    for n as long=lbound(pts) to ubound(pts)
        if incircle(pts(n).x,pts(n).y,10,mx,my) and mb=1 then
        mouse(n)
        end if
    next
 loop until len(inkey)
 timeendperiod
 sleep