Thanks for that info, Christoph!

Should we be naughtily sneaking a peek at the POV-Ray source code?

https://github.com/POV-Ray/povray/blob/master/source/backend/math/polysolv.cpp

The solve_cubic function in particular?

Post a reply to this message

"Nevado" <MYPUBLICNAME(a)sol.dk> wrote:

> I'll try and see if my high-school math is enough to handle Tor's code.

I didn't even try to follow the math - I just did a quick verbatim literal
translation from the c++ based SDL to JustBasic

See if you can plug some of this in where it needs to go.
I popped some comments in there and did the best I could on short notice, but I
figured that you "get" JustBASIC enough to know what to do with functions and
subroutines.  See what I tried to do in the SplineDraw.bas file in that respect.

I haven't tested this or looked for typos, but I figured this would get you
started.

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


Start = Times[0]
End   = Times[N]
Span  = End - Start

'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)


'#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) = -hh(i)/6*zz(i+1) + 2*zz(i)) +(yy(i+1) - yy(i))/hh(i)
next i

end sub ' or end function


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 Coeffs_X (Size, 3)
' // Copy kk to Coeffs_X
for i = 0 to N-1
 Coeffs_X (i, 0) = kk (i, 0)
 Coeffs_X (i, 1) = kk (i, 1)
 Coeffs_X (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 Coeffs_Y (Size, 3)

// Copy kk to Coeffs_Y
for i = 0 to N-1
 Coeffs_Y (i, 0) = kk (i, 0)
 Coeffs_Y (i, 1) = kk (i, 1)
 Coeffs_Y (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 Coeffs_Z (Size, 3)

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



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

  (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, Coeffs_X (i, 0), Coeffs_X (i, 1),
Coeffs_X (i, 2), X (i))

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

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

'#end // macro CalculateCoordinates
end sub

Color$ = "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     'replaces T.O.K's while loop
 T = Start + J/Intervals*Span
 ComponentX = 0
 ComponentY = 0
 ComponentZ = 0

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

 for  i = 0 to N-1
  if Times(i) <= T & T < Times(i+1) then
  ' you may have to check the proper BASIC syntax for this comparison logic
  sub CalculateCoordinates T, i
  '#declare Color = Colors[I];
  end if
 next i

 if Times(n) <= T then
  sub CalculateCoordinates T, N-1
  '#declare Color = (Colors[N-1]/2 + White);
 end if

 #grid, "size 5 ; color "; WORD$( Color$, n , "," ) ;" ; set "; XandY(Points,
1); " "; XandY(Points, 2)
 ' sphere {<ComponentX, ComponentY, ComponentZ>, Radius pigment { color Color }}
next J

Post a reply to this message

Error-checked it and fixed a bunch of comments, missing parentheses, and syntax
errors.
Added the basic stuff to get it to run in the compiler, but it doesn't do
anything, and I can't exit it without forcibly killing the process in the task
manager.

But it's a small improvement.  :)

Post a reply to this message

A bigger improvement would be if I actually attached the file.  :\

Post a reply to this message

Attachments:
Download 'cubicsplines.bas.txt' (5 KB)

"Bald Eagle" <cre### [at] netscapenet> wrote:
> A bigger improvement would be if I actually attached the file.  :\

I've now modified that JustBasic file so that it makes a 2D plot of the X(T),
Y(T) and Z(T) splines.

--
Tor Olav
http://subcube.com

Post a reply to this message

Attachments:
Download 'natural_cubic_splines_04.bas.txt' (6 KB)

"Tor Olav Kristensen" <tor### [at] TOBEREMOVEDgmailcom> wrote:
> "Bald Eagle" <cre### [at] netscapenet> wrote:
> > A bigger improvement would be if I actually attached the file.  :\
>
> I've now modified that JustBasic file so that it makes a 2D plot of the X(T),
> Y(T) and Z(T) splines.

T.O.K. = "Quick study".  Nice work.

Yeah, that's cool.
I'm sure it would be a lot cooler if I understood what to do with that...

Unfortunately, when I tried to:

for p = 0 to N
    #grid, "size 5"
    #grid, "down"
    #grid, "color pink"
    #grid, "set "; XX(p); " "; YY(p)
    #grid, "up"
next p

in that space to see where the original (x,y) coordinates would be, that didn't
work out so hot.

This must be some "T"-space, or "affine"-something-or-other that has to get
mapped back into "x,y"-space, or something like that.
Yeah.  That sounds good.

I'm getting the feeling that I'm like the short lady who needed my help getting
the butter off the back of the top shelf of the fridge in the supermarket....

Do you guys have crayons and sock-puppets?

Post a reply to this message

Sorry, went my own way with this, needing to build the code myself to understand
it. But I have something now that looks promising, using parametric functions
and solving with matrices.
I'm doing a little debugging on it now for the next couple of days.

Post a reply to this message

"Bald Eagle" <cre### [at] netscapenet> wrote:
> "Tor Olav Kristensen" <tor### [at] TOBEREMOVEDgmailcom> wrote:
> > "Bald Eagle" <cre### [at] netscapenet> wrote:
> > > A bigger improvement would be if I actually attached the file.  :\
> >
> > I've now modified that JustBasic file so that it makes a 2D plot of the X(T),
> > Y(T) and Z(T) splines.
>
> T.O.K. = "Quick study".  Nice work.
>
> Yeah, that's cool.
> I'm sure it would be a lot cooler if I understood what to do with that...
>
> Unfortunately, when I tried to:
>
> for p = 0 to N
>     #grid, "size 5"
>     #grid, "down"
>     #grid, "color pink"
>     #grid, "set "; XX(p); " "; YY(p)
>     #grid, "up"
> next p
>
> in that space to see where the original (x,y) coordinates would be, that didn't
> work out so hot.
>
> This must be some "T"-space, or "affine"-something-or-other that has to get
> mapped back into "x,y"-space, or something like that.
> Yeah.  That sounds good.
......

Hehe =)
Thank you

