## 3d without openGL

User projects written in or related to FreeBASIC.
BasicCoder2
Posts: 3566
Joined: Jan 01, 2009 7:03
Location: Australia

### Re: 3d without openGL

@dodicat,
Yes your 3d demos look great but I can't visualize or understand how your compact code works.

May I suggest you add keyboard control over the rotation and translation of the object? It is only then that it would be visually clear that the computations allow rotation around an arbitrary axis. In my wire frame rocket example I show the absolute world coordinates and also the rocket coordinates. The problem was to rotate around any one of the rocket's coordinates.
.
dafhi
Posts: 1345
Joined: Jun 04, 2005 9:51

### Re: 3d without openGL

@paul doe - Very intrigued and inspired by your offerings
dodicat
Posts: 6555
Joined: Jan 10, 2006 20:30
Location: Scotland

### Re: 3d without openGL

OK basiccoder2.
I use only the necessary for a Saturn type ring of stones.
Use the arrow keys.

Code: Select all

`Type V3    As Single x,y,z    As Uinteger colEnd TypeType 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 Angle3DEnd TypeType float As V3#define incircle(cx,cy,radius,x,y) (cx-x)*(cx-x) +(cy-y)*(cy-y)<= radius*radiusFunction 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 FunctionFunction 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,p.col)End Function Function perspective(p As V3,eyepoint As V3) As V3    Dim As Single   w=1+(p.z/eyepoint.z)    Return Type<V3>((p.x-eyepoint.x)/w+eyepoint.x,_    (p.y-eyepoint.y)/w+eyepoint.y,_    (p.z-eyepoint.z)/w+eyepoint.z,p.col)End Function  Sub QsortZ(array() As V3,begin As Long,Finish As Ulong)    Dim As Long i=begin,j=finish    Dim As V3 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 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 sleeptimeEnd Function        screen 19,32        const pi=4*atn(1)        redim as v3 a(0)                'set some points in a saturn type disc        for n as long=0 to 100000           var xp=rnd*800,yp=rnd*600           if incircle(400,300,150,xp,yp)=0 then              if incircle(400,300,250,xp,yp) then                  var u=ubound(a)                  redim preserve a(1 to u+1)                  a(ubound(a))=type<v3>(xp,yp,rnd*10-rnd*10,rnd*rgb(255,255,255))                  pset(xp,yp),a(ubound(a)).col              end if              end if          next          print "the points are on the x y plane with a small random z component"          print "flip them all 90 degrees around y axis and expand a bit"          print "Press a key"          sleep          cls                    redim as V3 rot(lbound(a) to ubound(a))'to hold all rotated points                    dim as Angle3D A3d=angle3D.construct(0,pi/2,0)   'flip all points by pi/2 on y axis                    for n as long=lbound(a) to ubound(a)         rot(n)=rotate(type(400,300,0),a(n),A3D,type(1.4,1.4,1.4)) '1.4 is a scale up         a(n)=rot(n)         pset(a(n).x,a(n).y),a(n).col     next     print "Thats 'em all flipped without perspctive"      print "Press a key"          sleep          dim as float ang          dim as long fps          'fix these two components at start          ang.z=pi/2 'tilt angle sideways          ang.y=-pi/7 'tilt angle back and forward      dim as string key              do              key=inkey              If key=Chr(255)+"K" Then ang.z-=.02     'left              If key=Chr(255)+"M" Then ang.z+=.02     'right              If key=Chr(255)+"P" Then ang.y-=.02     'down              If key=Chr(255)+"H" Then ang.y+=.02     'up                             ang.x+=.01  'the orbiting speed                             A3D=Angle3D.construct(ang.x,ang.y,ang.z)'get the six rotate components (sines, coses  ...)              screenlock              cls              draw string(50,50),"Framerate "&fps               draw string(50,150),"Use the arrow keys"              for n as long=lbound(a) to ubound(a)                rot(n)=rotate(type(400,300,0),a(n),A3D,type(1,1,1))                rot(n)=perspective(rot(n),type(400,300,1000))              next              qsortZ(rot(),lbound(rot),ubound(rot))              for n as long=lbound(rot) to ubound(rot)                pset(rot(n).x,rot(n).y),rot(n).col              next              screenunlock              sleep regulate(64,fps),1              loop until key=chr(27)                                `
BasicCoder2
Posts: 3566
Joined: Jan 01, 2009 7:03
Location: Australia

### Re: 3d without openGL

Demonstration without enough explanation for my programming limitations.
What I tried to do was modify your code to produce the three axes of the object (without the object). That is the axes were the object. I tried to make the z blue the x red and the y green but they appear to be random colors. I wanted to be able to rotate around any of those three axes.

Code: Select all

`Type V3    As Single x,y,z    As Uinteger colEnd TypeType 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 Angle3DEnd TypeType float As V3#define incircle(cx,cy,radius,x,y) (cx-x)*(cx-x) +(cy-y)*(cy-y)<= radius*radiusFunction 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 FunctionFunction 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,p.col)End Function Function perspective(p As V3,eyepoint As V3) As V3    Dim As Single   w=1+(p.z/eyepoint.z)    Return Type<V3>((p.x-eyepoint.x)/w+eyepoint.x,_    (p.y-eyepoint.y)/w+eyepoint.y,_    (p.z-eyepoint.z)/w+eyepoint.z,p.col)End Function  Sub QsortZ(array() As V3,begin As Long,Finish As Ulong)    Dim As Long i=begin,j=finish    Dim As V3 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 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 sleeptimeEnd Functionscreenres 1280,600,32const pi=4*atn(1)redim as v3 a(0)                'set some pointsfor n as single =-250 to 250    var xp=0,yp=0,zp = n    var u=ubound(a)    redim preserve a(1 to u+1)    a(ubound(a))=type<v3>(xp, yp, zp, rnd*rgb(0,0,255))    pset(xp,yp),a(ubound(a)).colnext nfor n as single =-250 to 250    var xp=0,yp=n,zp = 0    var u=ubound(a)    redim preserve a(1 to u+1)    a(ubound(a))=type<v3>(xp, yp, zp, rnd*rgb(255,0,0))    pset(xp,yp),a(ubound(a)).colnext nfor n as single =-250 to 250    var xp=n,yp=0,zp = 0     var u=ubound(a)    redim preserve a(1 to u+1)    a(ubound(a))=type<v3>(xp, yp, zp, rnd*rgb(0,255,0))    pset(xp,yp),a(ubound(a)).colnext n               'print "the points are on the x y plane with a small random z component"'print "flip them all 90 degrees around y axis and expand a bit"'print "Press a key"'sleep'cls          redim as V3 rot(lbound(a) to ubound(a))'to hold all rotated points          dim as Angle3D A3d=angle3D.construct(0,pi/2,0)   'flip all points by pi/2 on y axis          for n as long=lbound(a) to ubound(a)    rot(n)=rotate(type(400,300,0),a(n),A3D,type(0.5,0.5,0.5)) '1.4 is a scale up    a(n)=rot(n)    pset(a(n).x,a(n).y),a(n).colnext'print "Thats 'em all flipped without perspctive"'print "Press a key"     'sleepdim as float angdim as long fps'fix these two components at startang.z=pi/2 'tilt angle sidewaysang.y=-pi/7 'tilt angle back and forwarddim as string key    do    key=inkey    If key=Chr(255)+"K" Then ang.z-=.02     'left    If key=Chr(255)+"M" Then ang.z+=.02     'right    If key=Chr(255)+"P" Then ang.y-=.02     'down    If key=Chr(255)+"H" Then ang.y+=.02     'up    if key="z" then ang.x -=.02    if key="x" then ang.x +=.02                  A3D = Angle3D.construct(ang.x,ang.y,ang.z)'get the six rotate components (sines, coses  ...)        screenlock    cls    draw string(50,50),"Framerate "&fps    draw string(50,150),"Use the arrow keys"        for n as long=lbound(a) to ubound(a)        rot(n)=rotate(type(400,300,0),a(n),A3D,type(1,1,1))        rot(n)=perspective(rot(n),type(400,300,1000))    next        qsortZ(rot(),lbound(rot),ubound(rot))        for n as long=lbound(rot) to ubound(rot)        pset(rot(n).x,rot(n).y),rot(n).col    next        screenunlock    sleep regulate(64,fps),1loop until key=chr(27)`
paul doe
Posts: 1244
Joined: Jul 25, 2017 17:22
Location: Argentina

### Re: 3d without openGL

owen wrote:@ paul doe

This problem disappears if you use matrices

I couldn't agree with you more.
And not only that, using matrices avoids gimbal lock

in my example you can see it encounter gimbal lock by first rotating 90 degrees on the x axis, then 90 on y axis and then stop at 89 on the z.
when you get to 89 on the z axis and proceed to rotate to 90 look closely and notice a tiny glitch.

what a wonderful experience it was to try and avoid gimbal lock to no avail. well actually i kinda did it but what length i went to.

im glad that my example can serve the perfect example of what not to do. LOL

Exactly. Thanks for pointing it out. I can imagine that you had a lot of fun trying to avoid the lock with Euler angles X-D. Don't fret. We all used Euler, at one point of another. In physics is the same, you stick with an Euler integrator until the simulation craps out on you ('explodes' to infinity).

Fun fact: the Apollo 11 mission almost failed thanks to the cheap a\$\$es of the NASA, because they didn't wanted to install another gimbal in the landing module to cut costs. Look at the Wikipedia entry for Gimbal Lock, quite amusing LOL.

@dodicat: pretty nice!. Reminds me of Point Clouds. I had some code on that, will look for it.

@BasicCoder2: aha. I see it now. I think that you are confusing rotation about an arbitrary axis with rotation around the object's own axes. One is extrinsic (outside the object's coordinate system) and the other is intrinsic. I will now post some code that will help you clarify that (hopefully).
paul doe
Posts: 1244
Joined: Jul 25, 2017 17:22
Location: Argentina

### Re: 3d without openGL

dafhi wrote:@paul doe - Very intrigued and inspired by your offerings

Hi dafhi

Last edited by paul doe on Oct 11, 2017 11:33, edited 1 time in total.
paul doe
Posts: 1244
Joined: Jul 25, 2017 17:22
Location: Argentina

### Re: 3d without openGL

BasicCoder2 wrote:Demonstration without enough explanation for my programming limitations.
What I tried to do was modify your code to produce the three axes of the object (without the object). That is the axes were the object. I tried to make the z blue the x red and the y green but they appear to be random colors. I wanted to be able to rotate around any of those three axes.

BasicCoder2:

As I promised, here is some code that will help to understand how coordinate systems work. You will need five more files (besides the ones that I already posted). Check the next post.
Last edited by paul doe on Oct 11, 2017 11:35, edited 2 times in total.
paul doe
Posts: 1244
Joined: Jul 25, 2017 17:22
Location: Argentina

### Re: 3d without openGL

camera.bas

Code: Select all

