[solved]Math Guru II How to make a Normal little bit random FAST?

For other topics related to the FreeBASIC project or its community.
gothon
Posts: 224
Joined: Apr 11, 2011 22:22

Re: Math Guru II How to make a Normal little bit random FAST?

Postby gothon » Aug 18, 2016 1:25

D.J.Peters wrote:@dafhi @gothon first you think it's clever but in practice it looks different.

If you get any surface normal (it points in any direction)

1) how you would select a precalculated random vector ?
(as a minimum you have to get the right qudrant from sphere to chose one)


I suggested pre-computing a random low angle rotation matrix not a vector. A simple way to do this might be to make a small random rotation about the X, Y and Z axis, then multiplying all three rotations into one 3x3 matrix.

https://en.wikipedia.org/wiki/Rotation_ ... _rotations

Though there are probably other methods as-well since you might find a method giving you a more desirable distribution of rotations.

D.J.Peters wrote:2) more important how to rotate any new random vector does it fits the cone ?
(I think a hour about it but I don't get it)


Once you have a rotation matrix, you simply multiply it by your vector to apply the rotation.

Something like:
U.X = V.X * m_00 + V.Y * m_01 + V.Z * m_02
U.Y = V.X * m_10 + V.Y * m_11 + V.Z * m_12
U.Z = V.X * m_20 + V.Y * m_21 + V.Z * m_22

(9 multiplications and 6 additions)

D.J.Peters wrote:Show me the code if you like please


Yep, real code is best, but hold on as it will take me a bit of time to write it and debug it. :)
gothon
Posts: 224
Joined: Apr 11, 2011 22:22

Re: Math Guru II How to make a Normal little bit random FAST?

Postby gothon » Aug 18, 2016 5:44

Ok here is some code

Code: Select all

Type Vect3D
    As Single X, Y, Z
End Type

Type Mat3x3
    M(2, 2) As Single
End Type

'Multiply two 3x3 Matricies possably backwards? :P
Function MatMult3x3(A As Mat3x3, B As Mat3x3) As Mat3x3
    Dim C As Mat3x3
   
    For I As Integer = 0 To 2
        For J As Integer = 0 To 2
            For K As Integer = 0 To 2
                C.M(I, J) += A.M(I, K) * B.M(K, J)
            Next K
        Next J
    Next I
   
    Return C
End Function

'Fill parameters with a 2D rotation matrix
Sub RotMat2D(Angle As Double, ByRef M_00 As Single, ByRef M_01 As Single, ByRef M_10 As Single, ByRef M_11 As Single)
    M_00 = Cos(Angle): M_01 = -Sin(Angle) '<- Trig Functions
    M_10 = -M_01: M_11 = M_00
End Sub

Dim Normals() As Vect3D
ReDim Normals(5000000 - 1)

Dim Low_Angle_Rot(255, 9) As Mat3x3

' Precompute Rotation Matricies
For I As Integer = 0 To UBound(Low_Angle_Rot, 1)
    For J As Integer = 0 To UBound(Low_Angle_Rot, 2)
        Dim As Mat3x3 Rx, Ry, Rz, R_Tmp
        Rx.M(0, 0) = 1
        RotMat2D (J + 1) / 10 * (Rnd() - Rnd()), Rx.M(1, 1), Rx.M(1, 2), Rx.M(2, 1), Rx.M(2, 2)
        Ry.M(1, 1) = 1
        RotMat2D (J + 1) / 10 * (Rnd() - Rnd()), Ry.M(0, 0), Ry.M(0, 2), Ry.M(2, 0), Ry.M(2, 2)
        Rz.M(2, 2) = 1
        RotMat2D (J + 1) / 10 * (Rnd() - Rnd()), Rz.M(0, 0), Rz.M(0, 1), Rz.M(1, 0), Rz.M(1, 1)
        R_Tmp = MatMult3x3(Ry, Rx)
        Low_Angle_Rot(I, J) = MatMult3x3(Rz, R_Tmp)
    Next J
Next I

' Initialize normal vectors to something random but length 1
For I As Integer = 0 To UBound(Normals)
    Normals(I).X = Rnd * 2 - 1
    Normals(I).Y = Rnd * 2 - 1
    Normals(I).Z = Rnd * 2 - 1
    Var L = Sqr(Normals(I).X ^ 2 + Normals(I).Y ^ 2 + Normals(I).Z ^ 2)
    Normals(I).X /= L
    Normals(I).Y /= L
    Normals(I).Z /= L
Next I

#Define SMALL_RND (rnd()-rnd()) * surface_roughness
Var surface_roughness = 0.1

Var Start = Timer 'Randomize and Renormalize
For I As Integer = 0 To UBound(Normals)
    Normals(I).X += SMALL_RND
    Normals(I).Y += SMALL_RND
    Normals(I).Z += SMALL_RND
    Var L = Normals(I).X*Normals(I).X + Normals(I).Y*Normals(I).Y + Normals(I).Z*Normals(I).Z
    L=1/Sqr(L) : Normals(I).X*=L : Normals(I).Y*=L : Normals(I).Z*=L
Next I
Var T1 = Timer - Start

Dim As Double E = 0 'Determine how close to 1 the length is
For I As Integer = 0 To UBound(Normals)
    Var L = Normals(I).X*Normals(I).X + Normals(I).Y*Normals(I).Y + Normals(I).Z*Normals(I).Z
    E += Abs(1 - L)
Next I

Print
Print "Time to randomize and re-normalize vectors: " & T1
Print "Average Float Error: " & E / (UBound(Normals) + 1)

Start = Timer 'Multiply random rotation matricies
For I As Integer = 0 To UBound(Normals)
    Dim J As Integer = Int(Rnd * (UBound(Low_Angle_Rot) + 1))
    Dim N As Vect3D
    With Low_Angle_Rot(J, Int(surface_roughness* 10))
        N.X = Normals(I).X * .M(0,0) + Normals(I).Y * .M(0,1) + Normals(I).Z * .M(0,2)
        N.Y = Normals(I).X * .M(1,0) + Normals(I).Y * .M(1,1) + Normals(I).Z * .M(1,2)
        N.Z = Normals(I).X * .M(2,0) + Normals(I).Y * .M(2,1) + Normals(I).Z * .M(2,2)
    End With
    Normals(I).X = N.X
    Normals(I).Y = N.Y
    Normals(I).Z = N.Z
Next I
Var T2 = Timer - Start

E = 0 'Determine how close to 1 the length is
For I As Integer = 0 To UBound(Normals)
    Var L = Normals(I).X*Normals(I).X + Normals(I).Y*Normals(I).Y + Normals(I).Z*Normals(I).Z
    E += Abs(1 - L)
Next I

Print
Print "Time to rotate vectors: " & T2
Print "Average Float Error: " & E / (UBound(Normals) + 1)
Print
Print "Speedup factor: " & T1 / T2

Sleep

I get a speedup factor of around 2-3 on my computer, although that may also have something to do with there being only 1 Rnd() call in the loop instead of 6. I had to discretize surface_roughness to 10 levels so I could have different levels pre-computed. Though I just let it be 0.1 since there is no source of that information in this simple demonstration anyway.
D.J.Peters
Posts: 7838
Joined: May 28, 2005 3:28

Re: Math Guru II How to make a Normal little bit random FAST?

Postby D.J.Peters » Aug 18, 2016 6:21

@gothon I know how to use matrices I do every day and thank for your speed compare

