3D Geometry , basics

Post your FreeBASIC source, examples, tips and tricks here. Please don’t post code without including an explanation.
Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

After a number of years exploring other languages and numerous OS's I've returned
to update this code ; in particular I've included a polygon filling routine . This routine
was translated from C by another FreeBASIC user .
From what I'm able to determine thus far , the polygon filler is quite accurate , fast
and efficient.

Quite a few years ago I wrote similar code for openGL , this however required quite
a few external libraries ; whereas this doesn't and may even run in DOS.

Obviously there's a fair bit more to do , for now though this might suffice .

In the routine Persp , the variable d and the number 2 might be adjusted to tweek the
appearance of the graph .

Code: Select all

'
' -----------------------------------------------------------------------------
'
'   My graf 3d
'
'   MyGraf3_3d.bas
'
'
'    (c) copyright 2019 , sciwiseg@gmail.com ,
'
'             Edward.Q.Montague.  [ alias]
'
'
'
'
'
' -----------------------------------------------------------------------------
'
type point
         x as single
         y as single
         z as single
         u as single '  possible extension for special coord system
end type
'
const Pi=4*atn(1)

' replaces the defines above (single line Macro's in FB)
Const As Long   POINTS = 4 , POLYGONS = 20, SCR_W = 740, SCR_H = 680

'
dim as single x1,y1,z1,x2,y2,z2
dim as integer i,j,k
'
' ----------------------------------------------------------------------------
'
'
'                  Looking at a cube .
'
'            cc -1,1 _______<_______  1,1    start         z = -1
'                   |               |             back face.
'                   |               |
'                   v               ^
'                   |               |
'                   |_______________|
'                -1,-1        >        1,-1
'               
'
' -----------------------------------------------------------------------------
'
declare function rotx(q as point,angx as single) as point
declare function roty(q as point,angy as single) as point
declare function rotz(q as point,angz as single) as point
declare function tranx(q as point,movx as single) as point
declare function trany(q as point,movy as single) as point
declare function tranz(q as point,movz as single) as point
declare function persp(q as point,d as single) as point
'
declare function Trall( p1() as point,n as integer,edge() as integer, div as integer ,np as integer,ne as integer) as integer

declare sub gendata1(a() as long,n as integer,m as integer)
declare sub gendata2(a() as long,n as integer,m as integer)


declare sub genmatrix(gm() as single,n as integer,m as integer)
declare sub gmdata(gm() as single,n as integer,m as integer)
declare function f1(x as single , y as single)as single
declare sub trallg(gm() as single,dm() as single , n as integer,m as integer)
declare sub trand(dm() as single,a() as long,n as integer,m as integer)
declare Sub fill_polygon(a() As Long, ByVal c As ULong)
declare sub outline_polygon(a() As Long, ByVal c As ULong)

declare sub tqxyz(x as single , y as single ,n as integer , m as integer , byref q as point)
declare sub w2scrn(p as point,n as integer , m as integer,byref u as long,byref v as long)

'
' ============================================================================
'
dim as integer np , ne,n,m
'
restore  storeA
read np
'
dim as point p1(1 to np)
'
restore store1
for i =1 to np
   read p1(i).x
   read p1(i).y
   read p1(i).z
next i
'
'
restore  storeB
read ne
'
dim as integer edge(1 to ne,0 to 1)
'
restore store2
for i = 1 to ne
   read edge(i,0)
   read edge(i,1)
next i
'
' -----------------------------------------------------------------------------
'
'screen 12
'window (-1.5,-1.5)-(1.5,1.5)
'line (-1.4,-1.4)-(1.4,1.4),11,b
'
'cls
'k=Trall( p1() ,3,edge() , 32 ,np ,ne )
'sleep
'INITIALIZING GRAPHICS _________________________________________________
ScreenRes(SCR_W, SCR_H, 24)      'initialize graphics
'window(10,10)-(210,110)
Cls

n=520
m=520
'
dim a(0 to POINTS,0 to 1) as long
dim gm(0 to n , 0 to m,0 to 1) as single
dim dm(0 to n , 0 to m,0 to 1) as single
'screen 12
'window(0,0)-(SCR_W,SCR_H)


gendata2(a() ,n ,m )


sleep
end

genmatrix(gm() ,n ,m )
gmdata(gm() ,n ,m )
trallg(gm(),dm()  , n ,m )
trand(dm() ,a() ,n ,m )


sleep
end
'
' ===================================
'
'    number of vertices
'
storeA:
data 8 
'
'    vertex data , easier to keep track of
'  data when we use multiple data statements.
'
store1:
data  1,1,1
data  -1,1,1
data  -1,-1,1
data 1,-1,1
data  1,1,-1
data  -1,1,-1
data  -1,-1,-1
data  1,-1,-1
'
'   Number of edges.
'
storeB:
data 12
'
'  edge data
'
store2:
data 1,2
data 1,4
data 1,5
data 2,3
data 2,6
data 3,4
data 3,7
data 4,8
data 5,6
data 5,8
data 6,7
data 7,8
'
' -------------------------------------------------------------------------------
'
function rotx(q as point,angx as single) as point
'
'              Rotate around x axis .
'
static as point p
'
             p.x = q.x
             p.y= q.y*cos(angx)-sin(angx)*q.z
             p.z= q.z*cos(angx)+sin(angx)*q.y
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function roty(q as point,angy as single) as point
'
'              Rotate around y axis .
'
static as point p
'
            p.x = sin(angy)*q.z + cos(angy)*q.x
            p.y = q.y
            p.z = cos(angy)*q.z -sin(angy)*q.x
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function rotz(q as point,angz as single) as point
'
'              Rotate around z axis .
'
static as point p
'
            p.x = sin(angz)*q.y + cos(angz)*q.x
            p.y = cos(angz)*q.y-sin(angz)*q.x
            p.z = q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function tranx(q as point,movx as single) as point
'
'             Translate point along x axis
'
static as point p
'
              p.x=q.x + movx
              p.y=q.y
              p.z=q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function trany(q as point,movy as single) as point
'
'              Translate point along y axis .
'
static as point p
'
              p.x=q.x
              p.y=q.y + movy
              p.z=q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function tranz(q as point,movz as single) as point
'
'              Translate point along z axis .
'
static as point p
'
              p.x=q.x
              p.y=q.y
              p.z=q.z + movz
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function persp(q as point,d as single) as point
'
'              3d  perspective . 
'
'    Add 2.5 to the numerator when using any negative z value.
'
'     In this instance    -1 <= z  <= 1  ,  unit cube .
'
'    Therefore 2.5 is appropriate .
'
static as point p
'
     p.x = d*q.x/(q.z+2)
     p.y = d*q.y/(q.z+2)
     p.z = d
'
     return p
'
end function
'
' -----------------------------------------------------------------------------
'
function Trall( p1() as point,n as integer,edge() as integer, div as integer ,np as integer,ne as integer) as integer
'
'              Translate and rotate all vertices .
'     as an animation ,  for  n  cycles .
'
'    np number of points .
'    ne number of integers .
'
'
'   With  div number of angle divisions .
'
static as point p2(1 to np)
static as single theta,thi,x1,y1,z1,x2,y2,z2
static as integer i,j,k
static as integer i1,j1,k1
'
theta = Pi/div
'
for i=1 to n
  for j = 0 to div
  'cls
       thi = j*theta
   for k = 1 to np
     p2(k) = roty(p1(k),thi)
     p2(k)=persp(p2(k),0.8)
   next k   
 '  
    cls 
'
for i1 = 1 to ne
      j1 = edge(i1,0)
     k1 = edge(i1,1)
   x1 = p2(j1).x
   y1 = p2(j1).y
'   z1 = p2(j1).z   
   x2 = p2(k1).x
   y2 = p2(k1).y
'   z2 = p2(k1).z   
line(x1,y1)-(x2,y2),14
next i1     
'
sleep 100
  next j
next i
'
     return 0
'
end function
'
' --------------------------------------------------------------------------------
'
'
' ----------------------------------------------------------------
'
sub gendata1(a() as long,n as integer,m as integer)
'
'  Generate 3d data grid .
'
'  For perspective the range for x and y and z , is [-1,1]
'
'
'  Then a translation back to screen coordinates is required .  [-1,1] => [10,n] & [-1,1] => [10,m]  ?
'    x = (x + 1)*0.5*(n - 10) + 10
'    y = (y + 1)*0.5*(m - 10) + 10
'
'
'   w = f(x,y)
'   q.x = x
'   q.y = w
'   q.z = z
'
'   p = persp(q,d)
'
'   ( x', y',z , w') = persp( x , y , z , w )
'
'
static as integer i1,j1,k1,l1,qg,pg
static as single x,y,z,d
static as ulong colour,u,v
static as point p , q

d=0.8

'qg=10
'pg=10

pg=n/20
qg=m/10

colour=rgb(120,200,200)
for j1=10 to m step qg
    'z = 1 - 2*(j1-10)/(m-10) '  1  ->  -1 
    
  for i1=10 to n step pg
  
      k1=0
      x=i1
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
     ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
            
      k1=k1+1
      x=i1+pg
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
    ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1
      x=i1+pg
      y=j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )  
    ' f1(q.x,q.y) here ?
      p = persp(q ,d )   
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1     
      x = i1
      y = j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )   
    ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
      fill_polygon(a(), CULng(rnd*&hFFFFFF))
      outline_polygon(a() , colour)
 
   next i1
  '
next j1
'
end sub
'
' ----------------------------------------------------------------
'
sub gendata2(a() as long,n as integer,m as integer)
'
'  Generate 3d data grid .
'
'  For perspective the range for x and y and z , is [-1,1]
'
'
'  Then a translation back to screen coordinates is required .  [-1,1] => [10,n] & [-1,1] => [10,m]  ?
'    x = (x + 1)*0.5*(n - 10) + 10
'    y = (y + 1)*0.5*(m - 10) + 10
'
'
'   w = f(x,y)
'   q.x = x
'   q.y = w
'   q.z = z
'
'   p = persp(q,d)
'
'   ( x', y',z , w') = persp( x , y , z , w )
'
'
static as integer i1,j1,k1,l1,qg,pg
static as single x,y,z,d , theta
static as ulong colour,u,v,chrome
static as point p , q , s


theta = Pi/4

d=0.98
d=1.2

'qg=10
'pg=10

pg=n/50
qg=m/50

colour=rgb(120,200,200)
for j1=10 to m step qg
    'z = 1 - 2*(j1-10)/(m-10) '  1  ->  -1 
    
  for i1=10 to n step pg
  
      k1=0
      x=i1
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y= f1(q.x,q.y) 
      chrome = (q.y + 1)*32
      s=rotx(q ,theta )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
            
      k1=k1+1
      x=i1+pg
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y = f1(q.x,q.y) 
      s=rotx(q ,theta )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1
      x=i1+pg
      y=j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y = f1(q.x,q.y) 
      s=rotx(q ,theta )
      p = persp(s ,d )   
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1     
      x = i1
      y = j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )   
      q.y = f1(q.x,q.y) 
      s=rotx(q ,theta )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
     ' fill_polygon(a(), CULng(rnd*&hFFFFFF))
      
      fill_polygon(a(), rgb(0,0,0))
     outline_polygon(a() , colour)
 
   next i1
  '
next j1
'
end sub
'
' ----------------------------------------------------------------
'
sub genmatrix(gm() as single,n as integer,m as integer)
'
'   Generate a matrix with limits [x,-1,1][y,-1,1]
'
'
static as integer i,j
static as single x,y,z
'
for j=0 to m
    y=1-2*i/m
    z=y
    for i=0 to n
      x=-1+2*i/n
      gm(i,j,0)=x
      gm(i,j,1)=z
    next i
next j      
'
'
end sub
'
' -------------------------------------------------------------------
'
sub gmdata(gm() as single,n as integer,m as integer)
'
'  Generate data from function f1(x,y) .
'
'
static as integer i,j
static as single x,y,z
'
for j=0 to m
    for i=0 to n
      x=gm(i,j,0)
      y=gm(i,j,1)
      z=f1(x,y)
      gm(i,j,0) = z
    next i
next j   
'
end sub
'
' ---------------------------------------------------------------------
'
sub trallg(gm() as single,dm() as single , n as integer,m as integer)
'
'   translate , rotate , apply perspective to all of gm()
'
'
static as integer i,j
static as single x,y,z,d
static as point p,q
'
d=0.8
'
for j=0 to m
     q.z=1-2*j/m
    for i=0 to n
      q.x=-1 +2*i/n
      q.y=gm(i,j,0)
      p = persp(q ,d )       
     dm(i,j,0) = p.x
     dm(i,j,1) = p.y     
    next i
next j   
'
'
end sub
'
' -------------------------------------------------------------------
'
sub trand(dm() as single,a() as long,n as integer,m as integer)
'
'  translate dm() to a() , also translate to screen coordinates for fill_polygon() routine .
'
'
static as integer i1,j1,k1,l1,qg,pg
static as long u , v
static as single x,y,z
static as point p , q
static as ulong colour

qg=10
pg=10