Use Abscissa(XX(p)), Ordinate(YY(p)) instead, like this:

#grid, "set "; Abscissa(XX(p)); " "; Ordinate(YY(p))

To show Origo, use this:

#grid, "set "; Abscissa(0); " "; Ordinate(0)

http://en.wikipedia.org/wiki/Abscissa

I didn't want to confuse further by introducing more Xs and Ys.

I can explain more later about the scaling and "mapping" that these two
functions do.

--
Tor Olav
http://subcube.com

Post a reply to this message

I figured that there had to be a scaling to the screen relative to the range of
the point coordinates - I was just way too tired to "get it."

for p = 0 to N
    #grid, "size 5"
    #grid, "down"
    #grid, "color pink"
    #grid, "set "; Abscissa(XX(p), MinX, MaxX); " "; Ordinate(YY(p), MinY, MaxY)
    #grid, "up"
next p

Seems to work ok - except for the far right point, probably since the window
isn't full screen.  Multiplying ordinate and abscissa by 0.9 inside the function
scales that down a bit, and all 5 points are visible.

Now I just need to work out how to draw the right x,y curves in 2d and splice
out the parts between the endpoints, excluding the control points.  I'm guessing
studying the POV code that already graphs the curves might clue me in if I'm
alert enough to follow what's happening

Post a reply to this message

clipka <ano### [at] anonymousorg> wrote:

>      t
>
> is the time since the beginning of the segment.


Perhaps with your greater experience dealing with the structure of these things,
you could suggest how best I can "reverse engineer" a sphere_sweep to a spline.

It look like a regular spline has these t-values defined.

#declare MySpline =
spline {
  cubic_spline
    -.25, <0,0,-1>
    0.00, <1,0,0>
    0.25, <0,0,1>
    0.50, <-1,0,0>
    0.75, <0,0,-1>
    1.00, <1,0,0>
    1.25, <0,0,1>
    }

It seems that with a sphere sweep, the t-values must be constructed as part of
the spline calculation?  Can I access those t-values in some way to use my
sphere sweep data to generate a "regular" spline?
I'm looking to progress along a spline in a linear fashion - regularly spaced
intervals.

sphere_sweep {
 cubic_spline
 14,
 <-8.8750, 29.3750, 12>, 1.7500/2, // Control Point
 <-8.8750, 28.2500, 12>, 1.7500/2,
 <-8.8750, 27.1250, 12>, 1.2500/2,
 <-8.3750, 26.5000, 12>, 1.0000/2,
 <-8.7500, 25.6250, 12>, 0.7500/2,
 <-8.5313, 24.8438, 12>, 0.5625/2,
 <-8.0625, 24.1250, 12>, 0.5000/2,
 <-7.5000, 23.7500, 12>, 0.5000/2,
 <-6.7500, 23.6250, 12>, 0.5000/2,
 <-6.0000, 23.7500, 12>, 0.5000/2,
 <-5.4375, 24.0000, 12>, 0.5000/2,
 <-4.9375, 24.4375, 12>, 0.6250/2,
 <-4.4375, 25.3125, 12>, 1.1250/2,
 <-3.9375, 26.1875, 12>, 1.1250/2   //  Control Point
     tolerance 0.1
}

(What I'm trying to do at this point is "draw" between my sphere sweep and
another sphere sweep or spline to "fill-in" the space between the two.  Maybe
there's a better way that escapes me at the moment.)

Post a reply to this message