`#include once "platform.bi" '' platform specific definitions#include once "core.bi" '' core functions and defines#include once "math.bi" '' math functions and constants#include once "vec4.bi" '' vec4 type and operations#include once "mat4.bi" '' mat4 type and operations/'   this serves to define a projection plane      in the code, x and y coordinates are not used. It only serves to replace   the IRenderContext interface'/type rectangle   as double x   as double y   as double width   as double heightend type/'   simple camera class      responsibility: represents a pinhole-model camera. Can be used to navigate a scene FPS style      very useful for debugging purposes (especially geometry code)   mutability: mutable. Use with caution or solely for debugging. Just don't try to write Doom 5 with this.      10/10/2017: refactored the code because it's a mess.       '/class camera  public:      declare constructor( _         byval position as vec4, _         byval x as vec4, _         byval y as vec4, _         byval z as vec4, _         byval near as double = 1.0, _         byval far as double = 10000.0, _         byval projPlane as rectangle )                '' lots of crappy functions to play with         declare function getU() as vec4         declare function getV() as vec4         declare function getDir() as vec4         declare function getPos() as vec4                  declare sub move( byval offset as vec4 )         declare sub rotate( byval axis as vec4, byval angle as double )         declare function getRatioUV() as double         declare function getRatioVU() as double         declare sub setRatioUV( byval ratio as double )         declare sub setRatioVU( byval ratio as double)         declare function getZoomU() as double         declare function getZoomV() as double         declare sub setZoomU( byval a as double )         declare sub setZoomV( byval a as double )         declare sub zoom( byval a as double )         declare sub zoomU( byval a as double )         declare sub zoomV( byval a as double )         declare function getFOVV() as double         declare function getFOVU() as double         declare sub setFOVU( byval angle as double )         declare sub setFOVV( byval angle as double )         declare function getYaw() as double         declare function getPitch() as double         declare function getRoll() as double         declare sub setYaw( byval angle as double )         declare sub setPitch( byval angle as double )         declare sub setRoll( byval angle as double )         declare sub yawAbsolute( byval angle as double )         declare sub yaw( byval angle as double )         declare sub pitch( byval angle as double )         declare sub roll( byval angle as double )         declare function getCameraMatrix() as mat4         declare function getInverseCameraMatrix() as mat4         declare function transform( byval v as vec4 ) as vec4         declare function orthogonal( byval pnt as vec4, byval plane as rectangle, byref x as double, byref y as double, byref z as double ) as integer         declare function perspective( byval pnt as vec4, byval plane as rectangle, byref x as double, byref y as double, byref z as double ) as integer         declare function projectOnScreen( byval pnt as vec4, byval plane as rectangle, byref x as double, byref y as double, byref z as double ) as integer         declare function projectOnScreenOrtho( byval pnt as vec4, byval plane as rectangle, byref x as double, byref y as double, byref z as double ) as integer         declare function getScale() as double         declare sub setScale( byval dirLength as double )         declare sub scale( byval factor as double )         declare sub resetSkewU()         declare sub resetSkewV()         declare function getDistance( byval pnt as vec4 ) as double         declare sub setDistance( byval pnt as vec4, byval dist as double )         declare sub setLookDirection( byval newDir as vec4 )         declare sub lookAt( byval target as vec4 )         declare sub lookAt( byval target as vec4, byval up as vec4 )         declare sub setU( byval newU as vec4 )         declare sub setV( byval newV as vec4 )         declare sub setDir( byval newDir as vec4 )         declare sub setPos( byval newPos as vec4 )        declare function nearClip() as double        declare function farClip() as double        declare function projectionPlane() as rectangle          private:      declare sub generateMatrix()          '' members     as vec4 m_pos     as vec4 m_u     as vec4 m_v     as vec4 m_dir      as mat4 m_position           as mat4 m_camMatrix     as mat4 m_invCamMatrix        as double m_nearClip     as double m_farClip        as double m_fov     as rectangle m_projPlaneend classconstructor camera( _   byval position as vec4, _   byval x as vec4, _   byval y as vec4, _   byval z as vec4, _   byval near as double = 1.0, _   byval far as double = 10000.0, _   byval projPlane as rectangle )      m_pos = position   m_u = x   m_v = y   m_dir = z   m_nearClip = near   m_farClip = far   m_position = matrices.translation( position )   m_projPlane = projPlane      /'      adjust the aspect ratio of the camera to match the      projection plane resolution            sorry for this crap, I had to reimplement it this way to shorten      the code   '/   with m_projPlane      if( .width < .height ) then         setRatioUV( .width / .height )      else         setRatioVU( .height / .width )      end if   end with      generateMatrix()end constructorsub camera.generateMatrix()  '' generates the camera's rotation matrix  with m_camMatrix      .a = m_u.x : .b = m_v.x : .c = m_dir.x : .d = 0.0      .e = m_u.y : .f = m_v.y : .g = m_dir.y : .h = 0.0      .i = m_u.z : .j = m_v.z : .k = m_dir.z : .l = 0.0      .m = m_u.w : .n = m_v.w : .o = m_dir.w : .p = 1.0  end with      '' and calculates the inverse immediately   m_invCamMatrix = inverse( m_camMatrix )end subsub camera.move( byval offset as vec4 )   m_pos = m_pos + offset   end subsub camera.rotate( byval axis as vec4, byval angle as double )   '' rotates the camera   m_u = rotateAroundAxis( m_u, axis, angle )   m_v = rotateAroundAxis( m_v, axis, angle )   m_dir = rotateAroundAxis( m_dir, axis, angle )         generateMatrix()end subfunction camera.projectionPlane() as rectangle   return( m_projPlane )end functionfunction camera.getU() as vec4   return( m_u )end functionfunction camera.getV() as vec4   return( m_v )end functionfunction camera.getDir() as vec4   return( m_dir )end functionfunction camera.getPos() as vec4   return( m_pos )end functionfunction camera.getRatioUV() as double   return( m_u.length() / m_v.length() )end functionfunction camera.getRatioVU() as double   return( m_v.length() / m_u.length() )end functionsub camera.setRatioUV( byval ratio as double )   m_v.normalize()   m_v = m_v * m_u.length() / ratio      generateMatrix()end subsub camera.setRatioVU(byval ratio as double)   m_u.normalize()   m_u = m_u * m_v.length() / ratio      generateMatrix()end sub function camera.getZoomU() as double   return( m_dir.length() / m_u.length() )end functionfunction camera.getZoomV() as double   return( m_dir.length() / m_v.length() )end functionsub camera.setZoomU( byval a as double )   m_u = m_u / ( a / getZoomU() )      generateMatrix()end subsub camera.setZoomV( byval a as double )   m_v = m_v / ( a / getZoomV() )      generateMatrix()end subsub camera.zoom( byval a as double )   zoomU( a )   zoomV( a )end subsub camera.zoomU( byval a as double )   m_u = m_u / a      generateMatrix()end subsub camera.zoomV( byval a as double )   m_v = m_v / a      generateMatrix()end subfunction camera.getFOVV() as double   return( 2.0 * atan2( m_v.length(), m_dir.length() ) )end functionfunction camera.getFOVU() as double   return( 2.0 * atan2( m_u.Length(), m_dir.length() ) )end functionsub camera.setFOVU( byval angle as double )   setzoomU( 1.0 / tan( angle / 2.0 ) )end subsub camera.setFOVV( byval angle as double )   setzoomV( 1.0 / tan( angle / 2.0 ) )end subfunction camera.getYaw() as double   return( atan2( m_dir.x, m_dir.z ) )end functionfunction camera.getPitch() as double   return( atan2( m_dir.y, sqr( m_dir.x * m_dir.x + m_dir.z * m_dir.z ) ) )end functionfunction camera.getRoll() as double  return( vectorAngle( cross( vec4( 0.0, 1.0, 0.0 ), m_dir ), cross( m_v, m_dir ) ) )end functionsub camera.setYaw( byval angle as double )  '' to change yaw, you have to rotate around the "up" axis of the WORLD = the y axis  dim currentAngle as double = getyaw()    rotate( vec4( 0.0, 1.0, 0.0 ), -angle + currentAngle )end subsub camera.setPitch( byval angle as double )   dim currentAngle as double = getPitch()      rotate( m_u, angle - currentAngle )end subsub camera.setRoll( byval angle as double )   dim currentAngle as double = getRoll()      rotate( m_dir, angle - currentAngle )end subsub camera.yawAbsolute( byval angle as double )   rotate( vec4( 0.0, 1.0, 0.0 ), angle )end subsub camera.yaw( byval angle as double )   rotate( m_v, angle)end subsub camera.pitch(byval angle as double)   rotate( m_u, angle )end subsub camera.roll( byval angle as double )   rotate( m_dir, angle )end subfunction camera.getCameraMatrix() as mat4   return( m_camMatrix )end functionfunction camera.getInverseCameraMatrix() as mat4   return( m_invCamMatrix )end functionfunction camera.transform( byval v as vec4 ) as vec4   return( m_invCamMatrix * v )end functionfunction camera.orthogonal( byval pnt as vec4, byval plane as rectangle, byref x as double, byref y as double, byref z as double ) as integer   dim as mat4 PM   dim as double ratio      if( plane.width < plane.height ) then      ratio = plane.height   else      ratio = plane.width   end if      with PM      .a = ( 2 / ratio ) * 7 : .b = 0.0 : .c = 0.0 : .d = 0.0      .e = 0.0 : .f = -( 2 / ratio ) * 7 : .g = 0.0 : .h = 0.0      .i = 0.0 : .j = 0.0 : .k = ( -2.0 / ( m_farClip - m_nearClip ) ) * 7 : .l = ( ( m_farClip + m_nearClip ) / ( m_farClip - m_nearClip ) )      .m = 0.0 : .n = 0.0 : .o = 0.0 : .p = 1.0   end with      Dim as vec4 b   b = PM * pnt      dim as double px, py      px = ( b.x * plane.width ) + plane.width / 2   py = ( b.y * plane.height ) + plane.height / 2      if( px >= 0 and px <= plane.width - 1 and py >= 0 and py <= plane.height - 1 ) then       x = px       y = py       z = b.z          return( -1 )   else       x = px       y = py       z = b.z              return( 0 )   end ifend functionfunction camera.perspective( byval pnt as vec4, byval plane as rectangle, byref x as double, byref y as double, byref z as double ) as integer      '' the current matrix   dim as double fovX = radians( 90 )   dim as double fovY = radians( 90 )   dim as double fovTanX = 1 / tan( fovX / 2 )   dim as double fovTanY = 1 / tan( fovY / 2 )      '' projection matrix (assumes a right-handed coord system)  dim as mat4 PM    with PM      .a = fovTanX : .b = 0.0 : .c = 0.0 : .d = 0.0      .e = 0.0 : .f = fovTanY : .g = 0.0 : .h = 0.0      .i = 0.0 : .j = 0.0 : .k = -( ( farClip + nearClip ) / ( farClip - nearClip ) ) : .l = -( ( 2 * ( farClip * nearClip ) ) / ( farClip - nearClip ) )      .m = 0.0 : .n = 0.0 : .o = -1.0 : .p = 0.0  end with      '' project the point   dim as vec4 b   b = PM * pnt      dim as double px   dim as double py      px = ( b.x * plane.width ) / ( 1.0 * b.w ) + plane.width / 2   py = ( b.y * plane.height ) / ( 1.0 * b.w ) + plane.height / 2         if( px >= 0 and px <= plane.width - 1 and py >= 0 and py <= plane.height - 1 and b.z > m_nearClip ) then       x = px       y = py       z = b.z          return( -1 )   else       x = px       y = py       z = b.z          return( 0 )   end ifend functionfunction camera.projectOnScreen( byval pnt as vec4, byval plane as rectangle, byref x as double, byref y as double, byref z as double ) as integer   '' first transformation: translation   dim as mat4 TM = matrices.translation( -m_pos )   dim as vec4 a = vec4( TM * pnt )      '' second transformation: rotation   '' transform() multiplies the vector with the INVERSE camera matrix   '' to bring the point from camera space to world space   dim as vec4 b = vec4( transform( a ) )      '' Third transformation: projection   return( perspective( b, plane, x, y, z ) )end functionfunction camera.projectOnScreenOrtho( byval pnt as vec4, byval plane as rectangle, byref x as double, byref y as double, byref z as double ) as integer   dim as mat4 TM      '' first transformation: translation   with TM       .a = 1.0 : .b = 0.0 : .c = 0.0 : .d = -m_pos.x       .e = 0.0 : .f = 1.0 : .g = 0.0 : .h = -m_pos.y       .i = 0.0 : .j = 0.0 : .k = 1.0 : .l = -m_pos.z       .m = 0.0 : .n = 0.0 : .o = 0.0 : .p = 1.0   end with      dim as vec4 a = vec4( TM * pnt )      '' transform() multiplies the vector with the INVERSE camera matrix   '' to bring the point from camera space to world space   dim as vec4 b = vec4( transform( a ) )      '' third transformation: projection   return( orthogonal( b, plane, x, y, z ) )end functionfunction camera.getScale() as double   return( m_dir.length )end functionsub camera.setScale( byval dirLength as double )   scale( m_dir.length() / dirLength )end subsub camera.scale( byval factor as double )   m_dir = m_dir * factor   m_u = m_u * factor   m_v = m_v * factor      generateMatrix()end subsub camera.resetSkewU()   dim as double oldzoomU = getzoomU()   dim as double oldzoomV = getzoomV()      m_u = cross( m_dir, m_v )   m_v = cross( m_dir, -m_u )      setzoomU( oldzoomU )   setzoomV( oldzoomV )end subsub camera.resetSkewV()   dim as double oldzoomU = getzoomU()   dim as double oldzoomV = getzoomV()      m_v = cross( m_dir, m_u )   m_u = cross( m_dir, -m_v )      setzoomU( oldzoomU )   setzoomV( oldzoomV )end subsub camera.setLookDirection( byval newDir as vec4 )   dim as vec4 axis = vec4( ( cross( m_dir, newDir ) ) )      if axis.length() = 0 then      exit sub   end if      dim as double angle = vectorAngle( m_dir, newDir )      if angle <> 0 then      rotate( axis, angle )   end ifend subsub camera.lookAt( byval target as vec4, byval up as vec4 )   /'      preserve the lengths of the camera vectors            every time one messes with the camera vectors directly, it has to      remember to preserve their length, if you do not do this, the projection      gets distorted, as these vectors control different parameters of the      camera model         REMEMBER            u controls the horizontal FOV            v controls the vertical FOV            dir controls the focal zoom   '/   dim as double dl = m_dir.length()   dim as double ul = m_u.length()   dim as double vl = m_v.length()         m_dir = ( normalize( target - m_pos ) * dl )   m_U = ( normalize( cross( up, m_dir ) ) * ul )   m_V = ( normalize( cross( m_dir, m_u ) ) * vl )      generateMatrix()end subsub camera.lookAt( byval target as vec4 )   /'      preserve the lengths of the camera vectors            every time one messes with the camera vectors directly, it has to      remember to preserve their length, if you do not do this, the projection      gets distorted, as these vectors control different parameters of the      camera model         REMEMBER            u controls the horizontal FOV            v controls the vertical FOV            dir controls the focal zoom   '/   dim as double dl = m_dir.length()   dim as double ul = m_u.length()   dim as double vl = m_v.length()      '' compute the forward vector   dim as vec4 forward = vec4( normalize( target - m_pos ) )      '' compute temporal up vector based on the forward vector   '' watch out when look up/down at 90 degree   '' for example, forward vector is on the Y axis   dim as vec4 up      if( abs( forward.x ) < epsilon and abs( forward.z ) < epsilon ) then      if( forward.y > 0 ) then       '' forward vector is pointing +Y axis         up = vec4( 0.0, 0.0, -1.0 )      else         '' forward vector is pointing -Y axis         up = vec4( 0.0, 0.0, 1.0 )      end if   else      '' in general, up vector is straight up      up = vec4( 0.0, 1.0, 0.0 )   end if         '' compute the left vector   dim as vec4 leftV = vec4( normalize( cross( up, forward ) ) )      '' re-calculate the orthonormal up vector   up = normalize( cross( forward, leftV ) )      m_dir = forward * dl   m_u = leftV * ul   m_v = up * vl   generateMatrix()end subsub camera.setU( byval newU as vec4 )   m_u = newU      generateMatrix()end subsub camera.setV( byval newV as vec4 )   m_v = newV      generateMatrix()end subsub camera.setDir( byval newDir as vec4 )   m_dir = newDir      generateMatrix()end subsub camera.setPos( byval newPos as vec4 )   m_pos = newPosend subfunction camera.nearClip() as double   return( m_nearClip )end functionfunction camera.farClip() as double   return( m_farClip )end functionfunction camera.getDistance( byval pnt as vec4 ) as double   return( distance( m_pos, pnt ) )end functionpublic sub camera.setDistance( byval pnt as vec4, byval dist as double )  dim currentDist as double = distance( m_pos, pnt )  move( currentDist * normalize( pnt - m_pos ) )end sub`