colour=rgb(20,120,20)
for j1=10 to m step qg
    
  for i1=10 to n step pg
  
      k1=0
      p.x=dm(i1,j1,0)
      p.y=dm(i1,j1,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
      k1=k1+1
      p.x=dm(i1+pg,j1,0)
      p.y=dm(i1,j1,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
      
      k1=k1+1  
      p.x=dm(i1+pg,j1,0)
      p.y=dm(i1,j1+qg,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
       
      k1=k1+1
      p.x=dm(i1,j1,0)
      p.y=dm(i1,j1+qg,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

   fill_polygon(a(), CULng(rnd*&hFFFFFF))
   'outline_polygon(a() , colour)
 
   next i1
  '
next j1
'
'
end sub
'
' -------------------------------------------------------------------
'
function f1(x as single , y as single)as single
'
'       function to generate values upon [x,-1,1][y,-1,1]
'
'       [ ,-1,1] 
'
static as single r , z

        r=x*x+y*y
        if ( r > 0 ) then
             r = 5*sqr(r)
             z=-sin(r*Pi)/(r*Pi)
        else
             z=-1
        end if          

       return z

'
end function
'
' -----------------------------------------------------------------------------------------------------
'
Sub fill_polygon(a() As Long, ByVal c As ULong)
   'translation of a c snippet by Angad
   'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
   Dim As Long      i, j, k, dy, dx, x, y, temp
   Dim As Long      xi(0 to Ubound(a, 1))
   Dim As Single    slope(0 to Ubound(a, 1))
   
   'join first and last vertex
   a(Ubound(a, 1), 0) = a(0, 0)
   a(Ubound(a, 1), 1) = a(0, 1)

   For i = 0 To Ubound(a, 1) - 1

      dy = a(i+1, 1) - a(i, 1)
      dx = a(i+1, 0) - a(i, 0)

      If (dy = 0) Then slope(i) = 1.0
      If (dx = 0) Then slope(i) = 0.0

      If (dy <> 0) AndAlso (dx <> 0) Then slope(i) = dx / dy
   Next i

   For y = 0 to SCR_H - 1
      k = 0
      ' using FB's short-cut operators (which C doesn't have!)
      For i = 0 to Ubound(a, 1) - 1
         If (a(i, 1) <= y AndAlso a(i+1, 1) > y) OrElse _
             (a(i, 1) > y AndAlso a(i+1, 1) <= y) Then
            xi(k) = CLng(a(i, 0) + slope(i) * (y - a(i, 1)))
            k += 1
         End If
      Next i

      For j = 0 to k - 2
         'Arrange x-intersections in order
         For i = 0 To k - 2
            If (xi(i) > xi(i + 1)) Then
               temp = xi(i)
               xi(i) = xi(i + 1)
               xi(i + 1) = temp
            End If
         Next i
      Next j
      'line filling
      For i = 0 To k - 2 Step 2
         Line (xi(i), y)-(xi(i + 1) + 1, y), c
      Next i
   Next y
End Sub
'
' -----------------------------------------------------------------------------
' 
sub outline_polygon(a() As Long, ByVal c As ULong)
'
'  Draw an outtline for the polygon , in color c  .
'
  'translation of a c snippet by Angad
   'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
   Dim As Long      i, j,  x, y, u , v , temp
      
   'join first and last vertex
   a(Ubound(a, 1), 0) = a(0, 0)
   a(Ubound(a, 1), 1) = a(0, 1)

   For i = 0 To Ubound(a, 1) - 1
       x=a(i,0)
       y=a(i,1)
       u=a(i+1,0)
       v=a(i+1,1)
       line(x,y)-(u,v),c
   Next i


end sub
'
' ----------------------------------------------------------------------
'
sub w2scrn( p as point,n as integer , m as integer, byref u as long, byref v as long)
'
' input x,y,n,m
' output u , v
'
     p.x =  (p.x + 1)*0.5*(n - 10) + 10  ' [-1,1] -> [10,n]
     p.y =  (p.y + 1)*0.5*(m - 10) + 10  ' [-1,1] -> [10,m]
     u = clng(p.x)
     v = clng(p.y)

'
end sub
'
' -------------------------------------------------------------
'
sub tqxyz(x as single , y as single ,n as integer , m as integer , byref q as point)
'
'  translate  to  x ~ i1 , y ~ j1 to q.x , q.y , q.z
'

      q.z = 1 - 2*(y-10)/(m-10)
      q.x = -1 + 2*(x-10)/(n-10)
      q.y = -1 + 2*(y-10)/(m-10)

'
end sub







Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

Minor update , partial rotation of the complete 3d graph and a colour gradient as background.
May require different routine for perspective and an alternative to this painters algorithm .

Fascinating alternate approaches at games 3d FreeBASIC location.

Code: Select all


'
' -----------------------------------------------------------------------------
'
'   My graf 3d
'
'   MyGraf4_3d.bas
'
'
'    (c) copyright 2019 , sciwiseg@gmail.com ,
'
'             Edward.Q.Montague.  [ alias]
'
'
'
'
'
' -----------------------------------------------------------------------------
'
type point
         x as single
         y as single
         z as single
         u as single '  possible extension for special coord system
end type
'
const Pi=4*atn(1)

' replaces the defines above (single line Macro's in FB)
Const As Long   POINTS = 4 , POLYGONS = 20, SCR_W = 740, SCR_H = 680

'
dim as single x1,y1,z1,x2,y2,z2
dim as integer i,j,k
'
' ----------------------------------------------------------------------------
'
'
'                  Looking at a cube .
'
'            cc -1,1 _______<_______  1,1    start         z = -1
'                   |               |             back face.
'                   |               |
'                   v               ^
'                   |               |
'                   |_______________|
'                -1,-1        >        1,-1
'               
'
' -----------------------------------------------------------------------------
'
declare function rotx(q as point,angx as single) as point
declare function roty1(q as point,angy as single) as point
declare function roty(q as point,angy as single) as point


declare function rotz(q as point,angz as single) as point
declare function tranx(q as point,movx as single) as point
declare function trany(q as point,movy as single) as point
declare function tranz(q as point,movz as single) as point
declare function persp(q as point,d as single) as point
'
declare function Trall( p1() as point,n as integer,edge() as integer, div as integer ,np as integer,ne as integer) as integer

declare sub gendata1(a() as long,n as integer,m as integer)
declare sub gendata2(a() as long,n as integer,m as integer)


declare sub genmatrix(gm() as single,n as integer,m as integer)
declare sub gmdata(gm() as single,n as integer,m as integer)
declare function f1(x as single , y as single)as single
declare sub trallg(gm() as single,dm() as single , n as integer,m as integer)
declare sub trand(dm() as single,a() as long,n as integer,m as integer)
declare Sub fill_polygon(a() As Long, ByVal c As ULong)
declare sub outline_polygon(a() As Long, ByVal c As ULong)

declare sub tqxyz(x as single , y as single ,n as integer , m as integer , byref q as point)
declare sub w2scrn(p as point,n as integer , m as integer,byref u as long,byref v as long)

'
' ============================================================================
'
dim as integer np , ne,n,m
'
restore  storeA
read np
'
dim as point p1(1 to np)
'
restore store1
for i =1 to np
   read p1(i).x
   read p1(i).y
   read p1(i).z
next i
'
'
restore  storeB
read ne
'
dim as integer edge(1 to ne,0 to 1)
'
restore store2
for i = 1 to ne
   read edge(i,0)
   read edge(i,1)
next i
'
' -----------------------------------------------------------------------------
'
'screen 12
'window (-1.5,-1.5)-(1.5,1.5)
'line (-1.4,-1.4)-(1.4,1.4),11,b
'
'cls
'k=Trall( p1() ,3,edge() , 32 ,np ,ne )
'sleep
'INITIALIZING GRAPHICS _________________________________________________
ScreenRes(SCR_W, SCR_H, 24,2)      'initialize graphics
'window(10,10)-(210,110)
ScreenSet 1,0
Cls

n=520
m=520
'
dim a(0 to POINTS,0 to 1) as long
dim gm(0 to n , 0 to m,0 to 1) as single
dim dm(0 to n , 0 to m,0 to 1) as single
dim as ulong colour
'
gendata2(a() ,n ,m )
'
sleep
end

'genmatrix(gm() ,n ,m )
'gmdata(gm() ,n ,m )
'trallg(gm(),dm()  , n ,m )
'trand(dm() ,a() ,n ,m )


'sleep
'end
'
' ===================================
'
'    number of vertices
'
storeA:
data 8 
'
'    vertex data , easier to keep track of
'  data when we use multiple data statements.
'
store1:
data  1,1,1
data  -1,1,1
data  -1,-1,1
data 1,-1,1
data  1,1,-1
data  -1,1,-1
data  -1,-1,-1
data  1,-1,-1
'
'   Number of edges.
'
storeB:
data 12
'
'  edge data
'
store2:
data 1,2
data 1,4
data 1,5
data 2,3
data 2,6
data 3,4
data 3,7
data 4,8
data 5,6
data 5,8
data 6,7
data 7,8
'
' -------------------------------------------------------------------------------
'
function rotx(q as point,angx as single) as point
'
'              Rotate around x axis .
'
static as point p
'
             p.x = q.x
             p.y= q.y*cos(angx)-sin(angx)*q.z
             p.z= q.z*cos(angx)+sin(angx)*q.y
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function roty1(q as point,angy as single) as point
'
'              Rotate around y axis .
'
static as point p
'
            p.x = sin(angy)*q.z + cos(angy)*q.x
            p.y = q.y
            p.z = cos(angy)*q.z -sin(angy)*q.x
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function roty(q as point,angy as single) as point
'
'              Rotate around y axis .
'
static as point p
'
            p.x = sin(angy)*q.z + cos(angy)*q.x
            p.y = q.y
            p.z = cos(angy)*q.z -sin(angy)*q.x
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function rotz(q as point,angz as single) as point
'
'              Rotate around z axis .
'
static as point p
'
            p.x = sin(angz)*q.y + cos(angz)*q.x
            p.y = cos(angz)*q.y-sin(angz)*q.x
            p.z = q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function tranx(q as point,movx as single) as point
'
'             Translate point along x axis
'
static as point p
'
              p.x=q.x + movx
              p.y=q.y
              p.z=q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function trany(q as point,movy as single) as point
'
'              Translate point along y axis .
'
static as point p
'
              p.x=q.x
              p.y=q.y + movy
              p.z=q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function tranz(q as point,movz as single) as point
'
'              Translate point along z axis .
'
static as point p
'
              p.x=q.x
              p.y=q.y
              p.z=q.z + movz
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function persp(q as point,d as single) as point
'
'              3d  perspective . 
'
'    Add 2.5 to the numerator when using any negative z value.
'
'     In this instance    -1 <= z  <= 1  ,  unit cube .
'
'    Therefore 2.5 is appropriate .
'
static as point p
'
     p.x = d*q.x/(q.z+2.25)
     p.y = d*q.y/(q.z+2.25)
     p.z = d
'
     return p
'
end function
'
' -----------------------------------------------------------------------------
'
function Trall( p1() as point,n as integer,edge() as integer, div as integer ,np as integer,ne as integer) as integer
'
'              Translate and rotate all vertices .
'     as an animation ,  for  n  cycles .
'
'    np number of points .
'    ne number of integers .
'
'
'   With  div number of angle divisions .
'
static as point p2(1 to np)
static as single theta,thi,x1,y1,z1,x2,y2,z2
static as integer i,j,k
static as integer i1,j1,k1
'
theta = Pi/div
'
for i=1 to n
  for j = 0 to div
  'cls
       thi = j*theta
   for k = 1 to np
     p2(k) = roty(p1(k),thi)
     p2(k)=persp(p2(k),0.8)
   next k   
 '  
    cls 
'
for i1 = 1 to ne
      j1 = edge(i1,0)
     k1 = edge(i1,1)
   x1 = p2(j1).x
   y1 = p2(j1).y
'   z1 = p2(j1).z   
   x2 = p2(k1).x
   y2 = p2(k1).y
'   z2 = p2(k1).z   
line(x1,y1)-(x2,y2),14
next i1     
'
sleep 100
  next j
next i
'
     return 0
'
end function
'
' --------------------------------------------------------------------------------
'
'
' ----------------------------------------------------------------
'
sub gendata1(a() as long,n as integer,m as integer)
'
'  Generate 3d data grid .
'
'  For perspective the range for x and y and z , is [-1,1]
'
'
'  Then a translation back to screen coordinates is required .  [-1,1] => [10,n] & [-1,1] => [10,m]  ?
'    x = (x + 1)*0.5*(n - 10) + 10
'    y = (y + 1)*0.5*(m - 10) + 10
'
'
'   w = f(x,y)
'   q.x = x
'   q.y = w
'   q.z = z
'
'   p = persp(q,d)
'
'   ( x', y',z , w') = persp( x , y , z , w )
'
'
static as integer i1,j1,k1,l1,qg,pg
static as single x,y,z,d
static as ulong colour,u,v
static as point p , q

d=0.8

'qg=10
'pg=10

pg=n/20
qg=m/10

colour=rgb(120,200,200)
for j1=10 to m step qg
    'z = 1 - 2*(j1-10)/(m-10) '  1  ->  -1 
    
  for i1=10 to n step pg
  
      k1=0
      x=i1
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
     ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
            
      k1=k1+1
      x=i1+pg
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
    ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1
      x=i1+pg
      y=j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )  
    ' f1(q.x,q.y) here ?
      p = persp(q ,d )   
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1     
      x = i1
      y = j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )   
    ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
      fill_polygon(a(), CULng(rnd*&hFFFFFF))
      outline_polygon(a() , colour)
 
   next i1
  '
next j1
'
end sub
'
' ----------------------------------------------------------------
'
sub gendata2(a() as long,n as integer,m as integer)
'
'  Generate 3d data grid .
'
'  For perspective the range for x and y and z , is [-1,1]
'
'
'  Then a translation back to screen coordinates is required .  [-1,1] => [10,n] & [-1,1] => [10,m]  ?
'    x = (x + 1)*0.5*(n - 10) + 10
'    y = (y + 1)*0.5*(m - 10) + 10
'
'
'   w = f(x,y)
'   q.x = x
'   q.y = w
'   q.z = z
'
'   p = persp(q,d)
'
'   ( x', y',z , w') = persp( x , y , z , w )
'
'
static as integer i1,j1,k1,l1,qg,pg
static as single x,y,z,d , theta , thi
static as ulong colour,u,v,chrome
static as point p , q , s , t
static as integer i,j,k

theta = Pi/10
'theta=0
thi=Pi/6
'thi=0
d=0.98
d=1.2

'qg=10
'pg=10

pg=n/50
qg=m/50

colour=rgb(120,200,200)
for  theta=0 to Pi/2.5 step Pi/120 
'
cls
for j=0 to SCR_H
    k=int(200*j/SCR_H)
    colour=rgb(12,20,k)

   line(0,j)-(SCR_W,j),colour
next j  
'
colour=rgb(120,200,200)
'
for j1=10 to m step qg
    'z = 1 - 2*(j1-10)/(m-10) '  1  ->  -1 
    
  for i1=10 to n step pg
  
      k1=0
      x=i1
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y= f1(q.x,q.y) 
      chrome = (q.y + 1)*32
      
      t=roty(q ,theta )
      s=rotx(t ,thi )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
            
      k1=k1+1
      x=i1+pg
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y = f1(q.x,q.y) 
      t=roty(q ,theta )
      s=rotx(t ,thi )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1
      x=i1+pg
      y=j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y = f1(q.x,q.y) 
      t=roty(q ,theta )
      s=rotx(t ,thi )
      p = persp(s ,d )   
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1     
      x = i1
      y = j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )   
      q.y = f1(q.x,q.y) 
      t=roty(q ,theta)
      s=rotx(t ,thi )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
     ' fill_polygon(a(), CULng(rnd*&hFFFFFF))
      
      fill_polygon(a(), rgb(45,0,0))
     outline_polygon(a() , colour)

   next i1
  
   
  '
next j1
'
      Flip 1,0    'Copies our graph from page 1 to page 0


'
'sleep 10
'
next theta
'
end sub
'
' ----------------------------------------------------------------
'
sub genmatrix(gm() as single,n as integer,m as integer)
'
'   Generate a matrix with limits [x,-1,1][y,-1,1]
'
'
static as integer i,j
static as single x,y,z
'
for j=0 to m
    y=1-2*i/m
    z=y
    for i=0 to n
      x=-1+2*i/n
      gm(i,j,0)=x
      gm(i,j,1)=z
    next i
next j      
'
'
end sub
'
' -------------------------------------------------------------------
'
sub gmdata(gm() as single,n as integer,m as integer)
'
'  Generate data from function f1(x,y) .
'
'
static as integer i,j
static as single x,y,z
'
for j=0 to m
    for i=0 to n
      x=gm(i,j,0)
      y=gm(i,j,1)
      z=f1(x,y)
      gm(i,j,0) = z
    next i
next j   
'
end sub
'
' ---------------------------------------------------------------------
'
sub trallg(gm() as single,dm() as single , n as integer,m as integer)
'
'   translate , rotate , apply perspective to all of gm()
'
'
static as integer i,j
static as single x,y,z,d
static as point p,q
'
d=0.8
'
for j=0 to m
     q.z=1-2*j/m
    for i=0 to n
      q.x=-1 +2*i/n
      q.y=gm(i,j,0)
      p = persp(q ,d )       
     dm(i,j,0) = p.x
     dm(i,j,1) = p.y     
    next i
next j   
'
'
end sub
'
' -------------------------------------------------------------------
'
sub trand(dm() as single,a() as long,n as integer,m as integer)
'
'  translate dm() to a() , also translate to screen coordinates for fill_polygon() routine .
'
'
static as integer i1,j1,k1,l1,qg,pg
static as long u , v
static as single x,y,z
static as point p , q
static as ulong colour

qg=10
pg=10