Do you sure your random rotation forms a small cone around the original normal ?
(I will test if i'm back home)

Init 4 normals with for example an up vector 0,1,0 rotate it and plot it on screen (you can ignore z for this short test)
does it looks like a cone ?

Joshy
Image
Last edited by D.J.Peters on Sep 25, 2017 21:40, edited 1 time in total.
dodicat
Posts: 5988
Joined: Jan 10, 2006 20:30
Location: Scotland

Re: Math Guru II How to make a Normal little bit random FAST?

Postby dodicat » Aug 18, 2016 11:06

Try making the unit vectors in situ
(Excuse my Latin)
The time killer is getting a random + or -
I use doubles to check errors.

Code: Select all



screen 19,32
type pt
    as double x,y,z
end type


dim as long max=5000000

redim as pt p(1 to max)
redim as byte sign(1 to max+3)'get random signs
do
   
 dim as double t=timer

for n as long=1 to max+3
    sign(n)=IIf(Rnd>.5,1,-1)
    next

'dim as double t=timer

for n as long=1 to max 'create unit vectors on the fly
p(n)=type(sign(n)*rnd , _ 
          sign(n+1)*rnd*sqr(1-p(n).x*p(n).x), _
          sign(n+2)*sqr((1-(p(n).x*p(n).x+p(n).y*p(n).y))))
next n
print "Time taken ",timer-t,"  Press a key"


'get error residue
dim as single R
for n as long=1 to max
   r+=abs( 1-sqr(p(n).x^2+p(n).y^2+p(n).z^2))
    next n
print "ERROR ";r

dim as long ctrx,ctry,ctrz

for z as long=lbound(p) to ubound(p)
if (sqr(p(z).x^2+p(z).y^2+p(z).z^2) -1) >1e-5 then print z, sqr(p(z).x^2+p(z).y^2+p(z).z^2)'extra ckeck
if p(z).x>0 then ctrx+=1
if p(z).y>0 then ctry+=1
if p(z).z>0 then ctrz+=1
next z
print max/2, "Expected for each"
print ctrx,ctry,ctrz, "Got"
draw string(350,150),"Shade by .z"
draw string(350+250,150),"Bare bones"

Circle(400,300),100,4
For z As Long=LBound(p) To UBound(p)
   
   PSet(400+100*p(z).x,300+100*p(z).y),RGB(Abs(p(z).z)*255,Abs(p(z).z)*255,Abs(p(z).z)*25)
    PSet(400+250+100*p(z).x,300+100*p(z).y)
Next
print "press a key"
sleep
cls
loop until inkey=chr(27)
Tourist Trap
Posts: 2768
Joined: Jun 02, 2015 16:24

Re: Math Guru II How to make a Normal little bit random FAST?

Postby Tourist Trap » Aug 18, 2016 11:11

dodicat wrote:(Excuse my Latin)

Hi, dodicat,

I'm lost with your perspective. I don't know how to zoom/unzoom, or whatever (rotate view) with the perspective. I would need it to test some rodrigues' rotation (axis rotation). I have some problem obtaining the correct rotation, and to vizualise the things would help.

Code: Select all

'attempt of application of the Rodrigues rotation

'_________________________the dodi's perspective
type P3D
    as single   x,y,z
end type
function Perspective(P as P3D, Eyepoint as P3D) as P3D
    dim as single   w => 1 + ( P.z/Eyepoint.z )
    '
    return type<P3D>( (P.x-Eyepoint.x)/w + Eyepoint.x, _
                      (P.y-Eyepoint.y)/w + Eyepoint.y, _
                      (P.z-Eyepoint.z)/w + Eyepoint.z )
end function

'_________________________data
'origin h of the rotation axis (U)
dim as integer   hx = 295
dim as integer   hy = 400
dim as integer   hz = 160
dim as P3D h
    h.x = hx
    h.y = hy
    h.z = hz

'unit direction vector of axis (U)
dim as single   a   =   0
dim as single   b   =   100
dim as single   c   =   0
dim as single   norm    = sqr(a^2 + b^2 + c^2)
a = a/norm
b = b/norm
c = c/norm
dim as P3D u
    u.x = a
    u.y = b
    u.z = c

'point M to apply rotation on
dim as integer   mx   = 225
dim as integer   my   = 200
dim as integer   mz   = 160
dim as P3D M
    M.x = mx
    M.y = my
    M.z = mz

'transform
dim as P3D MM
dim as P3D MM2

'angle of the rotation
dim as double   axisRotationAngle

'the Rodrigues transform
function RodriguesRotation(H as P3D, U as P3D, M as P3D, angle as double) as P3D
    dim as P3D    result
    '
    'vector position of M relative to H
    dim as P3D G
        G.x = M.x - H.x
        G.y = M.y - H.y
        G.z = M.z - H.z
    '
    'K = 1 - cos angle
    dim as double   K = 1 - cos(angle)
    '
    dim as single aam = U.x*(U.x*G.x + U.y*G.y + U.z*G.z)
    dim as single bbm = U.y*(U.x*G.x + U.y*G.y + U.z*G.z)
    dim as single ccm = U.z*(U.x*G.x + U.y*G.y + U.z*G.z)
    '
    result.x = H.x + G.x*( U.x*U.x*K + cos(angle) ) + _
                     G.y*( U.x*U.y*K - U.z*sin(angle) ) + _
                     G.z*( U.x*U.z*K + U.y*sin(angle) )
    result.y = H.y + G.x*( U.x*U.y*K + U.z*sin(angle) ) + _
                     G.y*( U.y*U.y*K + cos(angle) ) + _
                     G.z*( U.y*U.z*K - U.x*sin(angle) )
    result.z = H.z + G.x*( U.x*U.z*K - U.y*sin(angle) ) + _
                     G.y*( U.y*U.z*K + U.x*sin(angle) ) + _
                     G.z*( U.z*U.z*K + cos(angle) )
   '
    return type<P3D>( result.x , result.y, result.z)
end function


'---------------------------------------_________
'---------------------------------------_________

screenRes 800, 550, 32

dim as P3D      pts(1 to 4)
dim as P3D      eyepoint

dim as double   angle
dim as integer  gmx,gmy,wheel
do
    angle += 0.04
    if angle>7 then angle = 0
   
    getmouse gmx,gmy,wheel
   
    eyepoint.x = (gmx - 100)
    eyepoint.y = 1.2*(gmy) + 20
    eyepoint.z = 455 - wheel*10
   
    pts(1) = type(50,400,0)
    pts(2) = type(540,400,0)
    pts(3) = type(295,400,160)
    pts(4) = type(295,200,160)
   
    MM      => RodriguesRotation(H, U, M, angle)        ''****
    MM2     => RodriguesRotation(H, U, M, 2*angle)      ''****

    screenlock
    cls
        circle (eyepoint.x, eyepoint.y - 20), 1
       
        for i as long=1 to 20
            var p1 => Perspective(pts(1), eyepoint)
            var p2 => Perspective(pts(2), eyepoint)
            pts(1).z += 20
            pts(2).z += 20
            line(p1.x, p1.y)-(p2.x, p2.y)
        next i
       
        var p3  => perspective(pts(3), eyepoint)
        var p4  => perspective(pts(4), eyepoint)
        var pH  => perspective(H, eyepoint)
        var pM => perspective(M, eyepoint)
        var pMM => perspective(MM, eyepoint)
        var pMM2 => perspective(MM2, eyepoint)
       
        circle (p3.x, p3.y), 4
        circle (p4.x, p4.y), 4
        circle (pM.x, pM.y), 4, rgb(200,0,200), , , , f
       
        circle (pMM.x, pMM.y), 4, rgb(100,200,100), , , , f
        MM = RodriguesRotation(H, U, M, angle + 0.01)
        MM2 = RodriguesRotation(H, U, M, angle + 0.01)
       
        circle (pMM.x, pMM.y), 2, rgb(100,200,100), , , , f
        circle (pMM2.x, pMM2.y), 2, rgb(200,200,100), , , , f
        MM = RodriguesRotation(H, U, M, angle + 0.02)
       
        pMM=perspective(MM,eyepoint)
        circle (pMM.x, pMM.y), 2, rgb(100,200,100), , , , f
       
        pMM=perspective(M,eyepoint)
        MM = RodriguesRotation(H, U, M, angle)

        circle (pH.x, pH.y), 4, rgb(100,100,200)
        line(p3.x,p3.y) - (p4.x,p4.y)
        line(int(pMM.x), int(pMM.y))-(int(pH.x), int(pH.y))
       
    screenunlock
    sleep 15
loop until len(inkey)

sleep


'(eof)
dafhi
Posts: 1258
Joined: Jun 04, 2005 9:51

Re: Math Guru II How to make a Normal little bit random FAST?

Postby dafhi » Aug 18, 2016 11:29

Finally got it done :-)

[update] translate() was not correct

Code: Select all

...
Last edited by dafhi on Aug 23, 2016 7:36, edited 3 times in total.
Tourist Trap
Posts: 2768
Joined: Jun 02, 2015 16:24

Re: Math Guru II How to make a Normal little bit random FAST?

Postby Tourist Trap » Aug 18, 2016 11:41

dafhi wrote:Finally got it done :-)

Very nice rendering.
dafhi
Posts: 1258
Joined: Jun 04, 2005 9:51

Re: Math Guru II How to make a Normal little bit random FAST?

Postby dafhi » Aug 18, 2016 11:51

Thanks man :-)
D.J.Peters
Posts: 7838
Joined: May 28, 2005 3:28

Re: Math Guru II How to make a Normal little bit random FAST?

Postby D.J.Peters » Aug 18, 2016 14:13

@dafhi
i'm back at home and tested your code looks ok
good job man

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

Re: Math Guru II How to make a Normal little bit random FAST?

Postby dodicat » Aug 18, 2016 15:08

Yea, nice Dafhi.
I cannot fill 5 million random unit vectors in any less than about 1/5 second on 32 bit.
64 bit is a tad faster
Tried the asm sqrt, but it is slower and only 32 bit.
randomize,2 seems the fastest rnd.
A lookup table for rnd, and use only once speeds it up a bit.

Code: Select all


randomize ,2
screen 19,32
type pt
    as double x,y,z
end type

dim as long max=5000000

redim as pt p(1 to max)
redim as byte sign(1 to max+3)'get random signs
redim as double rnds(1 to max)
'create a look up table

for n as long=1 to max
    rnds(n)= rnd-rnd
next
'===============================================
do
   
 dim as double t=timer

for n as long=1 to max+3
     sign(n)=IIf(Rnd(n)>.5,1,-1)
    next

''dim as double t=timer
'======================  CREATE THE UNIT VECTORS ==============
for n as long=1 to max ' on the fly
p(n)=type(sign(n)*rnds(n) , _    'use the look up table here only
          sign(n+1)*rnd*sqr(1-p(n).x*p(n).x), _
          sign(n+2)*sqr((1-(p(n).x*p(n).x+p(n).y*p(n).y))))
next n
print "Time taken ",timer-t


'get error residue
dim as double R
for n as long=1 to max
   r+=abs( 1-sqr(p(n).x^2+p(n).y^2+p(n).z^2))
    next n
print "ERROR ";r

dim as long ctrx,ctry,ctrz

for z as long=lbound(p) to ubound(p)
if (sqr(p(z).x^2+p(z).y^2+p(z).z^2) -1) >1e-5 then print z, sqr(p(z).x^2+p(z).y^2+p(z).z^2)'extra ckeck
if p(z).x>0 then ctrx+=1
if p(z).y>0 then ctry+=1
if p(z).z>0 then ctrz+=1
next z
print max/2, "Expected for each"
print ctrx,ctry,ctrz, "Got"
draw string(350,150),"Shade by .z"
draw string(350+250,150),"Bare bones"