object.bas

Code: Select all

`'' needed headers#include once "platform.bi"#include once "core.bi"#include once "math.bi"#include once "vec4.bi"#include once "mat4.bi"#include once "arrayList.bi"type edge   public:      declare constructor( byval v1 as vec4 ptr, byval v2 as vec4 ptr )            declare function vertex1() as vec4 ptr      declare function vertex2() as vec4 ptr         private:      as vec4 ptr m_v1      as vec4 ptr m_v2end typeconstructor edge( byval v1 as vec4 ptr, byval v2 as vec4 ptr )   m_v1 = v1   m_v2 = v2end constructorfunction edge.vertex1() as vec4 ptr   return( m_v1 )end functionfunction edge.vertex2() as vec4 ptr   return( m_v2 )end functionclass object3D   public:      declare constructor()            declare constructor( _         byval position as vec4 = vec4( 0.0, 0.0, 0.0 ), _         byval xAxis as vec4, _         byval yAxis as vec4, _         byval zAxis as vec4 )            declare destructor()            declare function transform( byval v as vec4 ) as vec4      declare function objectSpaceToWorldSpace( byval pnt as vec4 ) as vec4      declare function getObjectMatrix() as mat4      declare function getInverseObjectMatrix() as mat4      declare sub move( byval offset as vec4 )      declare sub rotate( byval axis as vec4, byval angle as double )      declare function getyaw() as double      declare function getPitch() as double      declare function getRoll() as double      declare sub setyaw( byval angle as double )      declare sub setPitch( byval angle as double )      declare sub setRoll( byval angle as double )      declare sub yawAbsolute( byval angle as double )      declare sub yaw( byval angle as double )      declare sub pitch( byval angle as double )      declare sub roll( byval angle as double )      declare function getDistance( byval pnt as vec4 ) as double      declare sub setDistance( byval pnt as vec4, byval dist as double )      declare sub setLookDirection( byval newDir as vec4 )      declare sub lookAt( byval lookAtMe as vec4 )      declare sub setU( byval newU as vec4 )      declare sub setV( byval newV as vec4 )      declare sub setDir( byval newDir as vec4 )      declare sub setPos( byval newPos as vec4 )      declare function getU() as vec4      declare function getV() as vec4      declare function getDir() as vec4      declare function getPos() as vec4            '' these two functions are wrong on so many levels that isn't even funny      '' made it this way to simplify a little      declare function vertices() as arrayList ptr      declare function edges() as arrayList ptr         private:      declare sub generateMatrix()            as vec4 m_pos = any      as vec4 m_u = any      as vec4 m_v = any      as vec4 m_dir = any            as mat4 m_objMatrix = any      as mat4 m_invObjMatrix = any            as arrayList m_vertices      as arrayList m_edgesend classconstructor object3d()   '' construct an object with the default axes   m_pos = vec4( 0.0, 0.0, 0.0 )   m_u = vec4( 1.0, 0.0, 0.0 )   m_v = vec4( 0.0, 1.0, 0.0 )   m_dir = vec4( 0.0, 0.0, 1.0 )      generateMatrix()end constructorconstructor object3d( _   byval position as vec4 = vec4( 0.0, 0.0, 0.0 ), _   byval xAxis as vec4, _   byval yAxis as vec4, _   byval zAxis as vec4 )      m_pos = position   m_u = xAxis   m_v = yAxis   m_dir = zAxis      generateMatrix()end constructordestructor object3d()   '' if there are vertices on the object, release them   dim as vec4 ptr v      if( m_vertices.count > 0 ) then      for i as integer = 0 to m_vertices.count - 1         v = m_vertices.get( i )         delete( v )      next   end if      '' same for the edges   dim as edge ptr e      if( m_edges.count > 0 ) then      for i as integer = 0 to m_edges.count - 1         e = m_edges.get( i )         delete( e )      next   end ifend destructorfunction object3d.vertices() as arrayList ptr   return( @m_vertices )end functionfunction object3d.edges() as arrayList ptr   return( @m_edges )end functionfunction object3d.getObjectMatrix() as mat4   return( m_objMatrix )end functionfunction object3d.getInverseObjectMatrix() as mat4   return( m_invObjMatrix )end functionsub object3d.generateMatrix()  '' generates the object rotation matrix  m_objMatrix = mat4( _      m_u.x , m_v.x , m_dir.x , 0.0, _      m_u.y , m_v.y , m_dir.y , 0.0, _      m_u.z , m_v.z , m_dir.z , 0.0, _      m_u.w , m_v.w , m_dir.w , 1.0 )      '' and compute the inverse immediately   m_invObjMatrix = inverse( m_objMatrix )end subsub object3d.move( byval offset as vec4 )   m_pos += offsetend subsub object3d.rotate( byval axis as vec4, byval angle as double )   '' rotates the object around an arbitrary axis   m_u = rotateAroundAxis( m_u, axis, angle)   m_v = rotateAroundAxis( m_v, axis, angle)   m_dir = rotateAroundAxis( m_dir, axis, angle)      generateMatrix()end subfunction object3d.getU() as vec4   return( m_u )end functionfunction object3d.getV() as vec4   return( m_v )end functionfunction object3d.getDir() as vec4   return( m_dir )end functionfunction object3d.getPos() as vec4   return( m_pos )end functionfunction object3d.getyaw() as double   return( atan2( m_dir.x, m_dir.z ) )end functionfunction object3d.getPitch() as double   return( atan2( m_dir.y, sqr( m_dir.x * m_dir.x + m_dir.z * m_dir.z ) ) )end functionfunction object3d.getRoll() as double  return( vectorAngle( cross( vec4( 0.0, 1.0, 0.0 ), m_dir ), cross( m_v, m_dir ) ) )end functionsub object3d.setyaw( byval angle as double )  '' to change yaw, you have to rotate around the "up" axis of the WORLD = the y axis  dim as double currentAngle = getyaw()  rotate( vec4( 0.0, 1.0, 0.0 ), -angle + currentAngle )end subsub object3d.setPitch( byval angle as double )  dim as double currentAngle = getPitch()  rotate( m_u, angle - currentAngle )end subsub object3d.setRoll( byval angle as double )   dim as double currentAngle = getRoll()   rotate( m_dir, angle - currentAngle )end subsub object3d.yawAbsolute( byval angle as double )   rotate( vec4( 0.0, 1.0, 0.0 ), angle )end subsub object3d.yaw( byval angle as double )   rotate( m_v, angle )end subsub object3d.pitch( byval angle as double )   rotate( m_u, angle )end subsub object3d.roll( byval angle as double )   rotate( m_dir, angle )end subfunction object3d.getDistance( byval pnt as vec4 ) as double   return( distance( m_pos, pnt ) )end functionsub object3d.setDistance( byval pnt as vec4, byval dist as double )   dim as double currentDist = distance( m_pos, pnt )   move( currentDist * normalize( pnt - pos ) )end subsub object3d.setLookDirection( byval newDir as vec4 )   dim axis as vec4 = vec4( ( cross( m_dir, newDir ) ) )      if axis.length() = 0 then      exit sub   end if      dim as double angle = vectorAngle( m_dir, newDir )      if( angle <> 0 ) then      rotate( axis, angle )   end ifend subsub object3d.lookAt( byval lookAtMe as vec4 )   dim as double dl = m_dir.length()   dim as double ul = m_u.length()   dim as double vl = m_v.length()      dim as vec4 w = vec4( lookAtMe - m_pos )         setDir( Normalize( w ) * dl )   setU( normalize( cross( vec4( 0.0, 1.0, 0.0 ), m_dir ) ) * ul )   setV( normalize( cross( m_dir, m_u ) ) * vl )end subsub object3d.setU( byval newU as vec4 )   m_u = newU   generateMatrix()end subsub object3d.setV( byval newV as vec4 )   m_v = newV   generateMatrix()end subsub object3d.setDir( byval newDir as vec4 )   m_dir = newDir   generateMatrix()end subsub object3d.setPos( byval newPos as vec4 )   m_pos = newPosend subfunction object3d.transform( byval v as vec4 ) as vec4   return( m_objMatrix * v )end functionfunction object3d.objectSpaceToWorldSpace( byval pnt as vec4 ) as vec4  dim as vec4 b = vec4( transform( pnt ) )  dim as vec4 a = vec4( b + m_pos )    return( a )end function`

UPDATE 2: 13/11/2017 - camera.bas: corrected some more bugs and added a few methods that I forgot to port. Sorry :'(
UPDATE 1: 10/11/2017 - camera.bas: killed a few bugses and sprinkled some more comments
Last edited by paul doe on Oct 13, 2017 10:27, edited 2 times in total.
paul doe
Posts: 1244
Joined: Jul 25, 2017 17:22
Location: Argentina

### Re: 3d without openGL

vec4.bi

Code: Select all