colour=rgb(20,120,20)
for j1=10 to m step qg
    
  for i1=10 to n step pg
  
      k1=0
      p.x=dm(i1,j1,0)
      p.y=dm(i1,j1,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
      k1=k1+1
      p.x=dm(i1+pg,j1,0)
      p.y=dm(i1,j1,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
      
      k1=k1+1  
      p.x=dm(i1+pg,j1,0)
      p.y=dm(i1,j1+qg,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
       
      k1=k1+1
      p.x=dm(i1,j1,0)
      p.y=dm(i1,j1+qg,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

   fill_polygon(a(), CULng(rnd*&hFFFFFF))
   'outline_polygon(a() , colour)
 
   next i1
  '
next j1
'
'
end sub
'
' -------------------------------------------------------------------
'
function f1(x as single , y as single)as single
'
'       function to generate values upon [x,-1,1][y,-1,1]
'
'       [ ,-1,1] 
'
static as single r , z

        r=x*x+y*y
        if ( r > 0 ) then
             r = 5*sqr(r)
             z=sin(r*Pi)*cos((x-0.5)*2*Pi)/(r*Pi)
        else
             z=-1
        end if          

       return z

'
end function
'
' -----------------------------------------------------------------------------------------------------
'
Sub fill_polygon(a() As Long, ByVal c As ULong)
   'translation of a c snippet by Angad
   'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
   Dim As Long      i, j, k, dy, dx, x, y, temp
   Dim As Long      xi(0 to Ubound(a, 1))
   Dim As Single    slope(0 to Ubound(a, 1))
   
   'join first and last vertex
   a(Ubound(a, 1), 0) = a(0, 0)
   a(Ubound(a, 1), 1) = a(0, 1)

   For i = 0 To Ubound(a, 1) - 1

      dy = a(i+1, 1) - a(i, 1)
      dx = a(i+1, 0) - a(i, 0)

      If (dy = 0) Then slope(i) = 1.0
      If (dx = 0) Then slope(i) = 0.0

      If (dy <> 0) AndAlso (dx <> 0) Then slope(i) = dx / dy
   Next i

   For y = 0 to SCR_H - 1
      k = 0
      ' using FB's short-cut operators (which C doesn't have!)
      For i = 0 to Ubound(a, 1) - 1
         If (a(i, 1) <= y AndAlso a(i+1, 1) > y) OrElse _
             (a(i, 1) > y AndAlso a(i+1, 1) <= y) Then
            xi(k) = CLng(a(i, 0) + slope(i) * (y - a(i, 1)))
            k += 1
         End If
      Next i

      For j = 0 to k - 2
         'Arrange x-intersections in order
         For i = 0 To k - 2
            If (xi(i) > xi(i + 1)) Then
               temp = xi(i)
               xi(i) = xi(i + 1)
               xi(i + 1) = temp
            End If
         Next i
      Next j
      'line filling
      For i = 0 To k - 2 Step 2
         Line (xi(i), y)-(xi(i + 1) + 1, y), c
      Next i
   Next y
End Sub
'
' -----------------------------------------------------------------------------
' 
sub outline_polygon(a() As Long, ByVal c As ULong)
'
'  Draw an outline for the polygon , in color c  .
'
  'translation of a translated  c snippet by Angad
   'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
   Dim As Long      i, j,  x, y, u , v , temp
      
   'join first and last vertex
   a(Ubound(a, 1), 0) = a(0, 0)
   a(Ubound(a, 1), 1) = a(0, 1)

   For i = 0 To Ubound(a, 1) - 1
       x=a(i,0)
       y=a(i,1)
       u=a(i+1,0)
       v=a(i+1,1)
       line(x,y)-(u,v),c
   Next i


end sub
'
' ----------------------------------------------------------------------
'
sub w2scrn( p as point,n as integer , m as integer, byref u as long, byref v as long)
'
' input x,y,n,m
' output u , v
'
     p.x =  (p.x + 1)*0.5*(n - 10) + 10  ' [-1,1] -> [10,n]
     p.y =  (p.y + 1)*0.5*(m - 10) + 10  ' [-1,1] -> [10,m]
     u = clng(p.x)
     v = clng(p.y)

'
end sub
'
' -------------------------------------------------------------
'
sub tqxyz(x as single , y as single ,n as integer , m as integer , byref q as point)
'
'  translate  to  x ~ i1 , y ~ j1 to q.x , q.y , q.z
'

      q.z = 1 - 2*(y-10)/(m-10)
      q.x = -1 + 2*(x-10)/(n-10)
      q.y = -1 + 2*(y-10)/(m-10)

'
end sub








Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

To partly round off the 3d graf examples , here's an animation .

The function f1(x,y,fp) might be whatever you're interested in .

I'm not introducing height colour at this stage .

I'm wondering if an array holding polygon vertices might work with OpenGL , anyway that's a bit
off topic .

Code: Select all


'
' -----------------------------------------------------------------------------
'
'   My graf 3d
'
'   MyGraf5_3d.bas
'
'
'    (c) copyright 2019 , sciwiseg@gmail.com ,
'
'             Edward.Q.Montague.  [ alias]
'
'
'
'
'
' -----------------------------------------------------------------------------
'
type point
         x as single
         y as single
         z as single
         u as single '  possible extension for special coord system
end type
'
const Pi=4*atn(1)

' replaces the defines above (single line Macro's in FB)
Const As Long   POINTS = 4 , POLYGONS = 20, SCR_W = 740, SCR_H = 680

'
dim as single x1,y1,z1,x2,y2,z2
dim as integer i,j,k
'
' ----------------------------------------------------------------------------
'
'
'                  Looking at a cube .
'
'            cc -1,1 _______<_______  1,1    start         z = -1
'                   |               |             back face.
'                   |               |
'                   v               ^
'                   |               |
'                   |_______________|
'                -1,-1        >        1,-1
'               
'
' -----------------------------------------------------------------------------
'
declare function rotx(q as point,angx as single) as point
declare function roty1(q as point,angy as single) as point
declare function roty(q as point,angy as single) as point


declare function rotz(q as point,angz as single) as point
declare function tranx(q as point,movx as single) as point
declare function trany(q as point,movy as single) as point
declare function tranz(q as point,movz as single) as point
declare function persp(q as point,d as single) as point
'
declare function Trall( p1() as point,n as integer,edge() as integer, div as integer ,np as integer,ne as integer) as integer

declare sub gendata1(a() as long,n as integer,m as integer)
declare sub gendata2(a() as long,n as integer,m as integer)


declare sub genmatrix(gm() as single,n as integer,m as integer)
declare sub gmdata(gm() as single,n as integer,m as integer)
declare function f1(x as single , y as single,fp as single)as single
declare sub trallg(gm() as single,dm() as single , n as integer,m as integer)
declare sub trand(dm() as single,a() as long,n as integer,m as integer)
declare Sub fill_polygon(a() As Long, ByVal c As ULong)
declare sub outline_polygon(a() As Long, ByVal c As ULong)

declare sub tqxyz(x as single , y as single ,n as integer , m as integer , byref q as point)
declare sub w2scrn(p as point,n as integer , m as integer,byref u as long,byref v as long)

'
' ============================================================================
'
dim as integer np , ne,n,m
'
restore  storeA
read np
'
dim as point p1(1 to np)
'
restore store1
for i =1 to np
   read p1(i).x
   read p1(i).y
   read p1(i).z
next i
'
'
restore  storeB
read ne
'
dim as integer edge(1 to ne,0 to 1)
'
restore store2
for i = 1 to ne
   read edge(i,0)
   read edge(i,1)
next i
'
' -----------------------------------------------------------------------------
'
'screen 12
'window (-1.5,-1.5)-(1.5,1.5)
'line (-1.4,-1.4)-(1.4,1.4),11,b
'
'cls
'k=Trall( p1() ,3,edge() , 32 ,np ,ne )
'sleep
'INITIALIZING GRAPHICS _________________________________________________
ScreenRes(SCR_W, SCR_H, 24,2)      'initialize graphics
'window(10,10)-(210,110)
ScreenSet 1,0
Cls

n=520
m=520
'
dim a(0 to POINTS,0 to 1) as long
dim gm(0 to n , 0 to m,0 to 1) as single
dim dm(0 to n , 0 to m,0 to 1) as single
dim as ulong colour
'
gendata2(a() ,n ,m )
'
sleep
end

'genmatrix(gm() ,n ,m )
'gmdata(gm() ,n ,m )
'trallg(gm(),dm()  , n ,m )
'trand(dm() ,a() ,n ,m )


'sleep
'end
'
' ===================================
'
'    number of vertices
'
storeA:
data 8 
'
'    vertex data , easier to keep track of
'  data when we use multiple data statements.
'
store1:
data  1,1,1
data  -1,1,1
data  -1,-1,1
data 1,-1,1
data  1,1,-1
data  -1,1,-1
data  -1,-1,-1
data  1,-1,-1
'
'   Number of edges.
'
storeB:
data 12
'
'  edge data
'
store2:
data 1,2
data 1,4
data 1,5
data 2,3
data 2,6
data 3,4
data 3,7
data 4,8
data 5,6
data 5,8
data 6,7
data 7,8
'
' -------------------------------------------------------------------------------
'
function rotx(q as point,angx as single) as point
'
'              Rotate around x axis .
'
static as point p
'
             p.x = q.x
             p.y= q.y*cos(angx)-sin(angx)*q.z
             p.z= q.z*cos(angx)+sin(angx)*q.y
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function roty1(q as point,angy as single) as point
'
'              Rotate around y axis .
'
static as point p
'
            p.x = sin(angy)*q.z + cos(angy)*q.x
            p.y = q.y
            p.z = cos(angy)*q.z -sin(angy)*q.x
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function roty(q as point,angy as single) as point
'
'              Rotate around y axis .
'
static as point p
'
            p.x = sin(angy)*q.z + cos(angy)*q.x
            p.y = q.y
            p.z = cos(angy)*q.z -sin(angy)*q.x
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function rotz(q as point,angz as single) as point
'
'              Rotate around z axis .
'
static as point p
'
            p.x = sin(angz)*q.y + cos(angz)*q.x
            p.y = cos(angz)*q.y-sin(angz)*q.x
            p.z = q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function tranx(q as point,movx as single) as point
'
'             Translate point along x axis
'
static as point p
'
              p.x=q.x + movx
              p.y=q.y
              p.z=q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function trany(q as point,movy as single) as point
'
'              Translate point along y axis .
'
static as point p
'
              p.x=q.x
              p.y=q.y + movy
              p.z=q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function tranz(q as point,movz as single) as point
'
'              Translate point along z axis .
'
static as point p
'
              p.x=q.x
              p.y=q.y
              p.z=q.z + movz
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function persp(q as point,d as single) as point
'
'              3d  perspective . 
'
'    Add 2.5 to the numerator when using any negative z value.
'
'     In this instance    -1 <= z  <= 1  ,  unit cube .
'
'    Therefore 2.5 is appropriate .
'
static as point p
'
     p.x = d*q.x/(q.z+2.25)
     p.y = d*q.y/(q.z+2.25)
     p.z = d
'
     return p
'
end function
'
' -----------------------------------------------------------------------------
'
function Trall( p1() as point,n as integer,edge() as integer, div as integer ,np as integer,ne as integer) as integer
'
'              Translate and rotate all vertices .
'     as an animation ,  for  n  cycles .
'
'    np number of points .
'    ne number of integers .
'
'
'   With  div number of angle divisions .
'
static as point p2(1 to np)
static as single theta,thi,x1,y1,z1,x2,y2,z2
static as integer i,j,k
static as integer i1,j1,k1
'
theta = Pi/div
'
for i=1 to n
  for j = 0 to div
  'cls
       thi = j*theta
   for k = 1 to np
     p2(k) = roty(p1(k),thi)
     p2(k)=persp(p2(k),0.8)
   next k   
 '  
    cls 
'
for i1 = 1 to ne
      j1 = edge(i1,0)
     k1 = edge(i1,1)
   x1 = p2(j1).x
   y1 = p2(j1).y
'   z1 = p2(j1).z   
   x2 = p2(k1).x
   y2 = p2(k1).y
'   z2 = p2(k1).z   
line(x1,y1)-(x2,y2),14
next i1     
'
sleep 100
  next j
next i
'
     return 0
'
end function
'
' --------------------------------------------------------------------------------
'
'
' ----------------------------------------------------------------
'
sub gendata1(a() as long,n as integer,m as integer)
'
'  Generate 3d data grid .
'
'  For perspective the range for x and y and z , is [-1,1]
'
'
'  Then a translation back to screen coordinates is required .  [-1,1] => [10,n] & [-1,1] => [10,m]  ?
'    x = (x + 1)*0.5*(n - 10) + 10
'    y = (y + 1)*0.5*(m - 10) + 10
'
'
'   w = f(x,y)
'   q.x = x
'   q.y = w
'   q.z = z
'
'   p = persp(q,d)
'
'   ( x', y',z , w') = persp( x , y , z , w )
'
'
static as integer i1,j1,k1,l1,qg,pg
static as single x,y,z,d
static as ulong colour,u,v
static as point p , q

d=0.8

'qg=10
'pg=10

pg=n/20
qg=m/10

colour=rgb(120,200,200)
for j1=10 to m step qg
    'z = 1 - 2*(j1-10)/(m-10) '  1  ->  -1 
    
  for i1=10 to n step pg
  
      k1=0
      x=i1
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
     ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
            
      k1=k1+1
      x=i1+pg
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
    ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1
      x=i1+pg
      y=j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )  
    ' f1(q.x,q.y) here ?
      p = persp(q ,d )   
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1     
      x = i1
      y = j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )   
    ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
      fill_polygon(a(), CULng(rnd*&hFFFFFF))
      outline_polygon(a() , colour)
 
   next i1
  '
next j1
'
end sub
'
' ----------------------------------------------------------------
'
sub gendata2(a() as long,n as integer,m as integer)
'
'  Generate 3d data grid .
'
'  For perspective the range for x and y and z , is [-1,1]
'
'
'  Then a translation back to screen coordinates is required .  [-1,1] => [10,n] & [-1,1] => [10,m]  ?
'    x = (x + 1)*0.5*(n - 10) + 10
'    y = (y + 1)*0.5*(m - 10) + 10
'
'
'   w = f(x,y)
'   q.x = x
'   q.y = w
'   q.z = z
'
'   p = persp(q,d)
'
'   ( x', y',z , w') = persp( x , y , z , w )
'
'
static as integer i1,j1,k1,l1,qg,pg
static as single x,y,z,d , theta , thi , fp
static as ulong colour,u,v,chrome
static as point p , q , s , t
static as integer i,j,k

theta = Pi/9
'theta=0
thi=Pi/7
'thi=0
d=0.98
d=1.2

'qg=10
'pg=10

pg=n/50
qg=m/50

colour=rgb(120,200,200)
for  fp=0.75 to 4 step 0.0041
'
'cls
'
'  Alternately put a previously generated or loaded background image .
'
for j=0 to SCR_H
    k=int(200*j/SCR_H)
    colour=rgb(12,20,k)

   line(0,j)-(SCR_W,j),colour
next j  
'
colour=rgb(120,200,200)
'
for j1=10 to m step qg
    'z = 1 - 2*(j1-10)/(m-10) '  1  ->  -1 
    
  for i1=10 to n step pg
  
      k1=0
      x=i1
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y= f1(q.x,q.y,fp) 
      chrome = (q.y + 1)*32
      
      t=roty(q ,theta )
      s=rotx(t ,thi )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
            
      k1=k1+1
      x=i1+pg
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y = f1(q.x,q.y,fp) 
      t=roty(q ,theta )
      s=rotx(t ,thi )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1
      x=i1+pg
      y=j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y = f1(q.x,q.y,fp) 
      t=roty(q ,theta )
      s=rotx(t ,thi )
      p = persp(s ,d )   
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1     
      x = i1
      y = j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )   
      q.y = f1(q.x,q.y,fp) 
      t=roty(q ,theta)
      s=rotx(t ,thi )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
     ' fill_polygon(a(), CULng(rnd*&hFFFFFF))
      
      fill_polygon(a(), rgb(45,0,0))
     outline_polygon(a() , colour)

   next i1
  
   
  '
next j1
'
      Flip 1,0    'Copies our graph from page 1 to page 0


'
'sleep 10
'
next fp
'
end sub
'
' ----------------------------------------------------------------
'
sub genmatrix(gm() as single,n as integer,m as integer)
'
'   Generate a matrix with limits [x,-1,1][y,-1,1]
'
'
static as integer i,j
static as single x,y,z
'
for j=0 to m
    y=1-2*i/m
    z=y
    for i=0 to n
      x=-1+2*i/n
      gm(i,j,0)=x
      gm(i,j,1)=z
    next i
next j      
'
'
end sub
'
' -------------------------------------------------------------------
'
sub gmdata(gm() as single,n as integer,m as integer)
'
'  Generate data from function f1(x,y) .
'
'
static as integer i,j
static as single x,y,z,fp
'
fp=1
'
for j=0 to m
    for i=0 to n
      x=gm(i,j,0)
      y=gm(i,j,1)
      z=f1(x,y,fp)
      gm(i,j,0) = z
    next i
next j   
'
end sub
'
' ---------------------------------------------------------------------
'
sub trallg(gm() as single,dm() as single , n as integer,m as integer)
'
'   translate , rotate , apply perspective to all of gm()
'
'
static as integer i,j
static as single x,y,z,d
static as point p,q
'
d=0.8
'
for j=0 to m
     q.z=1-2*j/m
    for i=0 to n
      q.x=-1 +2*i/n
      q.y=gm(i,j,0)
      p = persp(q ,d )       
     dm(i,j,0) = p.x
     dm(i,j,1) = p.y     
    next i
next j   
'
'
end sub
'
' -------------------------------------------------------------------
'
sub trand(dm() as single,a() as long,n as integer,m as integer)
'
'  translate dm() to a() , also translate to screen coordinates for fill_polygon() routine .
'
'
static as integer i1,j1,k1,l1,qg,pg
static as long u , v
static as single x,y,z
static as point p , q
static as ulong colour

qg=10
pg=10