Circle(400,300),100,4
For z As Long=LBound(p) To UBound(p)
   
   PSet(400+100*p(z).x,300+100*p(z).y),RGB(Abs(p(z).z)*255,Abs(p(z).z)*255,Abs(p(z).z)*25)
    PSet(400+250+100*p(z).x,300+100*p(z).y)
Next
print "press a key"
sleep
cls
loop until inkey=chr(27)
Tourist Trap
Posts: 2768
Joined: Jun 02, 2015 16:24

Re: Math Guru II How to make a Normal little bit random FAST?

Postby Tourist Trap » Aug 18, 2016 17:20

A little off-topic now, but even if can't be displayed in perspective, here is the Rodrigues rotation formula, both normal and for small angles. It seems well working. It may serve somewhere some day.

Code: Select all

#include "fbgfx.bi"

'useless here!
nameSpace cst
    const as double _pi => 4*atn(1)
    dim shared as double    presin(0 to 255)
    dim shared as double    precos(0 to 255)
    dim shared as double    cycleStep => 1*256/360
end nameSpace
sub PrecalcultaSinCos()
    dim as ubyte    currentStep
    dim as ubyte    furtherStep
    do
        currentStep=furtherStep
        cst.presin(currentStep) = sin(currentStep)
        cst.precos(currentStep) = cos(currentStep)
        furtherStep += cst.cycleStep
    loop until furtherStep<=currentStep
end sub

'starts here...
type P3D
    declare property PerspFactor() as double
    declare property PerspFactor(byval as double)
    declare property X2D() as integer
    declare property Y2D() as integer
    declare sub Rotate(byval AngleAroundX as double, _
                       byval AngleAroundY as double, _
                       byval AngleAroundZ as double)
    declare sub DrawPoint(byval Colour as ulong=rgb(200,175,130))
        as integer  _x
        as integer  _y
        as integer  _z
        as double   _perspFactor
    static as integer   _xPersp
    static as integer   _yPersp
    static as integer   _zPersp
end type
dim as integer   P3D._xPersp
dim as integer   P3D._yPersp
dim as integer   P3D._zPersp
property P3D.PerspFactor() as double
    '---->
    return THIS._perspFactor
end property
property P3D.PerspFactor(byval SetValue as double)
    THIS._perspFactor = SetValue
end property
property P3D.X2D() as integer
    '---->
    return ( THIS._x )
end property
property P3D.Y2D() as integer
    '---->
    return ( THIS._y )
end property
sub P3D.Rotate(byval AngleAroundXaxis as double, _
               byval AngleAroundYaxis as double, _
               byval AngleAroundZaxis as double)
    THIS._y = THIS._y*cos(AngleAroundXaxis) - THIS._z*sin(AngleAroundXaxis)
    THIS._z = THIS._y*sin(AngleAroundXaxis) + THIS._z*cos(AngleAroundXaxis)
    THIS._x = THIS._x*cos(AngleAroundYaxis) - THIS._z*sin(AngleAroundYaxis)
    THIS._z = THIS._x*sin(AngleAroundYaxis) + THIS._z*cos(AngleAroundYaxis)
    THIS._x = THIS._x*cos(AngleAroundZaxis) - THIS._y*sin(AngleAroundZaxis)
    THIS._y = THIS._x*sin(AngleAroundZaxis) + THIS._y*cos(AngleAroundZaxis)
end sub
sub P3D.DrawPoint(byval Colour as ulong=rgb(200,175,130))
    circle (THIS.X2D, THIS.Y2D), 8 - (6 - 6*THIS._z/600), Colour, , , , f
end sub
   

function RodriguesRotation(byref H as P3D, _
                           byref U as P3D, _
                           byref M as P3D, _
                           byval angle as double) _
                           as P3D
    '_the Rodrigues transform_
    '
    dim as P3D    result
    '
    'vector position of M relative to H
    dim as P3D G
        G._x = M._x - H._x
        G._y = M._y - H._y
        G._z = M._z - H._z
    '
    'K = 1 - cos angle
    dim as double   K = 1 - cos(angle)
    '
    result._x = H._x + G._x*( U._x*U._x*K + cos(angle) ) + _
                       G._y*( U._x*U._y*K - U._z*sin(angle) ) + _
                       G._z*( U._x*U._z*K + U._y*sin(angle) )
    result._y = H._y + G._x*( U._x*U._y*K + U._z*sin(angle) ) + _
                       G._y*( U._y*U._y*K + cos(angle) ) + _
                       G._z*( U._y*U._z*K - U._x*sin(angle) )
    result._z = H._z + G._x*( U._x*U._z*K - U._y*sin(angle) ) + _
                       G._y*( U._y*U._z*K + U._x*sin(angle) ) + _
                       G._z*( U._z*U._z*K + cos(angle) )
    '
    '---->
    return result
end function
function SmallRodriguesRotation(byref H as P3D, _
                                byref U as P3D, _
                                byref M as P3D, _
                                byval angle as double) _
                                as P3D
    '_the Rodrigues transform_
    '
    dim as P3D    result
    '
    'vector position of M relative to H
    dim as P3D G
        G._x = M._x - H._x
        G._y = M._y - H._y
        G._z = M._z - H._z
    '
    'K = 1 - cos angle
    dim as double   K = 0
    '
    result._x = H._x + G._x*( U._x*U._x*K + 1 ) + _
                       G._y*( U._x*U._y*K - U._z*(angle) ) + _
                       G._z*( U._x*U._z*K + U._y*(angle) )
    result._y = H._y + G._x*( U._x*U._y*K + U._z*(angle) ) + _
                       G._y*( U._y*U._y*K + 1 ) + _
                       G._z*( U._y*U._z*K - U._x*(angle) )
    result._z = H._z + G._x*( U._x*U._z*K - U._y*(angle) ) + _
                       G._y*( U._y*U._z*K + U._x*(angle) ) + _
                       G._z*( U._z*U._z*K + 1 )
    '
    '---->
    return result
end function




'- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - .
'.- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  -
randomize TIMER
screenRes 800, 540, 32, 2
screenSet 1, 0
color , rgb(45,58,72)
cls

PrecalcultaSinCos()

'_____________________________data____
'origin h of the rotation axis (U)
dim as single   hx = 400
dim as single   hy = 250
dim as single   hz = 0
dim as P3D h
    h._x = hx
    h._y = hy
    h._z = hz

'direction unit vector of the axis (U)
dim as single   a   =   0
dim as single   b   =   0
dim as single   c   =   100
dim as single   norm    = sqr(a^2 + b^2 + c^2)
dim as P3D abc
    abc._x = a
    abc._y = b
    abc._z = c
dim as P3D u
    u._x = a/norm
    u._y = b/norm
    u._z = c/norm

'point M to apply rotation on
dim as single shifter = 10

dim as single   mx   = 400 - shifter
dim as single   my   = 250
dim as single   mz   = 160
dim as P3D M
    m._x = mx
    m._y = my
    m._z = mz

'angle of the rotation
dim as double   axisRotationAngle
'_____________________________________


'- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - _
dim as P3D      randomSelection(0 to 260)
dim as double   randomAngle(0 to 260)
dim as integer  selectionIndex
dim as P3D      rotationResultPoint
dim as integer  gmX, gmY, gmWheel, gmBtn
dim as double   axisRotationDecorateAngle
do
    getMouse gmX, gmY, gmWheel, gmBtn
    line (0, 250)-(799, 250)
   
    m.DrawPoint(rgb(100,100,80))
    h.DrawPoint(rgb(140,210,180))
    rotationResultPoint = RodriguesRotation(h, u, m, axisRotationDecorateAngle)
    rotationResultPoint.DrawPoint()
    rotationResultPoint = RodriguesRotation(h, u, m, axisRotationDecorateAngle + 0.1)
    rotationResultPoint.DrawPoint()
    rotationResultPoint = RodriguesRotation(h, u, m, axisRotationDecorateAngle + 0.2)
    rotationResultPoint.DrawPoint()
    rotationResultPoint = RodriguesRotation(h, u, m, axisRotationDecorateAngle + 0.3)
    rotationResultPoint.DrawPoint()
    rotationResultPoint = RodriguesRotation(h, u, m, axisRotationDecorateAngle + 0.4)
    rotationResultPoint.DrawPoint(rgb(200,100,130))

    if shifter<2600 then
        randomAngle(selectionIndex) = (selectionIndex/200 + rnd())*cst._pi/8
        randomSelection(selectionIndex) = SmallRodriguesRotation(h, u, m, randomAngle(selectionIndex))
        for i as integer = 0 to selectionIndex
            randomSelection(i).DrawPoint(rgb(200,100,130)) 
        next i
    else
        for i as integer = 0 to 260
            randomSelection(i).DrawPoint(rgb(100,200,130)) 
        next i
    end if

    draw string (10,10), "current angle  --> "& left(str(axisRotationAngle), 4)
    draw string (10,20), "current radius --> "& str(shifter)
    draw string (10,30), "current index --> "& str(selectionIndex)
   
    screenSync()
    flip 1, 0
    cls

    if shifter<2600 then
        shifter += 10                   ''increase orbit radius
        selectionIndex += 1
        axisRotationAngle += cst._pi/8
        if axisRotationAngle>2*cst._pi then axisRotationAngle = 0
        m._x += shifter/100*(cos(axisRotationAngle) + sin(axisRotationAngle))
        m._y -= shifter/100*(sin(axisRotationAngle) - cos(axisRotationAngle))
        m._z = rnd()*180 + 10
    else
        m._x = mx - 120
        m._y = my
        m._z = mz
    end if
   
    axisRotationDecorateAngle += 0.01
    if axisRotationDecorateAngle>2*cst._pi then axisRotationDecorateAngle = 0
   
    sleep 15
loop until inkey()=chr(27)