`#ifndef __vec4__   #define __vec4__         #include once "math.bi"      /'      homogeneous 4 tuple vector type            it is meant to be used as an homogeneous 3D vector (like the ones      used by OpenGL or Direct3D)            conceptually, they are to be interpreted like this:            | x |      | y |      | z |      | 1 |            which is known as an homogeneous column vector            10/7/2017: fixed a typo in the + operator         now this was the original code. Take a very close look at it (the comment was removed         after the fix out of pure anger, as you can imagine):                        operator +( byref v as vec4, byref w as vec4 ) as vec4               '' addition operator (all rocket science stuff...)               return( vec4( v.x + w.x, v.y + w.y, v.x + w.z ) )            end operator                  goes to show you what happens when you are a #%\$@ imbecile and write code in a hurry.         STAY TUNED FOR THE NExT ExCITING TUTORIAL: how to write unbelievable \$%#@ 3D code!!!   '/      type vec4      public:         as double x         as double y         as double z         as double w               declare constructor( byval nx as double = 0.0, byval ny as double = 0.0, byval nz as double = 0.0, byval nw as double = 1.0 )         declare constructor( byref rhs as vec4 )         declare operator let( byref rhs as vec4 )                  declare operator cast() as string                  declare function length() as double         declare function squaredLength() as double         declare sub normalize()         declare function normal() as vec4         declare sub homogeneize()         declare function homogeneous() as vec4         declare function cross( byref v as vec4 ) as vec4         declare function dot( byref v as vec4 ) as double         declare function distance( byref v as vec4 ) as double   end type      constructor vec4( byval nx as double = 0.0, byval ny as double = 0.0, byval nz as double = 0.0, byval nw as double = 1.0 )      '' default constructor creates an homogeneous vector      x = nx      y = ny      z = nz      w = nw   end constructor      constructor vec4( byref rhs as vec4 )      '' copy constructor      x = rhs.x      y = rhs.y      z = rhs.z      w = rhs.w   end constructor      operator vec4.let( byref rhs as vec4 )      '' assignment constructor      x = rhs.x      y = rhs.y      z = rhs.z      w = rhs.w   end operator      operator vec4.cast() as string      '' human readable string representation (useful for debugging)      return( _         "| " & trim( str( x ) ) & " |" & chr( 13 ) & chr( 10 ) & _         "| " & trim( str( y ) ) & " |" & chr( 13 ) & chr( 10 ) & _         "| " & trim( str( z ) ) & " |" & chr( 13 ) & chr( 10 ) & _         "| " & trim( str( w ) ) & " |" & chr( 13 ) & chr( 10 ) )   end operator      function vec4.squaredLength() as double      /'         returns the squared length of this vector         useful when you just want to compare which one is bigger, as            this avoids having to compute a square root      '/      return( x * x + y * y + z * z )   end function      function vec4.length() as double      '' returns the length of this vector      return( sqr( x * x + y * y + z * z ) )   end function      sub vec4.normalize()      /'         normalizes the vector            note that the homogeneous coordinate (w) is not touched      '/      dim as double l = length()            if l > 0 then         x /= l : y /= l : z /= l      end if   end sub      function normalize( byref v as vec4 ) as vec4      '' for compatibility      return( v.normal() )   end function      function vec4.normal() as vec4      /'         returns this vector normalized but without altering itself         again the homogeneous coordinate is left alone      '/      dim as double l = length()      dim as vec4 v = vec4( this )            if l > 0 then         v.x /= l : v.y /= l : v.z /= l      end if            return( v )   end function      sub vec4.homogeneize()      /'         homogeneizes the vector            this is done by dividing the components by the homogeneous coordinate (w)      '/      x /= w : y /= w : z /= w : w /= w   end sub      function vec4.homogeneous() as vec4      /'         returns this vector homogeneized but without altering it      '/      return( vec4( x / w, y / w, z / w, w / w ) )   end function      function vec4.cross( byref v as vec4 ) as vec4      /'         returns the cross product (aka vectorial product) of this vector and         another vector v            note that there's no cross product defined for a 4 tuple, so            we simply use a standard 3d cross product, and make the w component 1      '/      return( vec4( _         v.y * z - v.z * y, _         v.z * x - v.x * z, _         v.x * y - v.y * x, _         1.0 ) )   end function      function cross( byref v as vec4, byref w as vec4 ) as vec4      '' returns the cross product between vectors v and w      return( v.cross( w ) )   end function      function vec4.dot( byref v as vec4 ) as double     '' returns the dot product (aka scalar product) of this vector and vector v      return( v.x * x + v.y * y + v.z * z )   end function      function dot( byref v as vec4, byref w as vec4 ) as double      '' returns the dot product between two vectors      return( v.dot( w ) )   end function      function vec4.distance( byref v as vec4 ) as double     /'        gets the distance of this vector with vector v           to calculate the distance, substract them and calculate the length of the resultant vector      '/      return( sqr( ( v.x - x ) ^ 2 + ( v.y - y ) ^ 2 + ( v.z - z ) ^ 2 ) )   end function      function distance( byref v as vec4, byref w as vec4 ) as double      return( sqr( ( v.x - w.x ) ^ 2 + ( v.y - w.y ) ^ 2 + ( v.z - w.z ) ^ 2 ) )   end function      operator -( byref v as vec4, byref w as vec4 ) as vec4     '' substraction operator     return( vec4( v.x - w.x, v.y - w.y, v.z - w.z ) )   end operator      operator -( byref v as vec4 ) as vec4      '' negation operator      return( vec4( -v.x, -v.y, -v.z ) )   end operator      operator +( byref v as vec4, byref w as vec4 ) as vec4      return( vec4( v.x + w.x, v.y + w.y, v.z + w.z ) )   end operator      operator *( byref v as vec4, byref a as double ) as vec4      /'         multiplication with a scalar            note that this is not a 'scalar product', but a multiplication with a number            vectors do not define multiplications per se, they define the dot product            and the cross product. To avoid confusion (and also correctness), the            multiplication operator is overloaded to a scaling of the vector      '/      return( vec4( v.x * a, v.y * a, v.z * a ) )   end operator      operator *( byref a as double, byref v as vec4 ) as vec4      '' same as above but with the parameters inversed      return( vec4( v.x * a, v.y * a, v.z * a ) )   end operator      operator /( byref v as vec4, byref a as double ) as vec4      '' division by a scalar. See note above on multiplying a vector      return( vec4( v.x / a, v.y / a, v.z / a ) )   end operator         function vectorAngle( byref v as vec4, byref w as vec4 ) as double      /'         returns the angle between two vectors using the dot product, in radians      note that the result of the dot product used here      should mathematically speaking, always give a result between -1 and 1. Due to imprecisions of      numerical calculations it might sometimes be a little bit outside this range however (especially      if you define Scalar to be float instead of double). If that happens, the acos function will give      an invalid result. So instead a protection was added that sets the value back to 1 or -1      (because, if the value became 1.0000023 for example, it was probably supposed to be 1 anyway).                 dotProduct( v, w ) = length( v ) * length( w ) * cos( angle )         angle = aCos( dotProduct / ( length( v ) * length( w ) ) )                              thus:                              angle = aCos( dot( normal( v ) * normal( w ) ) )      '/     dim as double angleCos = dot( v.normal(), w.normal() )          '' for acos, the value has to be between -1.0 and 1.0, but due to numerical imprecisions it     '' sometimes comes outside this range            angleCos = clamp( -1.0, 1.0, angleCos )        return( -acos( angleCos ) )   end function      function rotateAroundAxis( byref v as vec4, byref axis as vec4, byval angle as double ) as vec4    /'       rotate vector v around arbitrary axis for angle radians       it can only rotate around an axis through our object, to rotate around another axis:       first translate the object to the axis, then use this function, then translate back       in the new direction.      '/    if( ( v.x = 0 ) and ( v.y = 0 ) and ( v.z = 0 ) ) then       return vec4( 0.0, 0.0, 0.0 )    end if    dim nAxis as vec4 = vec4( axis.x, axis.y, axis.z )    nAxis.normalize()    '' calculate parameters of the rotation matrix    dim as double c = cos( angle )    dim as double s = sin( angle )    dim as double t = 1 - c    '' multiply w with rotation matrix    dim w as vec4        w.x = ( t * nAxis.x * nAxis.x + c ) * v.x _        + ( t * nAxis.x * nAxis.y + s * nAxis.z ) * v.y _        + ( t * nAxis.x * nAxis.z - s * nAxis.y ) * v.z    w.y = ( t * nAxis.x * nAxis.y - s * nAxis.z ) * v.x _        + ( t * nAxis.y * nAxis.y + c ) * v.y _        + ( t * nAxis.y * nAxis.z + s * nAxis.x ) * v.z    w.z = ( t * nAxis.x * nAxis.z + s * nAxis.y ) * v.x _        + ( t * nAxis.y * nAxis.z - s * nAxis.x ) * v.y _        + ( t * nAxis.z * nAxis.z + c ) * v.z            '' the vector has to retain its length, so it's normalized and      '' multiplied with the original length    w.normalize()    w = w * v.length()    return( w )   end function   #endif`

mat4.bi

Code: Select all