colour=rgb(20,120,20)
for j1=10 to m step qg
    
  for i1=10 to n step pg
  
      k1=0
      p.x=dm(i1,j1,0)
      p.y=dm(i1,j1,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
      k1=k1+1
      p.x=dm(i1+pg,j1,0)
      p.y=dm(i1,j1,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
      
      k1=k1+1  
      p.x=dm(i1+pg,j1,0)
      p.y=dm(i1,j1+qg,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
       
      k1=k1+1
      p.x=dm(i1,j1,0)
      p.y=dm(i1,j1+qg,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

   fill_polygon(a(), CULng(rnd*&hFFFFFF))
   'outline_polygon(a() , colour)
 
   next i1
  '
next j1
'
'
end sub
'
' -------------------------------------------------------------------
'
function f1(x as single , y as single,fp as single)as single
'
'       function to generate values upon [x,-1,1][y,-1,1]
'
'       [ ,-1,1] 
'
static as single r , z

        r=x*x+y*y
        if ( r > 0 ) then
             r = 5*sqr(r)*fp
             z=(sin(r*Pi)/(r*Pi))*cos((x-0.5)*2*Pi*fp)
        else
             z=-1
        end if          

       return z

'
end function
'
' -----------------------------------------------------------------------------------------------------
'
Sub fill_polygon(a() As Long, ByVal c As ULong)
   'translation of a c snippet by Angad
   'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
   Dim As Long      i, j, k, dy, dx, x, y, temp
   Dim As Long      xi(0 to Ubound(a, 1))
   Dim As Single    slope(0 to Ubound(a, 1))
   
   'join first and last vertex
   a(Ubound(a, 1), 0) = a(0, 0)
   a(Ubound(a, 1), 1) = a(0, 1)

   For i = 0 To Ubound(a, 1) - 1

      dy = a(i+1, 1) - a(i, 1)
      dx = a(i+1, 0) - a(i, 0)

      If (dy = 0) Then slope(i) = 1.0
      If (dx = 0) Then slope(i) = 0.0

      If (dy <> 0) AndAlso (dx <> 0) Then slope(i) = dx / dy
   Next i

   For y = 0 to SCR_H - 1
      k = 0
      ' using FB's short-cut operators (which C doesn't have!)
      For i = 0 to Ubound(a, 1) - 1
         If (a(i, 1) <= y AndAlso a(i+1, 1) > y) OrElse _
             (a(i, 1) > y AndAlso a(i+1, 1) <= y) Then
            xi(k) = CLng(a(i, 0) + slope(i) * (y - a(i, 1)))
            k += 1
         End If
      Next i

      For j = 0 to k - 2
         'Arrange x-intersections in order
         For i = 0 To k - 2
            If (xi(i) > xi(i + 1)) Then
               temp = xi(i)
               xi(i) = xi(i + 1)
               xi(i + 1) = temp
            End If
         Next i
      Next j
      'line filling
      For i = 0 To k - 2 Step 2
         Line (xi(i), y)-(xi(i + 1) + 1, y), c
      Next i
   Next y
End Sub
'
' -----------------------------------------------------------------------------
' 
sub outline_polygon(a() As Long, ByVal c As ULong)
'
'  Draw an outline for the polygon , in color c  .
'
  'translation of a translated  c snippet by Angad
   'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
   Dim As Long      i, j,  x, y, u , v , temp
      
   'join first and last vertex
   a(Ubound(a, 1), 0) = a(0, 0)
   a(Ubound(a, 1), 1) = a(0, 1)

   For i = 0 To Ubound(a, 1) - 1
       x=a(i,0)
       y=a(i,1)
       u=a(i+1,0)
       v=a(i+1,1)
       line(x,y)-(u,v),c
   Next i


end sub
'
' ----------------------------------------------------------------------
'
sub w2scrn( p as point,n as integer , m as integer, byref u as long, byref v as long)
'
' input x,y,n,m
' output u , v
'
     p.x =  (p.x + 1)*0.5*(n - 10) + 10  ' [-1,1] -> [10,n]
     p.y =  (p.y + 1)*0.5*(m - 10) + 10  ' [-1,1] -> [10,m]
     u = clng(p.x)
     v = clng(p.y)

'
end sub
'
' -------------------------------------------------------------
'
sub tqxyz(x as single , y as single ,n as integer , m as integer , byref q as point)
'
'  translate  to  x ~ i1 , y ~ j1 to q.x , q.y , q.z
'

      q.z = 1 - 2*(y-10)/(m-10)
      q.x = -1 + 2*(x-10)/(n-10)
      q.y = -1 + 2*(y-10)/(m-10)

'
end sub






dodicat
Posts: 7976
Joined: Jan 10, 2006 20:30
Location: Scotland

Re: 3D Geometry , basics

Post by dodicat »

Thank you Luxan.
Really nice.
I tried with other f1 functions.
I like the polygon fill (and outline), really fast and efficient.
I note you use a 2D array instead of a standard p.x,p.y, but it works well, and for p.x,p.y people, these points are easily converted to an array.
I had a similar thing a few years ago, but I used the disphelper lib for the parser.
But I have used a different parser here with the same code.

Code: Select all

 


'Surface plotter

#include "crt.bi"
Dim Shared e_input    As String  
Dim Shared e_tok      As String  
Dim Shared e_spelling As String  
Dim Shared e_error    As long 
Dim Pi As Double = 4 * Atn(1)

'==============  PARSER START  ==================================
'sinh, cosh, tanh from crt.bi
Function SEC(Byval x As Double) As Double
    SEC = 1 / Cos(x)
End Function

Function COSEC(Byval x As Double) As Double
    COSEC = 1 / Sin(x)
End Function

Function COT(Byval x As Double) As Double
    COT = 1 / Tan(x)
End Function

Function ARCSEC(Byval x As Double) As Double ''''''
    ARCSEC = Atn(x / Sqr(x * x - 1)) + Sgn((x) -1) * (2 * Atn(1))
End Function

Function ARCCOSEC(Byval x As Double) As Double
    ARCCOSEC = Atn(x / Sqr(x * x - 1)) + (Sgn(x) - 1) * (2 * Atn(1))
End Function

Function ARCCOT(Byval x As Double) As Double
    ARCCOT = Atn(x) + 2 * Atn(1)
End Function

Function sech(Byval x As Double) As Double
    sech = 2 / (Exp(x) + Exp(-x))
End Function

Function cosech(Byval x As Double) As Double
    cosech = 2 / (Exp(x) - Exp(-x))
End Function

Function coth(Byval x As Double) As Double
    coth = (Exp(x) + Exp(-x)) / (Exp(x) - Exp(-x))
End Function

Function arcsinh(Byval x As Double) As Double
    arcsinh = Log(x + Sqr(x * x + 1))
End Function

Function arccosh(Byval x As Double) As Double
    arccosh = Log(x + Sqr(x * x - 1))
End Function

Function arctanh(Byval x As Double) As Double
    arctanh = Log((1 + x) / (1 - x)) / 2 
End Function

Function arcsech(Byval x As Double) As Double
    arcsech = Log((Sqr(-x * x + 1) + 1) / x)
End Function

Function arccosech(Byval x As Double) As Double
    arccosech = Log((Sgn(x) * Sqr(x * x + 1) +1) / x)
End Function

Function arccoth(Byval x As Double) As Double
    arccoth = Log((x + 1) / (x - 1)) / 2
End Function

Function HAVERSINE(Byval x As Double) As Double
    HAVERSINE = (Sin(x/2))^2
End Function

function pie(byval x as double=1) as double
    return (4*atn(1))*x
    end function

Function e_function(Byref fun As String,Byval arg As Double) As Double
    Dim n As Double
    
    Select Case Lcase(fun)
    Case "abs": n = Abs(arg)
    Case "atn": n = Atn(arg)
    Case "cos": n = Cos(arg)
    Case "exp": n = Exp(arg)
    Case "ezp": n = Exp(arg)
    Case "fix": n = Fix(arg)
    Case "int": n = Int(arg)
    Case "log": n = Log(arg)
    Case "rnd": n = Rnd(arg)
    Case "sgn": n = Sgn(arg)
    Case "sin": n = Sin(arg)
    Case "sqr": n = Sqr(arg)
    Case "tan": n = Tan(arg)
    Case "haversine":n=haversine(arg)
    Case "cosec":n=cosec(arg)
    Case "sec":n=sec(arg)
    Case "cot": n=cot(arg)
    Case "asin":n=Asin(arg)
    Case "acos":n=Acos(arg)
    Case "atn":n=Atn(arg)
    Case "arcsec":n=arcsec(arg)
    Case "arccosec":n=arccosec(arg)
    Case "arccot":n=arccot(arg)
    Case "sinh":n=sinh(arg)
    Case "cosh":n=cosh(arg)
    Case "tanh":n=tanh(arg)
    Case "sech":n=sech(arg)
    Case "cosech":n=cosech(arg)
    Case "coth":n=coth(arg)
    Case "arcsinh":n=arcsinh(arg)
    Case "arccoth":n=arccoth(arg)
    Case "arctanh":n=arctanh(arg)
    Case "arcsech":n=arcsech(arg)
    Case "arccosech":n=arccosech(arg)
    Case "pi"      :n=pie(arg)
    Case Else
        If Not e_error Then
            Locate 1,1
            Print "UNDEFINED FUNCTION " + fun
            Print
            e_error = -1
        End If
    End Select
    e_function = n
End Function

Sub e_nxt()
    Dim is_keyword As long
    Dim c As String 
    e_tok = ""
    e_spelling = ""
    Do
        c = Left(e_input, 1)
        e_input = Mid(e_input, 2)
    Loop While c = " " Or c = Chr(9) Or c = Chr(13) Or c = Chr(10)
    
    Select Case Lcase(c)
    
    Case "0" To "9", "."
        e_tok = "num"
        Do
            e_spelling = e_spelling + c
            c = Left(e_input, 1)
            e_input = Mid(e_input, 2)
        Loop While (c >= "0" And c <= "9") Or c = "."
        e_input = c + e_input
        
    Case "a" To "z", "_"
        Dim As long is_id
        e_tok = "id"
        Do
            e_spelling = e_spelling + c
            c = Lcase(Left(e_input, 1))
            e_input = Mid(e_input, 2)
            is_id = (c >= "a" And c <= "z")
            is_id = is_id Or c = "_" Or (c >= "0" And c <= "9")
        Loop While is_id
        e_input = c + e_input
        is_keyword = -1
        Select Case Lcase(e_spelling)
        Case "and"
        Case "eqv"
        Case "imp"
        Case "mod"
        Case "not"
        Case "or"
        Case "xor"
        Case Else: is_keyword = 0
        End Select
        If is_keyword Then
            e_tok = Lcase(e_spelling)
        End If
        
    Case "<", ">"
        e_tok = c
        c = Left(e_input, 1)
        If c = "=" Or c = ">" Then
            e_tok = e_tok + c
            e_input = Mid(e_input, 2)
        End If
        
    Case Else
        e_tok = c
    End Select
    
    If e_spelling = "" Then
        e_spelling = e_tok
    End If
End Sub

Sub e_match (Byref token As String)
    If Not e_error And e_tok <> token Then
        Locate 1,1
        Print "EXPECTED " + token + ", got '" + e_spelling + "'"
        e_error = -1':end
    End If
    e_nxt()
End Sub

Function e_prs (Byval p As long) As Double
    Dim n   As Double  
    Dim fun As String  
    If e_tok = "num" Then
        n = Val(e_spelling)
        e_nxt()
    Elseif e_tok = "-" Then
        e_nxt()
        n = -e_prs(11)   ''   11 before 12  jhh
    Elseif e_tok = "not" Then
        e_nxt()
        n = Not e_prs(6) 
    Elseif e_tok = "(" Then
        e_nxt()
        n = e_prs(1)
        e_match(")")
    Elseif e_tok = "id" Then
        fun = e_spelling
        e_nxt()
        e_match("(")
        n = e_prs(1)
        e_match(")")
        n = e_function(fun, n)
    Else
        If Not e_error Then
            Locate 1,1
            Print "syntax error, at '" + e_spelling + "'"
            e_error = -1':end
        End If
    End If
    
    Do While Not e_error
        If     p <= 11 And e_tok = "^"   Then 
            e_nxt(): n = n ^ e_prs(12)
        Elseif p <= 10 And e_tok = "*"   Then 
            e_nxt(): n = n *   e_prs(11)
        Elseif p <= 10 And e_tok = "/"   Then 
            e_nxt(): n = n /   e_prs(11)
        Elseif p <= 9  And e_tok = "\"   Then 
            e_nxt(): n = n \   e_prs(10)
        Elseif p <= 8  And e_tok = "mod" Then 
            e_nxt(): n = n Mod e_prs(9)
        Elseif p <= 7  And e_tok = "+"   Then 
            e_nxt(): n = n +   e_prs(8)
        Elseif p <= 7  And e_tok = "-"   Then 
            e_nxt(): n = n -   e_prs(8)
        Elseif p <= 6  And e_tok = "="   Then 
            e_nxt(): n = n =   e_prs(7)
        Elseif p <= 6  And e_tok = "<"   Then 
            e_nxt(): n = n <   e_prs(7)
        Elseif p <= 6  And e_tok = ">"   Then 
            e_nxt(): n = n >   e_prs(7)
        Elseif p <= 6  And e_tok = "<>"  Then 
            e_nxt(): n = n <>  e_prs(7)
        Elseif p <= 6  And e_tok = "<="  Then 
            e_nxt(): n = n <=  e_prs(7)
        Elseif p <= 6  And e_tok = ">="  Then 
            e_nxt(): n = n >=  e_prs(7)
        Elseif p <= 5  And e_tok = "and" Then 
            e_nxt(): n = n And e_prs(6)
        Elseif p <= 4  And e_tok = "or"  Then 
            e_nxt(): n = n Or  e_prs(5)
        Elseif p <= 3  And e_tok = "xor" Then 
            e_nxt(): n = n Xor e_prs(4)
        Elseif p <= 2  And e_tok = "eqv" Then 
            e_nxt(): n = n Eqv e_prs(3)
        Elseif p <= 1  And e_tok = "imp" Then 
            e_nxt(): n = n Imp e_prs(2)
        Else
            Exit Do
        End If
    Loop
    e_prs = n
End Function

Function eval(Byref sp As String ) As Double
    Dim As Double value
    e_error = 0
    e_input = sp
    e_nxt()
    value = e_prs(1)
    If Not e_error Then Return value else e_error=0
    
End Function

Function FindAndReplace(Byref instring As String,Byref ReplaceThis As String,Byref WithThis As String) As String 
    Var lens1=Len(ReplaceThis),lens2=Len(WithThis)
    If lens1=lens2 Then lens1=0
    Dim As String s=instring 
    Dim As long position=Instr(s,ReplaceThis)
    While position>0
        If lens1 Then   
            s=Left(s,position-1) & WithThis & Mid(s,position+Lens1)
        Else
            Mid(s,position) = WithThis
        End If
        position=Instr(position+Lens2,s,ReplaceThis)
    Wend
    Function=s
End Function

Sub Setvariable(s As String,REPLACE_THIS As String,WITHTHIS As Double)' As String
    var WITH_THIS=Str(WITHTHIS)
    var position=Instr(s,REPLACE_THIS)
    While position>0
        s=Mid(s,1,position-1) & WITH_THIS & Mid(s,position+Len(REPLACE_THIS))
        position=Instr(position+Len(WITH_THIS),s,REPLACE_THIS)
    Wend
End Sub
#macro _window(topleftX,topleftY,bottomrightX,bottomrightY) 
minx=topleftX
maxx=bottomrightX
miny=bottomrightY
maxy=topleftY
#endmacro
#macro _axis(colour)
If Sgn(minx)<>Sgn(maxx) Then
    Line(0,(yres-(miny/(miny-maxy))*yres))-(xres,(yres-(miny/(miny-maxy))*yres)),colour 'x axis
    Endif
    If Sgn(miny)<>Sgn(maxy) Then
        Line(((minx/(minx-maxx))*xres),0)-(((minx/(minx-maxx))*xres),yres),colour 'y axis
        Endif
        Draw String(0,yres/2),Str(minx),colour
        Draw String(xres-8-8*(Len(Str(maxx))),yres/2),Str(maxx),colour
        Draw String(xres/2,0),Str(maxy),colour
        Draw String(xres/2,yres-16),Str(miny),colour
        #endmacro
        Sub inspect
            Dim As long mx,my,mw
            mx=70:my=230
            mw=2
            Dim As Ulong array(1 To 24641)
            Dim As long count
            For x As long=mx-00 To mx+600
                For y As long=my-20 To my+20
                    count=count+1
                    array(count)=Point(x,y)
                Next y
            Next x
            count=0
            'draw
            For x As long=mx-00 To mx+600
                For y As long=my-20 To my+20
                    count=count+1
                    var NewX=mw*(x-mx)+mx
                    var NewY=mw*(y-my)+my 
                    Line(NewX-mw/2,NewY-mw/2)-(NewX+mw/2,NewY+mw/2),array(count),BF
                Next y
            Next x
        End Sub
        'progress
        Type bar
            As long x,y,l,d,percent
            As Ulong col
        End Type
        Dim As long percentage
        #define progress(value,lower,upper) 100*(value-lower)/(upper-lower)
        
        Sub progressbar(b As bar)
            Line(b.x+1,b.y+1)-( (b.l*b.percent/100+b.x),b.y+b.d-1),b.col,bf
            Line(b.x,b.y)-(b.x+b.l,b.y+b.d),,b
        End Sub
        Dim  As bar b
        b=Type<bar>(100,300,600,20,0,Rgb(0,0,200))
        
        Dim Shared As Integer xres,yres
        Dim  As Single minx,maxx,miny,maxy,PLOT_GRADE=15000
        
        Screen 20,32
        Screeninfo xres,yres
        Type vector3d
            As Single x,y,z
        End Type
        #define vct type<vector3d>
        #define distance(p1,p2) sqr((p1.x-p2.x)*(p1.x-p2.x)+(p1.y-p2.y)*(p1.y-p2.y))
        #macro incircle(cx,cy,radius,x,y)
        (cx-x)*(cx-x) +(cy-y)*(cy-y)<= radius*radius
        #endmacro
        Dim Shared As vector3d eyepoint
        Dim Shared Rx(1 To 3,1 To 3) As Single
        Dim Shared Ry(1 To 3,1 To 3) As Single
        Dim Shared Rz(1 To 3,1 To 3) As Single
        Dim Shared pivot_vector(1 To 3) As Single
        Dim Shared new_pos(1 To 3) As Single
        Dim Shared temp1(1 To 3) As Single
        Dim Shared temp2(1 To 3) As Single
        Operator + (v1 As vector3d,v2 As vector3d) As vector3d
        Return Type<vector3d>(v1.x+v2.x,v1.y+v2.y,v1.z+v2.z)
        End Operator
        
        Operator -(v1 As vector3d,v2 As vector3d) As vector3d
        Return Type<vector3d>(v1.x-v2.x,v1.y-v2.y,v1.z-v2.z)
        End Operator
        
        Operator * (f As Single,v1 As vector3d) As vector3d
        Return Type<vector3d>(f*v1.x,f*v1.y,f*v1.z)
        End Operator
        Function r(first As Double, last As Double) As Double
            Function = Rnd * (last - first) + first
        End Function
        
        Function apply_perspective(p As vector3d) As vector3d
            Dim As Single   w=(p.z*(-1)/500+1)*.75'300
            Return vct((p.x-eyepoint.x)/w+eyepoint.x,(p.y-eyepoint.y)/w+eyepoint.y,(p.z-eyepoint.z)/w+eyepoint.z)
        End Function
        
        Sub framecounter
            Static As Double frame,fps
            frame=frame+1
            Static As Double t1,t2
            If frame>=fps Then
                t1 = Timer
                fps = frame/(t1-t2)
                Windowtitle "Frames per second = " & fps
                t2=Timer
                frame=0
            End If
        End Sub
        Function rotatepoint3d(Byval pivot As vector3d,_  
            Byval pt As vector3d,_  
            Byval angle As vector3d,_ 
            Byval dilator As Single=1) As vector3d  
            #macro mv(m1,v,ans)
            For i As long=1 To 3
                s=0
                For k As long = 1 To 3
                    s=s+m1(i,k)*v(k)
                Next k
                ans(i)=s
            Next i
            #endmacro
            #define cr 0.0174532925199433
            Dim angle_radians As vector3d=Type<vector3d>(cr*angle.x,cr*angle.y,cr*angle.z)
            
            Dim s As Single=Any
            pivot_vector(1)=(pt.x-pivot.x)*dilator
            pivot_vector(2)=(pt.y-pivot.y)*dilator
            pivot_vector(3)=(pt.z-pivot.z)*dilator
            
            'rotat1on matrices about the three axix
            Rx(1,1)=1:Rx(1,2)=0:Rx(1,3)=0
            Rx(2,1)=0:Rx(2,2)=Cos(angle_radians.x):Rx(2,3)=-Sin(angle_radians.x)
            Rx(3,1)=0:Rx(3,2)=Sin(angle_radians.x):Rx(3,3)=Cos(angle_radians.x)
            
            Ry(1,1)=Cos(angle_radians.y):Ry(1,2)=0:Ry(1,3)=Sin(angle_radians.y)
            Ry(2,1)=0:Ry(2,2)=1:Ry(2,3)=0
            Ry(3,1)=-Sin(angle_radians.y):Ry(3,2)=0:Ry(3,3)=Cos(angle_radians.y)
            
            Rz(1,1)=Cos(angle_radians.z):Rz(1,2)=-Sin(angle_radians.z):Rz(1,3)=0
            Rz(2,1)=Sin(angle_radians.z):Rz(2,2)=Cos(angle_radians.z):Rz(2,3)=0
            Rz(3,1)=0:Rz(3,2)=0:Rz(3,3)=1
            
            mv (Rx,pivot_vector,temp1)           
            mv (Ry,temp1,temp2)
            mv (Rz,temp2,new_pos)
            
            new_pos(1)=new_pos(1)+pivot.x
            new_pos(2)=new_pos(2)+pivot.y
            new_pos(3)=new_pos(3)+pivot.z
            
            Return Type<vector3d>(new_pos(1),new_pos(2),new_pos(3))
        End Function
        
        Sub blow(a() As vector3d,mag As Single)
            For z As long=1 To Ubound(a)-6
                a(z)=mag*a(z)
            Next z
        End Sub
        
        Sub translate(a() As vector3d,pt As vector3d)
            For z As long=1 To Ubound(a)
                a(z)=a(z)+vct(pt.x,pt.y,pt.z)
            Next z
        End Sub
        
        
        Function vertex(piv As vector3d,p1 As vector3d,ang As vector3d,dil As Single,col As Ulong) As Single
            var _temp1=rotatepoint3d(piv,p1,ang,dil)
            _temp1=apply_perspective(_temp1)
            Pset (_temp1.x,_temp1.y),col
            Return _temp1.z
        End Function
        
        Sub set_perspective(x As Single,y As Single,z As Single,minz As Single,maxz As Single)
            eyepoint=vct(x,y,z)
        End Sub
        #macro combsort(array,comp)
        Scope
            var size=Ub,switch=0,j=0
            Dim As Single void=size
            Do
                void=void/1.3: If void<1 Then void=1
                switch=0
                For i As long =1 To size-void
                    j=i+void
                    If comp(i)>comp(j) Then 
                        Swap array(i),array(j): switch=1
                        Swap comp(i),comp(j)
                        Swap col(i),col(j)
                    End If
                Next
            Loop Until  switch =0 And void=1
        End Scope
        #endmacro
        
        Redim Shared As vector3d e(0)
        Dim count As long
        Redim Shared As Ulong col(0)
        Dim As Single funct
        Dim As Single dist,scale=5
        _window(-scale,scale,scale,-scale)
        
        Dim As vector3d pt,cent=Type<vector3d>(0,0,0)
        
        Dim As String formula,worker
        
        Dim As String i,j
        Dim As long flag,ub
        Dim As long mx,my,mw,counter,mb
        dim as string blink
        dim as double t,k=1
        start:
       
        
        Do
             If Timer-t>.5 Then
        t=Timer
        k=-k
        If k=1 Then blink=" " Else blink="_"
            End If
            i=Inkey
            If Left(i,1)=Chr(08) Then j=Mid(j,1,Len(j)-1)
            Screenlock
            Cls 
            Draw String (.53*xres,.48*yres),"SCALERS"
            Circle(.53*xres,.53*yres),20,Rgb(100,100,100),,,,f
            Draw String(.525*xres,.525*yres),Chr(30)
            Circle(.58*xres,.53*yres),20,Rgb(100,100,100),,,,f
            Draw String(.575*xres,.525*yres),Chr(31)
            Getmouse mx,my,mw,mb
            If incircle(.53*xres,.53*yres,20,mx,my) Then   Circle(.53*xres,.53*yres),20,Rgb(255,255,255)
            If incircle(.58*xres,.53*yres,20,mx,my) Then   Circle(.58*xres,.53*yres),20,Rgb(255,255,255)
            
            If incircle(.53*xres,.53*yres,20,mx,my) And mb=1 Then 
                scale=scale+.1
                _window(-scale,scale,scale,-scale)
            End If
            If incircle(.58*xres,.53*yres,20,mx,my) And mb=1 Then 
                scale=scale-.1
                If scale<.1 Then scale=.1
                _window(-scale,scale,scale,-scale)
            End If
            
            Locate 3,3
            Print "  Example function sin(x)*cos(y)"
            Print "  You can also use d as a variable, which is distance from origin"
            Print "  E.G. (sin(2*d)/(2*d))*5"
            Print "  Use the scalers to adjust the X/Y plane"
            Print "  Enter a function in  x and y (or d) at the arrow "',formula
            print "  Enter end at the arrow to finish"
            If flag=1 Then
                Locate 10,12
                Print "Error in " & formula & " --please redo"
            End If
            Select Case Left(i,1)
            Case Chr(0) To Chr(254)
                If Left(i,1)<>Chr(08) Then
                    j=j+Left(i,1)
                End If
            End Select
            
            Locate 15,5
            Print "----> " & j +blink
            inspect
            _axis(Rgb(0,200,0))
            Screenunlock
            Sleep 1,1
        Loop Until i=Chr(13)
        j=Rtrim(j,Chr(13))
        formula=j
        formula=FindAndReplace(formula,"exp","ezp") 
        if formula="end" then end
        For x As Single=-scale To scale Step 2*scale/150
            For y As Single=-scale To scale Step 2*scale/150
                pt=Type<vector3d>(x,y)
                dist=distance(pt,cent)
                count=count+1
                Redim Preserve e(count)
                Redim Preserve col(count)
                col(count)=Rgb(255*(x+scale)/(2*scale),155*(y+scale)/(2*scale),50)
                worker=formula
                setvariable(worker,"x",x)
                setvariable(worker,"y",y)
                setvariable(worker,"d",dist)
                funct=eval(worker)
                e(count)=Type<vector3d>(x,y,funct)
            Next y
            percentage=progress(x,-scale,scale)
            b.percent=percentage
            progressbar(b)
        Next x
        ub=Ubound(e)
        Redim Preserve e(ub+6)
        'axis ends
        e(ub+1)=Type<vector3d>(-.5*(xres/2),0,0)
        e(ub+2)=Type<vector3d>(.5*(xres/2),0,0)
        e(ub+3)=Type<vector3d>(0,-.5*(yres/2),0)
        e(ub+4)=Type<vector3d>(0,.5*(yres/2),0)
        e(ub+5)=Type<vector3d>(0,0,-.5*(yres/2))
        e(ub+6)=Type<vector3d>(0,0,.5*(yres/2))
        Dim As Single dummy
        For z As long=1 To ub
            If e(z).z=0  Then dummy=dummy+1
        Next z
        If dummy=ub Then
            Beep
            Erase(e)
            Cls
            flag=1
            j=""
            Goto start
        End If
        flag=0
        blow(e(),20*5/scale)
        translate(e(),vct(xres/2,yres/2,0))
        set_perspective(xres/2,yres/2,0,-100,100)
        Dim As Single dilation
        Dim As vector3d piv,ang
        piv=eyepoint
        Dim As Single zeds(1 To Ub),_mw
        Dim As vector3d axis(6)
        Dim As Ulong colour
        Dim As Single startdilation=1
        dilation=startdilation
        Dim As Any Pointer im=imagecreate(xres,yres)
        Paint im,(0,0),Rgb(255,255,255)
        
        Do
            framecounter
            i=Inkey
            If i= Chr(255) + "K"  Then ang.y=ang.y+5
            If i= Chr(255) + "M"  Then ang.y=ang.y-5
            If i= Chr(255) + "P"  Then ang.x=ang.x-5
            If i= Chr(255) + "H"  Then ang.x=ang.x+5
            if i=" " then ang.x=0:ang.y=0
            Screenlock
            Cls
            Put(0,0),im
            Draw String (20,20), "Use up/down keys and mouse wheel, space to reset , esc to redo",Rgb(0,0,0)
            Draw String(20,50), "z= " & formula,Rgb(0,0,0)
            Draw String (.7*xres,.1*yres),"X Y plane = " & scale & " by " & scale,Rgb(0,0,0)
            Getmouse mx,my,mw
            _mw=mw/100
            combsort(e,zeds)
            dilation=startdilation+_mw
            counter=0
            
            'rotate axis
            For z As long=ub+1 To ub+6
                counter=counter+1
                axis(counter)=rotatepoint3d(piv,e(z),ang,dilation)
                axis(counter)=apply_perspective(axis(counter))
                
            Next z
            
            Line (axis(1).x,axis(1).y)-(axis(2).x,axis(2).y),Rgb(100,100,100)
            Draw String(axis(1).x,axis(1).y),"-X",Rgb(0,0,0)
            Draw String(axis(2).x,axis(2).y),"+X",Rgb(0,0,0)
            Line (axis(3).x,axis(3).y)-(axis(4).x,axis(4).y),Rgb(100,100,100)
            Draw String(axis(3).x,axis(3).y),"+Y",Rgb(0,0,0)
            Draw String(axis(4).x,axis(4).y),"-Y",Rgb(0,0,0)
            Line (axis(5).x,axis(5).y)-(axis(6).x,axis(6).y),Rgb(100,100,100)
            Draw String(axis(5).x,axis(5).y),"-Z",Rgb(0,0,0)
            Draw String(axis(6).x,axis(6).y),"+Z",Rgb(0,0,0)
            'rotate points and draw
            For z As long=1 To Ub
                zeds(z)=vertex(piv,e(z),ang,dilation,col(z))
            Next z
            Screenunlock
            Sleep 1,1
        Loop Until i=Chr(27)
        goto start
        imagedestroy im
        Sleep
       
         
Last edited by dodicat on Apr 14, 2019 15:13, edited 1 time in total.
Pitto
Posts: 122
Joined: Nov 19, 2012 19:58

Re: 3D Geometry , basics

Post by Pitto »

Hi Luxan,

it's nice to see an application of scanline algorithm I've translated in Freebasic from a C snippet some time ago (and Mr Swiss improved the translation).

I've used it in some 2d applications such FB LowPolyEditor, but seems to work fine also on 3d environment.

Good job so far

Best regards
Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

Thank you far your reply .
Sounds like you don't mind my using the scanline code.
Perhaps we might eventually make this part of the fbgfx library , or maybe a fbgfx3 library .

There is a limitation to what I'm doing with the painters algorithm , the graph can only be rotated
about the y axis by Pi/2 before the hidden surface feature becomes infeasible . This might be quite adequate for a fair few purposes .

Quite a few decades ago I'd seen a Borland C example where polygons were displayed and filled at a fast rate , most likely using the BGI library . This library wasn't really portable to BASIC.
Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

To dodicat .

I did a quick run of your program , it's quite good , especially your ability to parse a function to the
compiled program .

You're right about MyGraf5_3d.bas , and similar , running fast ; though just how fast is somewhat dependent upon your hardware.
I used a rather small increment in the outermost loop , as the whole loop was finishing too quickly
otherwise .

Some of the code in MyGraf5_3d.bas has yet to be utilized , especially the 2d arrays for storing data
or function values ; which might become quite large for small outer loop increments . That's if you were
to generate a single large array for all loops involved.
Presently the main array being utilized is

Const As Long POINTS = 4 , POLYGONS = 20, SCR_W = 740, SCR_H = 680

dim a(0 to POINTS,0 to 1) as long
.
Which , obviously , is quite small .

What type of function did you try for f1(x,y,fp) ; where fp is mostly for the purpose of animation .
Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

To Mr Swiss

Your participation is duly noted , you might want to post a comment .

To everyone

Now another minor update , in this instance I've introduced a palette to colour by height .
Also I've attempted to find the global maximum and minimum prior to graphing .
There's an ongoing difficulty with the positioning of the graph ; eventually I shall resolve this .

Code: Select all


'
' -----------------------------------------------------------------------------
'
'   My graf 3d
'
'   MyGraf6_3d.bas
'
'
'    (c) copyright 2019 , sciwiseg@gmail.com ,
'
'             Edward.Q.Montague.  [ alias]
'
'
'
'
'
' -----------------------------------------------------------------------------
'
type point
         x as single
         y as single
         z as single
         u as single '  possible extension for special coord system
end type
'
const Pi=4*atn(1)

' replaces the defines above (single line Macro's in FB)
Const As Long   POINTS = 4 , POLYGONS = 20, SCR_W = 740, SCR_H = 680

'
dim as single x1,y1,z1,x2,y2,z2
dim as integer i,j,k
'
' ----------------------------------------------------------------------------
'
'
'                  Looking at a cube .
'
'            cc -1,1 _______<_______  1,1    start         z = -1
'                   |               |             back face.
'                   |               |
'                   v               ^
'                   |               |
'                   |_______________|
'                -1,-1        >        1,-1
'               
'
' -----------------------------------------------------------------------------
'
declare function rotx(q as point,angx as single) as point
declare function roty1(q as point,angy as single) as point
declare function roty(q as point,angy as single) as point


declare function rotz(q as point,angz as single) as point
declare function tranx(q as point,movx as single) as point
declare function trany(q as point,movy as single) as point
declare function tranz(q as point,movz as single) as point
declare function persp(q as point,d as single) as point
'
declare function Trall( p1() as point,n as integer,edge() as integer, div as integer ,np as integer,ne as integer) as integer

declare sub gendata1(a() as long,n as integer,m as integer)
declare sub gendata2(a() as long,n as integer,m as integer)


declare sub genmatrix(gm() as single,n as integer,m as integer)
declare sub gmdata(gm() as single,n as integer,m as integer)
declare function f1(x as single , y as single,fp as single)as single
declare sub trallg(gm() as single,dm() as single , n as integer,m as integer)
declare sub trand(dm() as single,a() as long,n as integer,m as integer)
declare Sub fill_polygon(a() As Long, ByVal c As ULong)
declare sub outline_polygon(a() As Long, ByVal c As ULong)

declare sub tqxyz(x as single , y as single ,n as integer , m as integer , byref q as point)
declare sub w2scrn(p as point,n as integer , m as integer,byref u as long,byref v as long)

'
' ============================================================================
'
dim as integer np , ne,n,m
'
restore  storeA
read np
'
dim as point p1(1 to np)
'
restore store1
for i =1 to np
   read p1(i).x
   read p1(i).y
   read p1(i).z
next i
'
'
restore  storeB
read ne
'
dim as integer edge(1 to ne,0 to 1)
'
restore store2
for i = 1 to ne
   read edge(i,0)
   read edge(i,1)
next i
'
' -----------------------------------------------------------------------------
'
'screen 12
'window (-1.5,-1.5)-(1.5,1.5)
'line (-1.4,-1.4)-(1.4,1.4),11,b
'
'cls
'k=Trall( p1() ,3,edge() , 32 ,np ,ne )
'sleep
'INITIALIZING GRAPHICS _________________________________________________
ScreenRes(SCR_W, SCR_H, 24,2)      'initialize graphics
'window(10,10)-(210,110)
ScreenSet 1,0
Cls

n=520
m=520
'
dim a(0 to POINTS,0 to 1) as long
dim gm(0 to n , 0 to m,0 to 1) as single
dim dm(0 to n , 0 to m,0 to 1) as single
dim as ulong colour
'
gendata2(a() ,n ,m )
'
sleep
end

'genmatrix(gm() ,n ,m )
'gmdata(gm() ,n ,m )
'trallg(gm(),dm()  , n ,m )
'trand(dm() ,a() ,n ,m )


'sleep
'end
'
' ===================================
'
'    number of vertices
'
storeA:
data 8 
'
'    vertex data , easier to keep track of
'  data when we use multiple data statements.
'
store1:
data  1,1,1
data  -1,1,1
data  -1,-1,1
data 1,-1,1
data  1,1,-1
data  -1,1,-1
data  -1,-1,-1
data  1,-1,-1
'
'   Number of edges.
'
storeB:
data 12
'
'  edge data
'
store2:
data 1,2
data 1,4
data 1,5
data 2,3
data 2,6
data 3,4
data 3,7
data 4,8
data 5,6
data 5,8
data 6,7
data 7,8
'
' -------------------------------------------------------------------------------
'
function rotx(q as point,angx as single) as point
'
'              Rotate around x axis .
'
static as point p
'
             p.x = q.x
             p.y= q.y*cos(angx)-sin(angx)*q.z
             p.z= q.z*cos(angx)+sin(angx)*q.y
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function roty1(q as point,angy as single) as point
'
'              Rotate around y axis .
'
static as point p
'
            p.x = sin(angy)*q.z + cos(angy)*q.x
            p.y = q.y
            p.z = cos(angy)*q.z -sin(angy)*q.x
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function roty(q as point,angy as single) as point
'
'              Rotate around y axis .
'
static as point p
'
            p.x = sin(angy)*q.z + cos(angy)*q.x
            p.y = q.y
            p.z = cos(angy)*q.z -sin(angy)*q.x
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function rotz(q as point,angz as single) as point
'
'              Rotate around z axis .
'
static as point p
'
            p.x = sin(angz)*q.y + cos(angz)*q.x
            p.y = cos(angz)*q.y-sin(angz)*q.x
            p.z = q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function tranx(q as point,movx as single) as point
'
'             Translate point along x axis
'
static as point p
'
              p.x=q.x + movx
              p.y=q.y
              p.z=q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function trany(q as point,movy as single) as point
'
'              Translate point along y axis .
'
static as point p
'
              p.x=q.x
              p.y=q.y + movy
              p.z=q.z
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function tranz(q as point,movz as single) as point
'
'              Translate point along z axis .
'
static as point p
'
              p.x=q.x
              p.y=q.y
              p.z=q.z + movz
'
              return p
'
end function
'
' -----------------------------------------------------------------------------
'
function persp(q as point,d as single) as point
'
'              3d  perspective . 
'
'    Add 2.5 to the numerator when using any negative z value.
'
'     In this instance    -1 <= z  <= 1  ,  unit cube .
'
'    Therefore 2.5 is appropriate .
'
static as point p
'
     p.x = d*q.x/(q.z+2.25)
     p.y = d*q.y/(q.z+2.25)
     p.z = d
'
     return p
'
end function
'
' -----------------------------------------------------------------------------
'
function Trall( p1() as point,n as integer,edge() as integer, div as integer ,np as integer,ne as integer) as integer
'
'              Translate and rotate all vertices .
'     as an animation ,  for  n  cycles .
'
'    np number of points .
'    ne number of integers .
'
'
'   With  div number of angle divisions .
'
static as point p2(1 to np)
static as single theta,thi,x1,y1,z1,x2,y2,z2
static as integer i,j,k
static as integer i1,j1,k1
'
theta = Pi/div
'
for i=1 to n
  for j = 0 to div
  'cls
       thi = j*theta
   for k = 1 to np
     p2(k) = roty(p1(k),thi)
     p2(k)=persp(p2(k),0.8)
   next k   
 '  
    cls 
'
for i1 = 1 to ne
      j1 = edge(i1,0)
     k1 = edge(i1,1)
   x1 = p2(j1).x
   y1 = p2(j1).y
'   z1 = p2(j1).z   
   x2 = p2(k1).x
   y2 = p2(k1).y
'   z2 = p2(k1).z   
line(x1,y1)-(x2,y2),14
next i1     
'
sleep 100
  next j
next i
'
     return 0
'
end function
'
' --------------------------------------------------------------------------------
'
'
' ----------------------------------------------------------------
'
sub gendata1(a() as long,n as integer,m as integer)
'
'  Generate 3d data grid .
'
'  For perspective the range for x and y and z , is [-1,1]
'
'
'  Then a translation back to screen coordinates is required .  [-1,1] => [10,n] & [-1,1] => [10,m]  ?
'    x = (x + 1)*0.5*(n - 10) + 10
'    y = (y + 1)*0.5*(m - 10) + 10
'
'
'   w = f(x,y)
'   q.x = x
'   q.y = w
'   q.z = z
'
'   p = persp(q,d)
'
'   ( x', y',z , w') = persp( x , y , z , w )
'
'
static as integer i1,j1,k1,l1,qg,pg
static as single x,y,z,d
static as ulong colour,u,v
static as point p , q

d=0.8

'qg=10
'pg=10

pg=n/20
qg=m/10

colour=rgb(120,200,200)
for j1=10 to m step qg
    'z = 1 - 2*(j1-10)/(m-10) '  1  ->  -1 
    
  for i1=10 to n step pg
  
      k1=0
      x=i1
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
     ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
            
      k1=k1+1
      x=i1+pg
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
    ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1
      x=i1+pg
      y=j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )  
    ' f1(q.x,q.y) here ?
      p = persp(q ,d )   
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1     
      x = i1
      y = j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )   
    ' f1(q.x,q.y) here ?
      p = persp(q ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
      fill_polygon(a(), CULng(rnd*&hFFFFFF))
      outline_polygon(a() , colour)
 
   next i1
  '
next j1
'
end sub
'
' ----------------------------------------------------------------
'
sub gendata2(a() as long,n as integer,m as integer)
'
'  Generate 3d data grid .
'
'  For perspective the range for x and y and z , is [-1,1]
'
'
'  Then a translation back to screen coordinates is required .  [-1,1] => [10,n] & [-1,1] => [10,m]  ?
'    x = (x + 1)*0.5*(n - 10) + 10
'    y = (y + 1)*0.5*(m - 10) + 10
'
'
'   w = f(x,y)
'   q.x = x
'   q.y = w
'   q.z = z
'
'   p = persp(q,d)
'
'   ( x', y',z , w') = persp( x , y , z , w )
'
'
static as integer i1,j1,k1,l1,qg,pg
static as single x,y,z,d , theta , thi , fp ,max,min , Pd
static as ulong colour,u,v,chrome
static as point p , q , s , t
static as integer i,j,k
'
'
Dim pal(0 To 255) As Integer
'
j=64
for i=0 to 255 
    chrome=rgb(i,j,167)
    pal(i)=chrome
next i
'

'theta = Pi/9
theta=Pi/15
'theta=0
thi=Pi/7
'thi=0
d=0.98
d=1.2
d=2.4
d=2.8
d=3
'qg=10
'pg=10

pg=n/50
qg=m/50

pg=8
qg=8
'
' ................. max & min ..........................
'
max=-1000
min=1000
for  fp=0.75 to 4 step 0.0041
 for j1=10 to m step qg
  for i1=10 to n step pg
      x=i1
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
      y= f1(q.x,q.y,fp)  
    if ( y < min) then min=y
    if (y > max ) then max =y
  next i1
 next j1
next fp  
'
Pd=abs(max-min)
if (Pd=0) then Pd =1
'
' .....................................................
'
colour=rgb(120,200,200)
for  fp=0.75 to 4 step 0.0041
'
'cls
'
'  Alternately put a previously generated or loaded background image .
'
for j=0 to SCR_H
    k=int(200*j/SCR_H)
    colour=rgb(12,20,k)

   line(0,j)-(SCR_W,j),colour
next j 

'
colour=rgb(200,200,200)
'
for j1=10 to m step qg
    'z = 1 - 2*(j1-10)/(m-10) '  1  ->  -1 
    
  for i1=10 to n step pg
  
      k1=0
      x=i1
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y= f1(q.x,q.y,fp) 
     ' y=(exp(q.y))/(exp(1))
      y=(q.y-min)/Pd
      q.y=y-1
      'y=log(y*1000000+1)/log(1000001)
      k = int(y*255)
      
      t=roty(q ,theta )
      s=rotx(t ,thi )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
            
      k1=k1+1
      x=i1+pg
      y=j1
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y = f1(q.x,q.y,fp) 
      q.y=(q.y-min)/Pd-1
      
      t=roty(q ,theta )
      s=rotx(t ,thi )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1
      x=i1+pg
      y=j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )  
      q.y = f1(q.x,q.y,fp) 
      q.y=(q.y-min)/Pd-1
     
      t=roty(q ,theta )
      s=rotx(t ,thi )
      p = persp(s ,d )   
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

      k1=k1+1     
      x = i1
      y = j1+qg
      tqxyz(x  , y  ,n  , m  ,  q )   
      q.y = f1(q.x,q.y,fp) 
      q.y=(q.y-min)/Pd-1
      
      
      t=roty(q ,theta)
      s=rotx(t ,thi )
      p = persp(s ,d ) 
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
     ' fill_polygon(a(), CULng(rnd*&hFFFFFF))
      
     ' fill_polygon(a(), rgb(45,0,0))
     'chrome = rgb(chrome,20,20)
     fill_polygon(a(), pal(255-k))
     outline_polygon(a() , colour)

   next i1
  
   
  '