'.- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -  -
'- - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - .

getKey()

'(eof)
dafhi
Posts: 1258
Joined: Jun 04, 2005 9:51

Re: Math Guru II How to make a Normal little bit random FAST?

Postby dafhi » Aug 18, 2016 17:49

@dodicat - I create a set of random points. After that I choose randomly 5 mill times

Code: Select all

'  subs Cone_Precalc and translate() have all the details

type imagevars ' 2016 Jan 31
  as any ptr              im, pixels
  as integer              w,h,bpp,bypp,pitch,numpages,flags,rate, is_screen, pitchBy, wm,hm
  as single               midx,midy, midxm, midym, diagonal
  as string               driver_name
  declare sub             scr_inf
  Declare Sub             screen_init(wid As single=-1, hgt As single=-1, bpp as UInteger=32, numPages as integer=1, Flags as integer=0)
  Declare Function        create(wid As long = -1, hgt As long=0, color As Ulong=RGB(0,0,0)) As Any ptr
  declare sub             cls(col as ulong=&HFF000000)
  declare sub             checkers(pColor As ULong=RGBA(145,145,145,255), size As UInteger = 12)
  Declare Destructor
 private:
  Declare Sub             destroy
  declare sub             vars_common
end type
Destructor imagevars:  Destroy:  End Destructor
Sub imagevars.Destroy():  If im <> 0 Then ImageDestroy im: im = 0: endif:  End Sub
Sub imagevars.cls(col as ulong):  line (0,0)-(wm,hm),col,bf:  end sub
sub imagevars.vars_common
  wm=w-1: midx=w/2: midxm = wm/2
  hm=h-1: midy=h/2: midym = hm/2
  diagonal=sqr(w*w+h*h): if bypp<>0 then pitchBy=pitch\bypp
end sub
sub imagevars.scr_inf
  ScreenInfo w,h, bpp, bypp, pitch, rate, driver_name:  pixels = ScreenPtr
  vars_common
end sub
Sub imagevars.screen_init(_w As single, _h As single, _bpp as UInteger, _numpages as integer, _flags as integer)
  dim as integer ww,hh
  ScreenInfo ww,hh
  _w=abs(_w): if _w<=1 then _w*=ww: _flags or= 8 'borderless
  _h=abs(_h): if _h<=1 then _h*=hh: _flags or= 8
  numpages=_numpages: flags=_flags
  Destroy: ScreenRes _w,_h,_bpp,numpages,flags: scr_inf: is_screen = -1
  if numPages > 1 then screenset 0,1
End sub
Function imagevars.create(wid As long, hgt As long, col As ULong) As Any Ptr
  if hgt=0 then scr_inf: wid=w: hgt=h
  Destroy: im = ImageCreate( wid, hgt, col )
  ImageInfo im, w, h, bypp, pitch, pixels:  bpp = bypp * 8
  vars_common: is_screen = 0: Return im
End Function
Sub imagevars.checkers(pColor As ULong, size As UInteger)
  Dim As UInteger SizeDouble=size*2,SizeM=size-1
  For Y as integer = 0 To hm Step size
    For X as integer = -size * ((Y/SizeDouble)=Int(Y/SizeDouble)) To wm Step SizeDouble
      Line im,(X,Y)-(X+SizeM,Y+SizeM),pColor,BF
    Next
  Next
End Sub


#Macro Alpha256(ret,back, fore, am, a256)
  ret=((_
  (fore And &Hff00ff) * a256 + _
  (back And &Hff00ff) * am + &H800080) And &Hff00ff00 Or (_
  (fore And &H00ff00) * a256 + _
  (back And &H00ff00) * am + &H008000) And &H00ff0000) Shr 8
#EndMacro


type aadot
  as single               rad = 0.65, alpha = 1
  as ulong                outlined = 0
  as imagevars ptr        p
  as any ptr              pixels
  declare sub             render_target(pti as imagevars ptr)
  declare sub             draw(x as single=0, y as single=0, col as ulong=&HFFFFFFFF)
 private:
  as single               slope, slope_X2
end type
sub aadot.render_target(pti as imagevars ptr)
  if pti->is_screen then: pixels = ScreenPtr
  else: pixels = pti->pixels
  end if:  p = pti
end sub
sub aadot.draw(x as single, y as single, col as ulong)
  if p->h < 1 then exit sub
  slope = alpha*256: slope_X2 = slope*2
  dim as single           cone_hgt = rad * slope
  dim as integer          _X = x-rad, X2 = x+rad:  _X += _X * (_X < 0)
  dim as integer          _Y = y-rad, Y2 = y+rad:  _Y += _Y * (_Y < 0)
  dim as single           dxClip = slope * (x - _X), dy = slope * (y - _Y)
  X2 += (X2 - p->wm) * (X2 > p->wm)
  Y2 += (Y2 - p->hm) * (Y2 > p->hm)
  dim as integer          pwidm = X2 - _X
  dim as ulong ptr        pCorner = pixels:  pCorner += _Y * p->pitchBy + _X
  if outlined then
    for py as ulong ptr = pCorner to @pCorner[ (Y2-_Y)*p->pitchBy ] step p->pitchBy
      dim as single ySq = dy*dy, dx = dxClip
      for px as ulong ptr = py to py + pwidm
        dim as integer  alph = cone_hgt - sqr(dx*dx + ySq)
        if alph > slope then alph = slope_X2 - alph
        if alph > 0 then
          dim as integer  AlphM = 256 - alph
          Alpha256( *px, *px, col, AlphM, alph )
        end if
        dx -= slope
      next:  dy -= slope
    next
  else
    for py as ulong ptr = pCorner to @pCorner[ (Y2-_Y)*p->pitchBy ] step p->pitchBy
      dim as single ySq = dy*dy, dx = dxClip
      for px as ulong ptr = py to py + pwidm
        dim as integer  alph = cone_hgt - sqr(dx*dx + ySq)
        if alph > slope then:  alph = slope
        else:  alph += alph * (alph < 0)
        end if:  dim as integer  AlphM = 256 - alph
        Alpha256( *px, *px, col, AlphM, alph )
        dx -= slope
      next:  dy -= slope
    next
  end if
end sub


#ifndef Pi
Const                   TwoPi = 8 * atn(1)
Const                   pi    = 4 * Atn(1)
#EndIf

Type v3
  As double             x,y,z
  as double             rad
  as long               col
End Type

Type Axis3D
  As v3                 AxisX=(1,0,0), AxisY=(0,-1,0), AxisZ=(0,0,1)
  As double             x,y,z
End Type

#Macro RotC_(dsta,dstb,srca,srcb)
    temp = cosa_*srca - sina_*srcb
    dstb = cosa_*srcb + sina_*srca
    dsta = temp
#EndMacro
#Macro RotC(a_,dst,src,dota,dotb)
  Scope
    Dim As single a = (a_)
    Dim As double cosa_ = Cos(a), sina_ = Sin(a), temp
    RotC_( (dst.AxisX)dota, (dst.AxisX)dotb, (src.AxisX)dota, (src.AxisX)dotb )
    RotC_( (dst.AxisY)dota, (dst.AxisY)dotb, (src.AxisY)dota, (src.AxisY)dotb )
    RotC_( (dst.AxisZ)dota, (dst.AxisZ)dotb, (src.AxisZ)dota, (src.AxisZ)dotb )
  End Scope
#EndMacro

#Macro xRot(dst,a_)
  RotC(a_,dst,dst,.y,.z)
#EndMacro

#Macro yRot(dst,a_)
  RotC(a_,dst,dst,.z,.x)
#EndMacro

#Macro zRot(dst,a_)
  RotC(a_,dst,dst,.x,.y)
#EndMacro


type tMouse
  as integer        x,y,b, xp,yp,bp
  as single         scalar = 1
  declare function  buttons as integer
  declare property  dx as single
  declare property  dy as single
end type
property tmouse.dx as single
  return scalar * (x-xp)
end property
property tmouse.dy as single
  return scalar * (y-yp)
end property
function tmouse.buttons as integer
  bp = b: xp = x: yp = y
  getmouse x,y,,b
  return b
end function


sub render(buf as imagevars, back_img as imagevars, points() as v3, axis as Axis3D, zoom as single)
  screenlock
    put (0,0), back_img.im, pset
    draw string (0,0), "cone precalculated", rgb(255,255,255)
    draw string (0,10), "translation", rgb(0,255,0)
    dim as aadot dot
    dot.render_target @buf
    dot.alpha = 0.75
    dim as single z_point_scale = 0.011
    For p as v3 ptr = @Points(LBound(points)) to @Points(UBound(points))
      dim as single rz_ = (axis.z + axis.AxisZ.z*p->z + axis.AxisX.z*p->x + axis.AxisY.z*p->y)
      If rz_ > 0.1 Then
        Dim As Single rz = zoom/rz_
        Dim as single y = buf.midy - rz*(axis.y + axis.AxisY.y*p->y + axis.AxisZ.y*p->z + axis.AxisX.y*p->x)
        dim as single x = buf.midx + rz*(axis.x + axis.AxisX.x*p->x + axis.AxisY.x*p->y + axis.AxisZ.x*p->z)
        dot.rad = z_point_scale * rz * p->rad
        dot.draw x,y, p->col
      end if
    Next
  screenunlock
end sub


sub translate(byref ret as v3, vec as v3, precalc as v3)
      '.axisY.z = .axisX.y  ' Referenced from sub Cone_Precalc()
      '.axisZ.x = .axisX.y
     
      '.axisY.x = .axisX.z
      '.axisZ.y = .axisX.z
   
      '.axisY.y = .axisX.x
      '.axisZ.z = .axisX.x
  ret.x = vec.x * precalc.x + vec.y * precalc.z + vec.z * precalc.y
  ret.y = vec.y * precalc.x + vec.z * precalc.z + vec.x * precalc.y
  ret.z = vec.z * precalc.x + vec.x * precalc.z + vec.y * precalc.y
