nomainwin
open "Cubic Spline Solver POV-Ray : Tor Olav Kristensen, Christian Nevado, Bald Eagle" for graphics_nf_nsb as #grid
#grid, "trapclose [quit]"

Size = 5
N = Size - 1


'This is the weird way I came up to populate an array from a list
dim Times (Size)
Times$ = "-2, -1, 0, 0.5, 1"
for i = 0 to N
    Times (i) = val (WORD$ (Times$, i , "," ))
next i


dim Points (Size, 3)
Points$ = "-2, 0, 40, 0, 0, -2, 1, 0, 2, 0, 6, 2, 8, 2, 2"
for i = 0 to N
    Points (i,1) = val (WORD$ (Points$, 1+i*3 , "," ))
    Points (i,2) = val (WORD$ (Points$, 2+i*3 , "," ))
    Points (i,3) = val (WORD$ (Points$, 3+i*3 , "," ))
next i


StartTime = Times (0)
EndTime   = Times (N)
Span  = EndTime - StartTime

'T.O.K.'s way, using global variables

dim X (Size)
dim Y (Size)
dim Z (Size)


for i = 0 to N
X (i) = Points (i, 1)
Y (i) = Points (i, 2)
Z (i) = Points (i, 3)
next i

dim tt (Size)
dim yy (Size)

dim hh (Size)
dim uu (Size)
dim vv (Size)
dim zz (Size)

dim kk (Size, 3)





for i = 0 to N
    tt (i) = Times (i)
    yy (i) = X (i)
next i


'// Calculate coefficients for the X polynomials
call CalculateCoefficients hh, tt, kk, yy, uu, vv, zz

dim CoeffsX (Size, 3)
' // Copy kk to CoeffsX
for i = 0 to N-1
    CoeffsX (i, 0) = kk (i, 0)
    CoeffsX (i, 1) = kk (i, 1)
    CoeffsX (i, 2) = kk (i, 2)
next i


'// Calculate coefficients for the Y polynomials

'// Copy Y to yy
for i = 0 to N
    yy (i) = Y (i)
next i

call CalculateCoefficients hh, tt, kk, yy, uu, vv, zz


dim CoeffsY (Size, 3)

'// Copy kk to CoeffsY
for i = 0 to N-1
    CoeffsY (i, 0) = kk (i, 0)
    CoeffsY (i, 1) = kk (i, 1)
    CoeffsY (i, 2) = kk (i, 2)
next i



'// Calculate coefficients for the Z polynomials
'// Copy Z to yy
for i = 0 to N
    yy (i) = Z (i)
next i

call CalculateCoefficients hh, tt, kk, yy, uu, vv, zz


dim CoeffsZ (Size, 3)

'// Copy kk to CoeffsZ
for i = 0 to N-1
    CoeffsZ (i, 0) = kk (i, 0)
    CoeffsZ (i, 1) = kk (i, 1)
    CoeffsZ (i, 2) = kk (i, 2)
next i








Colors$ = "cyan, green, orange, yellow, blue"
' I haven't yet gotten a firm grasp of BASIC graphics - Nevado appears to be really familiar
' do a lookup using using the   WORD$( Color$, n , "," )   trick
' or convert this to an rgb value, etc.
'#declare Colors =
'  array[Size] {
'    color Cyan,
'    color Green,
'    color Orange,
'    color Yellow,
'    color Blue // Not used
'  }

Intervals = 400

Radius = 0.08

for J = 0 to Intervals - 1
    T = StartTime + J/Intervals*Span
    ComponentX = 0
    ComponentY = 0
    ComponentZ = 0

    if T < Times (0) then
        var = 0
        call CalculateCoordinates T, var
        '#declare Color = (Colors[0]/2 + White);
         Color$ = "black"
    end if

    for  i = 0 to N-1
        if Time (i) <= T and T < Times (i+1) then
        ' you may have to check the proper BASIC syntax for this comparison logic
        call CalculateCoordinates T, i
        Color$ = WORD$( Colors$, i , "," )
        end if
    next i

    if Times (n) <= T then
        var = N-1
        call CalculateCoordinates T, var
        Color$ = "black"
    end if

    #grid, "size 5 ; color "; Color$ ;" ; set "; ComponentX; " "; ComponentY
    ' sphere {<ComponentX, ComponentY, ComponentZ>, Radius pigment { color Color }}
next J

print #grid, "setfocus"
scan
#grid, "when characterInput [quit]"

[loop]
goto [loop]

[quit]
close #grid
stop

'#macro
function CubicPolynomial (T, A, B, C, D)

Coefficient = (D + T*(C + T*(B + T*A))) '// Equals (A*T^3 + B*T^2 + C*T + D)

'#end // macro CubicPolynomial
end function


'#macro
sub CalculateCoordinates T, I

    dT = T - Times (i)
        ComponentX = CubicPolynomial (dT, CoeffsX (i, 0), CoeffsX (i, 1), CoeffsX (i, 2), X (i))

    ComponentY = CubicPolynomial (dT, CoeffsY (i, 0), CoeffsY (i, 1), CoeffsY (i, 2), Y (i))

    ComponentZ = CubicPolynomial (dT, CoeffsZ (i, 0), CoeffsZ (i, 1), CoeffsZ (i, 2), Z (i))

'#end // macro CalculateCoordinates
end sub


'#macro CalculateCoefficients()
' Change this to a subroutine or a function and place at end of code
sub CalculateCoefficients hh, tt, kk, yy, uu, vv, zz

for i = 0 to N-1
    hh (i) = tt (i+1) - tt (i)
    kk (i, 1) = 6 * (yy (i+1) - yy (i)) / hh (i)
next i

for i = 2 to N-1
    uu (1) = 2 * (hh (0) + hh (1))
    vv (1) = kk (1, 1) - kk (0, 1) - hh (i-1) * vv (i-1) / uu (i-1)
next i

zz (n) = 0
for i = n-1 to 1 step -1
    zz (i) = (vv (i) - hh (i) * zz (i+1)) / uu (i)
next i
zz (0) = 0

for i = 0 to N-1
    kk (i, 0) = (zz (i+1) - zz (i)) / (6 * hh (i))
    kk (i, 1) = zz (i)/2
    kk (i, 2) = 0-hh (i) / 6 * (zz (i+1) + 2 * zz (i)) + (yy (i+1) - yy (i)) / hh (i)
next i

end sub ' or end function