`#ifndef __mat4__   #define __mat4__   #include once "core.bi"   #include once "platform.bi"   #include once "vec4.bi"   /'                             | a b c d |      4x4 matrix type      | e f g h |                            | i j k l |                            | m n o p |            9/30/2017: improved 4x4 matrix inverse calculation. With -gen gcc -O max         (the settings I always use) it is more than 60% faster than the previous         version. A byproduct of correcting the bug mentioned below.      9/29/2017: fixed determinant calculation (it was erroneously computed)         The funny thing was that this library is used in various applications,         including my 3D engine and various tools. When I was implementing a         feature, it kept doing weird things, and the bug was finally tracked to         an erroneous calculation of the determinant. The determinant was         correctly calculated in the 3x3   matrix code, but when I ported the code         to use OpenGL for rendering, it   didn't worked as intended. Goes to show         you that one is to be extra   careful with the math code, as it is the         foundation of the entire game engine.                  the calculations were cross checked with the help of this online resource:                        https://www.mathsisfun.com/algebra/matrix-calculator.html                  which has a very neat matrix calculator for various dimensions.   '/      type mat4      public:         as double a, b, c, d         as double e, f, g, h         as double i, j, k, l         as double m, n, o, p                  declare constructor( _            byval sa as double = 1.0, byval sb as double = 0.0, byval sc as double = 0.0, byval sd as double = 0.0, _            byval se as double = 0.0, byval sf as double = 1.0, byval sg as double = 0.0, byval sh as double = 0.0, _            byval si as double = 0.0, byval sj as double = 0.0, byval sk as double = 1.0, byval sl as double = 0.0, _            byval sm as double = 0.0, byval sn as double = 0.0, byval so as double = 0.0, byval sp as double = 1.0 _            )         declare constructor( byref NM as mat4 )         declare operator let( byref RHS as mat4 )         declare operator cast() as string                  declare function determinant() as double         declare sub transpose()         declare sub inverse() '' slower version         'declare sub inverse()         declare sub identity()   end type      constructor mat4( _            byval sa as double = 1.0, byval sb as double = 0.0, byval sc as double = 0.0, byval sd as double = 0.0, _            byval se as double = 0.0, byval sf as double = 1.0, byval sg as double = 0.0, byval sh as double = 0.0, _            byval si as double = 0.0, byval sj as double = 0.0, byval sk as double = 1.0, byval sl as double = 0.0, _            byval sm as double = 0.0, byval sn as double = 0.0, byval so as double = 0.0, byval sp as double = 1.0 _            )      /'         default constructor initializes the matrix to an identity, if no coefficients are specified.         this is far more useful than initializing it to all zeros      '/      a = sa : b = sb : c = sc : d = sd      e = se : f = sf : g = sg : h = sh      i = si : j = sj : k = sk : l = sl      m = sm : n = sn : o = so : p = sp   end constructor      constructor mat4( byref RHS as mat4 )      '' copy constructor      a = RHS.a : b = RHS.b : c = RHS.c : d = RHS.d      e = RHS.e : f = RHS.f : g = RHS.g : h = RHS.h      i = RHS.i : j = RHS.j : k = RHS.k : l = RHS.l      m = RHS.m : n = RHS.n : o = RHS.o : p = RHS.p   end constructor      operator mat4.let( byref RHS as mat4 )      '' assignment construction      a = RHS.a : b = RHS.b : c = RHS.c : d = RHS.d      e = RHS.e : f = RHS.f : g = RHS.g : h = RHS.h      i = RHS.i : j = RHS.j : k = RHS.k : l = RHS.l      m = RHS.m : n = RHS.n : o = RHS.o : p = RHS.p   end operator      operator mat4.cast() as string      '' the matrix in a human readable form (very useful for debugging purposes)      return( _         "| " & trim( str( a ) ) & " | " & trim( str( b ) ) & " | " & trim( str( c ) ) & " | " & trim( str( d ) ) & " |" & chr( 13 ) & chr( 10 ) & _         "| " & trim( str( e ) ) & " | " & trim( str( f ) ) & " | " & trim( str( g ) ) & " | " & trim( str( h ) ) & " |" & chr( 13 ) & chr( 10 ) & _         "| " & trim( str( i ) ) & " | " & trim( str( j ) ) & " | " & trim( str( k ) ) & " | " & trim( str( l ) ) & " |" & chr( 13 ) & chr( 10 ) & _         "| " & trim( str( m ) ) & " | " & trim( str( n ) ) & " | " & trim( str( o ) ) & " | " & trim( str( p ) ) & " |" & chr( 13 ) & chr( 10 ) _      )   end operator      function mat4.determinant() as double      /'         computes the determinant of the matrix using Laplace cofactor expansion                  the determinant of a 3x3 matrix is:            a * ( e * i - f * h ) - b * ( d * i - f * g ) + c * ( d * h - e * g )                  and a 4x4 matrix determinant is given by:            a *         b *         c *         d *            | f g h |   | e g h |   | e f h |   | e f g |            | j k l | - | i k l | + | i j l | - | i j k |            | n o p |   | m o p |   | m n p |   | m n o |                        where the '|' means the determinant of the inner 3x3 matrices. Note that the            cofactors are already factored in the calculation.                     the determinant is thus:            + ( a * (   f * ( k * p - l * o ) - g * ( j * p - l * n ) + h * ( j * o - k * n ) ) )            - ( b * ( e * ( k * p - l * o ) - g * ( i * p - l * m ) + h * ( i * o - k * m ) ) )            +   ( c * ( e * ( j * p - l * n ) - f * ( i * p - l * m ) + h * ( i * n - j * m ) ) )            - ( d * ( e * ( j * o - k * n ) - f * ( i * o - k * m ) + g * ( i * n - j * m ) ) )               '/      dim as double det = _            ( a * (   f * ( k * p - l * o ) - g * ( j * p - l * n ) + h * ( j * o - k * n ) ) ) _         - ( b * ( e * ( k * p - l * o ) - g * ( i * p - l * m ) + h * ( i * o - k * m ) ) ) _         +   ( c * ( e * ( j * p - l * n ) - f * ( i * p - l * m ) + h * ( i * n - j * m ) ) ) _         - ( d * ( e * ( j * o - k * n ) - f * ( i * o - k * m ) + g * ( i * n - j * m ) ) )      /'         this isn't matematically correct, just a programmer's dirty hack.         if the determinant of a matrix is 0, it means it has no inverse. In the code for         calculating the inverse, a division by the determinant is performed; and if it is          zero, a division by zero is performed on *every* element of the matrix, filling it         with positive or negative infinity values and rendering it useless. A matrix          without inverse is the matrix   itself, so setting the determinant value to 1          does the trick.      '/      if det = 0 then det = 1.0            return( det )   end function      sub mat4.transpose()      /'        transposes the matrix                   [ a b c d ]T    [ a e i m ]           [ e f g h ]  =  [ b f j n ]           [ i j k l ]     [ c g k o ]           [ m n o p ]     [ d h l p ]                why have it, if it is not used by the matrix code itself? Well, there is a        nice property of matrices, which has to do with rotations. If you can        be sure that the matrix contains only rotations, transposing it is the        same as taking its inverse, thus saving you *a lot* of computation      '/      this = mat4( a, e, i, m, b, f, j, n, c, g, k, o, d, h, l, p )   end sub      /'   sub mat4.inverse()      /'         computes the inverse of a 4x4 matrix                  this version is 60%+ faster and 400%+ uglier than the previous version.         it was made so by computing the determinant inside the method and         recycling as much   calculation as possible      '/      '' list of 2x2 determinants      dim as double kplo = k * p - l * o      dim as double jpln = j * p - l * n      dim as double jokn = j * o - k * n      dim as double iplm = i * p - l * m      dim as double iokm = i * o - k * m      dim as double injm = i * n - j * m      dim as double gpho = g * p - h * o      dim as double fphn = f * p - h * n      dim as double fogn = f * o - g * n      dim as double ephm = e * p - h * m      dim as double eogm = e * o - g * m      dim as double enfm = e * n - f * m      dim as double glhk = g * l - h * k      dim as double flhj = f * l - h * j      dim as double fkgj = f * k - g * j      dim as double elhi = e * l - h * i      dim as double ekgi = e * k - g * i      dim as double ejfi = e * j - f * i         '' list of 3x3 determinants      dim as double d1kplo = f * kplo      dim as double d1jpln = g * jpln      dim as double d1jokn = h * jokn      dim as double d2kplo = e * kplo      dim as double d2iplm = g * iplm      dim as double d2iokm = h * iokm      dim as double d3jpln = e * jpln      dim as double d3iplm = f * iplm      dim as double d3injm = h * injm      dim as double d4jokn = e * jokn      dim as double d4iokm = f * iokm      dim as double d4injm = g * injm      dim as double d5kplo = b * kplo      dim as double d5jpln = c * jpln      dim as double d5jokn = d * jokn      dim as double d6kplo = a * kplo      dim as double d6iplm = c * iplm      dim as double d6iokm = d * iokm      dim as double d7jpln = a * jpln      dim as double d7iplm = b * iplm      dim as double d7injm = d * injm      dim as double d8jokn = a * jokn      dim as double d8iokm = b * iokm      dim as double d8injm = c * injm      dim as double d9gpho = b * gpho      dim as double d9fphn = c * fphn      dim as double d9fogn = d * fogn      dim as double d10gpho = a * gpho      dim as double d10ephm = c * ephm      dim as double d10eogm = d * eogm      dim as double d11fphn = a * fphn      dim as double d11ephm = b * ephm      dim as double d11enfm = d * enfm      dim as double d12fogn = a * fogn      dim as double d12eogm = b * eogm      dim as double d12enfm = c * enfm      dim as double d13glhk = b * glhk      dim as double d13flhj = c * flhj      dim as double d13fkgj = d * fkgj      dim as double d14glhk = a * glhk      dim as double d14elhi = c * elhi      dim as double d14ekgi = d * ekgi      dim as double d15flhj = a * flhj      dim as double d15elhi = b * elhi      dim as double d15ejfi = d * ejfi      dim as double d16fkgj = a * fkgj      dim as double d16ekgi = b * ekgi      dim as double d16ejfi = c * ejfi            '' 4x4 determinant (inversed)      dim as double det = _            ( a * ( d1kplo - d1jpln + d1jokn ) _         - ( b * ( d2kplo - d2iplm + d2iokm ) ) _         + ( c * ( d3jpln - d3iplm + d3injm ) ) _         - ( d * ( d4jokn - d4iokm + d4injm ) ) )            '' if the determinant is 0, the matrix has no inverse      if det = 0 then exit sub            '' Multiplying with the reciprocal is slightly faster than dividing      dim as double invDet = 1.0 / det            '' minors      dim as double Ma = d1kplo - d1jpln + d1jokn      dim as double Mb = d2kplo - d2iplm + d2iokm      dim as double Mc = d3jpln - d3iplm + d3injm      dim as double Md = d4jokn - d4iokm + d4injm      dim as double Me = d5kplo - d5jpln + d5jokn      dim as double Mf = d6kplo - d6iplm + d6iokm      dim as double Mg = d7jpln - d7iplm + d7injm      dim as double Mh = d8jokn - d8iokm + d8injm      dim as double Mi = d9gpho - d9fphn + d9fogn      dim as double Mj = d10gpho - d10ephm + d10eogm      dim as double Mk = d11fphn - d11ephm + d11enfm      dim as double Ml = d12fogn - d12eogm + d12enfm      dim as double Mm = d13glhk - d13flhj + d13fkgj      dim as double Mn = d14glhk - d14elhi + d14ekgi      dim as double Mo = d15flhj - d15elhi + d15ejfi      dim as double Mp = d16fkgj - d16ekgi + d16ejfi            /'         adjugate (the adjugate is the transpose of the cofactored matrix of minors)                             Ma   -Me    Mi   -Mm         -Mb    Mf   -Mj    Mn          Mc   -Mg    Mk   -Mo         -Md    Mh   -Ml    Mp      '/      this = mat4( _          Ma * invDet, -Me * invDet,  Mi * invDet, -Mm * invDet, _         -Mb * invDet,  Mf * invDet, -Mj * invDet,  Mn * invDet, _          Mc * invDet, -Mg * invDet,  Mk * invDet, -Mo * invDet, _         -Md * invDet,  Mh * invDet, -Ml * invDet,  Mp * invDet )   end sub   '/   sub mat4.inverse()     /'        calculates the inverse of a 4x4 matrix                the minor, cofactor and adjugate matrices are all calculated and factored into the        code to make it shorter and more efficient     '/          '' multiplying by the reciprocal of the determinant is faster than dividing by it       dim as double invDet = 1 / determinant          this = mat4( _       ( f * k * p + g * l * n + h * j * o - f * l * o - g * j * p - h * k * n ) * invDet, _       ( b * l * o + c * j * p + d * k * n - b * k * p - c * l * n - d * j * o ) * invDet, _       ( b * g * p + c * h * n + d * f * o - b * h * o - c * f * p - d * g * n ) * invDet, _       ( b * h * k + c * f * l + d * g * j - b * g * l - c * h * j - d * f * k ) * invDet, _       ( e * l * o + g * i * p + h * k * m - e * k * p - g * l * m - h * i * o ) * invDet, _       ( a * k * p + c * l * m + d * i * o - a * l * o - c * i * p - d * k * m ) * invDet, _       ( a * h * o + c * e * p + d * g * m - a * g * p - c * h * m - d * e * o ) * invDet, _       ( a * g * l + c * h * i + d * e * k - a * h * k - c * e * l - d * g * i ) * invDet, _       ( e * j * p + f * l * m + h * i * n - e * l * n - f * i * p - h * j * m ) * invDet, _       ( a * l * n + b * i * p + d * j * m - a * j * p - b * l * m - d * i * n ) * invDet, _       ( a * f * p + b * h * m + d * e * n - a * h * n - b * e * p - d * f * m ) * invDet, _       ( a * h * j + b * e * l + d * f * i - a * f * l - b * h * i - d * e * j ) * invDet, _       ( e * k * n + f * i * o + g * j * m - e * j * o - f * k * m - g * i * n ) * invDet, _       ( a * j * o + b * k * m + c * i * n - a * k * n - b * i * o - c * j * m ) * invDet, _       ( a * g * n + b * e * o + c * f * m - a * f * o - b * g * m - c * e * n ) * invDet, _       ( a * f * k + b * g * i + c * e * j - a * g * j - b * e * k - c * f * i ) * invDet _      )      '' lovely, isn't it?   end sub      sub mat4.identity()      '' makes the matrix an identity matrix      a = 1.0 : b = 0.0 : c = 0.0 : d = 0.0      e = 0.0 : f = 1.0 : g = 0.0 : h = 0.0      i = 0.0 : j = 0.0 : k = 1.0 : l = 0.0      m = 0.0 : n = 0.0 : o = 0.0 : p = 1.0   end sub      operator *( byref A as mat4, byref B as mat4 ) as mat4      /'         multiply two 4x4 matrices                  remember that matrix multiplication is not commutative!            A * B != B * A      '/      return( mat4( _      A.a * B.a + A.b * B.e + A.c * B.i + A.d * B.m, _      A.a * B.b + A.b * B.f + A.c * B.j + A.d * B.n, _      A.a * B.c + A.b * B.g + A.c * B.k + A.d * B.o, _      A.a * B.d + A.b * B.h + A.c * B.l + A.d * B.p, _      A.e * B.a + A.f * B.e + A.g * B.i + A.h * B.m, _      A.e * B.b + A.f * B.f + A.g * B.j + A.h * B.n, _      A.e * B.c + A.f * B.g + A.g * B.k + A.h * B.o, _      A.e * B.d + A.f * B.h + A.g * B.l + A.h * B.p, _      A.i * B.a + A.j * B.e + A.k * B.i + A.l * B.m, _      A.i * B.b + A.j * B.f + A.k * B.j + A.l * B.n, _      A.i * B.c + A.j * B.g + A.k * B.k + A.l * B.o, _      A.i * B.d + A.j * B.h + A.k * B.l + A.l * B.p, _      A.m * B.a + A.n * B.e + A.o * B.i + A.p * B.m, _      A.m * B.b + A.n * B.f + A.o * B.j + A.p * B.n, _      A.m * B.c + A.n * B.g + A.o * B.k + A.p * B.o, _      A.m * B.d + A.n * B.h + A.o * B.l + A.p * B.p _      ) )      'return( mat4( _      '   A.a * B.a + A.b * B.e + A.c * B.i + A.d * B.m, _      '   A.a * B.b + A.b * B.f + A.c * B.j + A.d * B.n, _      '   A.a * B.c + A.b * B.g + A.c * B.k + A.d * B.o, _      '   A.a * B.d + A.b * B.h + A.c * B.l + A.d * B.p, _      '   A.e * B.a + A.f * B.e + A.g * B.i + A.h * B.m, _      '   A.e * B.b + A.f * B.f + A.g * B.j + A.h * B.n, _      '   A.e * B.c + A.f * B.g + A.g * B.k + A.h * B.o, _      '   A.e * B.d + A.f * B.h + A.g * B.l + A.h * B.p, _      '   A.i * B.a + A.j * B.e + A.k * B.i + A.l * B.m, _      '   A.i * B.b + A.j * B.f + A.k * B.j + A.l * B.n, _      '   A.i * B.c + A.j * B.g + A.k * B.k + A.l * B.o, _      '   A.i * B.d + A.j * B.h + A.k * B.l + A.l * B.p, _      '   A.m * B.a + A.n * B.e + A.o * B.i + A.p * B.m, _      '   A.m * B.b + A.n * B.f + A.o * B.j + A.p * B.n, _      '   A.m * B.c + A.n * B.g + A.o * B.k + A.p * B.o, _      '   A.m * B.d + A.n * B.h + A.o * B.l + A.p * B.p _      ') )   end operator      operator +( byref A as mat4, byref B as mat4 ) as mat4      '' adds two 4x4 matrices      return( mat4( _         A.a + B.a, A.b + B.b, A.c + B.c, A.d + B.d, _         A.e + B.e, A.f + B.f, A.g + B.g, A.h + B.h, _         A.i + B.i, A.j + B.j, A.k + B.k, A.l + B.l, _         A.m + B.m, A.n + B.n, A.o + B.o, A.p + B.p _      ) )   end operator      operator -( byref A as mat4, byref B as mat4 ) as mat4      '' substracts two 4x4 matrices      return( mat4( _         A.a - B.a, A.b - B.b, A.c - B.c, A.d - B.d, _         A.e - B.e, A.f - B.f, A.g - B.g, A.h - B.h, _         A.i - B.i, A.j - B.j, A.k - B.k, A.l - B.l, _         A.m - B.m, A.n - B.n, A.o - B.o, A.p - B.p _      ) )   end operator      operator -( byref A as mat4 ) as mat4      '' negates the matrix      return( mat4( _         -A.a, -A.b, -A.c, -A.d, _         -A.e, -A.f, -A.g, -A.h, _         -A.i, -A.j, -A.k, -A.l, _         -A.m, -A.n, -A.o, -A.p _      ) )   end operator      operator *( byref A as mat4, byref s as double ) as mat4      '' scalar multiplication      return( mat4( _         A.a * s, A.b * s, A.c * s, A.d * s, _         A.e * s, A.f * s, A.g * s, A.h * s, _         A.i * s, A.j * s, A.k * s, A.l * s, _         A.m * s, A.n * s, A.o * s, A.p * s _      ) )   end operator      operator *( byref s as double, byref A as mat4 ) as mat4      '' scalar multiplication      return( mat4( _         A.a * s, A.b * s, A.c * s, A.d * s, _         A.e * s, A.f * s, A.g * s, A.h * s, _         A.i * s, A.j * s, A.k * s, A.l * s, _         A.m * s, A.n * s, A.o * s, A.p * s _      ) )   end operator      operator /( byref s as double, byref A as mat4 ) as mat4      /'         scalar divide            the 'division' of a matrix by another matrix is not defined, the            equivalent operation is multiplying one matrix by the inverse of the other            on scalars though, it can be defined, mostly for convenience purposes      '/      return( mat4( _         A.a / s, A.b / s, A.c / s, A.d / s, _         A.e / s, A.f / s, A.g / s, A.h / s, _         A.i / s, A.j / s, A.k / s, A.l / s, _         A.m / s, A.n / s, A.o / s, A.p / s _      ) )   end operator      operator /( byref A as mat4, byref s as double ) as mat4      return( mat4( _         A.a / s, A.b / s, A.c / s, A.d / s, _         A.e / s, A.f / s, A.g / s, A.h / s, _         A.i / s, A.j / s, A.k / s, A.l / s, _         A.m / s, A.n / s, A.o / s, A.p / s _      ) )   end operator      operator *( byref v as vec4, byref A as mat4 ) as vec4      /'         multiply a vector with a row matrix, resulting in a row vector (like Direct3D)            a row vector looks like this:                        | x y z w |                        and is the format that Direct3D uses. What this means, code-wise, is that you            have to pre-multiply the vectors with the matrices, and some other stuff, like            transposing the matrices if you are using column vectors (as this library does)      '/      return( vec4( _          A.a * v.x + A.e * v.y + A.i * v.z + A.m * v.w, _         A.b * v.x + A.f * v.y + A.j * v.z + A.n * v.w, _         A.c * v.x + A.g * v.y + A.k * v.z + A.o * v.w, _         A.d * v.x + A.h * v.y + A.l * v.z + A.p * v.w _        ) )   end operator      operator *( byref A as mat4, byref v as vec4 ) as vec4      /'         multiply a vector with a column matrix, resulting in a column vector (like OpenGL)            a column vector looks like this                        | x |            | y |            | z |            | w |                        and is the format favored by OpenGL. In this library, column vectors are used, for            compatibility.      '/    'w.X = A.a * v.X + A.b * v.Y + A.c * v.Z + A.d * v.W    'w.Y = A.e * v.X + A.f * v.Y + A.g * v.Z + A.h * v.W    'w.Z = A.i * v.X + A.j * v.Y + A.k * v.Z + A.l * v.W    'w.W = A.m * v.X + A.n * v.Y + A.o * v.Z + A.p * v.W      return( vec4( _        A.a * v.X + A.b * v.Y + A.c * v.Z + A.d * v.W, _        A.e * v.X + A.f * v.Y + A.g * v.Z + A.h * v.W, _        A.i * v.X + A.j * v.Y + A.k * v.Z + A.l * v.W, _        A.m * v.X + A.n * v.Y + A.o * v.Z + A.p * v.W _         ) )   end operator      '' utility functions   function inverse( byref M as mat4 ) as mat4      '' returns the inverse of the provided matrix      dim as mat4 I = mat4( M )            I.inverse()            return( I )   end function         namespace matrices      function translation( byval t as vec4 ) as mat4         return mat4( _            1.0, 0.0, 0.0, t.x, _            0.0, 1.0, 0.0, t.y, _            0.0, 0.0, 1.0, t.z, _            0.0, 0.0, 0.0, 1.0 )       end function            function identity() as mat4         '' returns an identity matrix         return( mat4( _            1.0, 0.0, 0.0, 0.0, _            0.0, 1.0, 0.0, 0.0, _            0.0, 0.0, 1.0, 0.0, _            0.0, 0.0, 0.0, 1.0 ) )      end function   end namespace   #endif`