'  ret.x = vec.x * axisX.x + vec.y * axisY.x + vec.z * axisZ.x
'  ret.y = vec.y * axisY.y + vec.z * axisZ.y + vec.x * axisX.y
'  ret.z = vec.z * axisZ.z + vec.x * axisX.z + vec.y * axisY.z
End Sub

sub TranslateAll(a() as v3, surfnorm as v3 = (0, 0, -1), plots as integer = 5000000)
  dim as integer ub = ubound(a)\2, num_points = ub+1
  for i as integer = 0 to plots - 1
    var source = int(rnd * num_points)
    var dest = num_points+source
    translate a(dest), surfnorm, a(source)
    a(dest).rad = 1
    a(dest).col = rgb(0,255,0)
  next
End Sub

sub Cone_Precalc(a() as v3, num_points as integer=25, cone_half_angle as single = pi / 20)
  redim a(num_points-1)
  for i as integer = 0 to num_points-1
    var mincos = cos(cone_half_angle)
    var cosdelta = 1 - mincos
    with a(i)
   
      var cosa = mincos + rnd * cosdelta
      var sina = sqr(1 - cosa*cosa)
      var x_angle = rnd * twopi
     
      .x = cosa
      .y = sina * cos(x_angle)
      .z = sina * sin(x_angle)
     
      ' 1. imagined to X matrix
     
      '.axisX.x = cosa
      '.axisX.y = sina * cos(x_angle)
      '.axisX.z = sina * sin(x_angle)
     
      ' 2. perpendicular vectors
   
      '.axisY.z = .axisX.y
      '.axisZ.x = .axisX.y
     
      '.axisY.x = .axisX.z
      '.axisZ.y = .axisX.z
   
      '.axisY.y = .axisX.x
      '.axisZ.z = .axisX.x
     
      .rad = 1
      .col = rgb(255,255,255)
   
    end with
  next
end sub


Sub Main

  dim as imagevars   buf, img
  dim as Axis3D   Axis

  buf.screen_init 640, 480
  img.create      buf.w, buf.h, rgb(48,48,48)
  img.checkers rgb(92,92,92), 30

  dim as single   anim_fps = 26
  dim as single   phys_fps = 55
 
  Axis.z  = 3.0
  dim as single   zoom = 0.6 * buf.diagonal

  dim as single   rot_plane = pi / 9, rot_plane_inc = 0.0001
 
  dim as v3 points()
  cone_precalc points()
 
  redim preserve points( (ubound(points)+1) * 2 - 1)
 
  ' ~~~~~~ anim / physics / input ~~~~~~~
  dim as single   anim_f, phys_f, ianim = 1 / anim_fps, iphys = 1 / phys_fps

  dim as double   tNow = Timer, td, tp = tNow
  dim as double   tDemoExit = tNow+5, tMessageTime = 3

  dim as string   kstr
  dim as tMouse   mouse
 
 
  do
 
    if anim_f <= 0 then
      translateall points()
      render buf, img, points(), axis, zoom
      tMessageTime -= ianim
      anim_f += ianim
    end If

    sleep 15

    tp = tNow:  tNow = Timer:  td = tNow - tp
    anim_f -= td
   
    if mouse.buttons > 0 then
      mouse.scalar = Pi / buf.diagonal
      yrot(axis, -mouse.dx)
      xrot(axis, -mouse.dy)
    else
      dim as single rspeed = TwoPi / 8500
      phys_f += td
      while phys_f > 0
        zRot(axis, -rot_plane)
        yRot(axis, rspeed)
        zRot(axis, rot_plane)
        rot_plane += rot_plane_inc
        phys_f -= iphys
      wend
    end if
   
    kstr=left(inkey$,1)
   
  Loop until kstr = chr(27) 'esc
 
End Sub

Main
dafhi
Posts: 1258
Joined: Jun 04, 2005 9:51

Re: Math Guru II How to make a Normal little bit random FAST?

Postby dafhi » Aug 18, 2016 19:23

Code: Select all

' plots visualization with translation

type imagevars ' 2016 Jan 31
  as any ptr              im, pixels
  as integer              w,h,bpp,bypp,pitch,numpages,flags,rate, is_screen, pitchBy, wm,hm
  as single               midx,midy, midxm, midym, diagonal
  as string               driver_name
  declare sub             scr_inf
  Declare Sub             screen_init(wid As single=-1, hgt As single=-1, bpp as UInteger=32, numPages as integer=1, Flags as integer=0)
  Declare Function        create(wid As long = -1, hgt As long=0, color As Ulong=RGB(0,0,0)) As Any ptr
  declare sub             cls(col as ulong=&HFF000000)
  declare sub             checkers(pColor As ULong=RGBA(145,145,145,255), size As UInteger = 12)
  Declare Destructor
 private:
  Declare Sub             destroy
  declare sub             vars_common
end type
Destructor imagevars:  Destroy:  End Destructor
Sub imagevars.Destroy():  If im <> 0 Then ImageDestroy im: im = 0: endif:  End Sub
Sub imagevars.cls(col as ulong):  line (0,0)-(wm,hm),col,bf:  end sub
sub imagevars.vars_common
  wm=w-1: midx=w/2: midxm = wm/2
  hm=h-1: midy=h/2: midym = hm/2
  diagonal=sqr(w*w+h*h): if bypp<>0 then pitchBy=pitch\bypp
end sub
sub imagevars.scr_inf
  ScreenInfo w,h, bpp, bypp, pitch, rate, driver_name:  pixels = ScreenPtr
  vars_common
end sub
Sub imagevars.screen_init(_w As single, _h As single, _bpp as UInteger, _numpages as integer, _flags as integer)
  dim as integer ww,hh
  ScreenInfo ww,hh
  _w=abs(_w): if _w<=1 then _w*=ww: _flags or= 8 'borderless
  _h=abs(_h): if _h<=1 then _h*=hh: _flags or= 8
  numpages=_numpages: flags=_flags
  Destroy: ScreenRes _w,_h,_bpp,numpages,flags: scr_inf: is_screen = -1
  if numPages > 1 then screenset 0,1
End sub
Function imagevars.create(wid As long, hgt As long, col As ULong) As Any Ptr
  if hgt=0 then scr_inf: wid=w: hgt=h
  Destroy: im = ImageCreate( wid, hgt, col )
  ImageInfo im, w, h, bypp, pitch, pixels:  bpp = bypp * 8
  vars_common: is_screen = 0: Return im
End Function
Sub imagevars.checkers(pColor As ULong, size As UInteger)
  Dim As UInteger SizeDouble=size*2,SizeM=size-1
  For Y as integer = 0 To hm Step size
    For X as integer = -size * ((Y/SizeDouble)=Int(Y/SizeDouble)) To wm Step SizeDouble
      Line im,(X,Y)-(X+SizeM,Y+SizeM),pColor,BF
    Next
  Next
End Sub


#Macro Alpha256(ret,back, fore, am, a256)
  ret=((_
  (fore And &Hff00ff) * a256 + _
  (back And &Hff00ff) * am + &H800080) And &Hff00ff00 Or (_
  (fore And &H00ff00) * a256 + _
  (back And &H00ff00) * am + &H008000) And &H00ff0000) Shr 8
#EndMacro


type aadot
  as single               rad = 0.65, alpha = 1
  as ulong                outlined = 0
  as imagevars ptr        p
  as any ptr              pixels
  declare sub             render_target(pti as imagevars ptr)
  declare sub             draw(x as single=0, y as single=0, col as ulong=&HFFFFFFFF)
 private:
  as single               slope, slope_X2
end type
sub aadot.render_target(pti as imagevars ptr)
  if pti->is_screen then: pixels = ScreenPtr
  else: pixels = pti->pixels
  end if:  p = pti
end sub
sub aadot.draw(x as single, y as single, col as ulong)
  if p->h < 1 then exit sub
  slope = alpha*256: slope_X2 = slope*2
  dim as single           cone_hgt = rad * slope
  dim as integer          _X = x-rad, X2 = x+rad:  _X += _X * (_X < 0)
  dim as integer          _Y = y-rad, Y2 = y+rad:  _Y += _Y * (_Y < 0)
  dim as single           dxClip = slope * (x - _X), dy = slope * (y - _Y)
  X2 += (X2 - p->wm) * (X2 > p->wm)
  Y2 += (Y2 - p->hm) * (Y2 > p->hm)
  dim as integer          pwidm = X2 - _X
  dim as ulong ptr        pCorner = pixels:  pCorner += _Y * p->pitchBy + _X
  if outlined then
    for py as ulong ptr = pCorner to @pCorner[ (Y2-_Y)*p->pitchBy ] step p->pitchBy
      dim as single ySq = dy*dy, dx = dxClip
      for px as ulong ptr = py to py + pwidm
        dim as integer  alph = cone_hgt - sqr(dx*dx + ySq)
        if alph > slope then alph = slope_X2 - alph
        if alph > 0 then
          dim as integer  AlphM = 256 - alph
          Alpha256( *px, *px, col, AlphM, alph )
        end if
        dx -= slope
      next:  dy -= slope
    next
  else
    for py as ulong ptr = pCorner to @pCorner[ (Y2-_Y)*p->pitchBy ] step p->pitchBy
      dim as single ySq = dy*dy, dx = dxClip
      for px as ulong ptr = py to py + pwidm
        dim as integer  alph = cone_hgt - sqr(dx*dx + ySq)
        if alph > slope then:  alph = slope
        else:  alph += alph * (alph < 0)
        end if:  dim as integer  AlphM = 256 - alph
        Alpha256( *px, *px, col, AlphM, alph )
        dx -= slope
      next:  dy -= slope
    next
  end if