next j1
'
      Flip 1,0    'Copies our graph from page 1 to page 0


'
'sleep 10
'
next fp
'

color rgb(255,255,255)
locate 72,4 
print "  Done  "

      Flip 1,0    'Copies our graph from page 1 to page 0

'
end sub
'
' ----------------------------------------------------------------
'
sub genmatrix(gm() as single,n as integer,m as integer)
'
'   Generate a matrix with limits [x,-1,1][y,-1,1]
'
'
static as integer i,j
static as single x,y,z
'
for j=0 to m
    y=1-2*i/m
    z=y
    for i=0 to n
      x=-1+2*i/n
      gm(i,j,0)=x
      gm(i,j,1)=z
    next i
next j      
'
'
end sub
'
' -------------------------------------------------------------------
'
sub gmdata(gm() as single,n as integer,m as integer)
'
'  Generate data from function f1(x,y) .
'
'
static as integer i,j
static as single x,y,z,fp
'
fp=1
'
for j=0 to m
    for i=0 to n
      x=gm(i,j,0)
      y=gm(i,j,1)
      z=f1(x,y,fp)
      gm(i,j,0) = z
    next i
next j   
'
end sub
'
' ---------------------------------------------------------------------
'
sub trallg(gm() as single,dm() as single , n as integer,m as integer)
'
'   translate , rotate , apply perspective to all of gm()
'
'
static as integer i,j
static as single x,y,z,d
static as point p,q
'
d=0.8
'
for j=0 to m
     q.z=1-2*j/m
    for i=0 to n
      q.x=-1 +2*i/n
      q.y=gm(i,j,0)
      p = persp(q ,d )       
     dm(i,j,0) = p.x
     dm(i,j,1) = p.y     
    next i