utility.bi

Code: Select all

`/'   auxiliary (read: rushed) functions for performing some drawings on the   3D plane'/function clipSegment( byref cam as camera, byref s1 as vec4, byref s2 as vec4 ) as integer   /'      clip both vertices of a line if necessary            note that this 'clipping' is only done on the near plane, not on the scene            the sensible (and correct) way to do this would be to clip them against      the view frustum; but alas, in order to keep the code short and clear      the calculation of the frustum was removed from the camera class. Hence      this crappy function.            Clipping (especially against the near plane) is important, because if you      don't clip, when one of the line endpoints goes behind the near clipping      plane, it gets a negative coordinate (due to the projection), and the line      segment is no longer valid   '/   if( s1.z >= cam.nearClip() and s2.z >= cam.nearClip() ) then       ' Neither one is behind the camera, draw them both       return( -1 )   elseIf( s1.z < cam.nearClip() and s2.z >= cam.nearClip() ) then       ' First coordinate behind the camera, clip it       s1.x = s1.x + ((cam.nearClip() - s1.z) / (s2.z - s1.z)) * (s2.x - s1.x)       s1.y = s1.y + ((cam.nearClip() - s1.z) / (s2.z - s1.z)) * (s2.y - s1.y)       s1.z = cam.nearClip()          return( -1 )   elseIf( s1.z >= cam.nearClip() and s2.z < cam.nearClip() ) then       ' Second coordinate behind the camera, clip it       s2.x = s1.x + ((cam.nearClip() - s1.z) / (s2.z - s1.z)) * (s2.x - s1.x)       s2.y = s1.y + ((cam.nearClip() - s1.z) / (s2.z - s1.z)) * (s2.y - s1.y)       s2.z = cam.nearClip()          return( -1 )   else       ' Both coordinates behind the camera, don't draw       return( 0 )   end ifend functionsub drawLine3d( byref cam as camera, byref p1 As vec4, byref p2 As vec4, byval c as uint32 )    dim as vec4 pos1, pos2    pos1 = cam.transform( p1 - cam.getPos() )    pos2 = cam.transform( p2 - cam.getPos() )      dim as double x1, y1, x2, y2          dim as integer visible = clipSegment( cam, pos1, pos2 )    if( visible ) then      dim as double z1, z2      cam.perspective( pos1, cam.projectionPlane, x1, y1, z1 )      cam.perspective( pos2, cam.projectionPlane, x2, y2, z2 )      line( x1, y1 ) - ( x2, y2 ), c    end ifend sub`
paul doe
Posts: 1244
Joined: Jul 25, 2017 17:22
Location: Argentina

### Re: 3d without openGL

So, you will need to have, in the same folder:

platform.bi
core.bi
math.bi
vec4.bi
mat4.bi
arrayList.bi
utility.bi
camera.bas
object.bas

Once you've got all that assembled, copy and paste this demo code, in the same folder:

camera test.bas

Code: Select all

`#include once "fbgfx.bi"#include once "platform.bi" '' platform specific definitions#include once "core.bi" '' core functions and defines#include once "math.bi" '' math functions and constants#include once "vec4.bi" '' vec4 type and operations#include once "mat4.bi" '' mat4 type and operations#include once "arrayList.bi" '' convenience#include once "camera.bas" '' the camera class definition#include once "object.bas" '' the object definition#include once "utility.bi" '' to draw lines an' stuff directly in 3D space/'   3D playground      intended as a simple tutorial/framework to do 3D stuff without too much complication   it is a bare-bones implementation, actually. Mostly to test all the math/conventions    behind the representation of a 3D scene (using the term 'scene' veeeery loosely here)      conventions used      RIGHT HANDED coordinate system      positive rotations are COUNTER-CLOCKWISE      facets are defined in COUNTER-CLOCKWISE order      angles are in RADIANS - use radians( angle_in_degrees ) for convenience'/'' set a screen modeconst as integer scrW = 1200, scrH = 600screenRes( scrW, scrH, 32, , fb.gfx_alpha_primitives )windowTitle( fbVersion & " - 3D Playground" )/'   define a projection plane (in this case the screen itself)      the projection plane is where the image gets projected (duh)   in the pinhole camera model, the projection plane is actually the   so called 'near' plane. The near parameter of the camera is really   the distance of this plane to the origin of coordinates (in the   case of the camera, its position in the world coordinate system'/dim as rectangle projectionPlane = ( 0.0, 0.0, scrW, scrH )/'   instantiate a camera class   remember the parameters for the constructor:      position, in WORLD space   X axis of the camera's coordinate system   Y axis of the camera's coordinate system   Z axis of the camera's coordinate system. If you flip it, you will end with a      left-handed coordinate system. Not that it matters much, save for the fact that      all other methods treat the matrix as a right-handed one. If you use lookAt()      and suddenly the axes get flipped, this is the most probable reason   near clipping plane Z distance   far clipping plane Z distance   projection plane      a word of advice: the axes of the camera have to be of a certain length, for the   projection matrix depends on it. See the code below, where the U and V vectors   get scaled to match the relation of the projection plane. If you fool around with   them and suddenly the image looks like crap, this is the most probable reason.      to remember:      the U vector (the X one) controls the horizontal field of view      the V vector (the Y one) controls the vertical field of view      and the Z vector controls the zoom            of course, they must be perfectly perpendicular to each other, if not the resulting      image is sheared. Anyway, you don't have to mess with the constructor if you don't       like and/or understand what it does for the camera to do its job.      having said all that, play around! See how each function transforms the camera, and   the effect it has on the image. Look at the implementation to see how it's done, but   most important of all, have fun!'/dim as camera cam = camera( _   vec4( 0.0, 1.0, 3.0 ), _   vec4( 1.0, 0.0, 0.0 ), _   vec4( 0.0, 1.0, 0.0 ), _   vec4( 0.0, 0.0, -1.0 ), _   0.1, 10000.0, projectionPlane )/'   scene manager (yeah right)      here, it's just a list that contains the objects to be rendered, but in real code it   is an instance of the scene manager.'/dim as arrayList objects/'   let's create some objects, shall we? One, actually (yeah I'm cheap)      these are the vertices and edges to define a small paper plane that   you can manipulate (see the code in the main loop)   the objects ideally should be loaded from a file, but that will   complicate the implementation. If there's some interest, I will show   how you can load an obj file (or invent your own format if you wish)'/dim as object3d ptr obj = new object3d() '' it will be at the origin, at the startobj->vertices->add( new vec4( -0.1, 0.0, 0.0 ) )obj->vertices->add( new vec4(  0.1, 0.0, 0.0 ) )obj->vertices->add( new vec4(  0.0, 0.0, 0.3 ) )obj->vertices->add( new vec4(  0.0,  0.0, 0.0 ) )obj->vertices->add( new vec4(  0.0, -0.1, 0.0 ) )obj->vertices->add( new vec4(  0.0,  0.0, 0.3 ) )obj->edges->add( new edge( obj->vertices->get( 0 ), obj->vertices->get( 1 ) ) )obj->edges->add( new edge( obj->vertices->get( 1 ), obj->vertices->get( 2 ) ) )obj->edges->add( new edge( obj->vertices->get( 2 ), obj->vertices->get( 0 ) ) )obj->edges->add( new edge( obj->vertices->get( 3 ), obj->vertices->get( 4 ) ) )obj->edges->add( new edge( obj->vertices->get( 4 ), obj->vertices->get( 5 ) ) )obj->edges->add( new edge( obj->vertices->get( 5 ), obj->vertices->get( 3 ) ) )'' add the defined object to the sceneobjects.add( obj )'' move it a little so it doesn't lie exactly in the center of the worldobj->move( vec4( 1.0, 0.0, -1.0 ) )'' some variables used for interactiondim as integer mouseX, mouseY, mouseButton dim as integer oldMouseX, oldMouseY'' to make movement (somewhat) constantdim as double frameTime'' background color (change it if you don't like it)dim as uint32 backgroundColor = rgba( 8, 8, 8, 255 )color( , backgroundColor )'' (extremely crappy) main loopdo         '' grab an object to have fun   dim as object3d ptr o = objects.get( 0 )      '' update interaction   oldMouseX = mouseX   oldMouseY = mouseY      getMouse( mouseX, mouseY, , mouseButton )      /'      camera controls            click and draw the mouse to free look      w: forward      s: backward      a: strafe left      d: strafe right      q: up      e: down      space: look at the paper plane (hold to keep looking at it)            you see, they are very similar to that of a FPS   '/   if( multikey( fb.sc_e ) ) then      '' down      cam.move( -vec4( 0.0, 1.0, 0.0 ) * 20.0 * frameTime )   end if      if( multikey( fb.sc_q ) ) then      '' up      cam.move( vec4( 0.0, 1.0, 0.0 ) * 20.0 * frameTime )   end if   if( multikey( fb.sc_w ) ) then      '' forward      cam.move( normalize( cam.getDir() ) * 20.0 * frameTime )   end if      if( multikey( fb.sc_s ) ) then      '' backwards      cam.move( normalize( -cam.getDir() ) * 20.0 * frameTime )      end if   if( multikey( fb.sc_a ) ) then      '' left      cam.move( normalize( cam.getU() ) * 20.0 * frameTime )   end if   if( multikey( fb.sc_d ) ) then      '' right      cam.move( normalize( -cam.getU() ) * 20.0 * frameTime )   end if      if( multikey( fb.sc_space ) ) then      '' look at the object (in case we lost it)      cam.lookAt( o->getPos() )   end if   /'      paper plane controls               i: forward         k: backwards         j: turn left         l: turn right         u: turn up         o: turn down    '/   if( multikey( fb.sc_u ) ) then      o->rotate( normalize( o->getU() ), radians( 0.5 ) )      end if   if( multikey( fb.sc_o ) ) then      o->rotate( normalize( o->getU() ), radians( -0.5 ) )      end if   if( multikey( fb.sc_i ) ) then      o->move( normalize( o->getDir() ) * 10.0 * frameTime )   end if   if( multikey( fb.sc_k ) ) then      o->move( normalize( -( o->getDir() ) ) * 10.0 * frameTime )   end if   if( multikey( fb.sc_j ) ) then      o->rotate( normalize( o->getV() ), radians( 0.5 ) )   end if   if( multikey( fb.sc_l ) ) then      o->rotate( normalize( o->getV() ), radians( -0.5 ) )   end if      '' if the left mouse button is pressed, activate free look mode   if( mouseButton = 1 ) then      cam.rotate( vec4( 0.0, 1.0, 0.0 ), 5.0 * ( oldMouseX - mouseX ) / cam.projectionPlane.width )      cam.rotate( cam.getU, -5.0 * ( oldMouseY - mouseY ) / cam.projectionPlane.height )   end if      '' render the screen   frameTime = timer()   screenLock()            cls()            '' draw the floor to have a frame of reference      for zz as integer = -10 to 10 step 1         drawLine3d(   cam, vec4( -10, 0, zz ), vec4( 10, 0, zz ), rgba( 32, 32, 32, 255 ) )      next            for xx as integer = -10 to 10 step 1         drawLine3d( cam, vec4( xx, 0, -10 ), vec4( xx, 0, 10 ), rgba( 32, 32, 32, 255 ) )      next            /'         draw the absolute axes of the world                  x = red         y = green         z = blue      '/      drawLine3d( cam, vec4( 0.0, 0.0, 0.0 ), vec4( 3.0, 0.0, 0.0 ), rgba( 128, 0, 0, 255 ) )      drawLine3d( cam, vec4( 0.0, 0.0, 0.0 ), vec4( 0.0, 3.0, 0.0 ), rgba( 0, 128, 0, 255 ) )      drawLine3d( cam, vec4( 0.0, 0.0, 0.0 ), vec4( 0.0, 0.0, 3.0 ), rgba( 0, 0, 128, 255 ) )            '' renders the objects      for i as integer = 0 to objects.count - 1         dim as object3d ptr o                  o = objects.get( i )         '' render the edges         for j as integer = 0 to o->edges->count - 1            dim as edge ptr e = o->edges->get( j )                       drawLine3d( cam, o->objectSpaceToWorldSpace( *e->vertex1 ), o->objectSpaceToWorldSpace( *e->vertex2 ), rgba( 255, 255, 0, 255 ) )                  next      next   screenUnlock()      /'      it's important not to let this var to be negative, as it gets multiplied with      the camera vectors and could screw some calculations (movement, for example)      told you, the loop implementation is very crappy, but it seems like this is      esier for most people to grasp   '/   frameTime = abs( timer() - frameTime )      sleep( 1 )loop until multikey( fb.sc_escape ) '' loop until esc is pressed'' finally, if there's objects in the list, free themif( objects.count > 0 ) then   dim as object3d ptr o      for i as integer = 0 to objects.count - 1      o = objects.get( i )      delete( o )   nextend if`