end sub


#ifndef Pi
Const                   TwoPi = 8 * atn(1)
Const                   pi    = 4 * Atn(1)
#EndIf

Type v3
  As double             x,y,z
  as double             rad
  as long               col
End Type

Type Axis3D
  As v3                 AxisX=(1,0,0), AxisY=(0,-1,0), AxisZ=(0,0,1)
  As double             x,y,z
End Type

#Macro RotC_(dsta,dstb,srca,srcb)
    temp = cosa_*srca - sina_*srcb
    dstb = cosa_*srcb + sina_*srca
    dsta = temp
#EndMacro
#Macro RotC(a_,dst,src,dota,dotb)
  Scope
    Dim As single a = (a_)
    Dim As double cosa_ = Cos(a), sina_ = Sin(a), temp
    RotC_( (dst.AxisX)dota, (dst.AxisX)dotb, (src.AxisX)dota, (src.AxisX)dotb )
    RotC_( (dst.AxisY)dota, (dst.AxisY)dotb, (src.AxisY)dota, (src.AxisY)dotb )
    RotC_( (dst.AxisZ)dota, (dst.AxisZ)dotb, (src.AxisZ)dota, (src.AxisZ)dotb )
  End Scope
#EndMacro

#Macro xRot(dst,a_)
  RotC(a_,dst,dst,.y,.z)
#EndMacro

#Macro yRot(dst,a_)
  RotC(a_,dst,dst,.z,.x)
#EndMacro

#Macro zRot(dst,a_)
  RotC(a_,dst,dst,.x,.y)
#EndMacro


type tMouse
  as integer        x,y,b, xp,yp,bp
  as single         scalar = 1
  declare function  buttons as integer
  declare property  dx as single
  declare property  dy as single
end type
property tmouse.dx as single
  return scalar * (x-xp)
end property
property tmouse.dy as single
  return scalar * (y-yp)
end property
function tmouse.buttons as integer
  bp = b: xp = x: yp = y
  getmouse x,y,,b
  return b
end function


sub render(buf as imagevars, back_img as imagevars, points() as v3, axis as Axis3D, zoom as single, LB as integer = -1)
  screenlock
    put (0,0), back_img.im, pset
    draw string (0,0), "original"
    draw string (0,10), "translation", rgb(0,255,0)
    dim as aadot dot
    dot.render_target @buf
    dot.alpha = 0.75
    dim as single z_point_scale = 0.011
    if lb < 0 then lb = lbound(points)
    For p as v3 ptr = @Points(LB) to @Points(UBound(points))
      dim as single rz_ = (axis.z + axis.AxisZ.z*p->z + axis.AxisX.z*p->x + axis.AxisY.z*p->y)
      If rz_ > 0.1 Then
        Dim As Single rz = zoom/rz_
        Dim as single y = buf.midy - rz*(axis.y + axis.AxisY.y*p->y + axis.AxisZ.y*p->z + axis.AxisX.y*p->x)
        dim as single x = buf.midx + rz*(axis.x + axis.AxisX.x*p->x + axis.AxisY.x*p->y + axis.AxisZ.x*p->z)
        dot.rad = z_point_scale * rz * p->rad
        dot.draw x,y, p->col
      end if
    Next
  screenunlock
end sub


sub translate(byref ret as v3, vec as v3, precalc as v3)
  ret.x = vec.x * precalc.x + vec.y * precalc.z + vec.z * precalc.y
  ret.y = vec.y * precalc.x + vec.z * precalc.z + vec.x * precalc.y
  ret.z = vec.z * precalc.x + vec.x * precalc.z + vec.y * precalc.y
End Sub

sub TranslateAll(a() as v3, surfnorm as v3 = (0, 0, -1))
  dim as integer ub = ubound(a)\2, num_points = ub+1
  for i as integer = 0 to ub
    var dest = i+num_points
    translate a(dest), surfnorm, a(i)
    a(dest).rad = a(i).rad
    a(dest).col = rgb(0,255,0)
  next
End Sub

sub Cone_Precalc(a() as v3, num_points as integer=25, cone_half_angle as single = pi / 20)
  redim a(num_points-1)
  for i as integer = 0 to num_points-1
    var mincos = cos(cone_half_angle)
    var cosdelta = 1 - mincos
    with a(i)
      var cosa = mincos + rnd * cosdelta
      var sina = sqr(1 - cosa*cosa)
      var x_angle = rnd * twopi
      .x = cosa
      .y = sina * cos(x_angle)
      .z = sina * sin(x_angle)
      .rad = 1
      .col = rgb(255,255,255)
    end with
  next
end sub


Sub Main

  dim as imagevars   buf, img
  dim as Axis3D   Axis

  buf.screen_init 640, 480
  img.create      buf.w, buf.h, rgb(48,48,48)
  img.checkers rgb(92,92,92), 30

  dim as single   anim_fps = 26
  dim as single   phys_fps = 55
 
  Axis.z  = 3.0
  dim as single   zoom = 0.6 * buf.diagonal

  dim as single   rot_plane = pi / 9, rot_plane_inc = 0.0001
 
  dim as v3 points()
  cone_precalc points(), 250, pi
 
  redim preserve points( (ubound(points)+1) * 2 - 1)
  translateall points()
 
  ' ~~~~~~ anim / physics / input ~~~~~~~
  dim as single   anim_f, phys_f, ianim = 1 / anim_fps, iphys = 1 / phys_fps

  dim as double   tNow = Timer, td, tp = tNow
  dim as double   tDemoExit = tNow+5, tMessageTime = 3

  dim as string   kstr
  dim as tMouse   mouse
 
 
  do
 
    if anim_f <= 0 then
      var translations_lbound = 0'(ubound(points)+1)\2
      render buf, img, points(), axis, zoom, translations_lbound
      tMessageTime -= ianim
      anim_f += ianim
    end If

    sleep 15

    tp = tNow:  tNow = Timer:  td = tNow - tp
    anim_f -= td
   
    if mouse.buttons > 0 then
      mouse.scalar = Pi / buf.diagonal
      yrot(axis, -mouse.dx)
      xrot(axis, -mouse.dy)
    else
      dim as single rspeed = TwoPi / 8500
      phys_f += td
      while phys_f > 0
        zRot(axis, -rot_plane)
        yRot(axis, rspeed)
        zRot(axis, rot_plane)
        rot_plane += rot_plane_inc
        phys_f -= iphys
      wend
    end if
   
    kstr=left(inkey$,1)
   
  Loop until kstr = chr(27) 'esc
 
End Sub

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

Re: Math Guru II How to make a Normal little bit random FAST?

Postby dodicat » Aug 18, 2016 20:05

Thanks Dafhi.
I think you have about .3 secomds per 5000000 (from your code 2 up)
(within the TranslateAll)
Is that correct?
If I reduce a sphere to it's eight quadrants, I can speed mine up a bit.
By the way -gen gas is by far the fastest graphics.
Gcc is very slow.

Code: Select all



randomize ,2
screen 19,32
dim as any ptr row=screenptr
dim as integer pitch
dim as ulong ptr pixel
screeninfo ,,,,pitch

type pt
    as double x,y,z
end type
'=============  EIGHT SPHERICAL QUADRANTS ================

sub f1(a() as pt)
    for n as long=1 to ubound(a)/8
     a(n)=type(rnd, _   
               rnd*sqr(1-a(n).x*a(n).x), _
               sqr((1-(a(n).x*a(n).x+a(n).y*a(n).y))))
next n
end sub
sub f2(a() as pt)
    for n as long=ubound(a)/8 to 2*ubound(a)/8
       a(n)=type(-rnd , _   
                 rnd*sqr(1-a(n).x*a(n).x), _
                 sqr((1-(a(n).x*a(n).x+a(n).y*a(n).y))))
next n
end sub
sub f3(a() as pt)
    for n as long=2*ubound(a)/8 to 3*ubound(a)/8
       a(n)=type(-rnd , _   
                -rnd*sqr(1-a(n).x*a(n).x), _
                 sqr((1-(a(n).x*a(n).x+a(n).y*a(n).y))))
next n
end sub

sub f4(a() as pt)
    for n as long=3*ubound(a)/8 to 4*ubound(a)/8
       a(n)=type(rnd , _   
                -rnd*sqr(1-a(n).x*a(n).x), _
                 sqr((1-(a(n).x*a(n).x+a(n).y*a(n).y))))
next n
end sub

sub f5(a() as pt)
    for n as long=4*ubound(a)/8 to 5*ubound(a)/8
       a(n)=type(rnd , _   
                 rnd*sqr(1-a(n).x*a(n).x), _
                 -sqr((1-(a(n).x*a(n).x+a(n).y*a(n).y))))
next n
end sub

sub f6(a() as pt)
    for n as long=5*ubound(a)/8 to 6*ubound(a)/8
       a(n)=type(-rnd , _   
                 rnd*sqr(1-a(n).x*a(n).x), _
                 -sqr((1-(a(n).x*a(n).x+a(n).y*a(n).y))))
next n
end sub

sub f7(a() as pt)
    for n as long=6*ubound(a)/8 to 7*ubound(a)/8
       a(n)=type(-rnd , _   
                 -rnd*sqr(1-a(n).x*a(n).x), _
                 -sqr((1-(a(n).x*a(n).x+a(n).y*a(n).y))))
next n
end sub