next j   
'
'
end sub
'
' -------------------------------------------------------------------
'
sub trand(dm() as single,a() as long,n as integer,m as integer)
'
'  translate dm() to a() , also translate to screen coordinates for fill_polygon() routine .
'
'
static as integer i1,j1,k1,l1,qg,pg
static as long u , v
static as single x,y,z
static as point p , q
static as ulong colour

qg=10
pg=10



colour=rgb(20,120,20)
for j1=10 to m step qg
    
  for i1=10 to n step pg
  
      k1=0
      p.x=dm(i1,j1,0)
      p.y=dm(i1,j1,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
    
      k1=k1+1
      p.x=dm(i1+pg,j1,0)
      p.y=dm(i1,j1,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
      
      k1=k1+1  
      p.x=dm(i1+pg,j1,0)
      p.y=dm(i1,j1+qg,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v
       
      k1=k1+1
      p.x=dm(i1,j1,0)
      p.y=dm(i1,j1+qg,1)
      w2scrn(p ,n  , m ,u ,v ) 
      a(k1, 0) = u
      a(k1, 1) = v

   fill_polygon(a(), CULng(rnd*&hFFFFFF))
   'outline_polygon(a() , colour)
 
   next i1
  '
next j1
'
'
end sub
'
' -------------------------------------------------------------------
'
function f1(x as single , y as single,fp as single)as single
'
'       function to generate values upon [x,-1,1][y,-1,1]
'
'       [ ,-1,1] 
'
static as single r , z

        r=x*x+y*y
        if ( r > 0 ) then
             r = 2.5*sqr(r)*fp
             z=-(sin(r*Pi)/(r*Pi))+0.1*cos((x-0.5)*Pi*fp*0.5)*sin((y+1.5)*0.25*Pi*fp)
           '   z=2+0.1*cos((y-0.5)*Pi*fp)/(abs(x)+0.001)-sin((x+1.5)*0.5*Pi*fp)/(abs(y)+1.5)
            
             
        else
             z= 1
        end if          

       return z

'
end function
'
' -----------------------------------------------------------------------------------------------------
'
Sub fill_polygon(a() As Long, ByVal c As ULong)
   'translation of a c snippet by Angad
   'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
   Dim As Long      i, j, k, dy, dx, x, y, temp
   Dim As Long      xi(0 to Ubound(a, 1))
   Dim As Single    slope(0 to Ubound(a, 1))
   
   'join first and last vertex
   a(Ubound(a, 1), 0) = a(0, 0)
   a(Ubound(a, 1), 1) = a(0, 1)
'
 '  color c
'
   For i = 0 To Ubound(a, 1) - 1

      dy = a(i+1, 1) - a(i, 1)
      dx = a(i+1, 0) - a(i, 0)

      If (dy = 0) Then slope(i) = 1.0
      If (dx = 0) Then slope(i) = 0.0

      If (dy <> 0) AndAlso (dx <> 0) Then slope(i) = dx / dy
   Next i

   For y = 0 to SCR_H - 1
      k = 0
      ' using FB's short-cut operators (which C doesn't have!)
      For i = 0 to Ubound(a, 1) - 1
         If (a(i, 1) <= y AndAlso a(i+1, 1) > y) OrElse _
             (a(i, 1) > y AndAlso a(i+1, 1) <= y) Then
            xi(k) = CLng(a(i, 0) + slope(i) * (y - a(i, 1)))
            k += 1
         End If
      Next i

      For j = 0 to k - 2
         'Arrange x-intersections in order
         For i = 0 To k - 2
            If (xi(i) > xi(i + 1)) Then
               temp = xi(i)
               xi(i) = xi(i + 1)
               xi(i + 1) = temp
            End If
         Next i
      Next j
      'line filling
      For i = 0 To k - 2 Step 2
         Line (xi(i), y)-(xi(i + 1) + 1, y), c
      Next i
   Next y
End Sub
'
' -----------------------------------------------------------------------------
' 
sub outline_polygon(a() As Long, ByVal c As ULong)
'
'  Draw an outline for the polygon , in color c  .
'
  'translation of a translated  c snippet by Angad
   'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
   Dim As Long      i, j,  x, y, u , v , temp
      
   'join first and last vertex
   a(Ubound(a, 1), 0) = a(0, 0)
   a(Ubound(a, 1), 1) = a(0, 1)

   For i = 0 To Ubound(a, 1) - 1
       x=a(i,0)
       y=a(i,1)
       u=a(i+1,0)
       v=a(i+1,1)
       line(x,y)-(u,v),c
   Next i


end sub
'
' ----------------------------------------------------------------------
'
sub w2scrn( p as point,n as integer , m as integer, byref u as long, byref v as long)
'
' input x,y,n,m
' output u , v
'
     p.x =  (p.x + 1)*0.5*(n - 240) + 240  ' [-1,1] -> [10,n]
     p.y =  (p.y + 1)*0.5*(m - 100) + 100  ' [-1,1] -> [10,m]
     u = clng(p.x)
     v = clng(p.y)

'
end sub
'
' -------------------------------------------------------------
'
sub tqxyz(x as single , y as single ,n as integer , m as integer , byref q as point)
'
'  translate  to  x ~ i1 , y ~ j1 to q.x , q.y , q.z
'

      q.z = 1 - 2*(y-10)/(m-10)
      q.x = -1 + 2*(x-10)/(n-10)
      q.y = -1 + 2*(y-10)/(m-10)

'
end sub


dodicat
Posts: 7976
Joined: Jan 10, 2006 20:30
Location: Scotland

Re: 3D Geometry , basics

Post by dodicat »

Thank you Luxan.
This is non string parser, but experimenting with light shading with a given trig function.

Code: Select all



Type V3
    As Single x,y,z
End Type

Operator -(v1 As v3,v2 As v3) As v3 'v1-v2 
Return Type(v1.x-v2.x,v1.y-v2.y,v1.z-v2.z)
End Operator

Operator ^ (Byref v1 As v3,Byref v2 As v3) As v3 'cross product
Return Type(v1.y*v2.z-v2.y*v1.z,-(v1.x*v2.z-v2.x*v1.z),v1.x*v2.y-v2.x*v1.y)
End Operator

Type float As V3

Type box
    As v3 p(1 To 4)
    As Ulong c    'colour
    As Single z   '.z of a box
End Type

Type angle3D             'FLOATS for angles
    As Single sx,sy,sz
    As Single cx,cy,cz
    Declare Static Function construct(As Single,As Single,As Single) As Angle3D
End Type

Declare  Function InputFunction( As Double) As Double

Screenres 1024,768,32,2
Width 1024\8,768\16 'max dos font size
Color , Rgb(0,0,50)
'============ globals =============
Const pi=4*Atn(1)
Redim Shared As box b()
Redim Shared As box rot1()
Dim Shared As Angle3D A3d
Dim Shared As V3 CC       'grid centre
Dim Shared As Double df
Dim Shared As Integer xres,yres
Screeninfo xres,yres
'==================================
Function dot(v1 As v3,v2 As v3) As Single 
    Dim As Single d1=Sqr(v1.x*v1.x + v1.y*v1.y + v1.z*v1.z)
    Dim As Single d2=Sqr(v2.x*v2.x + v2.y*v2.y + v2.z*v2.z)
    Dim As Single v1x=v1.x/d1,v1y=v1.y/d1,v1z=v1.z/d1 'normalize
    Dim As Single v2x=v2.x/d2,v2y=v2.y/d2,v2z=v2.z/d2 'normalize
    Return v1x*v2x+v1y*v2y+v1z*v2z  'dot product
End Function

Function map(a As Single,b As Single,x As Single,c As Single,d As Single) As Single
    Return ((d)-(c))*((x)-(a))/((b)-(a))+(c)
End Function

Sub Qsortz(array() As box,begin As Long,Finish As Long)
    Dim As Long i=begin,j=finish
    Dim As box x =array(((I+J)\2))
    While I <= J
        While array(I).z > X .z:I+=1:Wend
            While array(J).z < X .z:J-=1:Wend
                If I<=J Then Swap array(I),array(J): I+=1:J-=1
            Wend
            If J >begin Then Qsortz(array(),begin,J)
            If I <Finish Then Qsortz(array(),I,Finish)
        End Sub
        
        Function Angle3D.construct(x As Single,y As Single,z As Single) As Angle3D
            Return   Type (Sin(x),Sin(y),Sin(z), _
            Cos(x),Cos(y),Cos(z))
        End Function
        
        Function Rotate(c As V3,p As V3,a As Angle3D,scale As float=Type(1,1,1)) As V3
            Dim As Single dx=p.x-c.x,dy=p.y-c.y,dz=p.z-c.z
            Return Type<V3>((scale.x)*((a.cy*a.cz)*dx+(-a.cx*a.sz+a.sx*a.sy*a.cz)*dy+(a.sx*a.sz+a.cx*a.sy*a.cz)*dz)+c.x,_
            (scale.y)*((a.cy*a.sz)*dx+(a.cx*a.cz+a.sx*a.sy*a.sz)*dy+(-a.sx*a.cz+a.cx*a.sy*a.sz)*dz)+c.y,_
            (scale.z)*((-a.sy)*dx+(a.sx*a.cy)*dy+(a.cx*a.cy)*dz)+c.z)
        End Function 
        
        Function perspective(p As V3,eyepoint As V3) As V3
            Dim As Single   w=1+(p.z/eyepoint.z)
            If w=0 Then w=1e-6
            Return Type<V3>((p.x-eyepoint.x)/w+eyepoint.x,_
            (p.y-eyepoint.y)/w+eyepoint.y,_
            (p.z-eyepoint.z)/w+eyepoint.z)
        End Function 
        
        Function setgrid(sx As Single,bx As Single,sy As Single,by As Single,st As Single,p() As box,fn As Function(x As Double) As Double) As v3
            #define U Ubound(p)
            #define dst(p1,p2) Sqr( (p1.x-p2.x)*(p1.x-p2.x) + (p1.y-p2.y)*(p1.y-p2.y) + (p1.z-p2.z)*(p1.z-p2.z) )   
            Redim p(0)
            Static As v3 centre
            Dim As Single cx,cy,ctr
            For y As Single=sy To by Step st
                For x As Single=sx To bx Step st
                    ctr+=1
                    cx+=x
                    cy+=y
                Next
            Next
            Static As Single q=15
            
            centre=Type(cx/ctr,cy/ctr)
            For y As Single=sy To by Step st
                For x As Single=sx To bx Step st
                    Redim Preserve p(1 To U+1)
                    p(u).p(1)=Type<v3>(x,y,      fn(dst( p(u).p(1),centre)))
                    p(u).p(2)=Type<v3>(x+st,y,   fn(dst( p(u).p(2),centre)))
                    p(u).p(3)=Type<v3>(x+st,y+st,fn(dst( p(u).p(3),centre)))
                    p(u).p(4)=Type<v3>(x,y+st,   fn(dst( p(u).p(4),centre)))
                    p(u).c=Rgb(x*q,x*q xor y*q,y*q ) 'pattern
                Next
            Next
            Return centre
        End Function
        
        Sub fill_polygon(a() As Long, Byval c As Ulong)
            'translation of a c snippet by Angad
            'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
            Dim As Long      i, j, k, dy, dx, x, y, temp
            Static As Long      xi(0 To Ubound(a, 1))
            Static As Single    slope(0 To Ubound(a, 1))
            
            'join first and last vertex
            a(Ubound(a, 1), 0) = a(0, 0)
            a(Ubound(a, 1), 1) = a(0, 1)
            
            For i = 0 To Ubound(a, 1) - 1
                
                dy = a(i+1, 1) - a(i, 1)
                dx = a(i+1, 0) - a(i, 0)
                
                If (dy = 0) Then slope(i) = 1.0
                If (dx = 0) Then slope(i) = 0.0
                
                If (dy <> 0) Andalso (dx <> 0) Then slope(i) = dx / dy
            Next i
            
            For y = 0 To yres - 1
                k = 0
                ' using FB's short-cut operators (which C doesn't have!)
                For i = 0 To Ubound(a, 1) - 1
                    If (a(i, 1) <= y Andalso a(i+1, 1) > y) Orelse _
                    (a(i, 1) > y Andalso a(i+1, 1) <= y) Then
                    xi(k) = Clng(a(i, 0) + slope(i) * (y - a(i, 1)))
                    k += 1
                End If
            Next i
            
            For j = 0 To k - 2
                'Arrange x-intersections in order
                For i = 0 To k - 2
                    If (xi(i) > xi(i + 1)) Then
                        'swap xi(i),xi(i+1)
                        temp = xi(i)
                        xi(i) = xi(i + 1)
                        xi(i + 1) = temp
                    End If
                Next i
            Next j
            'line filling
            For i = 0 To k - 2 Step 2
                Line (xi(i), y)-(xi(i + 1) + 1, y), c
            Next i
        Next y
    End Sub
    
    Sub convert(p() As v3,a() As Long) 'convert v3 to array of long
        Redim  a((Ubound(p)-Lbound(p))+1,1)
        Dim As Long ctr
        For n As Long=Lbound(p) To Ubound(p)
            a(ctr,0)=p(n).x
            a(ctr,1)=p(n).y
            ctr+=1
        Next n
    End Sub
    
    Sub drawboxes(b() As box)
        Static As Long a()
        For n As Long=Lbound(b) To Ubound(b)
            convert(b(n).p(),a())
            Var rd=Cast(Ubyte Ptr,@b(n).c)[2]
            Var gr=Cast(Ubyte Ptr,@b(n).c)[1]
            Var bl=Cast(Ubyte Ptr,@b(n).c)[0]
            Dim As v3 screencentre=(xres\2,yres\2)
            Var v1=b(n).p(2)-b(n).p(1)
            Var v2=b(n).p(3)-b(n).p(2)
            Var norm=v1^v2 'cross product
            Var dt=dot(norm,Type(1,.5,0))
            Var f=map(-1,1,dt,.2,1)
            fill_polygon(a(),Rgb(f*rd,f*gr,f*bl))
            
            #macro grid
            Line(b(n).p(1).x,b(n).p(1).y)-(b(n).p(2).x,b(n).p(2).y),Rgb(f*100,f*100,f*100)  
            Line(b(n).p(2).x,b(n).p(2).y)-(b(n).p(3).x,b(n).p(3).y),Rgb(f*100,f*100,f*100)  
            Line(b(n).p(3).x,b(n).p(3).y)-(b(n).p(4).x,b(n).p(4).y),Rgb(f*100,f*100,f*100)
            Line(b(n).p(4).x,b(n).p(4).y)-(b(n).p(1).x,b(n).p(1).y),Rgb(f*100,f*100,f*100)
            #endmacro
            'grid  'optional
        Next
    End Sub
    
    Function Regulate(Byval MyFps As Long,Byref fps As Long) As Long
        Static As Double timervalue,_lastsleeptime,t3,frames
        frames+=1
        If (Timer-t3)>=1 Then t3=Timer:fps=frames:frames=0
        Var sleeptime=_lastsleeptime+((1/myfps)-Timer+timervalue)*1000
        If sleeptime<1 Then sleeptime=1
        _lastsleeptime=sleeptime
        timervalue=Timer
        Return sleeptime
    End Function
    
    Sub setup(x1 As Single,x2 As Single,y1 As Single,y2 As Single,meshsize As Single)
        CC= setgrid(x1,x2,y1,y2,meshsize,b(),@InputFunction)'create grid, CC is the centre
        Redim rot1(Lbound(b) To Ubound(b))                   'working array
        A3d=angle3D.construct(0,-pi/2,0)  
        For n As Long=Lbound(b) To Ubound(b)
            For m As Long=1 To 4
                rot1(n).p(m)=rotate(CC,B(n).p(m),A3D,Type(10,10,10)) 'align boxes horizontally based
                rot1(n).c=B(n).c
                B(n).p(m)=rot1(n).p(m)
            Next m
        Next n
    End Sub
    
    Sub display()
        #define resetwheel(w,fl) fl=w
        #define wheel(w,f) w-f
        Screenset 1,0
        Dim As float ang=(0,-pi/7,pi/2)  'default
        Dim As Long fps
        Dim As String key
        Dim As Long mx,my,mw,mb,rflag
        Dim As Single sc=1
        
        Do
            df+=.03                   'push the wave out
            setup(485,515,385-50,415-50,.75) 'reset the grid 
            
            Getmouse mx,my,mw,mb
            If mb=2 Then 'reset
                ang.z=pi/2:ang.y=-pi/7
                resetwheel(mw,rflag) 
            End If
            mw=wheel(mw,rflag)
            If mx>0 Then sc=2+(mw/10)'scaler
            If sc<=0 Then sc=.1
            key=Inkey
            If key=Chr(255)+"K" Then ang.z-=.05     'left
            If key=Chr(255)+"M" Then ang.z+=.05     'right
            If key=Chr(255)+"P" Then ang.y-=.05     'down
            If key=Chr(255)+"H" Then ang.y+=.05     'up 
            
            ang.x+=.02  'the orbiting speed 
            
            A3D=Angle3D.construct(ang.x,ang.y,ang.z)      'set the rotate trigs
            For n As Long=Lbound(b) To Ubound(b)
                For m As Long=1 To 4
                    rot1(n).p(m) =rotate(CC,B(n).p(m),A3D,Type(sc,sc,sc))
                    rot1(n).p(m) =perspective(rot1(n).p(m),Type(cc.x,cc.y,400*sc))'eyepoint
                    If mb=1 Then rot1(n).p(m).x-=cc.x-mx: rot1(n).p(m).y-=cc.y-my'follow the mouse
                Next m
                rot1(n).z=(rot1(n).p(1).z+rot1(n).p(3).z)/2
            Next n
            
            qsortz(rot1(),Lbound(rot1), Ubound(rot1))
            
            Cls
            Draw String(50,50),"Framerate "&fps
            Draw String(50,150),"Use the arrow keys and wheel"
            Draw String(50,250),"Right mouse click to reset"
            drawboxes(rot1())
            Flip
            Sleep regulate(35,fps),1
        Loop Until key=Chr(27)
    End Sub
    Function InputFunction(x As Double) As Double
        Return (Sin(x/4-df)+Cos(4*(x/4-df))-Sin(6*(x/8-df)))
    End Function
    
    display()
    
    Sleep
    
    
    
    
    
      
The 32 bit compiler is MUCH better with gfx graphics than the 64 bit.
Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

Well that's a nice thing to wake up to .

The shading gives the impression of a light source , your hidden surface method is just what's
needed .

In this instance did you set the function in :

Code: Select all


    Function InputFunction(x As Double) As Double
        Return (Sin(x/4-df)+Cos(4*(x/4-df))-Sin(6*(x/8-df)))
    End Function

How possible is it to put your shading in a more height dependent method , as standard when viewing
scientific or mathematical graphs .
UEZ
Posts: 972
Joined: May 05, 2017 19:59
Location: Germany

Re: 3D Geometry , basics

Post by UEZ »

dodicat wrote:Thank you Luxan.
This is non string parser, but experimenting with light shading with a given trig function.

Code: Select all



Type V3
    As Single x,y,z
End Type

Operator -(v1 As v3,v2 As v3) As v3 'v1-v2 
Return Type(v1.x-v2.x,v1.y-v2.y,v1.z-v2.z)
End Operator

Operator ^ (Byref v1 As v3,Byref v2 As v3) As v3 'cross product
Return Type(v1.y*v2.z-v2.y*v1.z,-(v1.x*v2.z-v2.x*v1.z),v1.x*v2.y-v2.x*v1.y)
End Operator

Type float As V3

Type box
    As v3 p(1 To 4)
    As Ulong c    'colour
    As Single z   '.z of a box
End Type

Type angle3D             'FLOATS for angles
    As Single sx,sy,sz
    As Single cx,cy,cz
    Declare Static Function construct(As Single,As Single,As Single) As Angle3D
End Type

Declare  Function InputFunction( As Double) As Double

Screenres 1024,768,32,2
Width 1024\8,768\16 'max dos font size
Color , Rgb(0,0,50)
'============ globals =============
Const pi=4*Atn(1)
Redim Shared As box b()
Redim Shared As box rot1()
Dim Shared As Angle3D A3d
Dim Shared As V3 CC       'grid centre
Dim Shared As Double df
Dim Shared As Integer xres,yres
Screeninfo xres,yres
'==================================
Function dot(v1 As v3,v2 As v3) As Single 
    Dim As Single d1=Sqr(v1.x*v1.x + v1.y*v1.y + v1.z*v1.z)
    Dim As Single d2=Sqr(v2.x*v2.x + v2.y*v2.y + v2.z*v2.z)
    Dim As Single v1x=v1.x/d1,v1y=v1.y/d1,v1z=v1.z/d1 'normalize
    Dim As Single v2x=v2.x/d2,v2y=v2.y/d2,v2z=v2.z/d2 'normalize
    Return v1x*v2x+v1y*v2y+v1z*v2z  'dot product
End Function

Function map(a As Single,b As Single,x As Single,c As Single,d As Single) As Single
    Return ((d)-(c))*((x)-(a))/((b)-(a))+(c)
End Function

Sub Qsortz(array() As box,begin As Long,Finish As Long)
    Dim As Long i=begin,j=finish
    Dim As box x =array(((I+J)\2))
    While I <= J
        While array(I).z > X .z:I+=1:Wend
            While array(J).z < X .z:J-=1:Wend
                If I<=J Then Swap array(I),array(J): I+=1:J-=1
            Wend
            If J >begin Then Qsortz(array(),begin,J)
            If I <Finish Then Qsortz(array(),I,Finish)
        End Sub
        
        Function Angle3D.construct(x As Single,y As Single,z As Single) As Angle3D
            Return   Type (Sin(x),Sin(y),Sin(z), _
            Cos(x),Cos(y),Cos(z))
        End Function
        
        Function Rotate(c As V3,p As V3,a As Angle3D,scale As float=Type(1,1,1)) As V3
            Dim As Single dx=p.x-c.x,dy=p.y-c.y,dz=p.z-c.z
            Return Type<V3>((scale.x)*((a.cy*a.cz)*dx+(-a.cx*a.sz+a.sx*a.sy*a.cz)*dy+(a.sx*a.sz+a.cx*a.sy*a.cz)*dz)+c.x,_
            (scale.y)*((a.cy*a.sz)*dx+(a.cx*a.cz+a.sx*a.sy*a.sz)*dy+(-a.sx*a.cz+a.cx*a.sy*a.sz)*dz)+c.y,_
            (scale.z)*((-a.sy)*dx+(a.sx*a.cy)*dy+(a.cx*a.cy)*dz)+c.z)
        End Function 
        
        Function perspective(p As V3,eyepoint As V3) As V3
            Dim As Single   w=1+(p.z/eyepoint.z)
            If w=0 Then w=1e-6
            Return Type<V3>((p.x-eyepoint.x)/w+eyepoint.x,_
            (p.y-eyepoint.y)/w+eyepoint.y,_
            (p.z-eyepoint.z)/w+eyepoint.z)
        End Function 
        
        Function setgrid(sx As Single,bx As Single,sy As Single,by As Single,st As Single,p() As box,fn As Function(x As Double) As Double) As v3
            #define U Ubound(p)
            #define dst(p1,p2) Sqr( (p1.x-p2.x)*(p1.x-p2.x) + (p1.y-p2.y)*(p1.y-p2.y) + (p1.z-p2.z)*(p1.z-p2.z) )   
            Redim p(0)
            Static As v3 centre
            Dim As Single cx,cy,ctr
            For y As Single=sy To by Step st
                For x As Single=sx To bx Step st
                    ctr+=1
                    cx+=x
                    cy+=y
                Next
            Next
            Static As Single q=15
            
            centre=Type(cx/ctr,cy/ctr)
            For y As Single=sy To by Step st
                For x As Single=sx To bx Step st
                    Redim Preserve p(1 To U+1)
                    p(u).p(1)=Type<v3>(x,y,      fn(dst( p(u).p(1),centre)))
                    p(u).p(2)=Type<v3>(x+st,y,   fn(dst( p(u).p(2),centre)))
                    p(u).p(3)=Type<v3>(x+st,y+st,fn(dst( p(u).p(3),centre)))
                    p(u).p(4)=Type<v3>(x,y+st,   fn(dst( p(u).p(4),centre)))
                    p(u).c=Rgb(x*q,x*q xor y*q,y*q ) 'pattern
                Next
            Next
            Return centre
        End Function
        
        Sub fill_polygon(a() As Long, Byval c As Ulong)
            'translation of a c snippet by Angad
            'source of c code: http://code-heaven.blogspot.it/2009/10/simple-c-program-for-scan-line-polygon.html
            Dim As Long      i, j, k, dy, dx, x, y, temp
            Static As Long      xi(0 To Ubound(a, 1))
            Static As Single    slope(0 To Ubound(a, 1))
            
            'join first and last vertex
            a(Ubound(a, 1), 0) = a(0, 0)
            a(Ubound(a, 1), 1) = a(0, 1)
            
            For i = 0 To Ubound(a, 1) - 1
                
                dy = a(i+1, 1) - a(i, 1)
                dx = a(i+1, 0) - a(i, 0)
                
                If (dy = 0) Then slope(i) = 1.0
                If (dx = 0) Then slope(i) = 0.0
                
                If (dy <> 0) Andalso (dx <> 0) Then slope(i) = dx / dy
            Next i
            
            For y = 0 To yres - 1
                k = 0
                ' using FB's short-cut operators (which C doesn't have!)
                For i = 0 To Ubound(a, 1) - 1
                    If (a(i, 1) <= y Andalso a(i+1, 1) > y) Orelse _
                    (a(i, 1) > y Andalso a(i+1, 1) <= y) Then
                    xi(k) = Clng(a(i, 0) + slope(i) * (y - a(i, 1)))
                    k += 1
                End If
            Next i
            
            For j = 0 To k - 2
                'Arrange x-intersections in order
                For i = 0 To k - 2
                    If (xi(i) > xi(i + 1)) Then
                        'swap xi(i),xi(i+1)
                        temp = xi(i)
                        xi(i) = xi(i + 1)
                        xi(i + 1) = temp
                    End If
                Next i
            Next j
            'line filling
            For i = 0 To k - 2 Step 2
                Line (xi(i), y)-(xi(i + 1) + 1, y), c
            Next i
        Next y
    End Sub
    
    Sub convert(p() As v3,a() As Long) 'convert v3 to array of long
        Redim  a((Ubound(p)-Lbound(p))+1,1)
        Dim As Long ctr
        For n As Long=Lbound(p) To Ubound(p)
            a(ctr,0)=p(n).x
            a(ctr,1)=p(n).y
            ctr+=1
        Next n
    End Sub
    
    Sub drawboxes(b() As box)
        Static As Long a()
        For n As Long=Lbound(b) To Ubound(b)
            convert(b(n).p(),a())
            Var rd=Cast(Ubyte Ptr,@b(n).c)[2]
            Var gr=Cast(Ubyte Ptr,@b(n).c)[1]
            Var bl=Cast(Ubyte Ptr,@b(n).c)[0]
            Dim As v3 screencentre=(xres\2,yres\2)
            Var v1=b(n).p(2)-b(n).p(1)
            Var v2=b(n).p(3)-b(n).p(2)
            Var norm=v1^v2 'cross product
            Var dt=dot(norm,Type(1,.5,0))
            Var f=map(-1,1,dt,.2,1)
            fill_polygon(a(),Rgb(f*rd,f*gr,f*bl))
            
            #macro grid
            Line(b(n).p(1).x,b(n).p(1).y)-(b(n).p(2).x,b(n).p(2).y),Rgb(f*100,f*100,f*100)  
            Line(b(n).p(2).x,b(n).p(2).y)-(b(n).p(3).x,b(n).p(3).y),Rgb(f*100,f*100,f*100)  
            Line(b(n).p(3).x,b(n).p(3).y)-(b(n).p(4).x,b(n).p(4).y),Rgb(f*100,f*100,f*100)
            Line(b(n).p(4).x,b(n).p(4).y)-(b(n).p(1).x,b(n).p(1).y),Rgb(f*100,f*100,f*100)
            #endmacro
            'grid  'optional
        Next
    End Sub
    
    Function Regulate(Byval MyFps As Long,Byref fps As Long) As Long
        Static As Double timervalue,_lastsleeptime,t3,frames
        frames+=1
        If (Timer-t3)>=1 Then t3=Timer:fps=frames:frames=0
        Var sleeptime=_lastsleeptime+((1/myfps)-Timer+timervalue)*1000
        If sleeptime<1 Then sleeptime=1
        _lastsleeptime=sleeptime
        timervalue=Timer
        Return sleeptime
    End Function
    
    Sub setup(x1 As Single,x2 As Single,y1 As Single,y2 As Single,meshsize As Single)
        CC= setgrid(x1,x2,y1,y2,meshsize,b(),@InputFunction)'create grid, CC is the centre
        Redim rot1(Lbound(b) To Ubound(b))                   'working array
        A3d=angle3D.construct(0,-pi/2,0)  
        For n As Long=Lbound(b) To Ubound(b)
            For m As Long=1 To 4
                rot1(n).p(m)=rotate(CC,B(n).p(m),A3D,Type(10,10,10)) 'align boxes horizontally based
                rot1(n).c=B(n).c
                B(n).p(m)=rot1(n).p(m)
            Next m
        Next n
    End Sub
    
    Sub display()
        #define resetwheel(w,fl) fl=w
        #define wheel(w,f) w-f
        Screenset 1,0
        Dim As float ang=(0,-pi/7,pi/2)  'default
        Dim As Long fps
        Dim As String key
        Dim As Long mx,my,mw,mb,rflag
        Dim As Single sc=1
        
        Do
            df+=.03                   'push the wave out
            setup(485,515,385-50,415-50,.75) 'reset the grid 
            
            Getmouse mx,my,mw,mb
            If mb=2 Then 'reset
                ang.z=pi/2:ang.y=-pi/7
                resetwheel(mw,rflag) 
            End If
            mw=wheel(mw,rflag)
            If mx>0 Then sc=2+(mw/10)'scaler
            If sc<=0 Then sc=.1
            key=Inkey
            If key=Chr(255)+"K" Then ang.z-=.05     'left
            If key=Chr(255)+"M" Then ang.z+=.05     'right
            If key=Chr(255)+"P" Then ang.y-=.05     'down
            If key=Chr(255)+"H" Then ang.y+=.05     'up 
            
            ang.x+=.02  'the orbiting speed 
            
            A3D=Angle3D.construct(ang.x,ang.y,ang.z)      'set the rotate trigs
            For n As Long=Lbound(b) To Ubound(b)
                For m As Long=1 To 4
                    rot1(n).p(m) =rotate(CC,B(n).p(m),A3D,Type(sc,sc,sc))
                    rot1(n).p(m) =perspective(rot1(n).p(m),Type(cc.x,cc.y,400*sc))'eyepoint
                    If mb=1 Then rot1(n).p(m).x-=cc.x-mx: rot1(n).p(m).y-=cc.y-my'follow the mouse
                Next m
                rot1(n).z=(rot1(n).p(1).z+rot1(n).p(3).z)/2
            Next n
            
            qsortz(rot1(),Lbound(rot1), Ubound(rot1))
            
            Cls
            Draw String(50,50),"Framerate "&fps
            Draw String(50,150),"Use the arrow keys and wheel"
            Draw String(50,250),"Right mouse click to reset"
            drawboxes(rot1())
            Flip
            Sleep regulate(35,fps),1
        Loop Until key=Chr(27)
    End Sub
    Function InputFunction(x As Double) As Double
        Return (Sin(x/4-df)+Cos(4*(x/4-df))-Sin(6*(x/8-df)))
    End Function
    
    display()
    
    Sleep
    
    
    
    
    
      
The 32 bit compiler is MUCH better with gfx graphics than the 64 bit.
Very nice example which somehow looks familiar to me. ;-)

Well, the x64 compilation runs at least on my system much faster than the x86 counterpart when I use one of these settings:
-gen gcc -Wc -O2
-gen gcc -O 3
-gen gcc -Wc -Ofast -fpmode FAST -fpu SSE

The default compiler option either x86 or x64 are much slower.

On my older i5 4300U CPU it runs with one of the compile options above ~46 fps. With default compiler setting ~15 fps.

The gfx display with x64 is fluent and smooth.


@Luxan: thanks for your contribution. Might be useful for me... ^^
dodicat
Posts: 7976
Joined: Jan 10, 2006 20:30
Location: Scotland

Re: 3D Geometry , basics

Post by dodicat »

Yea UEZ, it is that sine wave thing coloured in.
Luxan
lines 186 to 189 for the shading.
First take two adjacent sides of a four sided polygon and cross product them to get a normal line to the face.
Then dot product this with Type(1,.5,0), which represents a light source point if you imagine a tiny grid at the origin, but you can make this any other vector line, perhaps give the z component a value.
z is +ve into your monitor.
The function is InputFunction()
If you try maybe sin(x)/x type of thing, then you should tend to maths singularity point yourself (at zero), because the polygons are sorted by the z value of their centre, and quicksort doesn't like -1.#IND due to division by zero
something like
if x=0 then x=1
Return 7*sin(x-df)/x
should do.


The useless 64 bit is probably only on this computer.
Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

Been running in RAM , save file difficulties , installed to other USB just setting up system again.
I'm now using the 1.05 version .

The next consideration , after alignment , was shading ; your explanation is very useful.
Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

To dodicat ,

I altered some code in the region you mentioned.
Your shading is quite apparent now and conveys a considerable amount of information , more perhaps
than a palette of distinct values.

Code: Select all


    Sub drawboxes(b() As box)
        Static As Long a()
        static as ulong colour
        colour=rgb(34 , 230  ,240)
         For n As Long=Lbound(b) To Ubound(b)
            convert(b(n).p(),a())
            Var rd=Cast(Ubyte Ptr,@b(n).c)[2] ' might contain checker pattern .
            Var gr=Cast(Ubyte Ptr,@b(n).c)[1]
            Var bl=Cast(Ubyte Ptr,@b(n).c)[0]
            Dim As v3 screencentre=(xres\2,yres\2)
            Var v1=b(n).p(2)-b(n).p(1)
            Var v2=b(n).p(3)-b(n).p(2)
            Var norm=v1^v2 'cross product
            Var dt=dot(norm,Type(1,.5,0))
            Var f=map(-1,1,dt,.2,1) ' shading function ?
            
          '  fill_polygon(a(),Rgb(f*rd,f*gr,f*bl))
            fill_polygon(a(),Rgb(f*32,f*142,f*200l))
        '   outline_polygon(a() , colour )
           
            #macro grid
            Line(b(n).p(1).x,b(n).p(1).y)-(b(n).p(2).x,b(n).p(2).y),Rgb(f*100,f*100,f*100) 
            Line(b(n).p(2).x,b(n).p(2).y)-(b(n).p(3).x,b(n).p(3).y),Rgb(f*100,f*100,f*100) 
            Line(b(n).p(3).x,b(n).p(3).y)-(b(n).p(4).x,b(n).p(4).y),Rgb(f*100,f*100,f*100)
            Line(b(n).p(4).x,b(n).p(4).y)-(b(n).p(1).x,b(n).p(1).y),Rgb(f*100,f*100,f*100)
            #endmacro
           grid  'optional
        Next
    End Sub


I also investigated the generation of palette values , based upon the human response to light .
Translated from fortran , originally from http://www.physics.sfasu.edu/astro/color.html .

Code: Select all

 sub gen_spectra(pal() as ulong,M as integer)
'
'   Generate spectral data .
' 
static as integer MAX , I
static as single GAMMA ,WL,R,G,B,HST
 '     
'       M=1024  '  I  limit
      redim pal(0 to M)   
 '     
       MAX=255
       GAMMA=.80
'
I=1
 do while (I<M)
'
'         WAVELENGTH = WL
'
            WL = 380. + csng(I * 400. / M)

            IF ((WL>=380.) AND (WL<=440.)) THEN 
              R = -1.*(WL-440.)/(440.-380.)
              G = 0.
              B = 1.
            ENDIF
            IF ((WL>=440.) AND (WL<=490.)) THEN
              R = 0.
              G = (WL-440.)/(490.-440.)
              B = 1.
            ENDIF
            IF ((WL>=490.) AND (WL<=510.)) THEN 
              R = 0.
              G = 1.
              B = -1.*(WL-510.)/(510.-490.)
            ENDIF
            IF ((WL>=510.) AND (WL<=580.)) THEN 
              R = (WL-510.)/(580.-510.)
              G = 1.
              B = 0.
            ENDIF
            IF ((WL>=580.) AND (WL<=645.)) THEN
              R = 1.
              G = -1.*(WL-645.)/(645.-580.)
              B = 0.
            ENDIF
            IF ((WL>=645.) AND (WL<=780.)) THEN
              R = 1.
              G = 0.
              B = 0.
            ENDIF
'
'      LET THE INTENSITY HST FALL OFF NEAR THE VISION LIMITS
'
         IF (WL>700.) THEN
            HST=.3+.7* (780.-WL)/(780.-700.)
         ELSEIF (WL<420.) THEN
            HST=.3+.7*(WL-380.)/(420.-380.)
         ELSE
            HST=1.
         ENDIF
'
'      GAMMA ADJUST AND WRITE  TO AN ARRAY
'
         R=(HST*R)^GAMMA
         G=(HST*G)^GAMMA
         B=(HST*B)^GAMMA
         pal(I-1)=rgb(255*R,255*G,255*B)
         I=I+1
      loop 
'
' 
 end sub



Luxan
Posts: 222
Joined: Feb 18, 2009 12:47
Location: New Zealand

Re: 3D Geometry , basics

Post by Luxan »

Is there something wrong with the rgba() alpha ?

Likely not, this code however doesn't display alpha adjusted
values.

Code: Select all

 
'      Hu_color2.bas

'       palette values , based upon the human response to light 
 
declare function pal_spectra(I as integer, M as integer, ap as ubyte) as ulong

'
' ----------------------------------------------------------------------
'
screen 14, 32    '  32 bit  color  depth  .
cls
view (20,20)-(620,320)

'  reccomend , M > = 400, to avoid artifacts .

dim as integer M,i,j
M = 400
window(-10,-10)-(M+10,110) 
line(-10,-10)-(M+10,110),11,b
line(0,0)-(M,100),rgb(255,255,255),b
sleep
dim as ulong pal(0 to M), pxc
'
dim as ubyte a
a = 32
'
 for j = 0 to 100
     a = cubyte(j)+32
   for i = 1 to M    
   pxc = pal_spectra(i, M , a )
   'pset(i,j),pxc
   line(i-1,j-1)-(i-1,j), pxc
 next i
next j

line(0,0)-(M,100),rgb(255,255,255),b

sleep

end
'
' ----------------------------------------------------------------------
'
a=32

for i=1 to M
 j=int(255*i/M)
 pxc = pal_spectra(i, M , a )
 line(i-1,0)-(i-1,101), pxc
next i

line(0,0)-(M,100),rgb(255,255,255),b

sleep

end
'
' ======================================================================
' 
 function pal_spectra(I as integer, M as integer, ap as ubyte) as ulong
 '
 '        Generate a particular palette value with alpha a
 '  
static as integer MAX 
static as single GAMMA ,WL,R,G,B,HST
static as ulong pxc
 '     
'       M=1024  '  I  limit
         
 '     
       MAX=255
       GAMMA=.80
'
'         WAVELENGTH = WL
'
            WL = 380. + csng(I * 400. / M)

            IF ((WL>=380.) AND (WL<=440.)) THEN 
              R = -1.*(WL-440.)/(440.-380.)
              G = 0.
              B = 1.
            ENDIF
            IF ((WL>=440.) AND (WL<=490.)) THEN
              R = 0.
              G = (WL-440.)/(490.-440.)
              B = 1.
            ENDIF
            IF ((WL>=490.) AND (WL<=510.)) THEN 
              R = 0.
              G = 1.
              B = -1.*(WL-510.)/(510.-490.)
            ENDIF
            IF ((WL>=510.) AND (WL<=580.)) THEN 
              R = (WL-510.)/(580.-510.)
              G = 1.
              B = 0.
            ENDIF
            IF ((WL>=580.) AND (WL<=645.)) THEN
              R = 1.
              G = -1.*(WL-645.)/(645.-580.)
              B = 0.
            ENDIF
            IF ((WL>=645.) AND (WL<=780.)) THEN
              R = 1.
              G = 0.
              B = 0.
            ENDIF
'
'      LET THE INTENSITY HST FALL OFF NEAR THE VISION LIMITS
'
         IF (WL>700.) THEN
            HST=.3+.7* (780.-WL)/(780.-700.)
         ELSEIF (WL<420.) THEN
            HST=.3+.7*(WL-380.)/(420.-380.)
         ELSE
            HST=1.
         ENDIF
'
'      GAMMA ADJUST AND WRITE  TO A ulong variable .
'
         R=(HST*R)^GAMMA
         G=(HST*G)^GAMMA
         B=(HST*B)^GAMMA
         pxc=RGBA(255*R,255*G,255*B,ap)
' 
        return pxc
' 
 end function
'
' ----------------------------------------------------------------------
' 

Post Reply