And there you have it, your own 3D playground. Read the comments in the code to see how it works. Should you have any questions, or need some explanation on what the code does, I'm all ears ;-)

Have fun!
Paul
BasicCoder2
Posts: 3566
Joined: Jan 01, 2009 7:03
Location: Australia

### Re: 3d without openGL

paul doe wrote:And there you have it, your own 3D playground. Read the comments in the code to see how it works. Should you have any questions, or need some explanation on what the code does, I'm all ears ;-)

You get an A+ on the code. It does exactly what I would have liked to have done.
It is unlikely I will ever understand it but I wonder to what extent I could use the code myself?
I imagine it could be extended to have a world of filled polygons to do something like the old F/A-18 Interceptor Amiga game.
.
dafhi
Posts: 1345
Joined: Jun 04, 2005 9:51

### Re: 3d without openGL

paul doe wrote:
dafhi wrote:@paul doe - Very intrigued and inspired by your offerings

Hi dafhi

It was obvious from your first demo 'main' that you've got a 'complete solution' - based on my own experience one typically needs a 'camera,' an object list, or 'world' or 'object manager' - something many of us haven't got into yet. Your next demo 'grid test,' combined with what it does, and how 'trivially' it does it is, well, unprecedented.

A lot of us are eager to step into 3d programming and you've dropped a deluxe-model scene composition lab on our doorstep
Last edited by dafhi on Oct 11, 2017 18:00, edited 1 time in total.
dodicat
Posts: 6555
Joined: Jan 10, 2006 20:30
Location: Scotland

### Re: 3d without openGL

Basiccoder2

In my rotate the first parameter is the fulcrum.
The second parameter is the point to rotate.
The third parameter is the angle which holds all six trigs necessary :
Sin(x),Sin(y),Sin(z)
Cos(x),Cos(y),Cos(z)

The fourth parameter is the scaler which is (1,1,1) default (no scaling)

You have set up your points around (0,0,0), so this is the fulcrum.
Also the eye point for perspective will be .x=0 and .y=0
I set the eye point .z to 600 to show a little more perspective.
In the final showing of points I add 640 to the x's and 300 to the y's so the shape is in the screen centre.
You don't really need the initial rotate by pi/2.
(This was to get the stones around Saturn into a Saturn disc.)
You have only orthogonal axis to show.
You can apply your .5 scale during the rotate in the loop.

Code: Select all

`Type V3    As Single x,y,z    As Uinteger colEnd TypeType 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 Angle3DEnd TypeType float As V3#define incircle(cx,cy,radius,x,y) (cx-x)*(cx-x) +(cy-y)*(cy-y)<= radius*radiusFunction 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 FunctionFunction 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,p.col)End Function Function perspective(p As V3,eyepoint As V3) As V3    Dim As Single   w=1+(p.z/eyepoint.z)    Return Type<V3>((p.x-eyepoint.x)/w+eyepoint.x,_    (p.y-eyepoint.y)/w+eyepoint.y,_    (p.z-eyepoint.z)/w+eyepoint.z,p.col)End Function  Sub QsortZ(array() As V3,begin As Long,Finish As Ulong)    Dim As Long i=begin,j=finish    Dim As V3 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 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 sleeptimeEnd Functionscreenres 1280,600,32const pi=4*atn(1)redim as v3 a(0)                'set some pointsfor n as single =-250 to 250    var xp=0,yp=0,zp = n    var u=ubound(a)    redim preserve a(1 to u+1)    a(ubound(a))=type<v3>(xp, yp, zp, rnd*rgb(0,0,255))    pset(xp,yp),a(ubound(a)).colnext nfor n as single =-250 to 250    var xp=0,yp=n,zp = 0    var u=ubound(a)    redim preserve a(1 to u+1)    a(ubound(a))=type<v3>(xp, yp, zp, rnd*rgb(255,0,0))    pset(xp,yp),a(ubound(a)).colnext nfor n as single =-250 to 250    var xp=n,yp=0,zp = 0     var u=ubound(a)    redim preserve a(1 to u+1)    a(ubound(a))=type<v3>(xp, yp, zp, rnd*rgb(0,255,0))    pset(xp,yp),a(ubound(a)).colnext n               'print "the points are on the x y plane with a small random z component"'print "flip them all 90 degrees around y axis and expand a bit"'print "Press a key"'sleep'cls          redim as V3 rot(lbound(a) to ubound(a))'to hold all rotated points          dim as Angle3D A3d'=angle3D.construct(0,pi/2,0)   'flip all points by pi/2 on y axis          'for n as long=lbound(a) to ubound(a)    'rot(n)=rotate(type(0,0,0),a(n),A3D,type(0.5,0.5,0.5)) '1.4 is a scale up    'a(n)=rot(n)    'pset(a(n).x,a(n).y),a(n).col'next'print "Thats 'em all flipped without perspctive"'print "Press a key"     'sleepdim as float angdim as long fps'fix these two components at startang.z=pi/2 'tilt angle sidewaysang.y=-pi/7 'tilt angle back and forwarddim as string key    do    key=inkey    If key=Chr(255)+"K" Then ang.z-=.02     'left    If key=Chr(255)+"M" Then ang.z+=.02     'right    If key=Chr(255)+"P" Then ang.y-=.02     'down    If key=Chr(255)+"H" Then ang.y+=.02     'up    if key="z" then ang.x -=.02    if key="x" then ang.x +=.02                  A3D = Angle3D.construct(ang.x,ang.y,ang.z)'get the six rotate components (sines, coses  ...)        screenlock    cls    draw string(50,50),"Framerate "&fps    draw string(50,150),"Use the arrow keys"        for n as long=lbound(a) to ubound(a)        rot(n)=rotate(type(0,0,0),a(n),A3D,type(.5,.5,.5))        rot(n)=perspective(rot(n),type(0,0,600))    next        qsortZ(rot(),lbound(rot),ubound(rot))        for n as long=lbound(rot) to ubound(rot)        pset(rot(n).x+640,rot(n).y+300),rot(n).col    next        screenunlock    sleep regulate(64,fps),1loop until key=chr(27) `

Others:
Matrices are too slow.
My rotate does the matrix * vector multiplication straight off in three lines.
For axial rotation Rodrigues method can also be used (I have posted many examples).
But it involves normalising and cross products which slows the whole process down.
Although it is easier in a way, you just stipulate an axis (2,-6,8) for example.
BasicCoder2
Posts: 3566
Joined: Jan 01, 2009 7:03
Location: Australia

### Re: 3d without openGL

dodicat wrote:For axial rotation Rodrigues method can also be used (I have posted many examples).

Maybe axial rotation was the term I was looking for? That the controls would rotate the points around the axes of the object not the axes of the world? I would not have recognized at the time which examples of yours does that. I have collected all your demos but of course none of them have user control or any tutorial explanations of how the equations were derived.

With paul doe's paper plane the keys U and O rotate it around one unmoving axis and the same with the J and L keys. In both cases one axis remains fixed while the other two rotate around it. He uses the I and K to move along the the third axis which is what an plane would do but I would also want to be able to rotate around that axis as well move along the direction of an axis.

By the way with your example I just noticed the rnd*rgb(255,0,0) should have been plain rgb(255,0,0) to get the three pure colors for each axis instead of random values meant for your Saturn ring demo.

.
Last edited by BasicCoder2 on Oct 11, 2017 18:07, edited 1 time in total.
dodicat
Posts: 6555
Joined: Jan 10, 2006 20:30
Location: Scotland

### Re: 3d without openGL

Hi basiccoder2.
Here is an example with a grid:
(I have drawn the axis of rotation)

Code: Select all