sub f8(a() as pt)
    for n as long=7*ubound(a)/8 to ubound(a)
       a(n)=type(rnd , _   
                 -rnd*sqr(1-a(n).x*a(n).x), _
                 -sqr((1-(a(n).x*a(n).x+a(n).y*a(n).y))))
next n
end sub

dim f(1 to 8)as sub(() as pt) ={@f1,@f2,@f3,@f4,@f5,@f6,@f7,@f8}
'=============================================================================
#macro ppset(_x,_y,colour)
pixel=row+pitch*(_y)+(_x) shl 2
*pixel=(colour)
#endmacro


dim as long max=5000000

redim as pt p(1 to max)

'===============================================
do
   
 dim as double t=timer
 'get the unit vectors
for n as long=1 to 8
  f(n)(p())
next n
'same as
'f1(p()):f2(p()):f3(p()):f4(p()):f5(p()):f6(p()):f7(p()):f8(p())

print "Time taken ",timer-t


'get error residue
dim as double R
for n as long=1 to max
   r+=abs( 1-sqr(p(n).x^2+p(n).y^2+p(n).z^2))
    next n
print "ERROR ";r

dim as long ctrx,ctry,ctrz

for z as long=lbound(p) to ubound(p)
if p(z).x>0 then ctrx+=1
if p(z).y>0 then ctry+=1
if p(z).z>0 then ctrz+=1
next z
print max/2, "Expected for each"
print ctrx,ctry,ctrz, "Got"
draw string(350,150),"Shade by .z"
draw string(350+250,150),"Bare bones"
draw string(0,100),"Vector p(10) = " +str(p(10).x)+" , "+str(p(10).y)+" , " +str(p(10).z)
screenlock

For z As Long=LBound(p) To UBound(p)
   PpSet(clng(400+100*p(z).x),clng(300+100*p(z).y),RGB(Abs(p(z).z)*255,Abs(p(z).z)*255,Abs(p(z).z)*25))
   PpSet(clng(400+250+100*p(z).x),clng(300+100*p(z).y),rgb(255,255,255))
Next
screenunlock
print "press a key"
sleep
cls
loop until inkey=chr(27)
Tourist Trap
Posts: 2768
Joined: Jun 02, 2015 16:24

Re: Math Guru II How to make a Normal little bit random FAST?

Postby Tourist Trap » Aug 18, 2016 20:24

Hi Dahfi, don't take it bad, but I wanted to really look at your nice work, but were enable to read it well without some "clean-up". Same code, just better fit my eyes.

And added some quick wheel zooming.

Code: Select all

'  subs Cone_Precalc and translate() have all the details
'[a bit refactored for comfort (of mine..)]

#ifndef _pi
    const as double _twoPi => 8*atn(1)
    const as double _pi    => 4*atn(1)
#endIf

#macro _ALPHA256(ret,back, fore, am, a256)
  ret=((_
          ( fore and &hFF00FF)*a256 + _
          ( back and &hFF00FF)*am + &h800080 ) and &hFF00FF00 or _
          (_
              ( fore and &H00FF00)*a256 + _
              ( back and &h00FF00)*am + &h008000 ) and &h00FF0000) shr 8
#endMacro

#macro _ROTC0(dsta, dstb, srca, srcb)
    temp = cosa_*srca - sina_*srcb
    dstb = cosa_*srcb + sina_*srca
    dsta = temp
#endMacro

#macro _ROTC(a_,dst, src, dota, dotb)
    scope
        dim as single   a       = ( a_ )
        dim as double   cosa_   = cos(a), sina_ = sin(a), temp
        _ROTC0( (dst._axisX)dota, (dst._axisX)dotb, (src._axisX)dota, (src._axisX)dotb )
        _ROTC0( (dst._axisY)dota, (dst._axisY)dotb, (src._axisY)dota, (src._axisY)dotb )
        _ROTC0( (dst._axisZ)dota, (dst._axisZ)dotb, (src._axisZ)dota, (src._axisZ)dotb )
    end scope
#endMacro

#macro _XROT(dst,a_)
    _ROTC(a_, dst, dst, ._y, ._z)
#endMacro

#macro _YROT(dst,a_)
    _ROTC(a_, dst, dst, ._z, ._x)
#endMacro

#macro _ZROT(dst,a_)
    _ROTC(a_, dst, dst, ._x, ._y)
#endMacro


type IMAGEVARS
' 2016 Jan 31
    declare destructor()
    declare sub ClearScreen(byval Col as ulong=&HFF000000)
    private:
    declare sub Destroy()
    declare sub VarsCommon()
    public:
    declare sub ScrInf()
    declare sub ScreenInit(byval Wid as single=-1, _
                           byval Hgt as single=-1, _
                           byval Bpp as uinteger=32, _
                           byval NumPages as integer=1, _
                           byval Flags as integer=0)
    declare function Create(byval Wid as long=-1, _
                            byval Hgt as long=0, _
                            byval Colour as ulong=rgb(0,0,0)) _
                            as any ptr
   
    declare sub Checkers(byval PColor as ulong=rgba(145,145,145,255), _
                         byval Size as uinteger=12)
        as any ptr      _im, _pixels
        as integer      _w, _h, _bpp, _
                        _bypp, _pitch, _numpages, _
                        _flags, _rate, _isScreen, _
                        _pitchBy, _wm, _hm
        as single       _midx, _midy, _
                        _midxm, _midym, _diagonal
        as string       _driverName
end type
destructor IMAGEVARS()
    Destroy()
end destructor
sub IMAGEVARS.Destroy()
    if THIS._im<>0 then
        imageDestroy(THIS._im)
        THIS._im = 0
    end if
end sub
sub IMAGEVARS.VarsCommon()
    THIS._wm    = THIS._w - 1
    THIS._midx  = THIS._w/2
    THIS._midxm = THIS._wm/2
    THIS._hm    = THIS._h - 1
    THIS._midy  = THIS._h/2
    THIS._midym = THIS._hm/2
    THIS._diagonal = sqr(THIS._w*THIS._w + THIS._h*THIS._h)
    if THIS._bypp<>0 then THIS._pitchBy = THIS._pitch\THIS._bypp
end sub
sub IMAGEVARS.ClearScreen(byval Col as ulong)
    line (0,0)-(THIS._wm, THIS._hm), Col, bf
end sub
sub IMAGEVARS.ScrInf()
    screenInfo THIS._w, THIS._h, THIS._bpp, _
               THIS._bypp, THIS._pitch, THIS._rate, _
               THIS._driverName
    THIS._pixels = screenPtr()
    THIS.VarsCommon()
end sub
sub IMAGEVARS.ScreenInit(byval Wid as single=-1, _
                         byval Hgt as single=-1, _
                         byval Bpp as uinteger=32, _
                         byval NumPages as integer=1, _
                         byval Flags as integer=0)
    dim as integer ww, hh
    screenInfo ww, hh
    '
    Wid = abs(Wid)
    if Wid<=1 then
        Wid *= ww
    end if
    Hgt = abs(Hgt)
    if Hgt<=1 then
        Hgt *=  hh
    end if
    THIS._w         = Wid
    THIS._h         = Hgt
    THIS._Bpp       = Bpp
    THIS._flags     or= 8
    THIS._numpages  = Numpages
    THIS._flags     = Flags
    '
    THIS.Destroy()
    '
    screenRes Wid, Hgt, Bpp, Numpages, Flags
    '
    THIS.ScrInf()
    THIS._isScreen = -1
    if NumPages> 1 then screenSet 0, 1
end sub
function IMAGEVARS.Create(byval Wid as long=-1, _
                          byval Hgt as long=0, _
                          byval Colour as ulong=rgb(0,0,0)) _
                          as any ptr
    if Hgt=0 then
        ScrInf()
        Wid = THIS._w
        Hgt = THIS._h
    end if
    '
    THIS.Destroy()
    '
    THIS._im => imageCreate( Wid, Hgt, Colour )
    imageInfo THIS._im, THIS._w, THIS._h, THIS._bypp, THIS._pitch, THIS._pixels
    THIS._bpp = THIS._bypp*8
    '
    THIS.VarsCommon()
    THIS._isScreen = 0
    '
    return THIS._im
end function
sub IMAGEVARS.Checkers(byval PColor as ulong, byval Size as uinteger)
    dim as uinteger sizeDouble  => Size*2
    dim as uinteger sizeM       => Size - 1
    for y as integer = 0 to THIS._hm step Size
        for x as integer = -Size*((y/sizeDouble)=int(y/sizeDouble)) to _
                            THIS._wm step _
                            sizeDouble
            line THIS._im, (x, y)-(x + sizeM, y + sizeM), PColor, bf
        next x
    next y
End Sub


type AADOT
    declare sub             RenderTarget(byval Pti as IMAGEVARS ptr)
    declare sub             DrawAadot(byval X as single=0, _
                                      byval Y as single=0, _
                                      byval C as ulong=&hFFFFFFFF)
        as single               _rad        = 0.65
        as single               _alpha      = 1
        as ulong                _outlined   = 0
        as IMAGEVARS ptr        _p
        as any ptr              _pixels
    private:
        as single               _slope
        as single               _slope_X2
end type
sub AADOT.RenderTarget(byval Pti as IMAGEVARS ptr)
    if pti->_isScreen then
        THIS._pixels = screenPtr()
    else
        THIS._pixels = Pti->_pixels
    end if
    '
    THIS._p = Pti