`Screen 19,32Const pi=4*Atn(1)Type V3    As Single x,y,z    As Ulong colEnd TypeType 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 Angle3DEnd TypeType float As V3Function 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 FunctionFunction 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,p.col)End Function Function perspective(p As V3,eyepoint As V3) As V3    Dim As Single   w=1+(p.z/eyepoint.z)    Return Type<V3>((p.x-eyepoint.x)/w+eyepoint.x,_    (p.y-eyepoint.y)/w+eyepoint.y,_    (p.z-eyepoint.z)/w+eyepoint.z,p.col)End Function  Sub QsortZ(array() As V3,begin As Long,Finish As Ulong)'NOT USED HERE.    Dim As Long i=begin,j=finish    Dim As V3 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 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 setgrid(sx As Long,bx As Long,sy As Long,by As Long,st As Long,p() As v3)            #define U Ubound(p)            Redim p(0)            For y As Long=sy To by Step st                Redim Preserve p(1 To U+1)                p(U)=Type(sx,y,0)                Redim Preserve p(1 To U+1)                p(U)=Type(bx,y,0)                Line(sx,y)-(bx,y)            Next            For x As Long=sx To bx Step st                Redim Preserve p(1 To U+1)                p(U)=Type(x,sy,0)                Redim Preserve p(1 To U+1)                p(U)=Type(x,by,0)                Line(x,sy)-(x,by)            Next        End Sub                Sub drawgrid(p() As V3)            Dim As Ulong clr=Rgb(0,100,255)            For n As Long=Lbound(p) To Ubound(p)-1                Line (p(n).x,p(n).y)-(p(n+1).x,p(n+1).y),clr                If n=Ubound(p)-3 Then clr=Rgb(0,200,0)'the axis green                n+=1            Next        End Sub                Redim As v3 a()'the main array of 3D points        Print "The initial grid"                setgrid(200,600,100,500,40,a())        'get dead centre of grid -- fulcrum        Dim As Single cx,cy,cz        For n As Long=Lbound(a) To Ubound(a)            cx+=a(n).x            cy+=a(n).y            cz+=a(n).z        Next        Dim As Long sz=(Ubound(a) - Lbound(a) +1)        cx=cx/sz        cy=cy/sz        cz=cz/sz        Print "Grid fulcrum:"        Print cx,cy,cz                'put in an axis (normal to the grid)        Redim Preserve a(Lbound(a) To (Ubound(a)+2))        a(Ubound(a)-1)=Type(cx,cy,0)        a(Ubound(a))  =Type(cx,cy,-200)        Print "Now rotate pi/2 round y axis"        Print "Press a key"        Sleep                Cls        Redim As V3 rot(Lbound(a) To Ubound(a))'to hold all rotated points                Dim As Angle3D A3d=angle3D.construct(0,pi/2,0)   'flip all points by pi/2 on y axis        For n As Long=Lbound(a) To Ubound(a)            rot(n)=rotate(Type(cx,cy,cz),a(n),A3D,Type(1.4,1.4,1.4)) '1.4 is a scale up            a(n)=rot(n)            Pset(a(n).x,a(n).y),a(n).col        Next        Print "Thats the grid rotated, with axis shown"        Print "Press a key"        drawgrid(a())        Sleep                                Dim As float ang        Dim As Long fps        'fix these two components at start        ang.z=pi/2 'tilt angle sideways        ang.y=-pi/7 'tilt angle back and forward        Dim As String key                Do            key=Inkey            If key=Chr(255)+"K" Then ang.z-=.02     'left            If key=Chr(255)+"M" Then ang.z+=.02     'right            If key=Chr(255)+"P" Then ang.y-=.02     'down            If key=Chr(255)+"H" Then ang.y+=.02     'up                         ang.x+=.01  'the orbiting speed                         A3D=Angle3D.construct(ang.x,ang.y,ang.z)'get the six rotate components (sines, coses  ...)            Screenlock            Cls            Draw String(50,50),"Framerate "&fps            Draw String(50,150),"Use the arrow keys"            For n As Long=Lbound(a) To Ubound(a)                rot(n)=rotate(Type(cx,cy,cz),a(n),A3D,Type(1,1,1))                rot(n)=perspective(rot(n),Type(cx,cy,1000))            Next            ''qsortZ(rot(),lbound(rot),ubound(rot)) 'definately do not sort            drawgrid(rot())                        Screenunlock            Sleep regulate(64,fps),1        Loop Until key=Chr(27)                        Sleep `

I am not quite sure what you seek.
I presume a shape and an axis to rotate the shape around.

My first shape was Saturn's rings.
My second shape was your x/y/z axis.
My third shape (here) is a grid, using only the end points.
Rodrigues rotation is more complex.
Here is an example, I did it a long while back .

Code: Select all

`'SetEnviron("fbgfx=GDI")#include "fbgfx.bi"Using fbType v3    As Single x,y,z    Declare Property length As Single    Declare Property unit As v3    Declare Function AxialRotate(As v3,As Single,As v3) As v3    Declare Function PointRotate(As v3,As v3,As v3=Type<v3>(1,1,1)) As v3    Declare Function perspective(eyepoint As v3) As v3    #define vct Type<v3>    #define dot *    #define cross ^End TypeType line3d    As v3 v1,v2    Declare property length As SingleEnd Type#define map(a,b,x,c,d) ((d)-(c))*((x)-(a))/((b)-(a))+(c)#define incircle(cx,cy,radius,x,y) (cx-x)*(cx-x) +(cy-y)*(cy-y)<= radius*radius'degrees to radians#define r  .01745329#define onscreen (mx>0) and (mx<xres) and (my>0) and (my<yres)Operator + (v1 As v3,v2 As v3) As v3Return Type<v3>(v1.x+v2.x,v1.y+v2.y,v1.z+v2.z)End OperatorOperator -(v1 As v3,v2 As v3) As v3Return Type<v3>(v1.x-v2.x,v1.y-v2.y,v1.z-v2.z)End OperatorOperator * (f As Single,v1 As v3) As v3 Return vct(f*v1.x,f*v1.y,f*v1.z)End OperatorOperator *(v1 As v3,f As Single) As v3Return f*v1End OperatorOperator * (v1 As v3,v2 As v3) As Single 'dot productReturn v1.x*v2.x+v1.y*v2.y+v1.z*v2.zEnd OperatorOperator ^ (v1 As v3,v2 As v3) As v3     'cross productReturn Type<v3>(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 OperatorProperty v3.length As SingleReturn Sqr(this.x*this.x+this.y*this.y+this.z*this.z)End Propertyproperty line3d.length As SingleDim As v3 p=this.v1-this.v2Return p.lengthEnd propertyProperty v3.unit As v3Dim n As Single=this.lengthIf n=0 Then n=1e-20Return vct(this.x/n,this.y/n,this.z/n)End PropertyFunction v3.AxialRotate(centre As v3,Angle As Single,axis As v3) As v3    Dim As v3 v=This-centre    Dim As v3 norm=axis.unit    Return (V*Cos(Angle)+(Norm cross V)*Sin(Angle)+Norm*(Norm dot V)*(1-Cos(Angle)))+centre End FunctionFunction v3.perspective(eyepoint As v3) As v3    Dim As Single   w=1+(this.z/eyepoint.z)    If w=0 Then w=1e-20    Return vct((this.x-eyepoint.x)/w+eyepoint.x,(this.y-eyepoint.y)/w+eyepoint.y,(this.z-eyepoint.z)/w+eyepoint.z)End FunctionFunction length(v As v3) As Single    Return Sqr(v.x*v.x+v.y*v.y+v.z*v.z) End FunctionFunction v3.PointRotate(_point As v3,Angle As v3,scale As v3=Type<v3>(1,1,1)) As v3    Dim As v3 p=vct(This-_point)    Dim As v3 rot,temp    Dim As Single s=Sin(angle.x),c=Cos(angle.x)    temp=vct((p.y)*C+(-p.z)*S,(p.z)*C+(p.y)*S)    rot.y=temp.x    s=Sin(angle.y):c=Cos(angle.y)    temp=vct((temp.y)*C+(-p.x)*S,(p.x)*C+(temp.y)*S)    rot.z=temp.x    s=Sin(angle.z):c=Cos(angle.z)    temp=vct((temp.y)*C+(-rot.y)*S,(rot.y)*C+(temp.y)*S)    rot.x=temp.x:rot.y=temp.y    Return vct((scale.x*rot.x+_point.x),(scale.y*rot.y+_point.y),(scale.z*rot.z+_point.z))End FunctionFunction segment_distance(l As line3d,p As v3,Byref ip As v3=vct(0,0,0)) As Single    Dim As Single linelength=length(l.v1-l.v2)    Dim As Single dist= length( (1/linelength)*((l.v1-l.v2) cross (p-l.v1)))    If length(p-l.v1) >= length(p-l.v2) Then        var leg=(p-l.v1)        var part=Sqr((length(leg))^2-dist^2)        var temp=part/linelength        'If temp>=1 Then temp=1:dist=length(p-l.v2)        ip=l.v1+(temp)*(l.v2-l.v1)        Return dist    Else        var leg=(p-l.v2)        var part=Sqr((length(leg))^2-dist^2)        var temp=part/linelength        'If temp>=1 Then temp=1:dist=length(p-l.v1)        ip=l.v2+(temp)*(l.v1-l.v2)        Return dist    End If    Return distEnd Function  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'=================================================================================Function keys(a As v3) As Integer    If Multikey(SC_LEFT) Then a.y=a.y+.005':return 1    If Multikey(SC_RIGHT) Then a.y=a.y-.005':return 1    If Multikey(SC_UP   ) Then a.x=a.x-.005':return 1    If Multikey(SC_DOWN ) Then a.x=a.x+.005':return 1    If Multikey(SC_Q) Then a.z=a.z-.005':return 1    If Multikey(SC_W) Then a.z=a.z+.005':return 1    If Multikey(SC_SPACE) Then Return 2    Return 1End Function#macro display()ans=pt.AxialRotate(centre,angle*r,normal_line)Dim As v3 ppvar dist=segment_distance(L1,ans,pp)ans=ans.perspective(vct(xres/2,yres/2,600))angle=angle+2If angle >=360 Then angle=0ScreenlockClsCircle (pt.x,pt.y),5,5,,,,fpp=pp.perspective(vct(xres/2,yres/2,600))Circle(pp.x,pp.y),3,6,,,,fDraw String(pt.x,pt.y),"Point ("&int(pt.x) &"," &int(pt.y) & "," &pt.z &")"For a As Single=1 To 360 Step 2    Var trace=pt.AxialRotate(centre,a*r,normal_line)    trace=trace.perspective(vct(xres/2,yres/2,600))    Pset(trace.x,trace.y)Next avar p=(-250*normal_line+centre).x,q=(-250*normal_line+centre).yL1.v1=vct(p,q,-250*(normal_line+centre).z)Line(p,q)-(centre.x,centre.y)Draw String(p,q),"Normal (" & Int(p) &"," &int(q) &"," &int(-250*(normal_line+centre).z) &")"p=(250*normal_line+centre).x:q=(250*normal_line+centre).yL1.v2=vct(p,q,250*(normal_line+centre).z)Line(p,q)-(centre.x,centre.y)Draw String(p,q),"Normal (" & Int(p) &"," &int(q) &"," &int(250*(normal_line+centre).z) &")"''Circle(centre.x,centre.y),5,3,,,,fVar rad=map(-300,300,ans.z,20,10)Circle(ans.x,ans.y),rad,,,,,fDraw String(10,10),"FPS " &fpsDraw String(10,30),"Use Up,Down,Left,Right,q and w to turn the axis"Draw String(10,50),"Up/down for x axis, Left/right for y axis, q/w for z axis"Draw String(10,70),"Space resets"Draw String(10,90),"Use the mouse to move the point"Draw String(10,110),"Distance of point from normal " &int(dist)Draw String(10,140),"NOTE - Positive Z is into the screen"ScreenunlockSleep snooze,1#endmacro#macro mouse()Dim As Integer x=mx,y=my,dx,dyWhile mb = 1    Display()    Getmouse mx,my,,mb    If onscreen Then        If mx<>x Or my<>y  Then            dx = mx - x            dy = my - y            x = mx            y = my            pt.x=x+dx            pt.y=y+dy        End If    End IfWend#endmacroDim As Integer xres,yresScreen 20Screeninfo xres,yresDim As v3 centre=Type<v3>(xres/2,yres/2,0),ptDim As v3 normal_lineDim As Single angleDim As v3 ans,Rangle,rot=vct(0,1,.5)+centreDim As v3 firstpoint=centre-vct(200,0,0)Dim As line3d L1pt=firstpointDim As long flag,fpsDim As Integer mx,my,mb,mwDo    Getmouse mx,my,mw,mb    var snooze=regulate(60,fps)    flag=keys(Rangle)        If flag=1 Then rot=rot.PointRotate(centre,Rangle):Rangle=vct(0,0,0):flag=0    If flag=2 Then rot=vct(0,1,.5)+centre:Rangle=vct(0,0,0):flag=0:pt=firstpoint     normal_line=vct(rot.x,rot.y,rot.z)-centre        display()        If incircle(pt.x,pt.y,10,mx,my) Then                mouse()            End IfLoop Until Multikey(SC_ESCAPE)Sleep   `

My equations for the ordinary rotations (not rodrigue) were derived from the 3D rotation matrix.
But Matrix multiplication by looping is far too slow, and the matrix is only 3 X 3 anyway, so it can be written straight down.
Thus the contorted looking return in the rotator.