end sub
sub AADOT.DrawAadot(byval X as single=0, _
                    byval Y as single=0, _
                    byval C as ulong=&hFFFFFFFF)
    if THIS._p->_h<1 then exit sub
    '
    THIS._slope       = THIS._alpha*256
    THIS._slope_X2    = THIS._slope*2
    '
    dim as single     coneHgt   = THIS._rad*THIS._slope
    dim as integer    _x        = X - THIS._rad
    dim as integer    x2        = X + THIS._rad 
        _x += _x*(_x<0)
    dim as integer    _y        = Y - THIS._rad
    dim as integer    y2        = Y + THIS._rad 
        _y += _y*(_y<0)
    dim as single     dxClip    = THIS._slope*(X - _x)
    dim as single     dy        = THIS._slope * (y - _Y)
        x2 += (x2 - THIS._p->_wm)*(x2>THIS._p->_wm)
        y2 += (y2 - THIS._p->_hm)*(y2>THIS._p->_hm)
    dim as integer    pwidm     = x2 - _X
    dim as ulong ptr  pCorner   = THIS._pixels
        pCorner += _y*THIS._p->_pitchBy + _x
    '
    if THIS._outlined then
        for py as ulong ptr = pCorner to _
                              @pCorner[ (y2-_y)*THIS._p->_pitchBy ] step _
                              THIS._p->_pitchBy
            dim as single ySq   = dy*dy
            dim as single dx    = dxClip
            for px as ulong ptr = py to py + pwidm
                dim as integer  alph = coneHgt - sqr(dx*dx + ySq)
                if alph>THIS._slope then alph = THIS._slope_X2 - alph
                if alph>0 then
                dim as integer  alphM = 256 - alph
                _ALPHA256( *px, *px, C, alphM, alph )
                end if
                '
                dx -= THIS._slope
            next px
            dy -= THIS._slope
        next py
    else
        for py as ulong ptr = pCorner to @pCorner[ (Y2-_Y)*THIS._p->_pitchBy ] step THIS._p->_pitchBy
            dim as single ySq = dy*dy, dx = dxClip
            for px as ulong ptr = py to py + pwidm
                dim as integer  alph = coneHgt - sqr(dx*dx + ySq)
                '
                if alph> THIS._slope then
                    alph = THIS._slope
                else
                    alph += alph*(alph<0)
                end if
                '
                dim as integer  alphM = 256 - alph
                _ALPHA256( *px, *px, C, alphM, alph )
                dx -= THIS._slope
            next px
            dy -= THIS._slope
        next py
    end if
end sub


type V3
    as double   _x, _y, _z
    as double   _rad
    as long     _c
end type


type AXIS3D
    as V3       _axisX = (1,0,0)
    as V3       _axisY = (0,-1,0)
    as V3       _axisZ = (0,0,1)
    as double   _x, _y, _z
end type


type TMOUSE
    declare property  Dx() as single
    declare property  Dy() as single
    declare function  Buttons() as integer
        as integer      _x, _y, _b
        as integer      _xp, _yp, _bp
        as single       _scalar = 1   
end type
property TMOUSE.dx as single
    '
    return THIS._scalar*(THIS._x - THIS._xp)
end property
property TMOUSE.dy as single
    '
    return THIS._scalar*(THIS._y - THIS._yp)
end property
function TMOUSE.buttons as integer
    THIS._bp = THIS._b
    THIS._xp = THIS._x
    THIS._yp = THIS._y
    '
    getMouse    THIS._x, THIS._y, , THIS._b
    '
    return THIS._b
end function


'. _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
'. + + + + + + + + + + + + + + + + + + + +

declare sub Render(byref Buf as IMAGEVARS, _
                   byref BackImg as IMAGEVARS, _
                   Points() as V3, _
                   byref Axis as AXIS3D, _
                   byval Zoom as single, _
                   byval Lb as integer=-1)
sub Render(byref Buf as IMAGEVARS, _
           byref BackImg as IMAGEVARS, _
           Points() as V3, _
           byref Axis as AXIS3D, _
           byval Zoom as single, _
           byval lB as integer=-1)
    screenLock
        put (0,0), BackImg._im, PSET
        draw string (0,0), "cone precalculated", rgb(255,255,255)
        draw string (0,10), "translation", rgb(0,255,0)
        '
        dim as AADOT    dot
                        dot.RenderTarget (@Buf)
                        dot._alpha = 0.75
        '
        dim as single z_point_scale = 0.011
        '
        if Lb<0 then Lb = lBound(points)
        for p as V3 ptr = @Points(Lb) to @Points( uBound(Points) )
            dim as single rz_ = ( Axis._z + _
                                  Axis._axisZ._z*p->_z + _
                                  Axis._axisX._z*p->_x + _
                                  Axis._axisY._z*p->_y )
            If rz_>0.1 Then
                dim as single rz    = Zoom/rz_
                dim as single y     = Buf._midy - _
                                      rz*(Axis._y + _
                                          Axis._axisY._y*p->_y + _
                                          Axis._axisZ._y*p->_z + _
                                          Axis._AxisX._y*p->_x)
                dim as single x     = Buf._midx + _
                                      rz*(Axis._x + _
                                          Axis._axisX._x*p->_x + _
                                          Axis._axisY._x*p->_y + _
                                          Axis._axisZ._x*p->_z)
                '
                dot._rad = z_point_scale*rz*p->_rad
                dot.DrawAadot(x, y, p->_c)
            end if
        Next p
    screenUnlock
end sub


declare sub Translate(byref Ret as V3, Vec as V3, Precalc as V3)
sub Translate(byref Ret as V3, Vec as V3, Precalc as V3)
    Ret._x = Vec._x*Precalc._x + Vec._y*Precalc._z + Vec._z*Precalc._y
    Ret._y = Vec._y*Precalc._x + Vec._z*Precalc._z + Vec._x*Precalc._y
    Ret._z = Vec._z*Precalc._x + Vec._x*Precalc._z + Vec._y*Precalc._y
end sub


declare sub TranslateAll(A() as V3, byref SurfNorm as V3=(0, 0, -1))
sub TranslateAll(A() as V3, byref SurfNorm as V3=(0, 0, -1))
    dim as integer ub           = uBound(A)\2
    dim as integer num_points   = ub+1
    for i as integer = 0 to ub
        var dest    = i+num_points
        '
        Translate ( A(dest), SurfNorm, A(i) )
        A(dest)._rad = A(i)._rad
        A(dest)._c = rgb(0,255,0)
    next i
end sub


declare sub ConePrecalc(A() as V3, _
                        byval NumPoints as integer=25, _
                        byval ConeHalfAngle as single=_pi/20)
sub ConePrecalc(A() as V3, _
                byval NumPoints as integer=25, _
                byval ConeHalfAngle as single=_pi/20)
    redim A(NumPoints - 1)
    for i as integer = 0 to NumPoints - 1
        var mincos      = cos(ConeHalfAngle)
        var cosdelta    = 1 - mincos
        var cosa        = mincos + rnd()*cosdelta
        var sina        = sqr(1 - cosa*cosa)
        var x_angle     = rnd()*_twoPi
        '
        with A(i)
            ._x     = cosa
            ._y     = sina*cos(x_angle)
            ._z     = sina*sin(x_angle)
            ._rad   = 1
            ._c     = rgb(255,255,255)
        end with
    next i
end sub


declare sub Main()
sub Main()
    dim as IMAGEVARS      buf, img
    dim as AXIS3D         axis
    '
    buf.ScreenInit(640, 480)
    img.Create      buf._w, buf._h, rgb(48,48,48)
    img.Checkers    rgb(92,92,92), 30
    '
    dim as single   anim_fps = 26
    dim as single   phys_fps = 55
    '
    axis._z  = 3.0
    dim as single   zoom            = 0.6*buf._diagonal
    dim as single   rot_plane       = _pi/9, _
                    rot_plane_inc   = .0001
    dim as V3       points()
    ConePrecalc( points(), 5000, _pi )
    '
    redim preserve points( (ubound(points) + 1)*2 - 1)
    TranslateAll( points() )
    '
    '~~~~~~ anim / physics / input ~~~~~~
    dim as single   anim_f, _
                    phys_f, _
                    ianim = 1/anim_fps, _
                    iphys = 1/phys_fps
    '
    dim as double   tNow = TIMER, _
                    td, _
                    tp = tNow
    dim as double   tDemoExit = tNow+5, _
                    tMessageTime = 3
    dim as string   kstr
    dim as tMouse   mouse
    '~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
    do
        'addition--------------------------
        dim as integer gmX, gmY, gmWheel
        getMouse gmX, gmY, gmWheel
        zoom = 0.6*buf._diagonal + gmWheel/(40 + log(abs(gmWheel - 100)))*buf._diagonal
        '----------------------------------
       
       
        if anim_f<=0 then
            var translations_lbound = 0 ''(ubound(points)+1)\2
            '
            Render(buf, img, points(), axis, zoom, translations_lbound)
            tMessageTime    -= ianim
            anim_f          += ianim
        end If
        '
        tp  =   tNow
                tNow = TIMER
        td  =   tNow - tp
        anim_f -= td
        '
        if mouse.Buttons()>0 then
            mouse._scalar = _pi/buf._diagonal
            _YROT(axis, -mouse.Dx)
            _XROT(axis, -mouse.Dy)
        else
            dim as single rspeed    = _twoPi/8500
            phys_f += td
            while phys_f>0
                _ZROT(axis, -rot_plane)
                _YROT(axis, rspeed)
                _ZROT(axis, rot_plane)
                rot_plane += rot_plane_inc
                phys_f -= iphys
            wend
        end if
        '
        kstr=left(inkey$,1)
        sleep 15
    loop until ( kstr=chr(27) )
    '
end sub




'. _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
'. _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _
Main()
'. _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _ _


'END_OF_SOURCE_FILE

Return to “Community Discussion”

Who is online

Users browsing this forum: No registered users and 1 guest